分子动力学在材料科学中的应用

分子动力学在材料科学中的应用

摘要：本文综述了几种常见条件下的分子动力学模拟方法以及分子动力学模拟的最新发展趋势.介绍用分子动力学模拟方法研究固休的休相结构,表面问题,界面问题以及薄膜形成过程等方面的研究成果。

关键词：分子动力学; 计算机模拟; 材料科学

1引言

分子动力学（Molecular Dyanmica,简称MD）用于计算以固体、液体、气体为模型的单个分子运动，它是探索各种现象本质和某些新规律的一种强有力的计算机模拟方法，具有沟通宏观特性与微观结构的作用，对于许多在理论分析和实验观察上难以理解的现象可以做出一定的解释［1］。MD方法不要求模型过分简化，可以基于分子（原子、离子）的排列和运动的模拟结果直接计算求和以实现宏观现象中的数值估算。可以直接模拟许多宏观现象，取得和实验相符合或可以比较的结果，还可以提供微观结构、运动以及它们和体系宏观性质之间关系的极其明确的图象［2］。MD以其不带近似、跟踪粒子轨迹、模拟结果准确［3］，而倍受研究者的关注，在物理、化学、材料、摩擦学等学科及纳米机械加工中得到广泛而成功的应用。本文主要评述MD方法在材料科学中的应用.

目前在材料微观结构的研究中，由于实验条件的限制，使得许多重要的微观结构的信息难以得到，如，对于由液态金属快速凝固的非晶转变过程，其微观结构的瞬时变化根本无法用实验仪器去测量。理论分析、实验测定及模拟计算已成为现代材料科学研究的3种主要方法[2]。20世纪90年代以来，由于计算机科学和技术的飞速发展，模拟计算的地位日渐突显。计算机模拟可以提供实验上尚无法获得或很难获得的信息。虽然计算机模拟不能完全取代实验，但可以用来指导实验的进行，从而促进理论和实践的发展，所以有必要对这一领域进行介绍。

2 分子动力学基本原理

分子动力学将连续介质看成由N个原子或分子组成的粒子系统，各粒子之间的作用力可以通过量子力学势能函数求导得出，忽略量子效应后，运用经典牛顿力学建立系统粒子运动数学模型，通过数值求解得到粒子在相空间的运动轨迹，然后由统计物理学原理得出该系统相应的宏观动态、静态特性。图1所示是MD

模拟过程。MD具体的做法是计算机上求运动方程的数值解。通过适当的格式对方程进行近似，使之适于在计算机求数值解。从使用连续变数和微分算符的描述过渡到使用离散变数和有限差分算符的描述，显然会有误差，误差的阶数取决于具体的近似机制，即所用的算法。模拟首先是规定初始条件。为了确定起见，可令初始位置在格子的格点上，而初始速度则从波尔兹曼分布得出。一个按上述办法建立的系统不会具有所要的能量，而且，很可能这个状态并不对应于一个平衡态。为了推动系统到达平衡，需要一个趋衡阶段。可以通过增加或从系统中移走能量，对运动方程向前积分若干时间步，使系统弛豫到平衡态。接着是物理量的计算阶段，沿着系统在相空间中的轨道计算一切令人感兴趣的量。

模拟中，MD采用周期边界条件和最小镜像原理，可以大幅度减少计算工作量［1］。周期边界条件是将一定数量的粒子N集中在一定的容积V中，这个容积V称为原胞，原胞周围的部分可以看作是原胞的复制，它们称作镜像细胞。这些镜像细胞的尺寸和形状与原胞完全相同，并且每个镜像细胞所包含的N个粒子是原胞中粒子的镜像，原胞在各个方向周期复制便形成了宏观物质样本。这样只需根据原胞周围的边界条件计算原胞内粒子的运动，幅度减少了工作量。

原子间作用势能模型的构造对于;MD法的应用至关重要。最简单的偶势模型只考虑两体作用，而与其它原子无关，在模拟中运算量小。20世纪80年代以来，各种经验或半经验的多体势模型迅速发展，特别是镶嵌原子法（EAM)［4］既克服了偶势的缺陷，又不会使计算量太大。

3分子模拟主要技术细节

3．1 分子间作用势

分子动力学模拟的首要条件就是要知道分子间的相互作用势，分子间作用势函数确定后，通过势函数对“求导即可得出分子间的作用力。所以必须对分：子间的作用势进行研究，目前存在多种势能模型。分子间势函数的发展经历了对势一多体势的过程。对势认为原子之间的相互作用是两两之间的作用，与其他原子的位置无关；而实际上，在多原子体系中，一个原子的位置不同，将影响空间一定范围内的电子云分布，从而影响其他原子之间的有效相互作用，故多原子体系的势函数需更准确地用多体势表示。

在分子动力学模拟的初期，人们经常采用的是对势。常用的对势模型包括Lennard—Jones势、Morse势、Born—Lande势及Johnson势。在对势模型中，系统能量表示为原子对相互作用能量的加和。对于简单的、高对称性的晶体，对势能够描述原子间的相互作用。对势模型的主要缺点是忽略了多体原子间的相互作用，更具体地说，它是忽略了由于原子局域环境的变化引起的原子间相互作用的变化，因此需要建立新的模型。

基于EAM势的势函数还有很多种。这些多体势大都用于金属的微观模拟。为了将EAM势推广到共价键材料，需要考虑电子云的非球形对称。于是，Baskes 等提出了修正型嵌入原子核法(MEAM)。经过修正，Baskes理论已能解决上述问题，但其确定参数的过程相当复杂，应用上仍受到很大的限制。Pasianot等试图在总能量中加入一修正项，以修正原EAM模型中2条基本假设所引起的能量差。这一模型成功地描述了Nb、Fe、Cr等元素，但应用起来很繁杂。张邦维等综合考虑以上EAM模型的优缺点后，提出了分析型EAM理论，成为一个普适分析型EAM模型。胡望宇等在张邦维分析型EAM模型的基础上进行了改进，提出了改进的分析型EAM模型(MAEAM)。此外，还有许多形式的多体势函数形式，如Jacobsen 等在等效介质原理(EMT)的基础上提出的另一种函数形式，由于其简单、有效，因此也得到了广泛的应用。势函数确定后，通过对势函数求导即可得出分子间的作用力。

3．2周期性边界条件

原子和分子体系的分子模拟的目的是提供一个宏观样本的物性信息。在具有自

由边界的三维N个粒子的体系中处于界面的分子数正比于N1脂。在使用有限的原子数来模拟实际体系中原子的运动时，必须考虑表面对体系中原子运动的影响。为避免这种影响，可以通过周期性边界条件来实现。将含有N个粒子的体积当作具有与其相同单元无限周期点阵的原始单元(如图1是一个二维的示意图，在这个二维图象中每个单胞被其他的8个单胞所包围；在三维方向上每个单胞就会被26个单胞所包围)，一给定粒子则与在此无限周期体系中的全部其他粒子相互作用。对有些模拟，在所有方向都用周期性边界条件是不合适的。比如在研究表面的分子吸附时，在与表面垂直方向上不能用周期性边界条件，而仅在平行于表面的2个方向需要应用周期性边界条件。

4 分子动力学在材料科学中的应用

4.1 分子动力学的适用范围

材料科学中计算模拟研究范围极为广泛，从埃量级的量子力学计算到连续介质层次的有限元或有限差分模型，可分为4个层次：电子、原子、显微组织和宏观层次（见图2）。;MD主要是原子尺度上研究体系中与时间和温度有关的性质的模拟方法。

最早将分子动力学方法用于材料研究中的是Vineyard 于1960年探讨材料辐射损伤的动力学规律。模拟结果给出了原子轨迹，这一工作使得过去对热力学性能的定性估计迈向对微观过程的定量研究.1964年Rahman用MD方法模拟液体氩，同时加进了周期性边界条件，结果他惊奇地发现可以用很少的粒子（864个）来反映真实系统的热力学性质。自此，凝聚态物质的分子动力学模拟成为可能，许多研究者纷纷投入这一研究工作。

最初应用是基于偶函数，如Lendard-Jones势函数和Morse势函数，模型简单，运算量小，而得以在材料科学中广泛应用。但由于其未考虑到体积相关项，在研究材料的弹性系数性质和预言金属的结合能及空位形成能时，难以获得准确的结果［5］。EAM多体势主要用于fcc型金属及其合金中，处理其结构、热力学、表面、缺陷及液态金属等问题，也应用于hcp及bcc型金属及合金，以及半导体Si。

一般，MD方法在中型机或微机上进行时，由于其内存和运算速度的限制，模拟研究只能限于500-1000个原子的小系统。因而模拟结果虽然也能揭示一些微观结构的特征和规律，但与实际的大系统情况有较大差异。在并行处理系统上对更大量的原子系统进行模拟研究［6］，其结果必然会接近于实际，从而对生产实践将会更有实际指导意义。

4.2 分子动力学的应用

4.2.1 金属的液态结构

在目前实验条件下，液态金属的结构及其变化尚很难精确测定。王鲁红等人

［7］采用F-S型多体势描述了8种hcp型金属的液态微观结构并与实验相比较，模拟结果表明，Mg、Co和Zn的势函数可以较好的描述其液态结构，Ti和Zr的势函数则不能；由Be和Ru的势函数描述的液态结构较为合理，Hf 则与一般液态结构特点不一致。李辉等人［8］采用EAM多体势模型，很好地描述了液态过渡金属Ni的结构变化特性。

4.2.2 薄膜形成过程

薄膜研究是当今科学研究的热点之一。目前在很多薄膜制备方法中，都应用了低能离子轰击技术，在这些方法中，低能离子/表面相互作用在控制薄膜的微观结构方面起着重要作用。由于离子/表面相互作用发生在时间间隔小于10-12s 内，因而特别适合于用MD方法对这一过程进行描述。薄腊研究是当今科学研究的热点之一目前在很多薄摸制备方法中,都应用了低能离子轰击技术,如离子束增强沉积,等离子休辅助化学汽相沉积、溅射沉积、离子辅助分子束处延等。Garrison[9],Kitabtake[9]和郊正明四等人分别用低能粒子轰击Si(ool)一2xl表面,由此可用分子动力学方法研究低能粒子对表面原子行为的影响。郊等人的研究表明,10ev,100ev粒子的轰击,一方面增强了表面原子形成二聚体(dimer)的能力,使表面二聚体键数增加,另一方面也使表面原子的排列更趋无序。

Ethier和Lewis[10]模拟了纯Si、Si

0.5Ge

0.5

和纯Ge在si(loo)-2xl再构表面

上用分子束处延(MBE)法生长膜的过程,其结果给出了薄膜质量与衬底温度之间的关系,即衬底温度较低时,形成的结构有序性较差;在高的衬底温度下,发生外延生长。对于纯Ge的外延生长,只有最初的三层结晶,以后便出现岛状结构,这在定性上与实验和理论结果相一致。

4.2.3界面研究

文献中大量报道了近十儿年来关于晶界的一些分子动力学模拟结果.目前有关界面的分子动力学模拟研究开展的不多.金属一金属界面的分子动力学模拟研究还有一些报道.李明研究了Ag/Ni界面处裂纹扩展行为.Ag/Ni和Cu/Ni界面在弯曲状态下的力学行为作了分子动力学模拟研究，结果给出的力学性能曲线与宏观规律相符合,并给出了进行界面模拟时计算元胞的原子数为1000个左右.同时证明了界面的存在对复合材料性能影响很大,界面结构不同复合材料的性能也不相同。Yang等人研究了Ni(100)面涂弋层的结构和动力学行为,并对嵌人原子法

(EAM)和对势结果作了比较.Rahrnan等人对上述问题作了进一步研究,给出了单层Ag在Ni(loo)面上的结构与温度之间的关系。即室温下,Ag在Ni衬底上前后滑动,距Ni上表面平均距离为 2.15埃,温度为1200K时,Ag在Ni上形成孔泡;1300K时,Ag在Ni上形成单晶.这些工作足以说明用分子动力学方法研究界面问题是可行的.采用分子动力学模拟方法进行界面研究是十分有意义的工作,尤其是目前正急于解决的金属一陶瓷界面问题。

4.2.4固体体相结构性能的研究

Parrinell和Rahman是最先采用分子动力学方法研究固体性能的,以此代替Milstein和Faber的静力学计算.其对象是Ni单晶在单轴压力下由面心立方结构向密排六方结构转变的过程,这是等温一等压分子动力学问题,所用势函数是Morse对势,其结果给出了应力一应变曲线,并与Milstein和Faber的汁算结果符合得相当好。Zhong[8]等人利用在嵌人原子法(Embedded一Atom一Method一EAM)的基础上发展起来的MBA势(Many一Body-Alloy)研究Pd一H系统的力学稳定性.所选取的势函数形式是Morse势.其中利用了Nose[15]和Rahman一Parrincllo分子动力学形式。结果表明氢脆的微观起因是氢饱和使某些区域的韧性和塑性增强,这与假定的氢增强局域塑性的机制相一致。

5 结论

从国内外分子动力学的模拟研究可以看出，随着计算机的普及及计算能力的提高，分子模拟方法的推广应用日益受到高度重视。模拟已逐渐从分子结构简单的体系扩展到分子结构复杂的体系，模拟粒子数已从一二百扩展到几千，并从单相扩展到多相、从均相扩展到非均相，位能函数也从简单扩展到复杂乃至要考虑量子效应。模拟方法也已经从Mc转至MD，以模拟真正的动力学过程。MD模拟方法也从平衡态模拟发展到外加场的非平衡态模拟。可以预料，分子模拟与材料的关系将更加紧密，它对材料科学工程的基础研究、工艺过程以及新产品的开发都会发挥更为明显的、不可替代的作用。

参考文献

[1］罗熙淳，梁迎春，董申. 分子动力学在纳米机械加工技术中的应用［J］中国机械工程，1999,10(60:692

[2］胡壮麒，王鲁红，刘轶电子和原子层次材料行为的计

算机模拟［J］材料研究学报，1099,12(1):1

[3] Nose S.J.Chem.Phys.1984.85:511

[4] Rice B D,Chcn J.Philos,Mag,Lctt.1990,61:139

[5］卢百平，张丰收，傅恒志材料科学中的分子动力学模拟

[6］编者西北工业大学2002年材料论坛论文集［C］+

西安：2002.1-4

[7］郭俊梅，邓德国，潘建生，等计算材料与材料设计［J］+

贵金属，1999,20(4):62

[8］董科军，刘让苏，郑采星，等液态金属凝固过程分子动力

学模拟的并行算法研究［J］+ 计算机科学与工程，2002,22(3):45

[9] Zhang W.Cai Y,Tomanck D. Phys.Rev.B,1992,46:8099

[10] Puska M J,Nicminen R M,Manninen M.Phys.Rev.B.1981,24：3037