机器学习势函数在Ni-Al合金缺陷物理中的精准应用与验证

1. 项目概述：当机器学习势函数遇上Ni-Al合金缺陷物理

在计算材料学的世界里，我们一直在寻找一个“圣杯”：一种既能保持量子力学级别的精度，又能处理成千上万个原子、模拟纳秒乃至微秒尺度动力学的工具。传统的第一性原理方法，如密度泛函理论，精度虽高，但计算成本使其局限于数百个原子、皮秒尺度的模拟。而经典的半经验势函数，如嵌入原子法，虽然计算飞快，但其精度和可迁移性常常令人担忧，尤其是在处理远离平衡态或复杂缺陷构型时。这个矛盾在像Ni-Al这样的关键工程合金体系中尤为突出，因为其高温强度、蠕变抗性等核心性能，与材料内部微小的点缺陷、位错、层错等微观结构息息相关。

近年来，机器学习势函数的出现，让我们看到了弥合这一鸿沟的曙光。它不像传统势函数那样依赖预设的物理模型，而是像一个“学霸”，通过“阅读”大量高精度的DFT计算结果，自己总结出原子间相互作用的复杂规律。其中，矩张量势作为一种基于系统性张量展开的机器学习势，因其严格的旋转、平移和排列不变性，以及在精度和效率间的良好平衡，受到了广泛关注。

那么，一个为Ni-Al合金精心训练的MTP，到底能不能担此大任？它能否准确捕捉空位、空位团簇、反位缺陷这些“微观破坏者”的形成与演化？又能否精确描绘堆垛层错这张决定材料塑性变形的“能量地图”？这正是我们本次工作的核心。我们不仅仅是在测试一个势函数，更是在探索机器学习方法在描述复杂材料体系缺陷物理时的能力边界、优势与潜在的陷阱。对于从事合金设计、核材料、高温结构材料研究的同行来说，一个可靠的势函数意味着我们可以以前所未有的精度和尺度，在计算机中“制备”并“测试”材料，大幅加速新材料的研发进程。

2. MTP模型构建与训练策略解析

2.1 矩张量势的核心思想与算法优化

要理解MTP为何能在精度和效率上取得平衡，我们需要深入其数学内核。传统MTP将每个原子i的局部环境用一组矩张量来描述，这些张量由邻居原子j的相对位置向量r_ij构成。MTP的核心是将原子能量表达为这些矩张量分量的标量函数，这些标量是通过将矩张量进行各种缩并（类似于张量点乘）得到的。然而，随着所描述相互作用的复杂度（即“级别” level）提高，独立分量的数量会组合爆炸式增长，导致计算极其昂贵。

我们本次工作的一个关键创新，在于引入了一种基于遗传算法的优化方案，对矩张量的缩并方式进行智能筛选和优化。简单来说，我们不再盲目地考虑所有可能的张量缩并方式，而是让算法自己去寻找那些对最终能量预测贡献最大、最有效的缩并组合。这个过程就像是为MTP这个“模型”设计最精简、最高效的“指令集”。

注意：这种优化并非简单的特征选择。它是在保证模型表达能力的数学完备性的前提下，剔除冗余的、线性相关的或贡献微小的基函数，从而在几乎不损失精度的情况下，将计算速度提升近一个数量级。这对于后续进行大规模分子动力学模拟至关重要。

具体操作上，我们使用了一个包含4-body交互（level=16）的基组。经过遗传算法优化后，独立矩张量分量的数量从传统方法的数千个减少到几百个。这意味着在每次力计算中，需要评估和组合的项数大大减少。实测表明，在相同的精度要求下，优化后的MTP力场调用速度比优化前快了近10倍。这为后续开展包含数千个原子、数纳秒时长的缺陷扩散模拟扫清了性能障碍。

2.2 训练集构建：覆盖相空间的艺术

一个机器学习模型的性能上限，很大程度上由其训练数据决定。对于旨在描述复杂缺陷的MTP，训练集必须尽可能广泛地覆盖Ni-Al合金在相图、成分、缺陷构型以及原子局部环境上的多样性。

我们的训练集主要来源于主动学习框架下的分子动力学采样。具体步骤如下：

1. 初始种子：我们首先基于Ni、Ni3Al以及一些中间成分的完美晶体结构，在高温（最高至2500 K）下运行短时间的AIMD模拟，采集原子轨迹。高温有助于原子获得足够的动能，探索更丰富的局部构型。
2. 主动学习迭代：将初始轨迹的构型输入MTP进行预测，同时用DFT计算这些构型的真实能量、力和应力。比较两者的差异，筛选出MTP预测误差（通常基于最大力误差或能量误差）超过预设阈值的构型。
3. 构型增强：将这些“难例”构型加入训练集，重新训练MTP。然后，用更新后的MTP开启新一轮MD采样，再次寻找预测不准的构型。如此循环，直至在验证集上的误差收敛到满意水平。
4. 针对性缺陷采样：除了无偏的MD采样，我们还手动构建了包含点缺陷（单空位、双空位、三空位、反位缺陷）和平面缺陷（不同剪切矢量的堆垛层错）的超晶胞模型，并用DFT计算其能量，将这些明确针对缺陷的数据点加入训练集。这对于准确预测缺陷形成能至关重要。

最终，我们的训练集包含了约15,000个构型，涵盖了从完美晶体到各种缺陷、从低温到高温、从富Ni到富Al的广阔相空间。验证集则包含约3,000个未参与训练的独立构型，用于客观评估模型的泛化能力。

2.3 精度验证：与DFT的终极对标

模型训练完成后，首要任务是与“金标准”DFT进行全面的对标。我们计算了训练集和验证集上能量、原子力、应力的均方根误差。

这个结果非常令人鼓舞。首先，验证集的误差普遍略低于训练集，说明模型没有过拟合，泛化性能良好。其次，关键指标原子力的RMSE达到了0.028 eV/Å的水平。这是什么概念？在分子动力学模拟中，原子力直接决定了原子的运动轨迹和体系的动力学性质。通常认为，力误差在0.05 eV/Å以下，就可以用于可靠的分子动力学模拟。我们的MTP显然满足了这个严苛的标准，为后续研究缺陷扩散等动力学过程奠定了坚实基础。

实操心得：在评估MTP精度时，不能只看能量误差。对于动力学模拟，力的精度更为关键。同时，应力的精度则关系到能否准确预测材料的弹性性质和力学响应。一个全面的误差评估报告是判断势函数是否“可用”于特定研究场景的前提。

3. 基础物性预测：从晶格常数到声子谱

在深入复杂的缺陷世界之前，我们必须确保MTP能准确复现材料的基本物理性质。这是检验一个势函数是否“物理”的试金石。我们系统计算了纯Ni和Ni3Al的多种基础物性，并与DFT、广泛使用的Mishin EAM势以及实验值进行了对比。

3.1 静态性质：晶格常数与弹性常数

我们首先关注的是平衡态性质。如表III所示，MTP预测的Ni和Ni3Al的晶格常数与DFT结果几乎完全一致（差异在0.01 Å以内）。与实验值相比，两者都略有低估。这主要源于两个因素：第一，实验数据是在室温下测量的，而我们的计算是在0 K下进行的，忽略了热膨胀效应；第二，也是更重要的，我们使用的DFT泛函（PBE）本身对Ni这类具有强磁性的过渡金属，在预测晶格常数时就存在系统性的轻微低估。MTP“继承”了DFT的精度，也“继承”了它的系统误差。

弹性性质是材料力学行为的基石。MTP在预测Ni和Ni3Al的弹性常数（C11, C12, C44）、体模量B、剪切模量G、杨氏模量E和泊松比ν方面，与DFT结果高度吻合。例如，对于Ni3Al，MTP预测的C11为236.50 GPa，与DFT的237.56 GPa几乎无差。相比之下，EAM势虽然在弹性常数上也表现不错（因为它本身就是针对这些性质拟合的），但在其他更复杂的性质上就开始显现差距。

3.2 能量-体积关系与声子谱

能量-体积曲线反映了材料在压缩或膨胀时的响应。如图4(a)(b)所示，在很宽的体积变化范围内（约±10%的晶格常数变化），MTP完美地复现了DFT计算出的能量-体积曲线。而EAM势在偏离平衡体积较远时，与DFT结果产生了明显的分歧。这说明MTP在描述原子间相互作用随距离变化的细节上更为精准。

声子谱是晶格动力学的直接体现，关系到材料的热学性质（如热容、热导率）和动力学稳定性。图4(c)(d)的声子色散关系对比显示，MTP预测的声子谱与DFT参考结果吻合得非常好，无论是声学支还是光学支的频率和形状。EAM势则在高频区出现了一定偏差。MTP的成功意味着它能够可靠地用于研究与温度相关的性质，如热膨胀、自由能计算，以及非谐晶格动力学效应。

注意事项：基础物性的准确复现是必须通过的“体检”。如果一个势函数连平衡晶格常数和弹性常数都预测不准，那么用它来研究缺陷、相变等远离平衡态的过程，其结果的可信度将大打折扣。在选用任何机器学习势或经验势之前，务必检查其在这些基础测试上的表现。

4. 点缺陷 energetics 的精准描绘

4.1 单空位与空位团簇

点缺陷，特别是空位，是材料中原子扩散的载体，对蠕变、辐照损伤、时效析出等过程有决定性影响。我们系统计算了Ni和Ni3Al中各种点缺陷的形成能和结合能。

单空位形成能：在纯Ni中，MTP预测的单空位形成能为1.425 eV，与DFT结果1.425 eV完全一致，也落在实验值（1.4-1.74 eV）的分散范围内。在Ni3Al中，MTP预测的Ni空位和Al空位形成能分别为1.588 eV和1.710 eV，表明形成Ni空位相对更容易。这与实验给出的平均形成能1.6±0.2 eV相符。EAM势虽然数值接近，但其预测的Ni3Al中Al空位形成能（1.544 eV）甚至低于Ni空位（1.773 eV），这与DFT和物理直觉（Al原子通常更易形成空位）相悖。

双空位结合能：双空位的结合能揭示了空位间的相互作用是吸引还是排斥。在纯Ni中，MTP准确预测了第一近邻（1NN）双空位具有微弱的吸引作用（结合能+0.015 eV），而第二、三近邻则变为排斥。在Ni3Al中，情况更为复杂：两个Ni空位在1NN处排斥，在2NN处吸引；而一个Ni空位和一个Al空位在1NN和2NN处都表现为排斥。MTP成功地捕捉到了所有这些细微的符号和趋势，而EAM势在部分情况下给出了错误的符号（如Ni3Al中2NN Ni空位对的结合能）。

三空位构型搜索：为了进一步挑战MTP，我们在一个3×3×3的超胞中对Ni和Ni3Al中所有可能的对称不等价三空位构型进行了 exhaustive search（分别有144和280种构型）。如图6所示，MTP不仅复现了DFT计算出的三空位形成能/结合能的整体分布趋势，还准确预测了能量最低（最稳定）和最高（最不稳定）的构型。例如，在Ni中，最稳定的三空位是在{111}面上形成三角形构型，次稳定的是沿<111>方向排成一条线。这些细节的吻合，证明了MTP在描述多体相互作用和复杂局部环境方面的强大能力。

4.2 反位缺陷与空位迁移

反位缺陷：在Ni3Al这类有序金属间化合物中，Ni原子占据Al位（NiAl）或Al原子占据Ni位（AlNi）的反位缺陷对其力学性能有重要影响。MTP预测的AlNi反位缺陷形成能为0.957 eV，NiAl为0.159 eV，与DFT结果（1.081 eV和0.033 eV）在量级和相对大小上一致。而EAM势严重高估了AlNi的形成能（1.436 eV），甚至错误地预测NiAl的形成能为负值（-0.168 eV），这在物理上是不合理的（意味着反位会自发大量形成）。

空位迁移势垒：我们采用爬坡弹性带方法计算了空位迁移的最小能量路径和势垒。如图5所示，MTP预测的迁移势垒曲线与DFT高度吻合。在纯Ni中，MTP给出的空位迁移势垒为1.000 eV，与DFT相同。EAM势则普遍高估了迁移势垒（Ni中为1.195 eV），尽管其Ni的数值更接近某些实验值（~1.24-1.38 eV），但这很大程度上是因为该EAM势在拟合时特意包含了迁移势垒作为目标。一个值得注意的发现是，DFT和MTP都预测Ni3Al中的空位迁移势垒低于纯Ni，这意味着空位在Ni3Al中可能扩散得更快。

深度解析：MTP在点缺陷描述上超越EAM的关键在于其“表达能力”。EAM势的函数形式相对简单，依赖于电子密度等近似，难以精确描述缺陷周围严重的晶格畸变和复杂的电荷重排。而MTP通过机器学习从DFT数据中直接学到了这种复杂的能量面，因此对于训练集充分覆盖的缺陷类型，其预测精度可以非常接近DFT。然而，这也对训练集的质量提出了极高要求——如果训练数据中缺乏某种缺陷的充分采样，MTP在该缺陷上的预测就可能出现较大偏差，正如我们在涉及Al空位的三空位构型中观察到的那样。

5. 空位团簇的基态构型与神奇尺寸

空位在辐照或变形条件下会聚集形成团簇，如空位环、堆垛层错四面体或空洞，这些团簇对材料的力学性能（如硬化、脆化）有深远影响。确定不同尺寸空位团簇的最稳定构型是一个极具挑战性的问题，因为可能的构型数量随空位数量指数增长，而DFT计算大尺寸团簇成本过高。

我们利用训练好的MTP，在一个包含2048个原子（8×8×8超胞）的Ni体系中，系统搜索了从单空位到24空位团簇的基态构型。主要方法是对空位环和空洞构型进行穷举搜索，对SFT及其他可能构型则采用分子动力学退火模拟来寻找能量最小值。

研究结果揭示了一些非常有趣的规律（图9，10）：

• 3空位：最稳定构型是空位环（三角形）。
• 4和8空位：最稳定构型是空洞。
• 5空位：形成一个6空位空洞构型，并在中心附近伴随一个自间隙原子。
• 6， 10， 15， 21空位：最稳定构型是完美的堆垛层错四面体构型。这些尺寸被称为“神奇尺寸”，因为它们能形成完整、无悬挂键的SFT结构，能量上特别有利。
• 其他尺寸（如7， 9， 11-14等）：其最稳定构型也是SFT，但并非完美，通常带有一个邻近的空位，可以看作是在某个神奇尺寸的SFT基础上添加或缺失一个空位。

这项工作清晰地展示了MTP在大规模体系构型搜索中的威力。它不仅能给出能量排序，还能直观地给出原子构型，帮助我们理解团簇的稳定化机制。例如，SFT的稳定性源于其通过层错面重构，有效地消除了大部分悬挂键，降低了体系能量。

6. 空位团簇的扩散动力学：一个意外的发现

静态能量分析之外，我们更关心这些缺陷在高温下的动态行为。我们在800 K下，对包含相同空位浓度（12个空位）但以不同团簇形式（单空位、双空位、三空位、四空位、六空位）存在的Ni体系进行了2纳秒的分子动力学模拟，并计算了每个体系的均方位移。

结果出人意料（图11）：

• 单空位：在800 K下几乎不移动，整个轨迹中仅观察到一次跳跃。
• 双空位和六空位团簇：保持结构稳定，未发生分解。
• 三空位和四空位团簇：不稳定！三空位在几次迁移后分解为一个双空位和一个单空位；四空位先分解为一个三空位和一个单空位，随后三空位进一步分解。因此，在模拟后期，双空位、三空位、四空位体系的MSD曲线斜率几乎相同，因为其扩散主力都变成了双空位和单空位。
• 扩散冠军：具有SFT结构的六空位团簇（6v-SFT）扩散最快！

这个发现颠覆了“缺陷越大越不移动”的简单认知。为了探究6v-SFT快速扩散的机理，我们进行了更细致的分析。发现6v-SFT存在两种中心对称的构型。在0 K下，这两种构型之间转换的势垒非常低，仅有约0.4 eV。在800 K的MD模拟中，我们观察到了这两种构型之间的频繁、快速转换。

通过滤波处理轨迹以消除热振动噪声，我们捕捉到了其独特的迁移机制（图12）。在构型转换的瞬间，SFT核心区域的一个原子会发生一次快速、长程的“跳跃”，跨越整个团簇核心区域，然后在几飞秒内返回。这个过程仿佛在团簇中心短暂地打开了一个“通道”，极大地促进了原子的协同运动。在一次持续约500 ps的亚稳态期间，这种快速原子迁移持续发生，导致整个包含117个原子的SFT团簇质心移动了约17 Å。一旦团簇恢复到稳定的6v-SFT构型，这种快速迁移便停止。

实操心得：这个案例生动地说明了动力学模拟的重要性。静态能量最低的构型不一定在动力学上最“懒惰”。复杂的多体协同运动机制可能带来意想不到的快速扩散行为。机器学习势函数使得在原子尺度上观测并理解这类复杂动力学过程成为可能，这是传统经验势难以企及的。

7. 平面缺陷：广义堆垛层错能面

堆垛层错能是理解材料塑性变形、位错分解和孪生行为的关键参数。我们计算了Ni和Ni3Al密排面（111）的广义堆垛层错能面。

7.1 计算方法与结果

采用剪切别名法，在六层厚的超胞中，对（111）面沿[1-10]和[11-2]方向施加不同的剪切位移，计算每个位移下的能量增量，即GSFE。对于DFT，我们在16×26的网格上计算；而对于MTP和EAM，我们使用了更密的120×120网格来获得平滑的能量面。

从图13的等高线图和表V的数据可以看出：

• MTP与DFT高度一致：无论是Ni还是Ni3Al，MTP预测的GSFE面在低能区域（如本征层错ISF、反相畴界APB附近）与DFT几乎重合。在高能区域，误差有所增大，这是因为无偏MD采样很难充分覆盖这些能量不利的构型，导致训练数据不足。
• 优于EAM：MTP的整体表现明显优于EAM势，特别是在Ni3Al的高能区域，EAM的预测偏差较大。
• 关键参数：MTP预测的Ni的本征层错能γ_ISF为115 mJ/m²，Ni3Al的复杂层错能γ_CSF为196 mJ/m²，反相畴界能γ_APB为235 mJ/m²，超点阵本征层错能γ_SISF为73 mJ/m²。这些值与DFT结果及分散的实验/文献数据基本相符。

7.2 位错分解机制

通过分析GSFE面上的最小能量路径，我们可以推断位错的分解方式：

• 在纯Ni中：一个全位错1/2[1-10]会分解为两个肖克利不全位错1/6[2-1-1]和1/6[1-21]，中间夹着一片本征层错。这是一个典型的面心立方金属的位错分解模式。
• 在Ni3Al中：由于存在有序结构，一个[1-10]超点阵位错首先分解为两个1/2[1-10]位错，中间被反相畴界隔开。然后，每个1/2[1-10]位错再进一步分解为两个肖克利不全位错，中间夹着复杂层错。

这些预测与文献中基于DFT和实验观察的位错理论完全一致，证明了MTP在描述平面缺陷和位错核心结构方面的可靠性。

注意事项：在计算Ni3Al的反相畴界能时，我们发现几何理想的APB构型在经过结构弛豫后，原子位置会发生微小偏移，能量进一步降低。因此，我们采用了弛豫后的能量最小值作为γ_APB，这导致我们的计算值（163 mJ/m²，DFT）低于一些文献报道的未弛豫或不同计算方法的数值。在比较不同工作的层错能时，需要特别注意计算方法和弛豫细节。

8. 极限测试：熔点预测与模型局限性

8.1 熔点预测

作为对MTP模型外推能力的终极测试，我们采用固液共存法预测了Ni和Ni3Al的熔点。我们构建了包含8000个原子、具有（001）固液界面的超胞，通过分子动力学模拟寻找两相共存的温度。

• Ni：MTP预测的熔点约为1597 K，低于实验值1728 K。
• Ni3Al：MTP预测的熔点为1630 K，与实验值1645 K吻合得很好。

需要强调的是，我们的训练集中并没有刻意包含任何熔体或液固界面的结构。MTP能够相对准确地预测熔点，尤其是对于成分和结构更复杂的Ni3Al，这充分展示了其出色的外推能力。对于纯Ni的偏差，可能源于PBE泛函本身对Ni的熔化性质预测存在系统误差，以及磁性效应的缺失。

8.2 模型局限性与挑战

尽管MTP在多数测试中表现优异，但我们也在工作中清晰地认识到了其当前的局限性：

1. 对能量接近的缺陷构型分辨力不足：当不同缺陷构型的形成能差异非常小（如小于0.03 eV/atom）时，MTP可能无法准确判断它们的相对稳定性顺序。这要求MLIP的精度需要达到meV/atom级别，而这通常需要更复杂的基组，会牺牲计算效率。
2. 训练集覆盖度决定预测上限：对于训练数据采样不足的区域，如涉及大量Al空位的复杂团簇，MTP的预测误差会增大。这凸显了通过主动学习持续扩展和优化训练集的重要性。
3. 继承DFT的系统误差：MTP的精度上限受限于其训练数据——DFT计算。例如，使用PBE泛函可能导致表面能、空位形成能的系统偏差，这些偏差会被MTP学习并继承。
4. 计算效率与精度的权衡：虽然我们的遗传算法优化大幅提升了效率，但相比于最简单的EAM势，MTP的计算成本仍然高出1-2个数量级。在需要极长时程或极大体系模拟时，需权衡精度与成本。

9. 总结与展望

通过这项系统性的研究，我们验证了基于遗传算法优化的矩张量势在描述Ni-Al合金体系，特别是其复杂缺陷物理方面的强大能力。它不仅能够高精度地复现DFT级别的各种基本物性（晶格、弹性、声子、熔点），更能准确地预测点缺陷、空位团簇的形成/结合能、迁移势垒，以及平面缺陷的层错能面。在动态行为上，它揭示了6v-SFT团簇独特的快速扩散机制，这是传统模拟方法难以发现的。

与广泛使用的Mishin EAM势相比，MTP展现出了全面的优势。EAM势在拟合目标属性上表现良好，但对于未在拟合目标中的复杂缺陷（如特定构型的空位团簇、层错能面细节），其误差可能很大。MTP则提供了一种更通用、更可靠的解决方案。

当然，没有完美的模型。当前MTP的挑战在于如何以可接受的计算成本，将精度推向meV/atom量级，以分辨能量极其接近的构型；以及如何通过更智能的主动学习策略，确保训练集能覆盖所有感兴趣的区域。

这项工作为开发面向Ni-Al合金的通用型机器学习势函数奠定了坚实基础。未来的方向很明确：将当前模型与主动学习框架深度结合，在模拟中动态发现并补充新的、难以预测的构型到训练集中，不断迭代优化。同时，将研究拓展到更广泛的成分范围、更多的缺陷类型（如位错、晶界），以及更极端的环境条件（如高温高压、辐照），最终构建一个真正“无所不能”的Ni-Al合金原子模拟利器。对于材料模拟实践者而言，这项研究提供了一个完整的范例——如何严谨地评估一个MLIP，理解其能力边界，并有效地将其应用于解决实际的科学问题。