集成电路制造过程中有一类重要的气相加工工艺，主要包括刻蚀、等离子体增强/化学气相沉积(plasma-enhanced/chemical vapor deposition,PE/CVD)、 物理气相沉积(physical vapor deposition,PVD)以及氧化扩散工艺等[1]，都是在具有高洁净度的真空腔室中完成的，而工艺的均匀性是衡量真空腔室性能的主要指标[2]。研究表明，采用喷淋头(showerhead)的真空腔室可以大大提高工艺的均匀性[3, 4, 5, 6]，喷淋头的性能也成为影响真空腔室性能的关键因素之一[7, 8]。喷淋头上通常布置有几千个微孔，每个微孔直径约为1 mm。在现有工艺操作压力下，喷淋头微孔内存在明显的稀薄气体效应。
直接模拟Monte Carlo (direct simulation Monte Carlo,DSMC) [9]方法对于稀薄气体流动模拟来说是一种有效的数值模拟方法[10]。从微观分子角度出发的DSMC方法在模拟过渡领域的气体流动问题上取得了巨大成功，它是通过引入真实物理模型来构造概率模型对Boltzmann方程碰撞项求解，研究表明其收敛于Boltzmann方程[11, 12]。DSMC方法本质上是基于气体稀薄和分子混沌假设。在这两个假设下，对于气体分子间的碰撞，DSMC方法只考察两个分子间的碰撞，也称为二元碰撞； 并在近似真实的碰撞时间Δt内，将分子的运动和碰撞解耦。在Δt时间步长内的运动是一个简单且直接的确定论步骤。两个分子碰撞完成后，通过随机分配方式，分配模拟分子碰撞后的速度和能量。就单个分子来说是随机的，但整个系统还是保持了动量和能量的守恒[13, 14]。
DSMC方法[9]是一种直接从气体流动的物理模型出发的方法，通过统计网格单元内模拟粒子的速度、内能等运动参量实现模拟气体流动的目的。粒子碰撞对的选择也是在网格内进行的。DSMC方法所用网格大体上可以分为两大类: 1) 以位置元方案为主的直角结构网格方法[15, 16, 17, 18]。该方法计算效率高，但不利于匹配复杂几何型面，边界处理上有误差[19]。2) 非结构网格方法。20世纪90年代初期，美国国家航空航天局(National Aeronautic and Space Administration，NASA) 发展了第一代非结构网格DSMC通用程序[10],模拟了三角翼的绕流。随后Wilmoth等对其性能进行了改善[20],计算了单级入轨飞船。同一时期,欧洲Laux等也开展了非结构网格DSMC方法的研究[21],取得了许多有意义的成果。王学德等研究了一类非结构网格的实现策略[22]。非结构网格由于其很好的适体性，使得边界处理非常容易，但由于非结构网格拓扑结构的高度无序性,需要采用较为繁琐的轨迹追踪方法来判断粒子所在的网格单元。与直角网格相比,计算效率至少低数倍乃至一个数量级[20]。
本文以Delaunay三角形为基本网格单元,研究了2维非结构网格DSMC算法的实现过程。首先，用矢量线段来描述计算域边界、粒子运动轨迹、三角形单元等，这样做的好处是： 在实现边界处理、粒子追踪等算法时只需少量的矢量和逻辑运算，提高了计算效率。其次，引入分区思想，将矩形辅助网格与三角形网格单元相耦合，缩小搜索区域，实现快速追踪模拟粒子，使算法的计算效率同使用直角网格的算法相接近。
1 DSMC程序主要算法DSMC方法是在Δt内将粒子的运动与粒子间的碰撞进行解耦处理，即将粒子的运动视为匀速直线运动，在其运动过程中没有和其他粒子碰撞，并计算模拟粒子在Δt间隔内到达的新位置，记录该粒子位置信息，得到模拟粒子所属的网格单元。然后，对该网格单元内的粒子进行碰撞计算。如果粒子与物面边界发生碰撞，则先计算与物面边界的碰撞，然后用碰撞后的速度计算粒子的运动。
本文自主编写的DSMC程序由8个模块构成，如图1所示。网格剖分采用Jonathan的开源 Triangle(Version 1.6)网格划分程序，利用 Tecplot360 实现了计算结果可视化，并行运算采用了openMP技术。
|
| 图 1 DSMC代码构成 |
粒子位置、边界顶点和三角形单元顶点信息存储为矢量点； 粒子运动轨迹、计算域边界、三角形单元边信息存储为若干矢量线段。粒子信息数据组织形式如图2所示。数组个数为网格单元总数，而网格单元内粒子是一个链表结构，链表节点储存单个粒子信息，包括粒子所在网格单元ID、 子网格序号、位置坐标、 3个方向速度和转动能与振动能等内能信息。
|
| 图 2 粒子信息数据组织形式 |
在每个Δt时间内，粒子按照一定速度到达新的位置，需要对粒子当前位置信息作出判断。首先，判断该粒子是否离开当前三角形网格单元，方法是沿着三角形网格单元的边按同一方向依次计算新的位置点与三角形网格单元顶点的叉积的符号，如果符号相同，则新的位置点位于该三角形网格单元内，否则在三角形网格单元外。如图3所示，这个算法[23]只用到了3次叉积计算，没有除法、三角函数和开根号等运算，因此计算效率、计算精度都很高，避免了浮点误差。如果该新的位置点不在当前三角形网格单元内，则区别下面两种情况: 1) 判断该粒子是否离开计算域，方法是依次计算粒子运动轨迹同边界线段的交点[24]，如图4所示。如果粒子运动轨迹同边界存在交点，则该粒子可能离开计算域。这时依据该边界线段物理属性进行相应的碰撞计算。2) 计算后如果粒子返回计算域，则与粒子未离开计算域情况一同计算该粒子位于计算域内哪个三角形网格单元中。如果每次都对整个计算域进行搜索的话，计算效率将非常低。本文所采用的方法是将整个计算域划分为一定密度的矩形区域，每个矩形区域信息只储存同当前矩形区域相关三角形网格单元的ID。在实际应用中，每个矩形区域平均只同 2~3个三角形网格单元相关。首先，确定粒子位于哪个矩形区域内，而后再依次判断同该矩形区域相关的三角形网格单元。该方法只比直角网格计算多一步简单的矢量叉积运算，计算效率同样很高，而且该分区思想可以非常容易应用并行计算。具体实现过程如图5所示。
|
| 图 3 网格单元内位置点的判断过程 |
|
| 图 4 粒子运动轨迹同物面边界碰撞的判断过程 |
|
| 图 5 粒子轨迹追踪过程 |
整个粒子轨迹追踪过程循环进行，当确定该粒子离开计算域或在计算域三角形网格单元内时，停止当前轨迹追踪计算，将新位置信息存入相应数据结构中。
1.3 物理模型在本文DSMC程序中，粒子之间碰撞采用变径硬球模型(variable hard sphere,VHS)计算，用Larsen-Borgnakke模型来计算粒子之间能量交换，无粘光滑壁面和对称界面采用镜面反射模型，粘性壁面则采用完全漫反射模型。空间网格平均尺度小于1/3分子平均自由程，时间间隔小于1/3粒子真实碰撞时间。
1.4 入口/出口边界条件标准的DSMC程序[14]主要用来计算高速或超高速气体流动问题。在计算域入口采用无穷远处来流的速度边界条件。喷淋头微孔内通常为亚音速流动，而且实际工况只能通过检测仪器获得压力和温度等宏观量，因此需要将压力和温度边界条件转换为传统的速度边界条件。Ikegawa等[25]和Nance等[26]提出在入口边界上采用“通量法”来计算粒子的分布速度，这种方法的缺点是程序处理上比较复杂,且容易发散； Liou和Fang[27]在“通量法”基础上在入口处采用统计平均的方法分布粒子，简化了程序，但收敛速度较慢。王沫然等[28]借用传统的计算流体动力学(computational fluid dynamics,CFD)特征化理论提出一种压力边界处理方法，提高了收敛速度。
在本文DSMC程序中，采用虚拟入口/出口网格单元的方法，具体做法是在计算域外侧将入口/出口网格单元镜向对称形成虚拟网格单元。以入口网格且单一种类气体分子为例，通过理想气体状态方程计算得到该虚拟网格单元的分子数密度ninlet，
| ${n_{{\rm{inlet}}}}{\rm{ = }}\frac{{{p_{{\rm{inlet}}}}}}{{K{T_{{\rm{inlet}}}}}}$ | (1) |
其中: 入口压力pinlet、 温度Tinlet通过实验得到； K为Boltzmann常数。
已知该虚拟网格单元的分子数密度和虚拟网格单元体积V，可计算得到该虚拟网格单元内粒子数,
| $N = {n_{{\rm{inlet}}}}V$ | (2) |
每个Δt间隔内统计得到离开该虚拟网格单元的粒子数为Nout，则进入该虚拟网格单元的粒子数Nin 与离开该虚拟网格单元的粒子数Nout相等，
| ${N_{{\rm{in}}}}{\rm{ = }}{N_{{\rm{out}}}}$ | (3) |
入口/出口网格单元宏观速度分量U为网格内所有粒子速度分量u的平均值μ′,
| $U = \mu ' = \frac{1}{{{N_p}}}\sum\limits_{p = 1}^{{N_p}} {{u_p}} $ | (4) |
进入虚拟网格单元粒子速度分量C(u,v,w)通过Maxwell速度分布应用取舍法[14]抽样得到:
| $M = \exp \left( {\frac{{ - {U^2}}}{{C_m^2}}} \right) \pm \sqrt \pi \frac{U}{{{C_m}}}\left[{1 \pm {\rm{erf}}\left( {\frac{U}{{{C_m}}}} \right)} \right]$ | (5) |
| ${f_{{\rm{in}}}}\left( C \right) = \frac{{2u}}{{\pi C_m^2M}}\exp \left\{ {\frac{{ - 1}}{{C_m^2}}\left[{{{\left( {u \mp U} \right)}^2} + {v^2} + {w^2}} \right]} \right\}$ | (6) |
其中: fin(C)为进入虚拟网格单元粒子速度分量C(u,v,w)的概率密度函数，${C_m} = \sqrt {2R{T_{{\rm{inlet}}}}/{M_s}} $为粒子最可几速度，R为气体常数，Ms为粒子摩尔数。式(5)和(6)在计算入口时取加号，计算出口时取减号。另外，两个宏观速度分量分别为v=0,w=0。
计算宏观速度时，为了减少因信息“涨落”带来的统计误差，采用了设置权重因子的方法，见式(7)。
| ${U^n} = \theta u' + \left( {1 - \theta } \right){U^{n - 1}}$ | (7) |
其中，θ为权重因子，本文取0.001。Un为当前时间步的宏观速度，Un-1为前一个时间步的宏观速度。
2 DSMC程序验证为验证本文DSMC 模拟程序的有效性，模拟计算了NASA Glenn中心的2维微机电系统(microelectromechanical system，MEMS)微喷管流动问题[29]。
MEMS微喷管喉部直径为0.3 mm，燃烧室温度为2 000 K,燃烧室压力为10 kPa,壁面温度为300 K,出口边界设为真空环境,模拟气体为氩气。碰撞模型为VHS。假设粒子在壁面发生漫反射，气体在壁面没有热量交换。速度云图和迹线图如图6所示。气体粒子流在MEMS喷管驻室内为亚音速，经过喷管平直段后被加速，粒子速度迅速跨越音速，在喷管出口处气体粒子速度为近3倍音速。DSMC模拟程序计算的沿喷管轴向速度与文[29]中算例的曲线对比如图7所示。本文DSMC模拟程序计算的结果与文[29]中算例给出的结果相比，平均偏差小于10%。
|
| 图 6 速度云图及粒子迹线 |
|
| 图 7 MEMS喷管内气流轴向速度分布 |
本文DSMC程序适用于任意边界计算域的稀薄气体动力学模型的模拟计算。下面应用该DSMC程序模拟计算了某种典型的化学气相沉积(chemical vapor deposition,CVD)反应器的喷淋头微孔孔径变化对流场分布均匀性的影响。
3.1 基本几何参数喷淋板上按一定几何形状布有几千个喷淋孔。考虑到DSMC方法对整个喷淋板建模会带来巨大的计算量，本文只对单个喷淋孔建模，这也是目前比较普遍的建模方法[27]。单个喷淋孔2维结构简图如图8所示。本文主要讨论了当喷淋孔半径r取0.5 mm、 0.75 mm和1.0 mm的3种情况，分别对应试验1、 试验2和试验3，喷淋孔其他结构参量不变。
|
| 图 8 单个喷淋孔2维结构简图 |
DSMC计算域包括喷淋孔及其出口下方3 mm区域，如图9所示。
|
| 图 9 DSMC计算域 |
采用的网格为非结构三角形网格，大约有 168 000 个。每个三角形单元又划分为3个子单元，气体粒子平均自由程小于平均网格单元尺寸，保证了粒子在每个网格单元内至少发生一次碰撞。模拟计算物理参数如表1所示。选用氩气作为试验测试气体，将喷淋孔的壁面假定为漫反射模型，粒子碰撞采用VHS模型。用入口条件初始化整个计算域，每个几何单元在初始化时平均分配有50个模拟粒子。由于轴对称几何建模过程中，近对称轴附近单元体积与远离轴单元体积差别比较大，在计算过程中会带来额外的计算开销，同时还会造成一定的统计涨落，因此对每个单元赋予不同的权值以保证每个单元在初始状态时具有相同数量的模拟粒子。将在入口处进入计算域的粒子数和出口处离开计算域的粒子数相等作为计算收敛判据。每个计算过程执行2 000 000次迭代，每执行20 000次迭代对流场宏观量进行统计采样一次。
| 物理参数 | 参数值 |
| 试验气体 | 氩气 |
| 入口压强/Pa | 200 |
| 入口温度/K | 370 |
| 出口压强/Pa | 40 |
| 出口温度/K | 300 |
| 壁面温度/K | 300 |
| 气体平均自由程(出口)/m | 4.9×10-5 |
| 时间步长/s | 1.0×10-9 |
| Knudsen数 | 1.25 |
| 网格单元数 | 168 000 |
| 粘性指数 | 0.81 |
| VHS模型粒子直径/m | 4.17×10-10 |
对于真空腔室工艺而言，最重要的是腔室内晶片上方的热场和流场分布的均匀性。为了讨论方便，借鉴膜厚均匀性的定义[31]，本文定义热场和流场的均匀性为平行晶片平面内任意点的温度、压强或速度与中心点温度、压强或速度的差异，
| $\sigma = \frac{{\left| {{c_i} - {c_r}} \right|}}{{{c_r}}}$ | (8) |
其中： ci为温度、速度或压强样本值； cr为归一化的标准值，本文cr是中心点位置相关物理参量。σ为均匀度或均匀度分布。由于喷淋头微孔计算模型具有轴对称特点，考虑温度、速度或压强分布时，本文只取气流在微孔出口下方3 mm径向线上的温度、速度和压强值，并进行了归一化处理。微孔半径为R。
3.3.1 喷淋头微孔径大小对速度场均匀性的影响气流速度分布云图如图10所示。气流在微孔平直段逐渐被加速，在出口处呈喷射状，速度达到最大，然后速度会有所降低。
|
| 图 10 喷淋头微孔气流速度云图 |
不同孔径的微孔径向速度均匀性分布如图11所示。由计算结果可以看出，3个试验的径向速度均匀性分布曲线趋势基本一致，径向速度均匀性平均偏差小于20%。越远离孔中心位置，径向速度均匀性偏差越大。减小微孔孔径，对其径向速度均匀性基本没有影响； 而增大微孔孔径，其径向速度均匀性会有所改善，但影响也不大。
|
| 图 11 径向速度均匀性分布 |
不同孔径的微孔径向压力场均匀性分布如图12所示。由计算结果可以看出，其径向压力均匀性分布曲线斜率差别较大。减小微孔孔径，其压力均匀性分布曲线斜率减小，均匀度变好； 而增大微孔孔径，其径向压力均匀性分布曲线的斜率增加较大，均匀度变差。
|
| 图 12 径向压力均匀性分布 |
不同孔径的微孔径向温度均匀性分布如图13所示。由计算结果可以看出，其径向温度均匀性平均偏差小于10%。减小微孔孔径，径向温度均匀性变差； 增大微孔孔径，径向温度均匀性变好。
|
| 图 13 径向温度均匀性分布 |
由以上讨论可以看出，微孔孔径变化，对径向速度、温度分布均匀性影响较小，而对径向压力分布均匀性影响较大。增大微孔孔径有利于提升径向速度分布和温度分布均匀性，而减小微孔孔径有利于提升径向压力分布均匀性。
4 结 论本文对基于2维Delaunay三角形非结构网格的DSMC算法进行了研究，提出的背景网格分区粒子追踪算法，具有较高的计算效率。编写了2维DSMC程序，可以模拟计算平面内任意复杂的计算域边界的稀薄气体流动问题。特别针对喷淋头微孔内气流亚音速流动特点，在入口/出口边界处理上采用了虚拟网格方法，成功地将压力、温度边界条件转化为DSMC标准边界条件。通过与经典文献算例对比，验证了该DSMC程序算法的有效性。最后，应用该DSMC程序研究了喷淋头微孔孔径变化对速度、温度和压力分布均匀性的影响。结果表明： 微孔孔径变化，对径向速度、温度分布均匀性影响较小，而对径向压力分布均匀性影响较大。
增大微孔孔径有利于提升径向速度分布和温度分布均匀性，而减小微孔孔径有利于提升径向压力分布均匀性。
对于3维非结构网格的DSMC算法研究以及喷淋头微孔结构优化设计等问题将在下一步工作中完善和提高。
| [1] | Quirk M, Serda J. 半导体制造技术 [M]. 韩郑生, 译. 北京: 电子工业出版社, 2004.Quirk M, Serda J. Semiconductor Manufacturing Technology [M]. HAN Zhengsheng, trans. Beijing: Publishing House of Electronics Industry, 2004. (in Chinese) |
| [2] | 谢欣峰. 高密度电浆用于十二吋晶圆0.22微米线宽铝内连线结构高选择比蚀刻研究 [D]. 台中: 逢甲大学, 2007.XIE Xinfeng. High Density Plasma Used in High Sensitivity Etching of Aluminum Inner Connection Structure with 0.22 Micrometer CD of 12 Inches Wafer [D]. Taichung: Feng Chia University, 2007. (in Chinese) |
| [3] | Chen I. Mass transfer analyses of the plasma deposition process [J]. Thin Solid Films, 1983, 101(1): 41-53. |
| [4] | Caquineau H, Despax B. Influence of the reactor design in the case of silicon nitride PECVD [J]. Chemical Engineering Science, 1997, 52(17): 2901-2914. |
| [5] | Nienhuis G J, Goedheer W. Modelling of a large scale reactor for plasma deposition of silicon [J]. Plasma Sources & Technology, 1999, 8(2): 295-298. |
| [6] | Sobbia R, Sansonnens L, Bondkowski J. Uniformity study in large-area showerhead reactors [J]. Journal of Vacuum Science & Technology, 2005, 23(4): 927-932. |
| [7] | Sansonnens L, Howling A A, Hollenstein C. A gas flow uniformity study in large-area showerhead reactors for RF plasma deposition [J]. Plasma Sources & Technology, 2000, 9(2): 205-209. |
| [8] | Howling A A, Legradic B, Chesaux M, et al. Plasma deposition in an ideal showerhead reactor: A two-dimensional analytical solution [J]. Plasma Sources & Technology, 2012, 21: 0150051-01-0150051-15. |
| [9] | Bird G A. Approach to translational equilibrium in a rigid sphere gas [J]. Physics of Fluids, 1963, 6(10): 1518-1519. |
| [10] | Celenligl M C, Moss J N. Application of the DSMC method to hypersonic flow about a delta wing [C]// International Symposium on Rarefied Gas Dynamics. Aachen, Germany, 1990. |
| [11] | Wagner W. A convergence proof for Bird's direct simulation Monte Carto method for the Boltzmann equation [J]. Journal of Statistical Physics, 1992, 66(3/4): 1011-1044. |
| [12] | Pulvircnti M, Wagner W, Zavelani M B. Convergence of particle schemes for the Boltzmann equation [J]. European Journal of Mechanics B: Fluids, 1994, 13(3): 339-351. |
| [13] | Bird G A. Molecular Gas Dynamics [M]. Oxford: Clarendon Press, 1976. |
| [14] | Bird G A. Molecular Gas Dynamics and the Direct Simulation of Gas Flows [M]. Oxford: Clarendon Press, 1994. |
| [15] | Bird G A. Application of the DSMC Method to the Full Shuttle Geometry [R]. AIAA 90-1692. Washington DC, 1990. |
| [16] | 樊菁, 彭世锍, 刘宏立, 等. DSMC 位置元方案中的表面元的程序标识及分子表面反射的确定论判据 [J]. 力学学报, 1999, 31(6): 671-676.FAN Jing, PENG Shiliu, LIU Hongli, et al. Program making of surface elements and the determination of molecular surface reflection in position element algorithm of DSMC method [J]. Chinese Journal of Theoretical and Applied Mechanics, 1999, 31(6): 671-676. (in Chinese) |
| [17] | 樊菁, 刘宏立, 沈青, 等. 直接统计模拟位置元方法中的分子表面反射确定论判据 [J]. 空气动力学学报, 2000, 18(2): 180-187.FAN Jing, LIU Hongli, SHEN Qing, et al. A molecular reflection deterministic criterion used in the position element algorithm of direct statistical simulation [J]. Acta Aerodynamica Sinica, 2000, 18(2): 180-187. (in Chinese) |
| [18] | 沈青, 樊菁, 胡振华, 等. 过渡领域三维绕流直接统计模拟中的一种新方案 [J]. 空气动力学学报, 1996, 14(3): 295-303.SHEN Qing, FAN Jing, HU Zhenghua, et al. A new version of position element space discretization in direct statistical simulation of three-dimensional flow in transitional regime [J]. Acta Aerodynamica Sinica, 1996, 14(3): 295-303. (in Chinese) |
| [19] | Stanton T K. Simple approximates formulas for back scattering of sound by spherical and elongated objects [J]. Journal of the Acoustical Society in America, 1989, 86 (4): 1499-1510. |
| [20] | Wilmoth R G, LeBeau G J, Carlson A B. DSMC Grid Methodologies for Computing Low Density, Hypersonic Flows about Reusable Launch Vehicles [R]. AIAA 96-1812. Washington DC, 1996. |
| [21] | Laux M, Fasoulas S, Messerschmid E W. Development of a DSMC Code on Planar Unstructured Grids with Automatic Grid Adaptation [R]. AIAA 95-2053. Washington DC, 1996. |
| [22] | 王学德, 伍贻兆, 夏健. 二维非结构网格DSMC方法实现及其应用 [J]. 南京航空航天大学学报, 2004, 36(6): 704-707.WANG Xuede, WU Yizhao, XIA Jian. Implementation of 2-D unstructured DSMC method and its application [J]. Journal of Nanjing University of Aeronautics & Astronautics, 2004, 36(6): 704-707. (in Chinese) |
| [23] | Berg M. Computational Geometry: Algorithms and Applications [M]. 北京: 世界图书出版公司, 2013.Berg M. Computational Geometry: Algorithms and Applications [M]. Beijing: World Publishing Corporation, 2013. |
| [24] | 刘汝佳, 黄亮. 算法艺术与信息学竞赛 [M]. 北京: 清华大学出版社, 2004.LIU Rujia, HUANG Liang. Algorithm Art and Informatics Competitions [M]. Beijing: Tsinghua University Press, 2004. (in Chinese) |
| [25] | Ikegawa M, Kobayashi J. Development of a rarefied flow simulator using the direct simulation Monte Carlo method [J]. JSME International Journal, 1990, 30: 463- 467. |
| [26] | Nance R P, Hash D B, Hassan H A. Role of Boundary Conditions in Monte Carlo Simulation of MEMS Devices [R]. AIAA 97-0375. Washington DC, 1997. |
| [27] | Liou W W, Fang Y. Implicit boundary conditions for direct simulation Monte Carlo method in MEMS flow predictions [J]. Computer Modeling in Engineering & Science, 2000, 4: 119-128. |
| [28] | 王沫然, 王金库, 李志信. DSMC方法的压力边界条件实现 [J]. 计算物理, 2004, 21(4): 316-320.WANG Moran, WANG Jinku, LI Zhixin. New implement of pressure boundary conditions for DSMC [J]. Chinese Journal of Computational Physics, 2004, 21(4): 316-320. (in Chinese) |
| [29] | Alexeenko A A, Levin D A, Fedosov D A, et al. Performance analysis of microthrusters based on coupled thermal-fluid modeling and simulation [J]. Journal of Propulsion and Power, 2005, 21(1): 95-101. |
| [30] | Hash D B, Mihopoulos T, Govindan T R, et al. Characterization of showerhead performance at low pressure [J]. Journal of Vacuum Science & Technology, 2000, 18(6): 2808-2813. |
| [31] | 麦克劳德. 光学薄膜技术 [M]. 周九林, 尹树百, 译. 北京: 国防工业出版社, 1974.Macleod H A. Optical Thin Film Technology [M]. ZHOU Jiulin, YIN Shubai, trans. Beijing: National Defence Industry Press, 1974.(in Chinese) |