第一性原理计算方法讲义

第一性原理计算方法

引言

前面讲述的有限元和有限差分等数值计算方法中，求解的过程中需要知道一些物理参量，如温度场方程中的热传导系数和浓度场方程中的扩散系数等，这些参量随着材料的不同而改变，需要通过实验或经验来确定，所以这些方法也叫做经验或者半经验方法。而第一性原理计算方法只需要知道几个基本的物理参量如电子质量、电子的电量、原子的质量、原子的核电荷数、布朗克常数、波尔半径等，而不需要知道那些经验或半经验的参数。第一性原理计算方法的理论基础是量子力学，即对体系薛定额方程的求解。

量子力学是反映微观粒子运动规律的理论。量子力学的出现，使得人们对于物质微观结构的认识日益深入。原则上，量子力学完全可以解释原子之间是如何相互作用从而构成固体的。量子力学在物理、化学、材料、生物以及许多现代技术中得到了广泛的应用。以量子力学为基础而发展起来的固体物理学，使人们搞清了“为什么物质有半导体、导体、绝缘体的区别”等一系列基本问题，引发了通讯技术和计算机技术的重大变革。目前，结合高速发展的计算机技术建立起来的计算材料科学已经在材料设计、物性研究方面发挥着越来越重要的作用。

但是固体是具有～1023数量级粒子的多粒子系统，具体应用量子理论时会导致物理方程过于复杂以至于无法求解，所以将量子理论应用于固体系统必须采用一些近似和简化。绝热近似（Born-Oppenheimei近似）将电子的运动和原子核的运动分开，从而将多粒子系统简化为多电子系统。Hartree-Fock近似将多电子问题简化为仅与以单电子波函数（分子轨道）为基本变量的单粒子问题。但是

其中波函数的行列式表示使得求解需要非常大的计算量；对于研究分子体系，他可以作为一个很好的出发点，但是不适于研究固态体系。1964年，Hohenberg 和Kohn 提出了严格的密度泛函理论（Density Functional Theory, DFT ）。它建立在非均匀电子气理论基础之上，以粒子数密度()r ρ作为基本变量。1965年，Kohn 和Sham 提出Kohn-Sham 方程将复杂的多电子问题及其对应的薛定谔方程转化为相对简单的单电子问题及单电子Kohn-Sham 方程。将精确的密度泛函理论应用到实际，需要对电子间的交换关联作用进行近似。局域密度近似（LDA ）、广义梯度近似（GGA ）等的提出，以及以密度泛函理论为基础的计算方法（赝势方法、全电子线形缀加平面波方法（FLAPW ）等）的提出，使得密度泛函理论在化学和固体物理中的电子结构计算取得了广泛的应用，从而使得固体材料的研究取得长足的进步。

第一性原理计算方法的应用

1、体系的能量。

进行第一性原理计算前，首先需要确定体系模型，即模型的晶胞和晶胞中原子的坐标。对于晶体具有周期对称性，具有三个基矢方向和基矢大小（晶格常数）。由于理论计算确定的平衡晶格常数和实验值有所差别，建立模型前需要确定平衡晶格常数。晶格常数的确定采用如下步骤：

通过改变三个基矢的大小，改变单胞的体积（81-119%）。通过第一性原理计算可以得到具有不同体积的模型的能量。通过拟合Murnaghan 方程，得到晶体的晶格常数以及单胞的能量：

'

0'000000''00()()11(1)B B V V V E V E V B B B V V ⎡⎤⎛⎫⎛⎫⎢⎥=+-+- ⎪ ⎪-⎝⎭⎝⎭⎢⎥⎣⎦

其中，0V 为基态平衡体积，00()E V 为基态下体系的结合能（相对于对应孤立原子能量）。V 为原胞体积，0B 为体模量，'0B 为体模量对压强的导数。如课件中图形所示，可以确定在一定体积下体系的能量达到极小值，即体系的基态能量，所对应的体积为体系的平衡体积，进而可以得到模型三个基矢的大小确定晶体的平衡晶格常数。这里需要指出的是不同的第一性原理计算方法给出的能量，代表的物理意义不同，但是本质上都可以反应体系的稳定性。如总能指构成体系的原子孤立时的能量减去原子成键放出的能量；结合能是以孤立原子的能量为零点，体系具有的总能，即原子构成晶体时放出的能量。

在上面求得的晶格常数的基础上，根据要研究的物理问题，确定体系中包含原子数目的多少，建立第一性原理计算模型。第一性原理计算的模型通常选取一个或几个单胞（超单胞）作为模型，选取的模型具有三个基矢方向，应保证沿着三个基矢方向平移可以构成无限大的晶体。

第一性原理计算输入的原子坐标有两种坐标形式，一种是笛卡尔坐标（Cartesian coordinates ），一种是分数坐标（fractional coordinates ）。如对于Ni3Al 高温合金，具有如图所示的晶体结构，铝原子位于立方体的顶点，镍原子位于立方体的面心位置。如果取一个单胞作为研究模型，则三个基矢1a ，2a ，3a 分别为（a , 0 , 0）（0, a , 0）（0, 0, a ），其中a 为体系的晶格常数。单胞中包含四个不等价原子：三个Ni 和一个Al 。采用笛卡尔坐标四个原子的坐标可表示为（0，0，0），（a /2, a /2, 0），（a /2, 0 , a /2），（0, a /2, a /2）。如果采用分数坐标表示为（0，0，0），（1/2, 1/2, 0），（1/2, 0 , 1/2），（0, 1/2, 1/2）。迪卡尔坐标(x,y,z)和分数坐标(a,b,c)之间关系为1a ·a+1a ·b+1a ·c=(x, y, z)，其中a,b,c 为一个原子的三个分数坐标，x,y,z 为该原子的笛卡尔坐标。图

中所示各点表示将晶格常数的大小a取不同值时得到的单胞体积作为横轴，而纵轴表示对应体积下将原子坐标输入进行第一性原理计算求得的体系能量。拟合后得到Ni3Al的平衡基矢大小以及体系能量。

对于研究合金中的掺杂问题，由于掺杂元素的量很少，所以建立的模型需要取多个单胞（超晶胞）。随着模型中原子数目的增加，第一性原理计算方法的计算量指数增加，对于掺杂量很低的情况，如0.1%，需要模型中至少取1000个原子来和实际相符合，这超出了第一性原理计算在目前的计算机上的计算能力（100原子左右），所以建立模型时需要考虑能够反映要研究的实际问题就可以。假设一个超单胞中只存在一个掺杂原子，这样相邻两个超晶胞中掺杂原子的间距为超晶胞的基矢大小。一般两个原子之间相隔三到四个原子层，原子之间的相互作用就可以认为非常小了。所以选取八个单胞构成的超单胞就可以基本反应掺杂量很低的掺杂问题了。

将建立好的模型，进行第一性原理计算可以得到体系的总量，对总能进行变换可以定义体系的内聚能、形成能以及择优占位能，进而可以对掺杂是否有利于形成，形成掺杂后对体系稳定的影响而进行分析，如课件中所列。内聚能：体系的总能减去所有原子孤立时的能量，即由于原子之间的相互作用而放出能量，从而内聚能为负值，越小表示形成的体系越稳定。通过和没有掺杂体系的内聚能相比较可以看出掺杂元素对体系稳定性的影响。形成能：体系的总能减去体系中各自元素对应的晶体中原子的能量。形成能可表示各种金属组成合金的能力。另外通过比较掺杂原子替代合金中不同元素原子时体系的形成能可以得到掺杂原子倾向位于合金的什么位置，这个差值就可以定义为择优占据能。这里需要强调，各种能量是根据要研究的具体问题来定义的，比如我们要就掺杂原子倾向于位于