摘要：第一性原理方法（first-principles method 或 ab initio method）是指一种用于数值求解薛定谔方程的方法。晶格动力学方法（lattice dynamics），是指通过原子间力常数（IFCs）计算声子色散关系和弛豫时间的方法。玻

第一性原理的 PBTE（Phonon Boltzmann Transport Equation，声子玻尔兹曼输运方程）其实就是下面三个计算方法的结合：

第一性原理方法（first-principles method 或 ab initio method）是指一种用于数值求解薛定谔方程的方法。

晶格动力学方法（lattice dynamics），是指通过原子间力常数（IFCs）计算声子色散关系和弛豫时间的方法。

玻尔兹曼输运方程（BTE）则是通过将声子色散关系和弛豫时间作为输入，计算体系总热导率的热传输模型。

第一性原理（first-principles）方法，又称从头算（ab initio）方法，指在量子力学框架下，通过最少的经验参数来求解材料中的电子结构、原子间相互作用等微观信息。它主要基于对薛定谔方程（或相应的有效方程）的数值解，从而得到材料的基态电子密度、能带结构、总能量以及其他相关物理性质。

在固体物理与材料科学中，第一性原理方法大多以 密度泛函理论（Density Functional Theory, DFT）为基础。由于完整的多体薛定谔方程难以直接求解，DFT 提供了一种在电子密度空间中处理多体问题的有效途径，被广泛应用于各种体系（如金属、半导体、绝缘体，以及分子、表面等）中。

密度泛函理论的核心思想源自于 Hohenberg-Kohn 定理和 Kohn-Sham 方程：

Hohenberg-Kohn 定理指出，体系基态的所有物理量都可以视为电子密度 的泛函；也就是说，只要知道体系的基态电子密度，就能唯一确定体系的哈密顿量和能量等信息。

Kohn-Sham 方程将真实多电子体系映射为一个关于独立电子的辅助体系，通过一个合适的交换-关联势（exchange-correlation potential）来近似地描述真实电子间的相互作用。这样，求解含有电子-电子相互作用的多体问题，就转化为求解一组类单电子方程：

其中：

为原子核对电子的库仑势；

为经典 Hartree 库仑势；

为交换-关联势，包含所有复杂的多体效应。

在实际计算中，人们需要选用适当的交换-关联泛函（如局域密度近似 LDA、广义梯度近似 GGA 以及更高阶杂化泛函等），并结合赝势（pseudopotential）或全电子方法（如 PAW、APW+lo）等来有效处理实空间和动量空间。

在研究热导时，关键的一步是获取材料的原子间力常数（interatomic force constants, IFCs）。这些力常数用于构造动力学矩阵并计算晶格动力学性质（如声子色散关系、声子寿命等）。从第一性原理出发，人们常用以下两种主流方法来获取 IFCs：

1.2.1 线性响应方法（Density Functional Perturbation Theory, DFPT）

在密度泛函微扰理论框架下，对晶体的周期性结构施加一个小的周期性扰动，求解线性化的 Kohn-Sham 方程，从而直接得到势能对原子位移的一阶、二阶甚至三阶响应。

DFPT 具备在声子波矢空间中直接计算各种力常数和声子模的优点，且无需在实空间中进行大规模有限位移计算。

通过 DFPT 可获得二阶（简谐）和三阶（非简谐）原子间力常数，为后续的声子色散和散射分析打下基础。

1.2.2 有限位移法（Finite Displacement Method, FDM）

在有限位移法中，通过在超胞（supercell）模型中对某个原子施加小的离散位移（例如 0.01 Å），计算体系总能量或力的变化，从而数值求解二阶或更高阶 IFCs。

需要构造足够大的超胞，以保证相互作用在边界处衰减到可以忽略；

需要对各种不同位移构型做多次能量/力计算，以提取完整的二阶或三阶 IFCs。

无论是 DFPT 还是有限位移法，都需要在 DFT 水平上做能量或力的计算。由于声子性质对计算精度相当敏感，需要在计算流程中仔细选择交换-关联泛函、平面波截断能（plane-wave cutoff）、-点采样密度等，以确保得到高精度的原子间力常数。

在简谐晶格动力学（harmonic lattice dynamics）中，通过二阶原子间力常数（IFCs）可以获得声子的色散关系。

一旦得到了色散关系，就可以计算每一个声子模 的比热容。

声子模 的群速度定义为频率对倒易空间坐标的梯度。

声子的弛豫时间 需要通过非简谐晶格动力学（anharmonic lattice dynamics）计算获得，该计算需要使用二阶（简谐）与高阶（非简谐）IFCs。

对于一个周期性晶体，若原子在平衡位置附近仅发生小幅度振动，则体系总势能 可以展开为泰勒级数：

其中：

：平衡态势能；

：第 个原子在 方向的位移；

：二阶原子间力常数；

：三阶原子间力常数；

：更高阶项。

由于在平衡态时对任意原子 都有 ，因此在势能展开中不会出现一次项。若在简谐近似下忽略所有高于二阶的项，则势能只保留二阶 IFCs，可以得到运动方程。

如果晶体中第 个原子是第 个晶胞中的第 个原子，而 对应于晶胞 中的第 个原子，那么根据牛顿第二定律，有：

其中：

：为同一晶胞中第 个原子的质量；

：在第 个晶胞中、第 个原子沿 方向的位移。

若假设解的形式为平面波：

将解代入运动方程，可得到如下形式的特征值方程：

其中， 是动力学矩阵，其表达式为：

求解特征值问题，就可以得到声子色散关系 ，以及相应的特征向量 。

在实际材料中，声子不会无限传播而无损耗，而是会因各种散射机制而出现有限的传播寿命。

在玻尔兹曼输运方程（BTE）框架下，这些散射过程常被归纳进一个总的散射率 或等效的弛豫时间。

在非金属晶体中，最主要的散射通常来自声子-声子相互作用，亦即晶格的非简谐振动效应。其理论根源在于势能展开中的三阶、四阶及更高阶项。

量子力学框架下，晶体哈密顿量可分为简谐项与非简谐项两部分。通常先对简谐部分做近似求解（将每个正常模视为量子简谐振子），再将非简谐相互作用作为微扰来处理。

三声子散射是最常见且主导的非简谐散射类型：一个声子可以与另一个声子“湮灭并产生”第三个声子（组合过程），或者一个声子“衰变”成两个声子（拆分过程）。

当声子波矢 , 和 以及相应模频率 ， 和 满足能量和动量守恒（含 Umklapp 过程）时，就会发生三声子散射。

用 表示声子模 （由波矢 和声子分支 组成）。

三声子过程的跃迁概率为：

其中：

表示动量守恒（含 Umklapp 矢量 的 Kronecker delta；

表示能量守恒的 Dirac delta；

分别对应组合和拆分两类过程；

为三阶非简谐相互作用矩阵元，取决于所涉及的三个声子的特征向量和三阶 IFCs。

当晶体中存在不同质量或不同类型的原子（如同位素、点缺陷等）时，会导致局部力常数或质量分布的微扰。

散射率通常用下述经验公式近似：

其中，。 是杂质原子 的浓度。平均质量为 。

根据傅里叶定律，材料的热导率 描述其传热能力，定义为：

其中， 是热流矢量， 是温度梯度。若材料具有各向异性，热导率 是一个张量。

若要从微观层面（如声子特性）来预测 ，通常需要使用玻尔兹曼输运方程（BTE）来描述声子群体在有限温度梯度下的非平衡分布。

若假设系统处于稳态，且温度梯度足够小，可将 BTE 线性化，得到下式：

其中 ，即偏离平衡态的部分。

为声子模 的分布函数，

表示平衡态（玻色-爱因斯坦）分布， 即偏离平衡态的微扰部分，

为该声子模的弛豫时间，

为声子模的群速度。

声子是玻色子，其平衡态分布服从 玻色-爱因斯坦统计：

其中， 是约化普朗克常数， 是玻尔兹曼常数， 是声子模的频率， 是温度。

由声子贡献的热流（忽略电子热导等其他贡献）可写为

其中， 是系统体积，可以通过单位胞体积 与波矢 的总数 得到：。

将上式与傅里叶定律 进行对比，可以得到热导率张量的表达式：

其中， 为单个声子模的体积比热容。

若材料近似各向同性，则可进一步简化为

单模弛豫时间近似（SMRTA）方法的假设是：在求解某一声子模 的偏离分布时，其他所有模都保持在平衡分布 。，即：

这会简化三声子散射项的计算，只需关注声子模 自身对平衡态的偏离。

单模弛豫时间近似下的弛豫时间：

其中，上标 0 表示此弛豫时间是基于 SMRTA 获得的第零近似解。

为了克服 SMRTA 的局限性，可采用 全迭代方法（Full Iterative Method）来求解 BTE。

在这种方法中，当计算某声子模 的偏离分布时，其他声子模 也允许偏离平衡态，因而多声子模之间的散射耦合需要自洽求解。

初始猜测

常以 SMRTA 给出的 作为第零次近似。

自洽迭代

将所有声子模的偏离量耦合起来写成一个方程组，通过迭代求解，直到收敛。这样能更精确地体现三声子 Normal 过程所导致的“声子再分配”效应。

由于全迭代方法包含更多相互作用细节，对高热导率材料或低温下 N 过程活跃的体系，往往能得到更准确的结果。不过，它对计算资源的需求也更高，算法实现更复杂。

计算声子热导率的一般流程包括以下几个关键步骤：

原子间力常数是进行晶格动力学和热输运计算的基础数据，可通过以下三种方式获得：

经典势函数：适用于结构简单、已有成熟势函数的材料；

机器学习势（ML potentials）：近年来快速发展，可兼顾精度和效率；

第一性原理方法（First-Principles/DFT）：最常用，能在无任何经验参数的前提下精确描述材料的基态性质。

准确提取 IFCs 并非易事，目前主要有两种方法：

有限位移法（Finite Displacement Method）：对超胞中原子施加微小扰动，结合差分法从获得的原子力中提取 IFCs。该方法适用于经典势、机器学习势以及第一性原理计算；

密度泛函微扰理论（DFPT）：一种基于线性响应的严格方法，能直接从第一性原理中计算 IFCs，但仅适用于 DFT 框架。

理论上，IFC 可存在于晶体中任意两个原子之间。但在实际计算中，通常需要人为设定截断半径，只保留一定距离内的相互作用。因此，需要测试不同截断范围对热导率结果的影响。

此外，由于数值噪声与截断带来的误差，从第一性原理获得的 IFCs 通常不完全满足 平移不变性（translational invariance）等晶体对称性。这种不变性对于热导率预测至关重要，必须对 IFCs 做微小修正，以强制满足相应对称性。

在获得二阶和三阶 IFCs 后，可使用非简谐晶格动力学方法结合玻尔兹曼输运方程（BTE）计算热导率。常用的求解方式包括：

单模弛豫时间近似（SMRTA）：假设除目标声子模外其他模处于平衡态，计算简单，适用于快速估算；

全迭代解法（Iterative Solution）：考虑所有声子模之间的相互影响，能更准确地反映三声子 N 过程的贡献，适用于高导热或低温系统。

该方法无需任何拟合参数，仅依赖材料的初始原子结构即可预测热导率，具有高度的通用性和准确性。大量研究表明，该方法预测结果与实验数据高度一致。

尽管该方法非常强大，但仍可能受到以下数值因素的影响：

三阶 IFC 的精度依赖于 DFT 本身的收敛性、计算泛函以及位移步长的选择；

IFC 截断距离的选取可能引入系统性误差；

第一布里渊区内数值积分精度（如 -点网格密度）也会影响最终热导率的准确性。

即便如此，该方法仍被广泛认为是目前预测晶格热导率最可靠的理论工具之一。它不仅可以输出总热导率，还可获得模分辨热导率贡献，并进一步用于研究界面热导、纳米结构热输运等复杂体系。

目前已有多个开源软件包可实现上述流程，包括：ShengBTE，phono3py等；

这些工具通常可与主流第一性原理软件包（如 VASP、Quantum ESPRESSO、ABINIT等）良好对接，实现从结构优化到 IFC 提取、再到热导率预测的自动化计算流程。

一些之前的相关推文：

Quantum ESPRESSO：能量自洽 scf 计算

VASP 计算基础与结构优化过程

QE计算：声子色散关系

石墨烯色散关系：理论与QE计算

QE+phono3py 热导率计算

参考文献：

1. Mingo, N.

et al. Ab Initio Thermal Transport. in

Length-Scale Dependent Phonon Interactions

2. Bao, H. et al.

ES Energy Environ.(2018).

来源：科学的十分