下载掌阅APP，畅读海量书库
立即打开
众所周知，量化计算在原子分子尺度上可归结为薛定谔方程（SchÖrdinger Equation），但薛定谔方程的求解却极其复杂和困难，密度泛函理论 [1，2] （Density Functional Theory-DFT）的建立和计算机技术的飞速发展使薛定谔方程的求解在很大程度上得到了解决。密度泛函理论于 20 世纪 60 年代，在 Thomas-Fermi [3，4] 理论的基础上发展而来，其最核心的部分是用粒子密度描述体系基态的物理性质来代替传统量子理论作为体系基本物理量的波函数，而粒子的密度只跟空间坐标有关，这将使复杂的 3 N 维波函数问题简化为 3 维的粒子密度问题，其最基本的近似是Born-Oppenheimer近似 [5，6] 。而粒子密度可以通过实验直接观测，使其应用起来简单且直观，具有相当诱人的应用前景。Born-Oppenheimer 近似又称绝热近似 [5，6] 。原子核的质量是电子质量的几千倍，因而使电子运动的速度远高于原子核的速度。电子处于高速运动中，而原子核只是在其平衡位置附近振动。当原子核间进行任一微小运动时，电子立即调整建立起与新位置的原子核力场相对应的运动状态，始终跟上原子核的运动。在讨论电子运动和结构时，离子的动能被忽略，可以不考虑核的运动，在计算总能时才需考虑；而考虑核的运动时，不需考虑电子在空间的具体分布，只考虑电子运动的平均结果即可。因而电子和核的运动可以分开考虑，多种粒子系问题可简化为多电子问题。绝热近似是固体量子理论中的基本近似，是能带理论及密度泛函理论中的基础性近似。一般将基于密度泛函理论的计算称为第一性原理（First-princinples）计算或从头算（Ab-initio），它是将多原子体系当做由电子和原子核构成的多粒子体系，在仅采用 m 0 、 e 、 ħ 、 c 、 k B 这 5 个常数而不引入任何经验参数情况下，对多粒子体系进行处理的方法，是量化计算的基础和核心。
多体理论是量子力学的核心问题，要确定含 N 个粒子的量子体系的性质 [7] ，原则上可以通过求解该体系的薛定谔方程
当体系的势场 U 与时间无关时，上面的薛定谔方程的解可以用分离变量法进行简化，即可以得到定态薛定谔方程
对于多粒子体系，上述方程仍然不能求解，我们必须借助一些物理模型和一些基本原理，如 Thomas-Fermi [3，4] 模型、Hohenberg-Kohn [1，2] 定理、Kohn-Sham方程及常用的分子轨道理论近似方法 [8] 等。
Thomas [3] 和 Fermi [4] 在 1927年提出了 Thomas-Fermi模型。其基本原理是：体系的动能可以通过电子密度来表达。他们考察理想的均匀电子气模型，把空间分割成足够小的立方体，在这些立方体中求解无限深势阱中粒子的方程，假设电子之间无相互作用，得到相应的能量和密度的表达式。把它们联系起来，简化后得到动能与粒子密度的关系式为
对于原子的情况，加上核吸引势和电子间的库仑势的作用，得到总能量与电子密度 ρ 的关系为
式中， Z 是核电荷数。该模型的表达式比较简单，但在实际计算中结果却不是很好，在原子的计算中并没有显示出比其他方法更突出的优越性，且在分子的计算中得不到原子间可能的成键。后来，很多学者在该模型中加入了各种修正项进行修改，如Dirac [9] 加入了交换作用项，该项对均匀电子气中的交换能进行处理。而 Slater [10，11] 近似认为，一个具有变化电子密度的体系的交换势可以近似地用一个局域密度依赖的项 n （ r ）/3来表示，这样的近似使得从头算可以计算真实的固体。
密度泛函理论的关键在于用电子密度分布，而不是把电子波函数分布作为试探函数，将总能量 E 表示为电子密度的泛函。这样的处理首先要在理论上确认，存在总能量对于电子密度分布的这样一个泛函。Hohenberg和 Kohn [1，2] 对非均匀电子气理论提出了两个定理。
定理 1 ：不计自旋的全同费密子系统的基态能量是粒子数密度函数的唯一泛函。
定理 2 ：能量泛函 E [ ρ ]在粒子数不变的条件下，对正确的粒子数密度函数 ρ （ r ）取极小值，并等于基态能量。
这里所处理的基态并非是简并态的，多电子体系哈密顿量分开写为动能部分、多电子系统相互作用部分和多电子系统之外的外场部分。
则 Hohenberg-Kohn定理证明体系总能存在对基态电子密度分布函数的泛函形式。
Kohn-Sham方程的具体形式为
式中，| Φ 是基态波函数， F [ ρ ]是与外场无关的部分，不管外场取什么形式， F [ ρ ]总是有共性的部分，因此 F [ ρ ]泛函的具体形式是整个密度泛函理论最关键的部分 [1，2] 。Hohenberg-Kohn 定理虽证明了总能量的确能通过求解最有利的基态电子密度分布函数而得到，但是总能量对于电子密度分布函数的具体泛函形式，以及如何才能利用以上泛函极值的性质求解总能量的问题，Hohenberg-Kohn定理并没有给出回答。而Kohn和Sham提出的Kohn-Sham方法最终将密度泛函理论和实际应用结合起来。Kohn-Sham 方程利用无相互作用粒子的哈密顿量，代替有相互作用粒子的哈密顿量，将粒子之间的相互作用全部归入交换关联泛函。
密度泛函理论整个框架中只有一个未知部分
即交换关联势的形式未知，实际应用中，密度泛函计算结果的精度取决于选取怎样的交换关联势。局域密度近似（LDA，Local Density Approximation） [12] 和广义梯度近似（GGA，Generalized Gradient Approximation） [13] 是两种常用的交换关联势的近似方法。LDA 是最实用、最有效的一种近似，它最早由 Slater 在 1951 年提出并应用 [10～11，14～15] 。该近似假定空间某点的交换关联能只与该点的电荷密度有关，且等于同密度的均匀电子气的交换关联能，可表示为
目前计算中最常用的交换关联势的局域密度近似是根据 Ceperley [16] 和 Alder [17] 用蒙特-卡罗方法计算均匀电子气的结果。
通过实际计算表明，使用LDA计算原子游离能、分子解离能的误差在10%～20%之内，对分子键长、晶体结构可准确到误差在 1%左右 [18] ，对大多数材料的计算具有相当的优势。但对于均匀电子气或者空间缓慢变化的、电子气相差太远的系统却不太适用。需更精确地考虑计入某处附近的电荷密度对交换关联能的影响，如要考虑密度的一级梯度对交换关联能的贡献。
可以取交换能为修正的 Becke泛函形式
[19，20]
（
，
β
是常数）
这称为广义梯度近似（GGA，Generalized Gradient Approximation） [13] 。计算中常用的 GGA方法的交换关联势有 Becke [19，20] 、Perdew-Wang 91 [21～23] 和 BLYP [24] 等。更进一步地，还可考虑密度的高阶梯度近似，称为 Meta-GGA或 Post-GGA；甚至可考虑非局域的交换关联作用，如 Vander Waals作用 [25，26] 。对这些问题，人们虽都有研究，但仍未找到一个足够精确的、适合的、简单的形式。
在 0K 时，完整晶体中粒子呈周期性规则排列，因而电子受到的外部势能也是周期性的，该周期与单位晶胞的长度 l 相同，即一个电子的表面势可表示为： V （ r ）= V （ r + l ），这是应用 Bloch定理的前提。Bloch定理利用晶体的周期性，将一个电子的无穷数量波函数的计算转化为简单晶体单位晶胞中电子数量的计算。通过 Bloch 定理，可能将无限晶体的波函数用 Bravais 晶格倒空间矢量的波函数来描述。波函数可表示为晶胞周期性和类波两部分的乘积
式中，第一项是类波部分，第二项是波函数的晶胞的周期性部分。可用有限数量的平面波来表示，平面波波向量是晶体的倒晶格向量
上式中， G 是倒晶格向量， G 对所有 l 定义为 G · l =2π m ， l 是晶体的晶格向量， m 取整数。因此每个电子波函数可写为平面波之和。
通过 Bloch 理论的运用，将电子无穷数量的问题简化为在周期晶胞的第一个Brillouin区内，用无穷倒空间向量 k 表示波函数的问题。
每个
k
点处的电子波函数可用离散的平面波基集来表示。理论上，这个傅里叶级数是无穷的，平面波的每个系数
c
i
（
k
+
G
）
都有一个动能
。通常，具有较小动能的平面波比那些具有很高动能的平面波更重要，而引入平面波能量截断值可将基集加以限制。这个能量截断将导致在计算体系总能量时出现误差，但通过加大能量截断值和增加集基的大小可以降低误差。能量截断值的大小取决于研究的体系，需要在研究之初通过测试来确定适当的能量截断值。
将电子波函数用平面波集基展开的另一优点是 Kohn-Sham 方程的形式非常简单。将式（2.15）代入 Kohn-Sham方程，可得
可以看到，动能的倒空间表示是对角化的，各个势能可用其傅里叶级数来描述。通常，Kohn-Sham 方程平面波展开式可通过 Hamiltonian 矩阵的对角化来求解，该矩阵元素 H k + G ， k + G ′ 由方程（2.17）中波形括号项给出，所以 Hamiltonian矩阵的大小是由能量截断值决定的。
通常，对一个由许多原子组成的固体，坐标空间根据波函数的不同可分成两部分，设存在某个截断半径 r c 。① r c 以内的核区域（芯区）。波函数由紧束缚的芯电子波函数组成，对周围其他原子是否存在不敏感，与近邻的原子的波函数相互作用很小。② r c 以外的电子波函数（称为价电子波函数），由于与周围其他原子相互作用而变化明显。而在原子的所有电子中，基本上只有价电子具有化学活性，相邻原子的存在和作用对芯电子的影响不大。因此，对于固体材料，从考虑原子之间相互作用角度来看，可将电子的波函数用两部分表示：在 r c 以外的价电子波函数仍然保留其真实波函数的形状；在 r c 以内的波函数用空间变化平缓的形状来替代，而得到的电子波函数称为赝波函数。要使得到的赝波函数成为原子的一个本征态，原子势需要同步改变成某种有效势，即赝势，对应的“赝势+赝波函数”系统称为赝原子。赝原子在描述真实原子自身性质时是不正确的，但它能近似正确地描述原子与原子之间的相互作用。近似正确性的好坏，取决于截断距离 r c 的大小， r c 越大，赝波函数越平缓，与真实波函数的差别越大，近似带来的误差越大； r c 越小，误差就越小。具体计算量的大小受截断半径 r c 的影响，但越小的 r c ，赝波函数振荡部分计入得就越多，需要的平面波展开基底就越多，计算量也将增大，因此高的精度与少的计算量两者总是矛盾的。与 LAPW、LMTO等精度较高的第一性原理计算方法相比，平面波赝势是计算量较小的方法，适用于计算精度要求不太高的体系。
在从头算赝势中，常用的有模守恒（Norm-Conserving）赝势 [27，28] 、超软赝势（Ultrasoft Pseudopotential） [29] 和 Trouiller-Martins [30，31] 赝势等。模守恒赝势的特点是，赝波函数和全电子波函数 r c 以内与真实波函数的形状和幅度都相同，即模守恒，在 r c 以内动能变化缓慢且没有较大的动能，该赝势能产生正确的电荷密度，适合做自洽计算。超软赝势的特色是让波函数变得更平滑，也就是所需的平面波基底函数更少，放松了对模守恒的要求，允许赝波函数在 r c 以内尽可能软，大大降低了截断半径。本书计算均采用 Vanderbilt [29] 型的超软赝势。Trouiller-Martins [30，31] 赝势的合理性在于构建更加平滑的赝势时使用了附加系数，赝势所对应的波函数不仅仅与真实势对应的波函数具有同样的能量本征值，而且在芯区以外与真实波函数的形状和幅度都相同。
在 Debye模型中，物质的非平衡 Gibbs函数 G *（ V ； P ， T ）可表示为 [32]
式中， E （ V ）是总能量， P 和 V 分别表示压强和体积， A vib 是 Helmholtz振动自由能。
考虑简谐近似 [32] ， A vib [33～35] 可以表示为
式中， Θ （ V ）是 Debye 温度， n 是分子式中原子的个数。对于各向同性固体来说，泊松比 σ =0.25 [36] ，Debye温度可以由下面公式得到
式中， M 是每个原胞中分子的质量， f （ σ ）的表达式由文献[34，37]给出， B s 是用来表示晶体压缩率的绝热体模量， B s 的表达式为 [31]
非平衡 Gibbs函数 G *（ V ； P ， T ）对体积求最小值，即
通过式（2.20）可以得到晶体的热状态方程 V （ P ， T ），等温弹性模量 B T 、热容 C v 、热膨胀系数 α 分别用如下形式表示。
弹性性质对于固体的力学性能、动力学性质有着重要的意义，而弹性常数不仅对地球物理、材料、化学等方面的研究非常重要，还跟热力学性质有关。弹性常数决定了晶体在外场作用下的响应，通过体积模量、剪切模量、杨氏模量、泊松比等性质来描述，这些性质在描述材料的强度等方面有很重要的作用。通过计算弹性常数还可以得到两个近邻原子平面之间的键的性质、结构的稳定性等。实验上可以通过超声波测量、中子散射、布里渊散射 [38，39] 等方法来测量弹性常数，而对于实验上不易测量的材料，可以通过理论方法来预测其力学稳定性 [39～48] 。
对于各向同性的、均匀的材料，应力和应变的响应关系可以通过 Hooke定律 [49] 来描述，表示为
式中的四阶张量 C ijkl 为弹性常数张量，共有 81 个张量元。晶体的弹性常数还可定义为能量对 Lagrangian无穷小应力张量 σ kl 的二阶导数，由 Voigt约定，应力张量可转化为
对应变张量为
由 Voigt约定简化后，Hooke [49] 定律的应力应变关系可表示为
可简写为
式中， c ij 就是本书要研究的弹性常数。 c ij 是一个 6×6 的矩阵，根据对称性，对于大多数晶体来说，弹性常数仅有 21个非零独立分量；一般地，对于立方晶体来说，这 21 个非零独立分量可以约化为 C 11 、 C 12 、 C 44 ，对于六角晶系的晶体的独立弹性常数为 C 11 、 C 12 、 C 13 、 C 33 、 C 44 。