热力学函数模拟计算：从声子态密度到相自由能的一条完整链路

Ni-Al金属间化合物的相稳定性，是高温合金设计里一个纠缠了几十年的老问题。Ni₃Al、NiAl、NiAl₃、Ni₂Al₃……相图上一排化合物，哪些在服役温度下真正热力学稳定、哪些只是动力学亚稳态，直接影响合金的强化相选择和热处理工艺窗口。实验热力学数据并非没有——但实验测的是有限温度点，DFT可以在整个温度轴上给出连续的G(T)曲线，这是模拟的独特价值。

这个项目的目标很明确：对Ni₃Al、NiAl两种最重要的金属间化合物，分别计算从0K到1500K的吉布斯自由能G(T)，并与CALPHAD数据库做交叉验证。

声子计算的物理入口

热力学函数的核心是配分函数，配分函数的核心是态密度。对晶体来说，声子态密度（PhDOS）决定了振动自由能：

Fvib(T)=∫0∞g(ω)[12ℏω+kBTln⁡(1−e−ℏω/kBT)]dωFvib​(T)=∫0∞​g(ω)[21​ℏω+kB​Tln(1−e−ℏω/kB​T)]dω

括号里的第一项是零点振动能，第二项是温度依赖的振动自由能。所以只要算出g(ω)，整个F(T)就跟着出来了。

g(ω)怎么来？DFT做声子计算的标准路径是：先用VASP做结构优化（ISIF=3全弛豫），确保残余力<0.001 eV/Å；然后用IBRION=8或Phonopy的有限位移法构建2×2×2超胞，逐原子做位移，用DFT算每个位移构型的Hellmann-Feynman力，组装成力常数矩阵，对角化得到声子频率和本征矢。

Ni₃Al在Phonopy框架下的声子色散没有虚频，说明L1₂结构在基态是动力学稳定的。NiAl（B2结构）同样全正频。这是后续热力学积分的基础前提——如果有虚频，要么是结构没优化到底，要么是相本身在0K下不稳定，F(T)也就无从说起。

准简谐近似的适用边界

标准简谐近似下，声子频率与体积无关，热膨胀不存在。准简谐近似（QHA）放开了这个约束：在每个温度T下，寻找使Helmholtz自由能F(V,T)最小的体积V_min(T)，热膨胀系数就从V_min(T)的温度导数自然导出。

操作上是对7个不同体积（平衡体积的±3%，步长1%）分别算全套声子态密度，然后对每个体积在任意温度下计算F(V,T) = E_0K(V) + F_vib(V,T)，最后对每个T做F-V曲线的多项式拟合找到最小值。

Ni₃Al在1200K时，准简谐修正使体积膨胀了约2.7%，对应的热膨胀系数约14.3×10⁻⁶ K⁻¹，与文献中Ni₃Al的实验值13.8×10⁻⁶ K⁻¹吻合得很好。QHA在熔点的50-70%以下通常比较可靠——更高的温度下非简谐效应（声子-声子散射导致的频率重整化）开始显著，准简谐描述力有不足。

从形成能到形成自由能

0K下的形成焓是起点：对Ni₃Al，ΔH_f(Ni₃Al) = E(Ni₃Al) – [3E(Ni) + E(Al)] = -0.44 eV/atom（fcc-Ni和fcc-Al为参考态）。这个值说明Ni₃Al的形成在能量上是有利的——这也是它为什么在自然界中自发出现的主要原因。

但0K的形成焓远远不够。在1200K的典型服役温度下，自由能修正至关重要。做完声子态密度后，G(Ni₃Al, 1200K) = E_0K + F_vib(1200K) + PV项。振动自由能贡献约-0.12 eV/atom（熵效应使自由能随温度下降），而Ni和Al单质也有各自的振动自由能下降。两者相减，ΔG_f在1200K约-0.36 eV/atom，仍然为负，但已经比0K的-0.44收窄了不少。

NiAl的趋势类似：0K形成焓-0.67 eV/atom，1200K时ΔG_f收窄到约-0.55 eV/atom。两种化合物的热力学稳定性排序在0K和1200K是一致的（NiAl > Ni₃Al），但稳定性的”裕度”随温度变化——这点在相图计算中直接影响相边界的位置。

与CALPHAD的对标

将DFT-QHA得到的G(T)曲线与TCNI热力学数据库的CALPHAD参数做对比，整体趋势吻合：Gibbs自由能差ΔG在0.05 eV/atom以内（约5 kJ/mol），对应在1200K时相边界位置的偏差约10K量级。对工程相图使用来说，这个精度是够用的。

两个数据源的分歧主要出现在高温段（>1000K），原因有三：CALPHAD的参数是用实验相边界反向拟合的，天然吸收了非简谐效应和电子熵贡献；DFT-QHA只考虑声子贡献，缺少电子热激发对自由能的修正；空位浓度在高温下的贡献也被忽略了。在1500K时这些修正项加起来可能达到0.03-0.05 eV/atom，解释了大部分偏差。

结论的边界

DFT声子计算+准简谐近似这条路，对金属间化合物的热力学函数模拟计算来说，是一条经过充分验证的方法链。0K到1000K的温度范围内精度可靠，可以直接用于相稳定性判据和合金设计。超过1000K时需要注意非简谐效应和电子熵贡献的缺失，如果要求高精度，需要引入从头算分子动力学（AIMD）做非简谐自由能修正。

对实际合金设计而言，这种计算的最大价值不在单个数字的精确性上，而在于可以系统扫描成分空间——三元、四元化合物，实验测不动，DFT可以。