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

分子动力学在材料科学中的应用摘要：分子动力学方法是进行物质原子或分子层次计算机模拟时所采用的一种基本方法。通过分子动力学模拟，可以给出原子尺度上材料及其演化过程细节的可能性，具有史无先例的准确性，使材料设计和性能的预测成为可能。本文介绍了分子动力学模拟的基本原理和算法。综述了用分子动力学模拟方法研究固休的相结构、表面问题、界面问题以及薄膜形成过程的方法。关键词：材料;分子动力学;计算机模拟；Abstract：MolecularDynamicsmethodisoneofthebasicmethodsusedtosimulatetheatomicormolecularlevelcomputersimulation.Bymoleculardynamicssimulation，wecangivethepossibilityofthedetailsoftheevolutionprocessofthematerialontheatomicscale，andithasthepossibilityofnoprecedent，sothatthematerialdesignandperformancepredictioncanbemade.Thispaperintroducesthebasicprincipleandalgorithmofmoleculardynamicssimulation，andsummarizesthemethodsofstudyingthephasestructure，surfaceproblem，interfaceproblemandtheformationprocessofthethinfilmbymoleculardynamicssimulation.Keywords：materials;moleculardynamics;computersimulation;分子动力学(MolecularDynamics，简称MD)用于计算以固体、液体、气体为模型的单个分子运动，它是探索各种现象本质和某些新规律的一种强有力的计算机模拟方法，具有沟通宏观特性与微观结构的作用，对于许多在理论分析和实验观察上难以理解的现象可以做出一定的解释。MD方法不要求模型过分简化,可以基于分子(原子、离子)的排列和运动的模拟结果直接计算求和以实现宏观现象中的数值估算。可以直接模拟许多宏观现象，取得和实验相符合或可以比较的结果，还可以提供微观结构、运动以及它们和体系宏观性质之间关系的极其明确的图象。MD以其不带近似、跟踪粒子轨迹、模拟结果准确而倍受研究者的关注，在物理、化学、材料、摩擦学等学科及纳米机械加工中得到广泛而成功的应用。本文主要讲述MD方法在材料科学中的应用及MD的最新研究进展。1、分子动力学基本原理分子动力学将连续介质看成由N个原子或分子组成的粒子系统，各粒子之间的作用力可以通过量子力学势能函数求导得出，忽略量子效应后，运用经典牛顿力学建立系统粒子运动数学模型，通过数值求解得到粒子在相空间的运动轨迹，然后由统计物理学原理得出该系统相应的宏观动态、静态特性，图1所示是MD模拟过程。MD具体的做法是计算机上求运动方程的数值解。通过适当的格式对方程进行近似，使之适于在计算机求数值解。从使用连续变数和微分算符的描述过渡到使用离散变数和有限差分算符的描述，显然会有误差，误差的阶数取决于具体的近似机制，即所用的算法。模拟首先是规定初始条件。为了确定起见可令初始位置在格子的格点上，而初始速度则从波尔兹曼分布得出。一个按上述办法建立的系统不会具有所要的能量，而且很可能这个状态并不对应于一个平衡态。为了推动系统到达平衡，需要一个趋衡阶段。可以通过增加或从系统中移走能量，对运动方程向前积分若干时间步，使系统弛豫到平衡态。接着是物理量的计算阶段，沿着系统在相空间中的轨道计算一切令人感兴趣的量。模拟中,MD采用周期边界条件和最小镜像原理,可以大幅度减少计算工作量。周期边界条件是将一定数量的粒子N集中在一定的容积V中，这个容积V称为原胞，原胞周围的部分可以看作是原胞的复制，它们称作镜像细胞。这些镜像细胞的尺寸和形状与原胞完全相同，并且每个镜像细胞所包含的N个粒子是原胞中粒子的镜像，原胞在各个方向周期复制便形成了宏观物质样本。这样只需根据原胞周围的边界条件计算原胞内粒子的运动幅度,减少了工作量。原子间作用势能模型的构造对于MD法的应用至关重要。最简单的偶势模型只考虑两体作用，而与其它原子无关，在模拟中运算量小。20世纪80年代以来，各种经验或半经验的多体势模型迅速发展，特别是镶嵌原子法(EAM)既克服了偶势的缺陷,又不会使计算量太大。2、分子动力学的几种基本方法2.1等能过程的分子动力学方法人们最初用分子动力学方法研究的对象是微正则系综，即休系的粒子数(N)休积(V)和能量(E)均守恒，且总线动量恒为零(P=0)。其出发点是用哈密顿方程来描述N个粒子的相互作用。在经典范畴里，哈密顿方程可以给出各种不同形式的运动方程，虽然这些方程在数学上是等价的，但在数值上却不等价,这里选择牛顿力学方程形式：式中的Firimi分别是第i个原子所受的力，第i个原子的位移和第i个原子的质量。对式(l)所示的二级微分方程进行一次积分，可求得粒子运动的速度;二次积分可求得粒子的位置(即质点的运动轨迹)。其初始条件要求给出初始位置及初始速度，初始位置决定势能对总能量的贡献，而初始速度则确定动能对总能量的贡献，求注解二级微分方程可利用Verlet算法(Summedform)，计算的成功与否依赖于初始位置和初始速度的设置。通常将所研究的休系置于理想晶休格上,根据玻尔兹曼分布来设置初始速度或有时直接将初始速度全部设为零。但是，一般情况下并不能准确地知道与给定能量相关的初始条件，因此按前述给定的初始位置和初始速度对体系进行调节，直至达到所需要的能量值。2.2等温过程的分子动力学方法在许多实际情形中，要研究的对象是等温系统而不是等能系统。这种情况相当于将弧立的等能系统放在一恒温的大热源中，使该系统与大热源之间可以交换能量从而保证温度恒定。这一恒温系统即为正则系统，其粒子数(N)、体积(V)、温度(T)均恒定，同时总线动量恒为零(P=0),温度恒定可视动能恒定。在模拟过程中保持温度恒定的办法是对每一步的速度进行重新标度，即在速度前均乘上一因子β，使体系的动能保持为定值。该因子可由下式确定：式中,Kβ为玻尔兹曼常数;Tref为参考温度;mivi分别为元素i的原子质量和原子速度。Haile和Gupta又将上述方法做了改进，将动能恒定这一限制加到运动方程中去，作为运动方程的一特殊情形，他们提出了重新标度动量的方法。2.3等压过程的分子动力学方法如果对孤立体系(即N、E、F均守恒的休系)施加压强，则体积V就是与压强有关的变量,这时体系不再是孤立的。如果体系与外界是绝热的，那么粒子数N和压强P都是守恒量，但总的内能是不守恒的，此时的守恒量是焓H，因此这时的体系实际上是等温-等压休系(即NPH)。为使体积涨落变为可能，引进新的变量PeV（Pe为外加压力,V为体系的休积)，这是一个与时间相关的新的变量，将它看作势能，在拉格朗日方程中增加这一项，因此动力学方程变为：等压分子动力学模拟方法首先由Anderson给出的。随后由Pairnelo、Ranamn将等压分子动力学推广至允许计算元胞的形状改变，并将其应用到固体结构变化中，此后,Heycs用与Andesron相似的方法做等压分子动力学模拟(用真实变量代替无量纲变量)。2.4等温一等压过程的分子动力学方法在等压过程中，如果系统与外部不是绝热的，则此即为等温一等压系统，即NPT系统。该系统的分子动力学模拟算法与NPH系统的算法基本相同。在1、2节中我们已经明确了等温过程是通过调节每一步的速度来达到动能恒定从而温度恒定，将这一思想用于等压过程，即调节等压过程的每一步的速度来实现等温一等压的，对该过程的另一种处理方式是由Andesron提出的，用MD和MC桥接的办法来模拟等温一等压系统，该方法通过随机碰撞来改变粒子速度，碰撞后粒子速度的分布再用于产生正则系统。2.5其它分子动力学方法此外,Car和Parriello将密度函数理论与分子动力学方法统一起来，形成CP分子动力学。该方法除具有前述的分子动力学的优势外，还可以计算固体材料的基态电子结构。由于计算量过大，目前其应用还受到一定条件的限制。但随着计算机的发展，CP分子动力学将会越来越为更多研究者所采用。还有，Hoover和Evans等人从非平衡态出发得出约束分子动力学方法，又称为HLNE法,在这一方法中作用力方程式中增加了一项：αPi，其中,qi为坐标;Pi为动量;α为一系数，其值为：方程虽不再是正则形式，但在坐标空间可以产生正则分布。由此看来，利用分子动力学方法所能研究的问题非常广泛。目前，已有许多人利用这一方法做了大量的研究工作。这里大致归纳为四类问题作如下简要综述。3、分子动力学的应用3.1金属的液态结构在目前实验条件下，液态金属的结构及其变化尚很难精确测定。王鲁红等人采用F-S型多体势描述了8种hcp型金属的液态微观结构并与实验相比较,模拟结果表明Mg、Co和Zn的势函数可以较好的描述其液态结构，Ti和Zr的势函数则不能;由Be和Ru的势函数描述的液态结构较为合理,Hf则与一般液态结构特点不一致。李辉等人采用EAM多体势模型，很好地描述了液态过渡金属Ni的结构变化特性。3.2高压相稳定性预测高压试验比较困难，但在近自由电子近似下，用赝势理论可以预测高压相稳定性。金朝晖研究了高压下Mg的相变，在常压Mg以1×1012和6×1012K/s的速度冷却时获得hcp结构，冷却速度为1.2×1013K/s时，得到非晶结构。在高压(45GPa)下，Mg以8×1012和5×1013K/s的速度冷却后为bcc结构,冷却速度为1×1014K/s时可得非晶结构。3.3热力学参数的计算刘洪波利用EAM模型几等体积MD方法模拟计算Cu-Ni合金的Gibbs自由能，结果表明，计算的合金混合Gibbs自由能与实验值较一致。3.4薄膜形成过程薄膜研究是当今科学研究的热点之一。目前在很多薄膜制备方法中，都应用了低能离子轰击技术，在这些方法中，低能离子/表面相互作用在控制薄膜的微观结构方面起着重要作用。由于离子/表面相互作用发生在时间间隔小于10-12s内，因而特别适合于用MD方法对这一过程进行描述。3.5金属界面研究目前有关界面的MD研究的不多。罗旋曾对Ag/Ni和Cu/Ni界面在弯曲状态下的力学行为作了研究，结果给出的力学性能曲线与宏观规律相符合，这说明用MD方法研究界面问题是可行的。采用MD方法进行界面研究是十分有意义的工作，尤其是目前正急需解决的金属-陶瓷界面问题。3.6凝固及玻璃化转变的模拟由于实验室急冷技术所能达到的冷却速度一般小于107K/s，对于许多合金和金属不能达到其玻璃转化所需的冷却速度，对玻璃态结构的研究受到了限制。用EAM框架下的有效配对势模拟了Pd的一级相变及玻璃转化过程，得到的热力学参数可以同实验室进行比较，还得到了对应于结晶和玻璃转化的冷却速度。又如，常规的快冷技术冷却速度一般107K/s，无法使Ni3Al非晶化。王鲁红用计算方法知道在什么情况下能非晶化，弥补了实验的不足。另外MD法还在固液相变、固态相变、表面能计算、晶体生长等几个方面得到应用。4结论综上所述，分子动力学模拟方法在材料科学领域已越来越活跃，它已成为材料科学工作者的一个强有力的工具。有些化学工作者也尝试用MD方法来解决物质的结构与性质等问题。最近，第一性分子动力学方法正在兴起，它对经典分子动力学方法产生很大的挑战。以往的经典分子动力学方法只能给出原子(离子)的位置，以及晶体的结构，但第一性分子动力学却能同时给出晶体结构及基态性质。随着计算机的发展，计算量不再成为人们最感困难的问题时，第一性分子动力学模拟必将越来越为材料科学工作者所重视。参考文献[1]罗熙淳,梁迎春,董申；分子动力学在纳米机械加工技术中的应用[J].中国机械工程,1999,1