材料加工计算机模拟与应用

《材料加工计算机模拟与

应用》

学号：

流水号：

姓名：

完成日期：

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

摘要从分子水平研究材料的开发和设计是21世纪材料科学与工程的一个重要方向。综述了分子模拟中的Monte Carlo分子模拟和分子动力学模拟2种方法及其在材料科学中的应用，概述了目前国内外分子模拟的研究进展，评述了当前使用的数学模型的特点及局限性，展望了分子模拟技术的发展方向。

关键词分子模拟分子动力学法(MD) 蒙特卡罗法(Mc) 材料科学

引言

从分子水平研究材料的开发和设计，无疑是21世纪材料科学与工程的一个重要方向[1]。许多材料新技术的出现，如膜分离技术、界面技术、纳米技术等，都涉及到复杂的分子结构及其伴生的复杂现象和性能，在这些方面，除了准确的物性数据外，更要对各种复杂现象的机理有深刻的理解。因此材料科学工作者造韧需要从分子水平来研究系统的微观结构及其宏观热力学性质、传递性质。但是，目前在材料微观结构的研究中，由于实验条件的限制，使得许多重要的微观结构的信息难以得到，如，对于由液态金属快速凝固的非晶转变过程，其微观结构的瞬时变化根本无法用实验仪器去测量。近年来，利用计算机模拟技术研究材料日益成为人们感兴趣的课题。理论分析、实验测定及模拟计算已成为现代材料科学研究的3种主要方法[2]。20世纪90年代以来，由于计算机科学和技术的飞速发展，模拟计算的地位日渐突显。计算机模拟可以提供实验上尚无法获得或很难获得的信息。虽然计算机模拟不能完全取代实验，但可以用来指导实验的进行，从而促进理论和实践的发展，所以有必要对这一领域进行介绍。本文拟从分子模拟的种类及其基本原理、国内外最近研究进展状况作一简要的介绍。

1、分子模拟基本方法

分子模拟基本上有2种方法。一是Monte Carlo法(简称MC)，二是分子动力学法(简称MD)。Monte Carlo法和分子动力学法是利用统计力学的计算方法，用计算机模拟宏观体系的微观结构和运动，在此基础上用数值运算统计求和的方法。这2类方法不要求模型过于简化，可以基于分子(原子、离子)的排列和运动的模拟结果直接计算求和以实现宏观现象中的数值估算。它们一方面可以直接模拟许多宏观现象，取得和实验相符合或者可以比较的结果；另一方面可以提供微观结构、运动以及它们和体系宏观性质间的关系及其明确的图像，从中产生一些新概念或者新理论。

1．1 MC法基本原理

MC法最初是在原子能研究中为了计算中子散射过程而发展起来的，因而也可以用来计算某些动力学问题。基本模拟过程是在一定系统条件下，将系统内粒子进行随机的位移、转动，或粒子在两相间转移位置。根据给定的分子势能函数，进行粒子间内能的加和。采用Metropolis取样方法[3]，生成一系列体系的微观粒子随机构型，从而逐渐趋近于平衡时的Boltzmann分布。由于MC模拟的粒子位移是虚拟的，不代表粒子的真正运动历程，所以该法不能用于传递性质的模拟。通常采用MC方法所得到的粒子瞬时分布很接近实际情况，而模拟得到的粒子的运动方式却与实际情况有差异，这是因为Mc方法是采用随机数来描述粒子的运动，并使其符合Bolzmann分布的。因此，用MC方法研究物系平衡性质是可靠的，而用它研究动力学性质就必须谨慎。

1．2 MD法基本原理

分子动力学方法(MD)是另一种主要的计算机模拟方法，目前在材料科学、物理、化学等学科的各个领域得到广泛应用。基本模拟过程是在一定系及已知分子势能函数条件下，从计算分子间作用力着手，求解牛顿运动方程，得到体系中各分子微观状态随时间的真正变化，再将粒子的位置和动量组成的微观状态对时间平均，即可求出体系的压力、能量、粘度等宏观性质以及组成粒子的空间分布等微观结构。系统的初始位形和初始速度可以通过实验数据、或理论模型、或两者的结合来决定。如果被模拟的系统具有初始密度分布和温度分布，而没有固定的晶格结构，则每个原子的位置可以从初始密度分布，用舍选法或Metropolis等

方法得到。每个原子的初速度，则可以从初始温度分布下的Maxwell- Boltzmann 分布来随机选取。Maxwell- Boltz- mann分布可以用0～1之间均匀分布的随机数发生器的输出通过简单的变换而得到[5]。该方法既可计算体系的平衡性质，也可计算体系的各种动力学性质。在MD模拟中，又可根据是否对体系加一外力场而区分平衡MD模拟(EMD)和非平衡MD．模拟(NEMD)。外加力场的NEMD模拟方法目前尚处于初级阶段。对于非平衡系统，其分子动力学模拟的过程包括初始条件和边界条件的确定、牛顿方程的有限差分求解和作为时间函数的感兴趣量的提取。对于平衡系统，其分子动力学模拟的过程与非平衡系统的差别在于感兴趣量及边界条件与时间无关。从原理上说，利用MD法研究动力学过程更为合适，所以MC方法可以对系统平衡态的某些问题进行模拟研究，而对于快速凝固系统处于非平衡态的变化过程则必须采用MD方法。

2、分子模拟主要技术细节

2．1 分子间作用势

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

在分子动力学模拟的初期，人们经常采用的是对势。常用的对势模型包括Lennard-Jones势、Morse势、Born-Lande势及Johnson势。Lennard-Jones(L-J)势是最老的原子间作用势之一,是由J．E．Lennard-Jones[6]于1924年提出来L-J 势基本形式。在对势模型中，系统能量表示为原子对相互作用能量的加和。对于简单的、高对称性的晶体，对势能够描述原子间的相互作用。然而用对势模型计算立方晶体的弹性常数(Cauchy关系)，而在实际金属中很少有能够满足Cauchy 关系的。对势模型的主要缺点是忽略了多体原子间的相互作用，更具体地说，它是忽略了由于原子局域环境的变化引起的原子间相互作用的变化，因此需要建立新的模型。

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

模型的基础上进行了改进，提出了改进的分析型EAM模型(MAEAM)。此外，还有许多形式的多体势函数形式，如Jacobsen等在等效介质原理(EMT)的基础上