（材料学专业论文）氮化碳微观结构及团簇沉积的分子动力学模拟.pdf

江苏大学硕士学位论文 摘要 1 3 c 3 n 4 是首先在理论上预测的，结构类似于b s i 3 n 4 ，硬度有可能超过金刚 石的一种非天然共价化合物。在原子、分子尺度上研究其微观结构变化及团簇沉 积c n x 薄膜的动力学过程对探讨薄膜生长的物理机制，改进薄膜质量是至关重要 的。 本论文对晶态、液态和非晶态c 3 n 4 的微观结构进行了分子动力学模拟。模拟 结果显示：在c 3 n 4 晶体中，c 和n 原子的配位数分别为3 9 2 和3 0 7 ，c 为四面体 配位；液态c 3 n 4 中，c 和n 原子的配位数分别为3 6 1 和2 5 2 ，键角和键长的分布 范围明显变宽，结构处于无序状态；非晶c 3 n 4 的结构中c n 键的键角和键长分 布与晶态c 3 n 4 基本一致，说明构成非晶c 3 n 4 与晶态c 3 n 4 的结构基元相同。非晶 c 3 n 4 中c 和n 的配位数分别为3 4 6 和2 2 6 ，比c 3 n 4 晶体下降，表明结构中石墨 型构型增加。计算得到1 3 c 3 n 4 的体弹性模量为4 3 6 1 g p a ，与l i u 采用可变晶格模 型分子动力学c s m d ) 从头计算法计算得到的体弹性模量4 3 7 g p a 基本一致。 用分子动力学方法和第一性原理方法模拟了单个c 、n 原子以及c 2 、n 2 、c n 团簇在硬质合金表面的沉积过程，为小原子团簇与表面的作用搭建了两种不同的 吸附模型，即物理吸附模型和化学吸附模型，基体温度分别取为2 9 8 k 、5 0 0 k 和 7 0 0 k 。模拟中，计算了不同基体温度下团簇与表面的相互作用能，并对团簇在衬 底表面的吸附结构、吸附能以及成键情况进行了讨论；还研究了原子( 团簇) 在 s i ( 1 0 0 ) 表面和高速钢表面的吸附机理，以探讨不同衬底材料对c n 。薄膜生长的 影响。模拟结果显示：c 、n 原子与三种晶面都成化学吸附，但c 、n 在高速钢表面 的吸附最稳定，硬质合金表面对n 原子的选择性吸附要明显强于c 原子，这对硬 质合金表面的氮化和提高薄膜中氮的含量极为有利：w c ( 1 0 0 ) 表面由于c 2 的吸附 具有金属性，c n 吸附在w c 0 0 0 ) 表面上时在w 和c n 之间存在比较强的离子键。 在硬质合金、单晶硅和高速钢表面，模拟出了对生长c n x 薄膜有利的衬底温度分 别为5 0 0 k 、7 0 0 k 和2 9 8 k 。 关键词：c 3 n 4 ，计算机模拟，分子动力学，团簇，薄膜，态密度 江苏大学硕士学位论文 a b s t r a c t 1 3 。p h a s e c a r b o nn i t r i d ei saa r t i f i c i a lc o v a l e n c ec o m p o u n d l i ua n dc o h e n i n d i c a t e dt h a ti t sm i c r o s t r u c t u r ei ss i m i l a r 、析t 1 1 1 3 p h a s es i l i c o nn i t r i d ea n di t sh a r d n e s s m a ye x c e e di nd i a m o n d t od i s c u s sp h y s i c a lm e c h a n i s mo ft h i nf i l mg r o w t ha n d i m p r o v et h et h i n f i l mq u a l i t y ，i ti sc r u c i a lt os t u d yt h e1 3 - p h a s ec a r b o nn i t r i d e s m i c r o s t r u c t u r e sa n dd y n a m i c so f c a r b o nn i t r i d et h i nf i l m s d i s c o v e rm o d d ei nm a t e r i a l ss t u d i oh a sb e e na p p l i e dt os t u d ym i c r o s t r u c t u r e so f c r y s t a l 1 i q u i da n da m o r p h o u sp h a s e so f c a r b o nn i t r i d e t h es i m u l a t e dr e s u l t si n d i c a t et h a t t h ec o o r d i n a t i o nn u m b e ro fci nt h es t r u c t u r eo fc 3 n 4i s3 9 2a n dt h a to fni s3 0 7 t h e c o o r d i n a t i o np o l y h e d r o no fci st e t r a h e d r o n i nt h el i q u i ds t a t eo fc 3 n 4 ，t h ec o o r d i n a t i o n n u m b e ro fca n dni s3 61a n d2 5 2r e s p e c t i v e l y t h ed i s t r i b u t i o no fb o n da n g l ea n db o n d l e n g t hw i d e n sr e m a r k a b l ya n dt h es t r u c t u r ei sd i s o r d e r e df o rl i q u i dc 3 n 4 乃ed i s t r i b u t i o n o fb o n da n g l ea n db o n dl e n g t hi na m o r p h o u sp h a s ec 3 ；n 4i ss i m i l a rt ot h a ti nc r y s t a lc 3 n 4 i ti n d i c a t e st h a tt h ec o o r d i n a t i o np o l y h e d r o ni nc r y s t a lc 3 i n 4i ss a m ea st h a ti na m o r p h o u s p h a s ec 3 n 4 t h ec o o r d i n a t i o nn u m b e ro fca n dn i na m o r p h o u sp h a s ec 3 n 4i s3 4 6a n d 2 2 6r e s p e c t i v e l y , w h i c hi sl e s st h a nt h a ti nc r y s t a ls t r u c t u r e t h ed e p o s i t i o nc o u r s e so fc ，na t o ma n dc 2 ，n 2 ，c nc l u s t e r si nw c ( 1 0 0 ) s u r f a c e a r es i m u l a t e db ym o l e c u l a rd y n a m i c sa n dt h ef i r s t - p r i n c i p l e sp s e u d o p o t e n t i a lt o t a le n e r g y m e t h o d s p h y s i c a la n dc h e m i c a la d s o r p t i o nm o d e la r eb u i l d e dt os t u a yt h ec l u s t e r s ( a t o m s ) i n t e r a c t i o n sw i t hu n d e r l a ys u r f a c ei n2 9 8 k ，5 0 0 k ，7 0 0 k t h ei n t e r a c t i o ne n e r g yi s c a l c u l a t e da n d c l u s t e r s a d s o r p t i o nm e c h a n i s mi nw e ( 10 0 ) s u r f a c ei ss t u d i e d t od i s c u s st h e t h i nf i l m sg r o w t hi nd i f f e r e n tu n d e r l a ys u r f a c e ，t h ea d s o r p t i o ns t r u c t u r e 、a d s o r p t i o ne n e r g y a n db o n d i n gb e t w e e na t o m sa r ea l s os t u d y e di ns i ( 10 0 ) a n dh i g h - s p e e ds t e e ls u r f a c e t h e r e s u l t ss h o wca n dna t o mc a nb ea t 乜a c m dc h e m i c a l l yb yw c ( 10 0 ) ，s i ( 10 0 ) a n dh s s ( 1l0 ) $ l l r f a c e ，b u tca n dna t o ma d s o r b i n go nh s s ( 11 0 ) i ss t a b l e s t na t o ma d s o r b i n go n w c ( 10 0 ) s u r f a c 圮i sb e t t e rt h a nc a t o m t h ew e ( 10 0 ) s u r f a c ea d s o r b e db yc 2c l u s t e rs h o w s m e t a lp r o p e r t y t h e r ei sai n t e n s ei o m cb o n db e t w e e nwa n dc nw h e nc nc l u s t e ra d s o r b o nw c ( 11 0 ) s u r f a c e i ti sa v a i l a b l et oc n xt h i nf i l m sg r o w t ho nw c ( 11 0 ) s u b s t r a t ei n5 0 0 k , o ns i ( 11 0 ) s u b s t r a t ei n7 0 0 ka n do nh s s ( 1l o ) s u b s t r a t ei n2 9 8 t c k e yw o r d s ：c a r b o nn i t r i d e ，c o m p u t e rs i m u l a t i o n , m o l e c u l a rd y n a m i c s ，c l u s t e r s ， t h i nf i l m ，d e n s i t yo fs t a t e i i 学位论文版权使用授权书 本学位论文作者完全了解学校有关保留、使用学位论文的规定，同 意学校保留并向国家有关部门或机构送交论文的复印件和电子版，允许 论文被查阅和借阅。本人授权江苏大学可以将本学位论文的全部内容或 部分内容编入有关数据库进行检索，可以采用影印、缩印或扫描等复制 手段保存和汇编本学位论文。 保密亩，在年解密后适用本授权书。 本学位论文属于 不保密口。 学位论文作者签名：复j 寸互参 加年f 月 叮时 ( 2 1 5 ) 式中，盯为分子直径。只有当原子间距离小于o 2 时才有力的作用，当距离 大于o 1 2 时原子之间没有相互作用，这时粒子的运动轨迹为直线。这种模型形式 简单，数学处理容易，对有效短程斥力的相互作用给出了有效的近似。 ( 2 ) 连续对势 l e n n a r d j o n e s ( l j ) 势【2 5 】 口( 俨4 气 ( 马坨一( 翌) s 】 (216)vi 啄俨4 咄宁1 2 等) 6 】 q 式中，毛表示相互作用力的强度的参数，气是表示原子大小的参数。 b o r n m a y e r - h a g g i n s ( b m h ) 势 w = 孚坪唧 半】一軎一争 泣 m o r s e 势【2 6 1 l ( 勺) = 彳p 之口吩一引一2 e 一口，i 一局】 ( 2 1 8 ) j o h n s o n 势【2 7 】 西 ，( ) = 一4 ( ，；，一e ) 3 + q 勺一乜 ( 2 1 9 ) l j 势是为描述惰性气体分子之间相互作用力而建立的，因此它表达的作用力 较弱，描述的材料的行为也就比较柔韧：b m h 势是用来表征碱金属卤化物、碱土 金属卤化物等离子结合型物质的原子间相互作用的；m o r s e 势与j o h n s o n 势经常用 8 江苏大学硕士学位论文 来描述金属固体。 2 2 - 2 多体势 2 0 世纪8 0 年代初提出了可弥补对势在金属中使用的缺点的多体势函数。d a w 和b a s k e s 2 8 】基于密度泛函近似( d f a ) 的理论提出了嵌入原子模型( e m b e d e d a t o m m o t l l o d ，简称e a m ) 及相应的各向同性n 体e a m 势，应用于f e c 金属、b c c 金 属、以及几乎是d 壳层完全填满的过渡金属，获得了成功；m w f i n n i s 和 j e s i n c l a i r 2 9 】基于t i g h t b i n d i n gm o d e l 的二次矩近似提出了经验多体势，应用于 b c c 金属和过渡金属得到了良好的结果。r a j o h n s o n 【3 0 】开发了对f c c 金属只考虑 最近邻相互作用的e a m 势的解析函数，但是，此模型对所有金属给出了具有相同 的各向异性比的结果。另外，k w j a c o b s e n 3 1 】等人立足于更基本的理论思想提出 了有效介质理论( e f f e c t i v em e d i u mt h e o r y ) 模型。所有这些势函数可全部用下式 来表达： ( 1 ) e t o l = 去( 勺) + f ( ) ( 2 2 0 ) t = l = l i = 1 式中，e ( r , j ) 是对势，n j 是表示原子在位置f 的局域电子密度，f 相当于局域电子 密度的能量。 ( 2 ) s t i l l i n g e r - w a b e r 势 3 2 】 对于s i 和g e 等半导体，其键合强度依赖于周围原子的配置。对此类原子仅 考虑二体力是不够的，必须引入三体力及其以上的所谓多体力，1 9 8 5 年，s t i l l i n g e r 和w a b e r 提出了s w 3 2 】势，其表达式如下： ( 。) 2 萎占f ( 軎) + i 莓 i t f g ( 詈，罟，詈) ( 2 2 ) i i l ， i u uu 式中第一项为二体力，第二项为三体力。 ( 3 ) t e r s o f f 判3 3 】 t e r s o f f 势是键合强度依赖于周围原子配置的势，其基本函数为莫氏函数，可 依据键合强度与配位数的关系来构造，此函数由下式给出： ( 乃) = f c ( ，；，) 【以e x p ( 一九) 一气e x p ( - a , j r , # ) ( 2 2 2 ) i 勺 式中，乃是原子f 和原子的键合( i _ ，键) 距离，z ( 0 ) 是相互作用截断函数。 9 江苏大学硕士学位论文 2 3 分子动力学模拟的初始条件 分子动力学模拟计算是一个初值问题，选择合理的初始构型和初始速度将使 得模拟体系快速驰豫到所要求的平衡状态。系统的初始位置和速度可以通过实验 数据，或理论模型，或者两者的结合来决定。如果被模拟的系统具有初始密度分 布n ( r ) 和温度分布t ( 0 ，而没有固定的晶格结构，则每个原子的位置可以从初始密 度分布n ( r ) ，用舍入法或m e t r o p o l i s 等方法得到。每个原子的初速度，则可以从初 始温度下的m a x w e l l b o l t z m a n 来随机选取；如果系统有固定的晶格结构，则把n 个粒子的初始位置放在晶格的格点上，然后加以扰动。 2 4 周期性边界条件 由于计算机的运算能力有限，模拟计算系统的粒子数不可能很大，这就会导 致模拟系统粒子数少于真实系统的所谓“尺寸效应”问题。为了选取尽可能多的粒 子数而又不至于使计算工作量过分庞大，在统计物理中，对于平衡态分子动力学 模拟计算引用三维周期性边界条件，在这个具有自由边界的三维n 个粒子的体系 中处于界面的分子数正比于n 。m 。 图2 1周期性边界条件示意图 f i g 2 1p e r i o d i cb o u n d a r yc o n d i t i o n s 在使用有限的原子数来模拟实际体系中原子的运动时，就必须考虑表面对体 系中原子运动的影响。为避免这种影响，可以通过周期性边界条件来实现。将含有 n 个粒子的体积当作具有与其相同单元的无限周期点阵的原始单元( 如图2 1 是一 1 0 江苏大学硕士学位论文 个二维的示意图，在这个二维图象中每个单胞被其他的8 个单胞所包围；三维方 向每个单胞就会被2 6 个单胞包围) ，一给定粒子则与在此无限周期体系中的全部 其他粒子相互作用 3 4 1 。 根据周期性边界条件理论，模拟是在一个盒子( b o x ) 中进行的，原子从盒子 的一边溢出就会从与这个边相对的那个边进入盒子中，这样所有的原子无时无刻 都是处于这个盒子中的，也就是说这个盒子不存在自由表面。 理论上说，任何对称结构的几何图形都可以作为周期性边界条件，但为了简 单和方便起见，一般选择立方盒子作为周期性边界条件，此外，还可以通过关闭 或打开在x ，y ，z 三个方向上的周期性边界条件来模拟三维块状体系；模拟二维 薄片体系时，则只打开其中任意两个方向的周期性边界条件；要模拟一维长条体 系，可以通过只打开其中任意一个方向的周期性边界条件来实现。 2 5 运动方程积分 分子动力学模拟方法是属于统计物理学中的一种计算方法，在m d 方法处理 过程中方程组的建立是通过对物理体系的微观数学描述给出的。在这个微观的物 理体系中，每个分子都各自服从经典的牛顿力学，每个分子运动的内窠力学是用 理论上的哈密顿量或者拉格朗日量来描述，也可以直接用牛顿方程来描述。从上 面讨论的方程形式可以看出，分子动力学方法的基本方程是线性的或非线性的二 阶常微分方程，对此人们研究了多种方法。 1 9 6 7 年，y e r l e t 3 5 1 提出了一种算法。其公式为： ，o + a 0 2 厂( f ) 一，o 一f ) + ( a t ) 2 u f ( t ) ( 2 2 3 ) 聊 1 ，“) r ( t + a t ) - r ( t - a t ) ( 2 2 4 ) 2 a t 式( 2 2 3 ) 右边第一、二项是( 出) o 的量级项，而第三项是( 址) 2 的量级 项。这样若将大的绝对值项和小的绝对值项混合，则随着重复计算次数的增多， 其误差发生积累。 为了防止v e r l e t 算法中的误差积累，h o n e y c u t t 提出了所谓蛙跃( l e a p - - - f r o g ) 算法【3 6 】。有关公式如下： 江苏大学硕士学位论文 川+ 争m 一a t + 半 ( 2 2 5 ) r ( t + 百a t ) - 一r ( t ) _ v ( r + 等 ( 2 2 6 ) 址 、 2 ，0 + a t ) ：，( f ) + f d o + i a t ) ( 2 2 7 ) 式( 2 2 5 ) 和式( 2 2 7 ) 中，速度和位置矢量的变化依赖于步长等。正如以上 二 所述，如果使用v e r l e t 法或蛙跃法，则可以求解关于n e v 系综的线性常微分方程， 但不能求解除此之外的n t v 、n h p 和n t p 等系综给出的非线性常微分方程。 s w o p e 提出的v e l o c i t y - - v e r l e t 算法【3 刀可以同时给出位置、速度与加速度， 并且不牺牲精度。这种算法的优点是给出了显式速度项，且计算量适中，目前应 用较广泛。 g e a r 算法【3 8 1 的基本要点是先用泰勒展开，由 r i n 彳导至l j r i n + l ，随后计算力， 求得修正后的位置和速度。其预测步长为： r f ( t + a t ) = a r , ( t ) ( 2 2 8 ) 校正步长为： ( f + a t ) ：足o + 址) + ! 会 厂 r o + f ) 卜足o + f ) ( 2 2 9 ) 式( 2 2 8 ) 和式( 2 2 9 ) 中，址为时间步长，4 为由泰勒展开系数组成的矩阵， a 为根据精度和计算稳定性要求所确定的系数。2 6平衡系综的控制方法 2 6 1 调压技术 分子动力学模拟通常要求体系的压力保持不变，而体积可以随温度的变化而 改变，也就是说要进行等压模拟。压力通常比一些量如总能量波动范围要大的多。 这是因为压力与维里系数有关，而维里系数是位置与势能函数导数的乘积，它随 维里系数的变化比内能快的多。 、 宏观系统通过改变它的体积来保持压力为常值，体积变化量与等温压缩系数 有关，易压缩的物质体积波动较大；反之，在常体积模拟中，压缩系数越小的物 质，压力的波动越大。在等压模拟中，可以通过改变模拟原胞的三个方向或一个 方向的尺寸来实现体积的变化，常用的有以下三种调压技术： 1 2 江苏大学硕士学位论文 ( 1 ) b e r e n d s e n 方法【3 9 】 这种方法是假想把系统与一“压浴”相耦合，取模拟原胞体积的标度因子为旯， 而原子坐标的标度因子为兄m ，于是： 兄= l + r 盟( p 一) ( 2 3 0 ) 。p 式中，咋是耦合参数，忍耐是“压浴”的压力，p ( t ) 是f 时刻的真实压力。新的位 置由r j = 2 v 3 r , 给出，此式可以用于各向同性或各向异性。应用于各向异性方法，可 使得原胞的三个方向的尺寸能独立地变化。 ( 2 ) a n d e r s o n 方法【1 7 】 在压力耦合的系统方法中，一个额外的自由度，即原胞的体积，首先被 a n d e r s o n 引入到系统中。这个自由度( 相当于一个活塞作用在此系统上) 的动能 是q ( d v d t ) 2 2 ，这里q 是活塞的质量”，活塞的势能为p v ，其中p 是期望的 压力，v 是系统的体积，质量小的活塞引起原胞的快速波动，而质量大的则相反。 在模拟过程中，体积可以变化，而平均体积则由内部的压力与期望值达到平衡时 决定。 2 6 2 调温技术 在一个分子动力学模拟中通常需要对温度进行控制，因为开始进行模拟时， 按初始温度的麦克斯韦尔一波尔兹曼分布对每个原子赋予了一个初速度，也就是 说体系获得了初始动能。在模拟体系经过驰豫达到平衡态的过程中，系统中多余 的势能会转化为动能，从而使体系温度升高，这就需要在运行过程中对温度进行 调解。温度是一个热力学参量，体系的有效温度是由动能来决定的： k = 三翻v i l 2 一竿 3 ， 式中，厂为体系的自由度，为玻尔兹曼常数。 对一个由n 个粒子组成的非周期体系，自由度厂= 3 一6 ，所减去的6 个自由 度为质心的平移和旋转，则： r ： 兰茎 ( 2 3 2 ) ( 3 n 一6 ) k b 而对周期体系，质心旋转自由度不能忽略，自由度f = 3 n 一3 ，则： 江苏大学硕士学位论文 丁： 三丝 ( 2 3 3 ) ( 3 n 一3 ) k b 目前对温度的控制有4 种方式，分别如下所示： ( 1 ) 速度标度【4 0 】 改变温度最简单和最直观的方法是直接对速度进行标度。如果，时刻的温度是 t ( t ) ，速度乘以因子兄后，温度的变化为： 扯- 1 ) 矾刮急 3 4 ， 控制温度最简单的方法是在每一步乘上标度因子兄。这里，z r e q 是期望的参考 温度。 ( 2 ) b e r e n d s e n 热浴【3 9 】 这种保持温度的方法是b e r e n d s e n 在1 9 8 4 年提出的，它是假想系统和一个温 度为期望值的虚拟热浴相接触，对速度每一步进行标度，以保持温度的变化率与 热浴和系统的温差( 瓦枷一丁( f ) ) 成比例，每一步的温度变化为： 缸= 里( k 一荆) ( 2 3 5 ) 因此，速度的标度因子为： 冉+ 孚 焉_ 1 】 3 6 ， 当耦合参数f 等于步长o t 时，这个算法就等同于简单的速度标度法。f 的合 适值大约为o 4 p s ，若步长仃f 为l f s ，则给出里0 0 0 2 5 。这种方法的优点是它允 许系统在期望的温度值上下波动。 ( 3 ) g a u s s i a n 热浴 这种方法的基本思想是在运动方程中加上一项“摩擦力 项，这个力又与系 统的速度联系在一起，运动方程如下： 可m , d 2 r , = 互一彘v j ( 2 3 7 ) 这里，磊= ( _ 吒) ( 唯毗) ，以保持总动能守恒。为了把温度调到期望 的温度下，需要将初始的温度设定为期望值。 1 4 江苏大学硕士学位论文 ( 4 ) n o s e h o o v e r 热浴 1 9 4 2 】 此种方法和g a u s s i a n 方法采用相同的运动方程，但f 值给出的方式不同： 警= 一石1 。i t - 1 】 ( 2 3 8 ) 当钿寸0 0 时，孝= 0 ，就相当于绝热系统，而当一所，f = 彘，就等价 于g a u s s i a n 方法。t n h 的取值要根据实际模拟情况来选取。 2 7 时间步长( t i m es t e p ) 的选取 在m d 模拟中，时间步长的选取是很重要的，它在某种程度上决定了模拟所 需的时间，以及计算误差。一般来讲，为了减少误差，步长应取得小一点，但是 太小的话，系统模拟进程会很慢，这需要积累一定的经验。目前看来，对于一些 研究比较透彻的分子，如氩，水和甲烷等分子，时间步长已经有了公认的合理值， 通常为飞秒级。 2 8 模拟结果的分析方法 由分子动力学基本方程出发进行直接计算的输出结果是关于原子的位置坐标 和速度。但是，在材料科学中，不仅要很好地知道原子水平上的三维结构，而且 还必须计算出表征宏观物性( 热力学性质、光学性质等) 的物理量值。 物理参数值可根据位置坐标和速度通过统计处理求出。在统计物理中，可利 用系综微观量的统计平均值来计算物性参数值，即： ( 彳) 棚= 击j 弦洲( ) p ( p ，) ( 2 3 9 ) 式中，e r i s 是系综( e n s e m b l e ) 的英文缩写，p ( p ，) 是分布函数( 几率密度) 。几 率密度因所考察的系综不同而不同。 在分子动力学方法中，人们使用了所谓的时间平均等于系综平均的各态历经 ( e r g o d i c ) 假设，即： ( 彳) 鲫= ( 彳) 腑 ( 2 4 0 ) 在热力学统计物理学中，各态历经假设没有严格地证明，它的正确性是由它 的推论被实验结果所证实而被肯定。 江苏大学硕士学位论文 ( 彳) 咖。l i m 。1 - - ，a ( p ，) d r = 砌 ( 2 4 1 ) o b s 。o b s 时，= 文 泣4 2 ， 一矗 2 8 2 配位数( c o o r d i n a t i o nn u m b e r ) = 2 r 4 风g ( ，) d r ( 2 4 4 ) 对双体分布函数从到k 进行积分( r o 指g ( ，) 曲线_ k 第- - 个非零点所对应的 ，值，指g ( ，) 第一峰所对应的，值) 可得到最近邻配位数n ，式( 2 4 4 ) 是采 用对称的方法计算配位数，其前提是假定双体分布函数的第一峰对称。 1 6 江苏大学硕士学位论文 2 8 3 角分布函数 角分布函数恰好能弥补双体分布函数的特点，因而对于确定一个体系的微观 结构也非常重要【4 6 4 8 1 。在以任意一个粒子为中心、半径为，的球面上，对位于球面上 的任意两个原子，它们与中心粒子间夹角为秒的几率是a ( c o s o ，r ) s i n 8 ，那么定义角分布函 数a ( c o s 0 ，) 为： 2 9力场 彳( c o s 9 , r ) 2 瓦1 荟n 蔷。( 。( 一坝c o s 秒一等) ( 2 4 5 ) c = o ( r 吨) o ( ，一) ( 2 4 6 ) k = lf ( j 任何模拟的关键部分就是力场的选择。力场是原子尺度上的一种势能场，由 一套势函数与一套力常数构成，是用来精确计算分子结构和能量的计算方法，它 通过原子核的位置来计算体系的能量。力场近似描述了原子核运动上的势能超曲 面，它常被用来调整特有体系。因此，力场的选择取决于所要研究的结构类型。 用量子力学和相对论对一个分子进行全面精确描述，这是一个强有力的挑战。 理论上讲，相对量子理论可以完美描述物质的所有性质，但是，由于分子体系是 微观粒子多体问题，应用量子力学处理问题时必须进行适当近似，根据 b o r n o p p e n h e i m e r 近似，有如下s c h r s d in g e r 方程： h f ( r ，) = 胖( r ，) ( 2 4 7 ) 这里，日是h a m i l t o n 算符，甲是波函数，e 通常称为势能面，仅是原子核的 函数，甲( r ，- ) 是原子核坐标尺和电子坐标，的函数。 式( 2 4 7 ) 描述了势能表面的原子核运动，求解此式的努力是极其巨大的，因 此，经常对势能表面作经验处理，称之为力场。由于原子核是相对重的对象，量 子力学可忽略不计，因此公式( 2 4 7 ) 可以用牛顿运动方程代替： 一= d v ：了m d 2 r ( 2 4 8 ) 一一= l 4 j d r d t | 对势能表面作经验处理后的公式( 2 4 8 ) 称为分子动力学。 与一个力场结合的结构坐标就建立了一个能量表达式( 或目标函数) ，这个能 1 7 江苏大学硕士学位论文 量表达式描述了作为原子坐标函数的特定结构的势能表面： = + k 朋+ e 佣一6 删 ( 2 4 9 ) 式中，玩耐为系统总能量，下标c r o s s t e r m 表示由近邻原子所引起的键或角度扭转， 一。训为非键作用能。 2 10m a t e r i a l ss t u d i o 软件 本课题所采用的模拟软件是m a t e r i a l ss t u d i o ( m s ) 4 9 】一它是美国a c c e l r y s 公司在2 0 0 0 年推出的新一代分子模拟软件，是专门为材料、化学领域研究者开发 的一款可在p c 上运行的模拟软件。能帮助研究人员解决当今化学及材料工业中的 许多重要问题，使得任何研究人员可以轻易获得与世界一流研究机构相一致的材 料模拟能力。可以建立三维结构模型，对各种晶体，无定型以及高分子材料的性 质及相关过程进行深入的研究。m s 可以进行复杂的动力学模拟和先进的量子力学 计算，并得到切实可靠的数据，它采用了m i c r o s o f t 的标准用户界面，允许用户 通过各种控制面版直接对计算参数和计算结果进行设置和分析。其中的m a t e r i a l s v i s u a l i z e r 是m s 的核心模块，它提供了搭建分子、晶体及高分子材料结构模型所 需要的所有工具，运用这些工具可以操作、观察及分析结构模型；计算和显示重 要的结构参数；处理图表、表格或文本等形式的数据。 d i s c o v e r 模块【5 0 1 是m s 软件的计算引擎。它使用多种分子力学和动力学方法， 以仔细推导的力场为基础，可准确地计算出最低能量构型，分子体系的结构和动 力学轨迹等；它解释分子结构和行为的关系，提供对分子相互作用的了解，预测 液体、固体和气体的关键性质。 c a s t e p 5 l 】模块是由剑桥大学凝聚态理论研究组开发的一套先进的量子力学 程序，可进行化学和材料科学方面的研究。c a s t e p 使用基于密度泛函理论( d f t ) 的总能量一平面波赝势方法，在材料科学的数学模型中，离子势被赝势( 即只作 用于系统价电子的有效势) 所替代。电子波函数通过一平面波基组扩展，电子一 电子相互作用中的交换和相关效应通过局域密度近似( l d a ) 或广义梯度近似 ( g g a ) 得以包括。结合赝势和平面波基组的应用，使对体系中所有原子上的作 用力的计算变得极为容易，这使得对分子、固体、表面及界面的离子构型的有效 优化成为可能。基于总能量赝势方法，c a s t e p 根据系统中原子的类型和数目， 1 8 江苏大学硕士学位论文 即可预测出包括晶格常数、几何结构、弹性常数、能带、态密度、电荷密度及光 学性质在内的各种性质。 综上所述，可见m a t e r i a l ss t u d i o 软件满足本课题的要求。所有计算工作都是 利用m s 软件平台在i n t e lpi vp c 机上完成。 2 11分子动力学模拟在材料科学中的应用 1 9 9 6 年，c h u b a c i 5 2 】等人利用自开发的动力学程序t r i m 模拟了c 、n 离子束 在s i ( 1 0 0 ) 面沉积的物理过程，模拟中离子束入射能量为1 0 k e v ，c 离子束入射量 为1 6 x1 0 1 7 i o n c m 2 ，n 离子束入射量为2 0 x1 0 1 7 i o n c m 2 ，c n 比率为o 8 ，模拟 结果与c - n 膜的x p s 数据显示的c n 比率相当一致，还对c n 比率为0 6 和1 0 时的成膜过程进行了模拟。1 9 9 8 年，f r a u e n h e i m 删等人使用基于第一性原理的紧 束缚分子动力学方法研究了各种c n 体系的行为，范围包括在化学气相沉积金刚 石里的掺杂行为和n 在化学气相沉积金刚石生长上的影响以及超硬c 3 n 4 晶体可能 结构的形成等。2 0 0 0 年，吕华威 5 4 】采用s t i l l i n g e r - w e b e r 二体和三体作用势，模拟 了衬底温度为1 0 0 0 k 时团簇沉积c n 薄膜的过程。分子动力学研究显示，团簇初 始能量和衬底温度对团簇在衬底表面的扩散迁移能力以及所沉积薄膜的结构产生 巨大的影响。 郭玉宝【5 5 】等人使用m a t e r i a l ss t u d i o 程序包中的c v f f 力场，对甘氨酸分子在 单壁纳米管中的吸附和扩散行为进行了分子动力学模拟，并对甘氨酸分子在纳米 碳管中的构象和能量进行了优化。模拟计算结果表明：甘氨酸在纳米碳管中的构 象发生了伸缩和扭转，这种构象的改变将会导致氨基酸生物性能的改变；纳米碳 管对氨基酸分子具有较强的吸附作用，其中纳米碳管和氨基酸分子之间的一兀 相互作用增加了纳米碳管对氨基酸的吸附能。t s a j 【5 6 】采用基于第一性原理的分子 动力学方法，研究了g a n 表面气态g a 和n 原子的吸附特征。模拟结果显示：引 入的g a 原子相对不活动，而引入的n 原子却很活泼；g a 在g a n 表面的吸附系数 要比n 原子高很多。 h uy a n h o n g 【5 刀等人研究了在金刚石表面分子束沉积薄膜的成核和生长过程， 考虑了c 2 h 4 团簇、c l o h l 6 分子和c 2 0 这三种入射物质作为对比，随机入射角度分 别为0 0 、1 5 0 、4 5 0 和6 0 0 ，入射能量为4 0 0 e v 。结果显示：c 2 0 最易形成无定形薄 1 9 江苏大学硕士学位论文 膜，而分子团簇c 2 h 4 效果最差；尽管三种入射物的化学键不同，但在入射能量为 4 0 0 e v 时它们沉积的薄膜的结构却很相似。r e s e n d e s s 】等人使用分子动力学方法 模拟了c u l 3 团簇在c u ( 0 1 0 ) 表面的沉积过程。模拟结果显示：元胞尺寸对模拟 结果的影响很小，即使用小尺寸模拟也能获得可靠的结果和对沉积现象的定量描 述；还计算了8 0 0 k 下团簇中心的扩散系数，结果表明团簇中心的扩散系数比单个 粒子下降了1 1 0 。s e u n g h y e o bl e e ( 5 9 1 等人使用t e r s o f f 势能函数研究了无定形碳 薄膜的原子结合特征。模拟结果表明：入射能量在5 0 7 0 e v 之间时，获得了高致 密、s p 3 杂化和有残留应力的薄膜，这与实验观察到的t a c 薄膜相一致；在距基板 2 1 a 时原子的沉积数量最大；由碳原子碰撞而产生的局部热量能有效地分散在晶 格上时，就能获得高硬度和致密的无定形碳薄膜。这些工作的不断深入开展必将 对沉积高质量的c n x 薄膜有着深远的意义。 2 12 本文研究目的和意义 薄膜技术是材料表面改性、新材料制备和研究的主要技术手段，磁控溅射镀 膜是其中被广泛使用且有效的方法之一。在溅射镀膜中，成膜粒子到达衬底后结 合成膜的过程及其机理一直是人们所感兴趣的研究内容，这一领域的理论和实验 研究都非常活跃。本论文结合课题组课题磁控溅射法制备c n x 薄蒯6 0 铆】的研究， 采用分子动力学方法和密度泛函理论，对b c 3 n 4 的微观结构、原子及原子团簇与 基体表面的相互作用作了如下分析： 1 研究b c 3 n 4 的微观结构变化情况 在微观领域材料的硬度由体积模量决定，而体积模量又依赖于材料中化学键 的特性，正是由于化学键的强度和压缩系数在材料的抗形变中起着决定性的作用。 因此，通过分子动力学研究晶态、液态和非晶态c 3 n 4 的键长、键角和原子配位数 变化情况，对c n x 薄膜结构分析与检测具有一定的理论意义。 2 团簇在基体表面吸附机理的研究 原子和分子团簇，简称团簇( c l u s t e r s ) 是由几个乃至上千个原子、分子或离 子通过一定的物理或化学结合力组成的相对稳定的微观和亚微观聚集体。一方面， 人们可以利用团簇实现薄膜制备，如离化团簇沉积【6 2 】。实验上已发现可以利用低 能团簇制备出致密、光滑、膜基结合力优异的薄膜材料；另一方面，在薄膜制备 2 0 江苏大学硕士学位论文 过程中，薄膜的生长经常受到原子团簇的影响，甚至影响合成薄膜的质量，如溅 射过程中载能粒子轰击靶材表面时，溅射出的各种原子和团簇等。因此，研究团 簇在基体表面的吸附机理无论是对利用低能团簇制备优质薄膜，还是控制团簇对 薄膜的影响，均具有重要意义 6 3 - - 6 6 】。 2 l 江苏走学硕士学位论文 第三章氮化碳微观结构的分子动力学模拟 1 3 c 3 n 。作为新- - 4 超硬材料和功能材料将日盏受到重视。但由于氯化碳是一 种新材料，目前的技术尚不能合成单相的样品，其结构研究受到很大制约，特别 是对于非晶态氯化碳的结构，有许多问题有待研究。而分子动力学计算机模拟是 研究复杂凝聚态系统的有力工具，不仅能得到原了的运动细节，还可以提供实验 上尤法获得或很难获得的信息，它是联系物质微观信息和宏观性质的一种计算方 法，对c n 。薄膜结构分析与检测具有重要意义。 3 1 动力学模拟基本过程 实验采用d i s c o v e r 模块对晶态、液态和非晶态c n 4 的微观结构进行分子动力 学模拟。模拟中我们以p s i 3 n 4 的结构为基础，建立目c n 4 的晶胞结构，如图 31 所示。并调用d i s c o v e r 模块对其进行结构优化， 随后以优化后的构型为基础， 在x 、y 和z 方向分别扩展2 倍、2 倍和5 倍，得到一个超品胞，如图3 2 所示， 它包含1 4 2 2 5 - 2 8 0 个粒子。以下的模拟工作部在这个超晶胞e 进行。 在对0c ，n 4 的模拟巾采用等温分子动力学力法，使用剧期性边界条件，应用 v e r l e tv e l o c i t y 方祛求解牛顿运动方程组，超晶饱中原子的初始速度服从麦克斯韦 尔玻尔珐曼分布，时间步长为20 f s ，初始位置由o k 时的c 3 n 4 晶体结构构型确定， 体系温度控制采用a n d e r s o n 温控法。 圈318 c 3 n a 的晶胞结构 f i g 31u n l tc e l l o fb c 3 n 圈3 2bc 、n 。的超晶胞 f i g3 2 s u p e r c e l lo f sa n 江苏大学硕士学位论文 3 2 模拟结果与分析 3 2 1 晶体c 。n 结构分析 晶态模拟以o k 时的构型为基础，在3 0 0 k 驰豫1 0 ，0 0 0 步达到平衡，每隔1 0 0 0 步输出一次体系的构型，并记录每一步粒子的瞬态位置和速度，以计算体系的径 向分布函数、键角分布和其它性质，从而研究c 成键机理【6 8 】。模拟结果如下： - 4 2 0 0 0 4 2 2 0 0 - j b 2 0 0 - 4 3 4 0 0 o5 0 0 01 0 0 0 01 o o s i m u l a t i o nt i m e f s