（物理电子学专业论文）异质外延生长、表面弛豫和吸附等问题的分子动力学方法模拟研究.pdf

中文摘要 本文在简要介绍了分子动力学模拟的历史、发展、应用领域、在物理体系模 拟计算过程中的特点，以及分子动力学模拟的基本思想与基本方法后，主要介绍 运用分子动力学模拟方法，对异质外延、表面弛豫、表面吸附，以及碳纳米管与 硅样品的相互作用进行的研究。主要内容如下： 1 用分子动力学方法，模拟研究了异质薄膜的外延生长，包括有机分子的 范德瓦尔斯外延和金属之间的外延。 ， f 外延生长异质薄膜通常要求材料之间晶格匹配，然而当利用范德瓦尔斯作用 外逛生长时，晶格匹配要求显著降低。实验上已经得出了在胁s 基板上外延生 长c 。薄膜的结果，本文用分子动力学计算机模拟方法对胁是基板上外延生长 c 。薄膜进行了研究，证实范德瓦尔斯外延可以克服较大的晶格失配问题。类似 的模拟了c ，。m o s ，、c ，m o s ：等体系的外延生长，这些都属范德瓦尔斯外延。 金属之间的异质外延与金属之间的作用强弱关系很大，本文对铜、银、金之间的 外延生长进行了模拟研究，结果获得了外延形式的膜层结构，但三种金属之间存 在细微差别，这可以从作用强弱和失配度两方面去解释。1 2 模拟研究了s i ( i 0 0 ) - 2 i 、c ( 1 0 0 ) 一2 l 表面以下各层的弛豫，还对金 刚石表面吸附氢进行了动力学模拟研究。 无论是清洁表面，还是吸附表面，都存在与体材料不一样的结构。晶体表面 结构的研究，对于改进外延生长的质量和人工控制各种生长表面结构无疑都是十 分重要的。本文以研究体系总势能取极小值为出发点，通过牛顿一拉夫逊最优化 方法对s i ( 1 0 0 ) 一2 l 、金刚石c ( 1 0 0 ) - 2 l 体系作了模拟研究，得到了表面以 下五层的弛豫量，与文献上用弹性能取极小值得出的弛豫量相符。对金刚石表面 吸附氢的动力学模拟，得出了在金刚石成膜过程中表面吸附氢的稳定结构。) 3 模拟研究了碳纳米管作为s t m 针尖与s i 样品的作用。 ， f 由于碳纳米管很好的弹性力学性质以及独特的电学性质，将其作为s t m 的探 针，是一个重要的应用，有必要对碳纳米管对s i 样品表面的作用作一定的研究。 本文模拟研究了针尖与样品在各种作用距离下的作用，得出了针尖一样品作用平 衡距离；并且在常规操作距离下，研究了针尖对样品本身形貌的影响，主要表现 为样品原子受针尖吸引作用而改变位置。对针尖直接接触样品时的相互作用的模 拟研究，得出了与国外文献相似的结果。针尖在接触到样品s i 时，破坏了s i 表面的结构，从而在针尖拉离s 1 表面的时候，8 1 原子对被针尖带离3 1 表面。l 关键词：分子动力学稹拟，范德瓦尔斯作用，异质峁延，表击弛豫 表面吸附，碳纳米管，能量最优化 。， 。1 a b s t r a c t i n t h i st h e s i s ，t h eh i s t o r y , a p p l i c a t i o na r e a ，f i m d a m e n t a lt h e o r ya n dp r o c e s s ，a n d s p e c i a l t yo fm o l e c u l a rd y n a m i c s ( m d ) s i m u l a t i o na r ei n t r o d u c e da tf i r s t t h e n ，t h e s t u d yo fh e t e r o e p i t a x y , s u r f a c er e c o n s t r u c t i o na n da d s o r p t i o n ，a n dt h e i n t e r a c t i o n b e t w e e nc a r b o nn a n o t u b ea n ds i l i c o na r ed i s c r i b e db ym dm e t h o ds i m u l a t i o n t h e m a i ns i m u l a t e dr e s u l t sa r es h o w na sf o l l o w s ： 1 u s i n gm ds i m u l a t i o n ，h e t e r o e p i t a x yi n c l u d i n ge p i t a x i a lg r o w t hp r o c e e d i n g w i t hv a nd e rw a a l si n t e r a c t i o n sa n dh e t e r o e p i t a x yb e t w e e nm e t a l sh a v eb e e ns t u d i e d o n eo ft h em a j o ro b s t a c l e st or e a l i z i n ge p i t a x i a lh e t e r o s t r u c t u r ei st h el a t t i c e m a t c h i n g b e t w e e nc o n s t i t u e n tm a t e r i a l s h o w e v e r w h e nt h ee p i t a x i a lg r o w t h p r o c e e d sw i t hv a i ld e rw a a l si n t e r a c t i o n s ，t h el a t t i c e m a t c h i n gc o n d i t i o ni sd r a s t i c a l l y r e l a x e d i ti sr e p o r t e dt h a tt h ee x p e r i m e n to fc 6 0g r o w i n go nm o s 2s u b s t r a t eh a s b e e ns u c c e s s f u l l yc a r r i e do u t i nt h i st h e s i s ，t h em ds i m u l a t i o nm e t h o dw a sa p p l i e dt o t h i sk i n do fh e t e r o g e n e o u se p i t a x yt op r o v et h a tv a nd e rw a a l se p i t a x yw a ss u i t a b l e f o rt h eh i g h l yl a t t i c e m i s m a t c h i n gs y s t e m a l s o ，c 7 0 m o s 2a n dc 5 7 m o s 2s y s t e m s w e r es i m u l a t e d ，w h i c ha r ea l s oh i g h l yl a t t i c e m i s m a t c h i n gs y s t e m sp r o c e e d i n gv a n d e rw a a l si n t e r a c t i o n s i n t e r a c t i o np l a y sa l li m p o r t a n tr o l ei nh e t e r o e p i t a x yb e t w e e n m e t a l s ，c u ，a ga n da uw a ss e l e c t e dt op e r f o r mt h es i m u l a t i o ns t u d ya n dg o tf c c ， s t r u c t u r eo ft h ef i l m t h ed i f f e r e n tr e s u l t sa m o n gf i l m so ft h e s et h r e em e t a l sc o u l db e e x p l a i n e db yi n t e r a c t i o na n dm i s - m a t c h i n g 2 t h es u r f a c er e l a x a t i o no fs i ( 1 0 0 ) - 2 1a n dd i a m o n dc ( 1 0 0 ) - 2 1h a v eb e e n s t u d i e d ，a n dt h ea d s o r p t i o no fh y d r o g e no nd i a m o n d ( 111 ) s u r f a c eh a v eb e e n s i m u l a t e db ym dm e t h o d t h es t r u c t u r eo fs u r f a c e ，w h e t h e rc l e a no ra d s o r b e do n e ，w a sd i f f e r e n tf r o mt h a t o fb u l km a t e r i a l s t h es t u d yo ft h es u r f a c es t r u c t u r ei si m p o r t a n tf o re p i t a x i a lg r o w t h a n ds u r f a c eo fa l lk i n d so ft h i nf i l m s i nt h et h e s i s ，r e l a x a t i o nb e l o wt h es u r f a c eo f s i ( 1 0 0 ) 一2 1 a n d c ( 1 0 0 ) 一2 1 w e r es t u d i e d b yn e w t o n r a p h s o ne n e r g y m i n i m i z a t i o nm e t h o d m ds i m u l a t i o no fc o m p l e t eh y d r o g e nc o v e r a g ed i a m o n d ( 111 ) s u r f a c ec o n f i r m e dt h et h e o r yt h a tt h es u r f a c eh a sb u l k 。t e r m i n a t e d ( 1 i ) s t r u c t u r e w i t hd a n g l i n gb o n d st e r m i n a t e db yh y d r o g e na t o m s ，c ( i o o ) 一( 1 1 ) h 3 t h ei n t e r a c t i o nb e t w e e ns i l i c o na n dc a r b o nn a n o t u b ea st h et i po fs t mw a s s i m u l a t e d d u et oi t sg o o de l a s t o m e c h a n i c a lp r o p e r t i e sa n du n i q u ee l e c t r i c a l p r o p e r t i e s ， c a r b o nn a n o t u b eh a sb e e nu s e da st h et i po fs t ms ot h es t u d yo fi n t e r a c t i o nb e t w e e n c a r b o nn a n o t u b ea n ds i l i c o ns a m p l ei sn e c e s s a r y i nt h el a t e rp a r to ft h i st h e s i s ，t h e i n t e r a c t i o nb e t w e e nt h et i pa n ds i l i c o ns a m p l ew e r es t u d i e dw h e nt h e yw e r ea ta l l k i n d so fd i s t a n c e s ，a n dt h er e s u l to ft h et i p s a m p l eb a l a n c ed i s t a n c ew a sg o tf r o mt h e s i m u l a t i o n w h e nt h et i pa n ds a m p l ea r ea tn o r m a ld i s t a n c e ，c h a n g e st oa t o mp o s i t i o n o ft h es a m p l ew e r es t u d i e d f i n a l l y , am o d e lw a sm a d et os t u d yt h ei n t e r a c t i o nw h e n c a r b o nn a n o t u b et i pb a r e l yt o u c ht h es u r f a c e ，w i t ht h er e s u l ts a m ea st h a to ff o r e i g n r e s e a r c ht e a m w h e nt h et i pt o u c h e dt h es a m p l e ，i td e s t r o y e dt h es i l i c o ns u r f a c e s t r u c t u r em a k i n gs o n ! es i s ib o n d sr u p t u r e d s ow h e nt h et i pp u l l e db a c kf r o mt h e s u r f a c e ，ad i m m e ro fs i l i c o na t o m sw a sp i c k e du pb yt h et i p k e y w o r d s ：m o l e c u l a rd y n a m i c ss i m u l a t i o n ，v a nd e rw a a l si n t e r a c t i o n ， h e t e r o e p i t a x y , s u r f a c er e l a x a i o n ，a d s o r p t i o n ，c a r b o nn a n o t u b e s ， e n e r g ym i n i m i z a t i o n 第一章前言 梗旦大擎顽士论支 第一章前言 1 分子动力学计算机模拟方法概述 十九世纪中叶以前，物理学基本上是- f 实验科学，因而，实验物理学成了 物理学的基本分支。二十世纪初，随着量子力学和相对论的诞生，理论方法和数 学物理方法日益重要，使得理论物理学成为- - i 独立分支。理论物理研究中，在 给出合适的简化物理模型后，以主要研究手段是解析演绎推导方法为其基本特 点。但是随着科学技术的高速发展，使所研究的体系越来越复杂，传统的解析推 导方法已不敷应用；研究对象的复杂性和向极端条件( 极高压、极高温、极低温 等) 发展，实验方法也将因为设备的昂贵、条件的危险性而变得十分困难，二十 世纪中叶计算机的发明与计算技术的突飞猛进，大容量快速计算机的出现为物理 研究提供了有效的手段，以计算机为必不可少的工具，在合理的物理模型下，以 对物理问题进行数值研究和数值模拟为主要研究方法的“计算物理学”应运而生， 成为物理学的第三大分支 1 。计算物理有的称为模拟物理或计算机模拟物理。 对物理体系进行计算机模拟，有以下几点优势： ( 1 ) 计算机模拟方法研究一个具体问题时，并不需要像解析推导研究方法那 样，对具体问题的物理模型做某种必要的近似，所以计算机模拟适宜于研究 更接近实际的复杂问题。 ( 2 ) 计算机模拟方法可对模型和实验进行比较，从而为评估一个模型的正确 与否提供了一种手段 2 。 ( 3 ) 计算机模拟方法可以给出某些实验上无法或者难以测量的量，从而沟通 理论和实验。 ( 4 ) 计算机模拟可以用图形显示，因而具有形象直观的优点。 计算机模拟的方法主要分成两大类：一类是随机性模拟方法，例如蒙特卡罗 方法( m c ) ，实现了g i b b s 的统计力学途径；第二类是确定性方法，例如分子动力 学方法( m d ) ，实现了b o l t z l l l g r l n 的统计力学途径。分子动力学方法是通过数值解 n 个分子的运动方程得到这些分子每个时刻的坐标与动量，即相空间中的轨迹， 苎二主堕主 堡! 兰! 堑三生生 从而用统计方法计算出体系的宏观性质。 分子动力学方法从1 9 5 7 年首次运用到物理学领域至今，已经有将近半个世 纪的时间。二十世纪五十年代和六十年代，有四篇代表性论文，可以作为分子动 力学发展初期的里程碑： ( 1 ) 1 9 5 7 年，a l d e r 发表了第一篇报道分子动力学的论文，文章采用刚球模 型描述气体相互作用，粒子通过瞬间碰撞相互作用，两次碰撞之间粒子直线运动。 程序在i b m7 0 4 机器上运行。 ( 2 ) 1 9 6 0 年，g i b s o n 发表第一篇基于时间积分的连续势能的m d 计算的论文， 研究由辐射损伤引起的缺陷，模拟了一个5 0 0 原子的体系，同样是在i b m7 0 4 机器上运行，每步要执行1 分钟。即使从今天看来，模拟结果也是很好。 ( 3 ) 1 9 6 4 年，r a h m a n 在研究液体a r 时采用l e n n a r d j o n e s 势，模拟了一个 8 6 4 个原子的体系，程序在c d c3 6 0 0 计算机上运行。并且r a h m a n 的代码较为通 用，被后来的很多分子动力学程序引用。 ( 4 ) 1 9 6 7 年，l o u dv e r l e t 同样采用l e n n a r d j o n e s 势研究液体a r ，并且 计算出对相关函数，用来研究结构并验证液态理论，模拟中所采用的时间积分算 法后来被称为“v e r l e t 算法”。 半个世纪以来，分子动力学已经成为一种有力的模拟研究的方法，取得了 系列重要的成果。特别是近年来，由于计算机技术的迅速发展，使得分子动力学 方法在表面界面问题研究、薄膜生长、超晶格材料、纳米材料、材料力学等方面 得到了广泛的应用 3 5 。下面仅列举部分应用： ( 1 ) 表面研究：二十世纪八十年代以来，由于微分析实验手段的发展，表面 物理兴起，对于理解表面的特殊现象( 例如表面再构、表面熔化、表面扩散、表 面吸附等) ，分子动力学提供了有力的研究手段。 ( 2 ) 晶体与薄膜生长过程的研究：成膜动力学研究可以获得薄膜生长过程中 的薄膜结构、成膜机理、实验条件对成膜的影响等。 ( 3 ) 团簇研究：从几个到几千个原子组成的团簇，是分子系统和固体之间的 桥梁，由于其有限的体积，物理性质会与固体材料不一样。特别是金属团簇经常 可以作为化学反应的催化剂，是研究的热点。 ( 4 ) 生物有机分子：蛋白质动力学研究，对于在原子水平上理解蛋白质性质 以及分子设计有很大帮助。 2 墨二苎堕兰 堡：! 兰望主丝圭 ( 5 ) 液体：包括液体的相图研究、液体的疏运特性研究。 2 本论文工作简述 本论文首先简要介绍了分子动力学模拟在物理体系模拟计算过程中的特点， 以及分子动力学模拟的基本思想与基本方法后，着重介绍运用分子动力学模拟方 法，对异质外延、表面弛豫、表面吸附，以及碳纳米管作为s t m 针尖与硅样品的 相互作用进行的研究。主要内容简述如下： 一用分子动力学方法，模拟研究了有机分子的范德瓦尔斯外延和金属之间 的异质外延。 外延生长异质薄膜通常要求材料之问品格匹配，然而实验上已经得出了大的 晶格失配下坳& 基板上外延生长c 。薄膜的结果，这是利用范德瓦尔斯作用进行 的外延生长。本文用分子动力学计算机模拟方法对坳& 基板上外延生长c 。薄膜 进行了研究，证实范德瓦尔斯外延可以克服较大的晶格失配问题。类似 c ，。m o s ，、c ，m o s ，等体系的外延生长，都属范德瓦尔斯外延，本文也对此 进行了模拟研究。金属之间的异质外延匹配要求较高，而且成膜方式与金属之间 的作用强弱关系很大，本文对铜、银、金之间的外延生长进行了模拟研究，结果 获得了外延形式的膜层结构，但各种外延膜之间存在细微差别，这可以从作用强 弱和失配度两方面去解释。 二对s i ( 1 0 0 ) 一2 l 表面的再构和以下几层原子的弛豫进行了模拟研究， 还对金刚石c ( 1 0 0 ) - 2 1 表面弛豫以及金刚石薄膜表面吸附氢进行了动力学模 拟研究。 晶体表面存在与体材料不一样的结构，其结构的研究，对于改进外延生长的 质量和人工控制各种生长表面结构都是非常重要的。自从1 9 5 9 年报道了s i ( 1 0 0 ) 表面2 l 再构以后，有过很多的研究，有人运用单纯形法对表面弹性能作最小 化计算得出表面以下几层的弛豫量。本文以研究体系总势能取极小值为出发点， 通过牛顿一拉夫逊最优化方法对s i ( 1 0 0 ) 一2 xl 、金刚石c ( 1 0 0 ) 一2 l 体系作了模 拟研究，得到了表面以下五层的弛豫量。对金刚石表面吸附氢的动力学模拟，说 明了氢吸附在金刚石成膜过程中起着重要的作用，可以防止金刚石薄膜石墨化。 第一章前言 摸g 戈擎顿士话支 三模拟研究了碳纳米管作为s t m 针尖与硅样品的作用。 碳纳米管作为s t m 的探针是一个重要的应用，有必要对碳纳米管对s i 样品 表面的作用作一定的研究。本文模拟研究了针尖与样品在各种作用距离下的作 用，得出了针尖一样品作用平衡距离；并且在常规操作距离下，研究了针尖对样 品本身形貌的影响，并且给出了样品表面层原子受针尖作用位置的改变值。对针 尖直接接触样品时的相互作用的模拟研究，得出了与国外文献相似的结果。针尖 在接触到样品s i 时，使得一些s i _ s i 键断裂，再经过新的再构，靠近针尖的s i 原子对与针尖c 形成s i c 键，从而在针尖拉离s i 表面的时候，s i 原子对被针 尖带离s i 表面，实现针尖对样品的原子操作。 4 苎三主坌三塑塑兰兰苎i 堕丝查兰 堡! 兰! _ _ 生_ ! ：兰! 生 第二章分子动力学计算机模拟方法 1 分子动力学模拟的基本方法 分子动力学方法的出发点是物理体系的确定的微观描述。这种描述可以是哈 密顿描述或拉格朗日描述，也可以是直接用牛顿运动方程表示的描述。分子动力 学方法是通过运动方程来计算体系的性质，结果得到的既有体系的静态特性，也 有动态特性。 分子动力学方法通过直接数值求解构成体系的n 个分子的运动方程，得到这 些分子每个时刻的坐标与动量，即相空间中的轨迹，从而用统计方法计算出多体 体系的宏观性质。 在分子动力学中，研究体系中的每个分子都遵循牛顿第二定律： 即 ( 2 1 1 ) ( 2 1 2 ) 式( 2 1 1 ) 中，f 是分子i 受到的其它分子总的作用力，它由分子所处的势能函数 v ( r ) 决定。口，是分子f 的加速度。 方程( 2 1 2 ) 有不少数值解法，如v e r l e t 算法、f r o gl e a p 算法、龙格一库塔法、 欧拉法等。 在本论文的分子动力学模拟方法中，首先根据具体的研究体系，给出体系中 每个分子的初始位置，和由设定的体系平衡时的温度，给出每个分子的初始速度 值，即相空间中的n 个分子的初始点。然后，在合适的时间步长下，对体系中每 个分子按运动方程( 2 1 2 ) 随时间作动力学运动模拟，使体系逐步趋于平衡，直 至由体系的总能量决定的温度达到设定的体系平衡时的温度( 目标温度) 为止， 此时研究体系结构稳定，该平衡状态下的体系结构即为所求的结构。 r 一2 一出 一 m 聊 = 耻 塑打 d 一 苎三兰坌三型垄兰兰苎垫坚塑查兰堡! ! ! 堡兰丝一墨一 2 势能模型一原子间作用势 从f 2 1 2 ) 式可以看出分子动力学方法模拟物理过程的关键部分是原子f f i 相互 作用势的选取。考虑量子力学和相对论效应的分子运动的完全数学描述是非常艰 难的任务，通常分子动力学模拟研究中势场的选取是基于经验数据的，这些经验 数据包含了相对论效应和量子效应。 最简单也是最常用的一种相互作用势是l e m l a r d j o n e s 作用势，在很多固体 和液体中描述原子与原子相互作用时，和经验数据符合的较好。l e n n a r d j o n e s 作用势包含两个部分，一个是吸引项，另一个是排斥项： y ( r ) = 4 e ( t 7 ) ”一( 马6 】 rr ( 2 2 1 ) 式中，r 为两原子间距离，f 处于势阱最深处的v ( r ) 值，仃为v ( r ) 等于零时的r 值。由此可以看出l e n n a r d j o n e s 作用势是一种两体作用势，实际体系中只考虑 两体作用是不够的，必须考虑多体作用 6 。 多体作用势考虑粒子位置、两粒子相互作用、三体作用等，如下式表达： 矿( r ) = k ( t ) + ( r i ，r j ) + v 3 ( r , ，r j ，r k ) + ( 2 2 2 ) lj j i 1 7 1k j i 式中第一项代表体系处于外场中，后面的各项分别表示两体作用、三体作用 实际体系中一般三体作用以后的各项都可以忽略。 对于所要研究的有机分子、无机晶体、金属等体系，作用势的一般表达式可 以分五项： v = v 6 册4 + 矿嗍咖+ v + v ”删6 + v 打( 2 2 3 ) 总的势能v “包含键能v b “、角能矿删。、扭矩能v ”、非成键能y ”6 、氢键能 矿”。对于具体的体系会选用各自不同的、合适的表达式。 一般力场中的几种典型的势能项如下图： 6 第二章分子动力拳计算机模拟方法 2 旦克擎顽士话盘 酋g * 矾粕r 转翻s w - j i l l 降n 1 2 v 神翻w a 甜s 圜 绀阿一蚓 ”* 凇 圳m 警 q 线虹盼抽喇k 咖 织赫姜气k ， 卜- 啪细i 腑妇1 0 i i n g 小茹渗 噜j 。硒一瓤- e 棚口粕嗍 ，、 啄扛；铷，一略 一獬 蚰肿辩龋i 螺 烈 。k 丢女。一鼓 秘0 p 吖删蒋a 嗡毫尉鲰i 。o + 一砷 囊 t 豫d t 埘叫抖黼删辞啪辫斟对礴悄- n 轴静舶t 图2 ，1 几种典型的势能项 3 运动方程的数值解法 分子动力学通过数值解微分方程来计算粒子的运动轨迹 7 。体系中的粒子 在符合牛顿方程的规律下运动。经典的牛顿运动方程如下式： _ d r ( t ) ：了d 2 r ( t ) ：盟( 2 3 1 ) 出d t 2m 一一 其中，是作用于该粒子的力，m 是该粒子质量，v 是该粒子在，时刻，位置的速 度。时刻f a t 的位置由t a y l o r 展开法可得到： ，( t + a t ) _ ，掣，+ 圭学驴+ d ( ，3 ) ( 2 3 2 ) 将r o + f ) 和，o 一f ) 相加，可得到t + 出时刻的位置( 由f 一缸时刻和f 时刻的位 置值给出，且忽略高次项) ： r ( r + f ) ：一r ( f 一f ) + 2 r ( f ) + ：史2 f 2( 2 - 3 3 ) m 苎三主坌王兰生兰生兰塑壁苎查鲞 堡：查兰壁l 兰丝 这就是v e r l e t 算法 8 。v e r l e t 算法中，速度v ( ，) 由下式得出： v ( f ) = r ( t + a 1 t ) - 万r ( 一t - a t ) 另外一种是l e a p f r o g 方法 9 一l o ，表达式如下 帕+ 尹a t 叫，一争等， r 。+ r ) = r o ) + v ( f + 一a 2 t ，a t 而t 时刻的速度v ( t 1 由下式得出： ( 2 3 4 ) ( 2 3 5 ) ( 2 3 6 ) 砸，：翠叫睁三等， b ，刀 如果要求精确度更高一些，可以采用a b m 4 ( a d a m s b a s h f o r t h m o u l t o n f o u r t ho d e r ) 方法，这是一种预测一校正法，精度达到时间步长的五次方。这种 方法要利用前三步的结果，所以不能自动开始计算前三步的结果由 r u n g e - k u t t a 方法产生。由于存储了更多信息，所以需要更多的计算机内存。a b m 4 算法如下： 预测量为： 啪埘) 甜( r ) + 尝( 5 5 v 一5 9 v ( 卜岔) + 3 7 v ( 卜2 出) - 9 v ( 卜3 a t ) ) + o ( a t 5 ) ( 2 3 8 ) 由预测量0 + ，) 计算出速度v ( h f ) ，然后得到校正量： r o + r ) = ，o ) + 9 v ( f + f ) + 1 9 v ( f ) 一5 v o 一f ) + v o 一2 f ) ) + o ( a t 5 ) ( 2 3 9 ) 4 实际模拟过程中的几个关键的技术问题 分子动力学模拟过程中，根据具体研究体系以及所需要的物理结果、精度等 的不同，要注意以下一些细节问题。 苎三兰坌三茎垄兰兰墨垫堡塑查鲞 ! 兰三兰! l 二堡坐 4 1 边界处理 模拟体系的边界如何处理? 一种可能是什么都不用做，因为体系在此终结， 边界的原子比体内原子近邻原子少，这一般是研究团簇时碰到的情况。对于其它 情况，不管模型有多大，其所包含的原子数n 都远远小于宏观物质的1 0 2 3 的数量 级，这使得模拟体系中表面原子与所有原子的比率比实际体系大得多。 解决的办法是用周期性边界条件，此时粒子处于一个盒子中，此盒子朝三个 方向周期性的平移，可以让体系达到无穷。这个盒子称为分子动力学模拟的元胞。 于是元胞中一个位置为r 的粒予，其实代表了一组无穷多的粒子，这些粒子的共 同特征是： ，= ，0 + f a + m b + 珩。 ( 2 4 1 ) 式中，l ，t o ，n 是从一。到+ 一的整数，a ，b ，c 是元胞的三个边的矢量。于是 当元胞中一个粒子的移动，就代表了一组粒子的移动，而在模拟的程序中，只需 要计算一个粒子，就可以得到一组粒子的位置和速度。元胞的选取应该合适，不 能太大，否则粒子数太多计算时间太长；也不能太小，否则误差太大容易给出错 误结果。 实际上元胞的边界并不是真正的边界。计算势能时，不光要考虑元胞中的粒 子与粒子的作用，还要考虑元胞中的粒子与相邻元胞中的粒子的作用。而且，元 胞中的粒子有可能跑出元胞，当然这不会是体系的粒子数变化，因为一个粒子从 一个面跑出去，必然会有另一个粒子从相对的面跑进元胞，这就能使粒子数守恒。 4 2 势能截断和长程近似 势能计算中最耗费内存的项是非成键作用项，因为对于一个n 原子的体系， 存在n 木( n - 1 ) 2 这样的作用。所以为了节省计算时间，当两个原子的距离超过某 个值r ，( 比如l o a ) 时，往往忽略不计这种非成键作用。特别是对于范德瓦尔斯 作用，当两个原子距离1 0 a 时，作用力几乎为零。但是，如果粒子带有静电荷， 比如体系包含带正负电的离子，非成键项中包含了静电作用，这种作用能反比于 距离r ，是一种长程作用，所以截断距离r 应该大一些。 苎三兰坌三黧查兰! 墨塾堡塑查鲎 垒! ! ! 一塑! 坠 简单的势能截断会产生这样的问题：势能在r 。处不连续。为了解决这个问题 可以把势能平移一下，让势能在r 。处为零，以l j 势为例： 附) ： 纵力疋 ( 。r c ) 但4 2 ) ( r r 。) 、 以上处理虽然能使势能连续，但是势能的一阶导数依然不连续，也就是力不 连续。对于精度要求高的体系，需要力连续。对此可以做如下所示的处理： e 10 o r s ( o 一+ - - 、j 、 、lr c r i 厂 | - 一 、v ( r ) s ( r )沙 f 一一 图2 2 势能截断的校正 s ( r ) 称为转换函数，须满足： 1 ，在势能变化剧烈处保持恒定( 例如等于1 0 ) ； 2 ，在截断距离r 。处光滑的趋于零； 3 ，距离大于r ，时为零。 转变后的势能为v ( r ) s ( r ) ，在r 。处， l ，v ( r ) s ( r ) = v ( r 。) - s ( r 。) = 0 ，势能连续； 2 ，( y ( r ) s ( r ) ) = v + ( r 。) s ( r 。) + 矿( r 。) s 。( r 。) = 0 ，力连续。 苎三兰坌! 兰垄兰茎墨塾堡垫查兰 堡三兰! 二二生j 盟 4 3 粒子的初始位置、速度以及时间步长的选取 分子动力学模拟的初始位置要合理的给出，这样容易让体系较快的趋于平 衡。可以对体系在温度为零时的势能用能量最小化算法，给出合理的结构，也就 给出了粒子的初始位置。 粒子的初始速度可以由体系的温度控制。粒子的速度符合 m a x w e l 卜b o l t z m s n n 分布： m ， f ( v 肛( 刍儿面4 扩咖 r 2 4 3 ) 上式表示，质量为m 的粒子在温度为t 时，速度为v 的几率为f ( v ) 。 时间步长的选择要合适。一方面，时间步长应该尽可能的大，这样可以在一 定的模拟步数内得到大的时间范围的体系状态；另一方面，由于截断了t a y l o r 高次项，所以时间步长应尽可能的小，才能得到足够精确的粒子的位置和速度。 如果采取精确到3 次项的算法，可以取体系原子振动的最高频率作为时间步长， 对于所模拟的体系，时间步长取a t ：l f s 。 5 粒子的运动状态的确定以及体系宏观性质的统计 获得体系平衡时的热力学性质是分子动力学的目的之一。一个宏观变量a 沿 着计算轨迹的值为a ( t ) ，作时间平均就可以得到a 的热力学值： 4 = l i m ，+ 。彳1t j 4 ( f ) 衍 ( 2 5 1 ) 宏观变量4 可以是粒子位置和动量的任何函数。通过时间平均，可以得到体系的 一阶性质( f i r s t o r d e r p r o p e r t y ) ，例如内能、动能、压力等。 体系的瞬态动能由下式决定： k ( f ) = 去m 。【v ，o ) 】2 ( 2 5 2 ) 苎三兰坌兰丝生兰盐墨塾堡丝主兰 堡! 查! 生主丝 由动能可根据能量均分定理得出温度： ，：三生( 2 5 3 ) 3 n k 8 式( 2 5 - 3 ) 中，k 为体系温度，k 。为玻耳兹曼常数。 位能一般可以通过对关联函数表示： u = 2 矽f “( ，) g ( ，) r 2 d r ( 2 5 4 ) 式( 2 5 ，4 ) 中，p 是体系的密度，为体系的总粒子数，g ( r ) 为对关联函数，g ( r ) d r 表示在离，= o 的粒子，处单位体积元发现另一粒子的几率。若n ( r ) 为r 到r 十甜之 间的平均粒子数，则 g ( r ) = 万v 丽n ( r ) ( 2 5 5 ) 式( 2 5 5 ) 中，v 为体积。对关联函数g ( r ) 的计算可以通过下式 时，= 古( 姜，煮盯。，) ( 2 5 6 ) 苎三兰坌! 兰垄兰苎苎壑堂垫主兰 堡! 查! 堡l 兰坠 6 分子动力学方法模拟流程 本论文分子动力学计算使用m s i ( m o l e c u l a rs i m u l a t i o n si n c ) 公司的m a t e r i a l s s t u d i o2 0 软件。首先使用v i s u a l i z e r 模块搭建所要研究的体系结构，并使用 d i s 。o v e r 模块选择合适的力场建立体系的能量表达式，再运用最速下降法、牛顿 法等对体系作能量最优化处理，作为分子动力学模拟的初始结构。然后设定体系 的温度、时间步长、终止步数、截断半径等，开始进行动力学计算，最后可以对 所得的结构进行分析。图2 3 是模拟计算中的主要流程图。 图2 3 分子动力学方法模拟主要流程图 第二章分子动力肇计算机模拟方法 槎g 土擎顽 论支 7 本章小结 本章首先简单介绍了分子动力学模拟的历史、发展、应用领域，以及在物理 体系模拟计算过程中的特点，具体介绍了在实际模拟计算中要主要的细节问题， 以及不同体系的势能表达形式，给出了联系微观和宏观的物理量的表达形式，并 且以框图的形式给出了本论文分子动力学模拟的步骤。 苎三兰墨些竺墨生：! 兰圭苎堡壑 堡! ! 兰兰一生! 蔓l 1 引言 第三章薄膜的异质外延生长模拟 异质外延生长是指不相同的材料相互之间外延生长，能够形成自然界中没有 的人工结构材料。异质外延生长已经成为薄膜技术中的一个越来越重要的方面。 但是实验上只能获得有限的一些外延薄膜，原因是材料间的晶格需匹配，才容易 进行外延生长。晶格匹配的要求，严重限制了外延材料的组合。在半导体材料如 g a a s 和s i 中，晶格匹配限制尤其严格，因为共价键的长度和方向不容易改变( 见 图3 1 a ) 。离子晶体材料中，晶格匹配的条件会放松一些【1 1 - 1 2 。这时晶格匹配 条件通过一个厚度约几纳米的过渡层来实现( 见图3 1 b ) 。如果异质外延生长是 通过范德瓦尔斯力作用而结合的话，晶格匹配条件可以明显放宽( 见图3 1 c ) ， 典型的例子是在种表面没有悬挂键的基体表面上生长另一种层状材料，这种类 型的外延生长被称为范德瓦尔斯外延0 3 ，1 5 。 ( a ) ii ii iiii 1 1 毒 寸一 中1 一 呻争一 一范德瓦尔斯沟 ( c ) 图3 1 膜层间的界面 ( a ) 为共价键( b ) 为离子键( c ) 为范德瓦尔斯作用 过瀵罄 ( b ) 苎三主苎堕竺墨堕! ! 苎兰苎堡垫 堡! ! 丝兰兰二! 衬底和薄膜之间的失配度，被定义为： 厂：坠型 ( 3 1 1 ) “j 其中，d 。、d ，分别为衬底和薄膜的晶格常数 1 4 1 。 2m o s ：基板上外延生长c 。薄膜的分子动力学模拟 2 1 实验结果和理论意义 在m o s 2 的基板上外延生长g 。薄膜是晶格失配情况下异质外延生长的一个 典型例子，实验上已经由m a s a t o s h is a k u r a i 等取得了成功 1 6 】。由于c 6 0 的球状 分子内部各键已经饱和，g 。分子与分子之间是依靠范德瓦尔斯力结合在一起的， 据报道其晶体是面心立方结构，分子间距( 晶格常数) 为1 0 0 2 埃【1 7 1 。这一长 度远远大于一般基板材料的晶格长度，很难找到与其晶格匹配的基板，然而尽管 m o s ：与c 6 0 分子间晶格失配度大约为3 ，通过r h e e d 检测证明c 6 0 充分利用范 德瓦尔斯外延生长的优点，成功的生长在m o s 2 上面。经过仔细分析r h e e d 图 案，m o s 2 上面生长的是具有h c p 结构的层状c 矗。测量薄膜中c 岛分子间距为 1 0 0 5a ，与c 0 晶体的分子间距相当接近。 本文应用分子动力学模拟方法来研究m o s ：基板上外延生长c 。薄膜的过程， 从而验证了在范德瓦尔斯作用下确实能够在晶格失配条件下进行异质外延的薄 膜生长。 2 2 物理模型的建立和模拟计算的基本步骤 2 2 1 物理模型的建立 本文的主要工作是模拟c 6 。在m o s 2 基板上成膜时的合理结构与分布。主要 分两个步骤实现： 1 系统准备过程。 ( 1 ) 分子动力学计算初始位置的确定： 墨三兰墨垦竺墨堕盐垩兰苎堡丝 堡! 查! 堑生笙生 用能量最小化方法对c 。随机分布在m o s ：表面的体系进行优化，得到系 统的一个近似合理结构，作为分子动力学计算初始位置值； ( 2 ) 分子动力学计算初始速度的确定： 确定目标温度，即系统平衡时的温度。根据此温度，给出一组符合麦克 斯韦分布的速度值，赋给系统各原子作为分子动力学计算的初始速度值。 2 系统平衡过程。 从上面的系统准备过程，可知在分子动力学模拟开始的瞬问，由系统各个原 子的动能总和计算出的温度等于最终系统平衡时的温度，但是此时系统的结构显 然不是此温度下的最终平衡结构，系统按照牛顿第二定律随时间进行动力学运 动。系统在趋向平衡过程中温度会出现大的波动，甚至偏离平衡温度，但是最终 达到平衡时，温度会在目标温度处上下作微小波动，保持稳定。此时c 。m o s ：体 系结构稳定。 2 2 2 势函数的确定 本文研究对象c 。m o s ：体系的经验势y ( r ) 的表达式采用下面的形式 1 8 ： y ( r ) = 【五2 ( 6 - b 。) 2 十k 3 ( 6 - b 。) 3 + k 。( 6 一b o ) 4 6 型 + h 2 ( 口- 0 0 ) 2 + h j ( 曰一o o ) 3 + h 4 ( 目- 0 0 ) 4 】 口 型 + 巧【1 - c o s ( 一升) 】+ 【1 一c o s ( 2 一趟) 】十( 1 - c o s ( 3 砂一对) 】 型 + k ；x 2 + 。( 6 - b 。) ( 6 一6 ；) 十( 曰- 0 0 ) ( 0 一或) jbb p 口 丝盟盟 + 只，( 6 - b 。) ( 占- 0