1. 项目概述为什么势函数是辐射损伤模拟的“命门”搞材料模拟尤其是辐射损伤这块的同行估计都深有体会分子动力学模拟的结果准不准八成得看你的原子间势函数靠不靠谱。这玩意儿说白了就是一套描述原子之间怎么“推拉拽扯”的数学公式是模拟的物理基础。你用它来算晶格振动、扩散这些近平衡态的性质可能感觉还行但一旦涉及到高能粒子轰击、原子被撞飞、产生大量缺陷这种极端非平衡过程很多势函数就开始“露怯”了。最近几年机器学习势函数火得不行号称能逼近第一性原理计算的精度又比直接做第一性原理分子动力学快好几个数量级听起来简直是辐射损伤模拟的“梦中情势”。但问题是它真的那么神吗在不同场景下表现是否一致这就是我们这次要深扒的核心。我们这次就拿金属镍开刀。镍是核反应堆结构材料、聚变堆面向等离子体材料的关键组分理解它在辐射下的行为至关重要。手头有几套势函数两套基于不同赝势Ni和Ni pv训练的表格化高斯近似势也就是tabGAP这是机器学习势的代表还有两套经典的嵌入原子法势分别是Stoller和Bonny势。我们的目标很明确就是把这些势函数扔到辐射损伤模拟的“高压锅”里煮一煮看看它们在阈值位移能、单次级联碰撞、以及最要命的重叠级联累积损伤这些关键测试中到底表现如何。结果发现水还挺深不是所有机器学习势都能在所有环节“通吃”短程排斥相互作用这个细节成了决定模拟结果走向的关键胜负手。2. 核心思路从静态测试到动态“实战”的评估体系评估一个势函数在辐射损伤模拟中的适用性不能只看它算出来的晶格常数、弹性模量准不准。那只是“体检”的常规项目。真正的考验在于它能否准确描述从两个原子无限接近时的剧烈排斥短程到正常晶格距离的平衡相互作用中程再到缺陷形成、迁移等涉及能量变化的复杂过程。我们的评估逻辑是层层递进的从简单到复杂从静态到动态构建一个完整的“压力测试”链条。2.1 评估维度的逻辑拆解首先最基础的检查是看势函数对原子间相互作用的描述是否“健康”。我们画了二聚体曲线就是计算两个镍原子在不同距离下的总能量。这能直观看出短程排斥部分的“硬度”。如果这里就和第一性原理计算的结果差很远那后续的高能碰撞模拟基本可以不用看了结果肯定跑偏。紧接着我们做了准静态拖拽计算模拟一个原子在晶格中沿着特定方向被缓慢拖开时周围晶格原子的弛豫和能量变化。这个测试比二聚体曲线更进了一步它包含了晶格环境的影响能反映势函数对点缺陷形成能的预测趋势。通过了这些“静态”或“准静态”测试才能进入真正的动态模拟环节。第一步是计算阈值位移能。这是辐射损伤领域的一个核心参数指的是将一个晶格原子永久撞离其位置所需的最小动能。我们不是只算几个高对称方向就完事而是系统性地计算了不同入射方向上的阈值位移能生成一个三维的“能量球面”。这能全面评估势函数对各向异性效应的预测能力。然后就是模拟真实的辐射事件——初级离位级联。用一个初始动能为几千电子伏特的初级撞出原子去轰击晶格用分子动力学模拟整个级联过程能量沉积、热峰形成、冷却、以及最终稳定下来的缺陷空位和间隙原子对即Frenkel对数量。这一步看的是势函数在极端非平衡、高能条件下的瞬时响应能力。最后也是最接近真实辐照环境的测试重叠级联模拟。在现实中材料是持续受到辐照的新的碰撞级联会不断在已经存在损伤的区域产生。我们模拟了成千上万次级联事件在同一个模拟盒子中依次发生观察缺陷浓度是如何随着辐照剂量累积的以及最终会形成什么样的缺陷结构比如位错环、层错四面体等。这个测试考察的是势函数在长时间尺度、复杂缺陷相互作用下的表现是对其综合预测能力的终极考验。2.2 势函数选型背后的考量为什么选这四套势函数这里头有对比的巧思。两套tabGAP势函数Ni和Ni pv的核心区别在于训练它们所用的第一性原理数据采用了不同的赝势来处理镍的电子结构。Ni pv赝势 explicitly 处理了镍的3p半芯电子而标准的Ni赝势没有。这个看似技术性的选择会直接影响对电子云排布尤其是内层电子在原子非常接近时的描述从而显著影响短程排斥力。用它们俩对比就是想剥离出“机器学习框架相同但底层量子力学输入不同”所带来的影响。而Stoller和Bonny这两套EAM势是经过多年发展、在辐射损伤领域被广泛使用和引用的经验势。它们代表了在机器学习势流行之前我们所能获得的、经过一定校验的“传统智慧”。尤其是Bonny势它在描述镍及其合金的辐照肿胀等方面有一定声誉。将这组传统势函数与新的机器学习势进行对比不仅能判断机器学习势是否真的超越了前辈更能揭示不同势函数在物理描述上的根本差异点在哪里。注意势函数的评估一定要有参照系。单纯说“我的势函数预测缺陷浓度是X”没有意义必须和可靠的实验数据、更高级别的第一性原理计算结果或者其他经过广泛验证的势函数进行对比。本次工作中我们尽可能地将模拟结果与已有的实验阈值位移能数据、以及第一性原理分子动力学的结果进行交叉验证。3. 关键测试一短程相互作用与阈值位移能阈值位移能是辐射损伤研究的“敲门砖”它直接决定了有多少初始动能会转化为稳定的原子离位。如果势函数连这个基本参数都预测不准后续的级联模拟就像建立在沙滩上的城堡。3.1 二聚体曲线与准静态拖拽短程排斥的“体检报告”我们先从最基础的二聚体曲线看起。把两个镍原子从很远拉到几乎重叠计算总能量变化。理想情况下这条曲线在短程部分应该与第一性原理计算的结果高度吻合。我们的结果显示两套tabGAP势函数在短程部分出现了肉眼可见的分歧。基于Ni pv赝势训练的tabGAP势其排斥“墙”更“硬”即随着距离减小能量上升得更陡峭而基于标准Ni赝势的tabGAP势则相对“软”一些。这个差异在二聚体距离小于0.5埃时变得非常明显。相比之下Stoller的EAM势与Ni pv tabGAP的短程行为更接近而Bonny势则展现出截然不同的、更“硬”的排斥特征。准静态拖拽测试进一步放大了这种差异。当我们模拟一个原子沿100、110、111等高对称方向被缓慢拖离时需要克服的能垒高度直接反映了势函数对晶格恢复力的描述。从提供的能量差图可以清晰看到不同势函数给出的能量曲线形状和峰值差异显著。例如在111方向某些势函数预测的拖拽能峰可高达120 eV以上而另一些则低于80 eV。这直接预示着它们在阈值位移能上会有巨大差别。实操心得不要只看势函数论文里报告的平衡态性质拟合误差。一定要亲自画一下关键方向的二聚体曲线和拖拽曲线与可靠的DFT或更高精度的计算结果对比。短程部分的微小差异在级联模拟中会被指数级放大。很多现成的势函数文件可能未充分优化这一区域需要使用者保持警惕。3.2 阈值位移能的方向各向异性与平均值我们系统计算了数千个不同入射方向上的阈值位移能绘制成三维极坐标图。图像非常直观地展示了镍中TDE的强烈各向异性沿着密排方向如110通常需要更低的能量就能造成离位而沿着开放方向如111则需要更高的能量。这是由晶格通道效应决定的——原子在密排方向运动时更容易找到“缝隙”穿过而不引起巨大的晶格畸变。具体到数值不同势函数的预测可谓五花八门。对于100方向实验值大约在38 eV左右。我们的tabGAP势预测值在21-23 eV明显偏低Stoller势预测27 eVBonny势高达49 eV。在110方向实验值约21 eVtabGAP势在19-25 eV区间Stoller势23 eVBonny势又是49 eV。最夸张的是111方向实验值认为大于60 eVtabGAP势预测29-39 eVStoller势35 eV而Bonny势达到了惊人的69 eV其平均值更是接近100 eV远超其他势函数和实验认知。这些差异说明了什么首先Bonny势在几乎所有方向都给出了异常高的TDE这可能源于其势函数形式中人为加强了原子间的束缚力以提高某些平衡态性质如空位形成能的拟合精度但却牺牲了非平衡碰撞过程的准确性。其次两套tabGAP势虽然基于相同的机器学习框架但由于训练数据源赝势不同在TDE预测上也出现了不可忽视的差异尤其是在111方向。这凸显了训练数据的量子力学计算精度对最终势函数在极端条件下性能的决定性影响。表1不同势函数预测的镍阈值位移能对比单位eV势函数类型100方向110方向111方向整体平均值实验/参考值~38~2160~33-69Ni tabGAP21192942.9Ni pv tabGAP23253946.5Stoller EAM27233547.7Bonny EAM49496999.74. 关键测试二单次初级离位级联模拟阈值位移能是“门槛”而单次级联模拟则是看“进门之后”会发生什么。我们模拟了从1 keV到10 keV不同初始动能的初级撞出原子引发的级联过程。4.1 缺陷产生数量与PKA能量的关系模拟50次取平均后我们绘制了稳定Frenkel对数量随PKA能量变化的曲线。一个有趣的发现是缺陷产生数量并不与平均阈值位移能简单相关。例如平均TDE最低的Ni tabGAP势在几乎所有PKA能量下预测产生的空位数都是最高的。而平均TDE最高的Bonny势产生的缺陷数量却并非最低在某些能量区间甚至超过了Stoller势。这打破了“TDE越低产生缺陷越多”的简单直觉。原因在于缺陷最终数量不仅取决于初始离位是否容易发生TDE更取决于级联发展过程中热峰的能量耗散、原子混合、以及后续的空位-间隙原子复合效率。一个“较软”的势函数如Ni tabGAP其原子间排斥力较弱可能导致级联区域原子的碰撞截面更大能量沉积更分散形成更弥散、非球形的热峰。这种热峰冷却较慢原子有更多时间扩散和重组反而可能抑制了近距离的空位-间隙原子复合导致最终保留下来的缺陷更多。4.2 级联时空演化与热峰特征通过分析5 keV PKA下缺陷数量随时间演化的曲线我们获得了更动态的图景。所有势函数都展示了典型的级联发展三阶段快速损伤产生期~0.1 ps内、热峰冷却与复合期~0.1-1 ps、以及缺陷稳定期1 ps。然而不同势函数在峰值损伤热峰最大规模和最终稳定缺陷数量上的排序并不一致。Ni tabGAP势的峰值损伤反而是最低的但最终缺陷残留最多。这印证了上面的推测它的级联更“软”、更弥散热峰温度可能相对较低但范围广在冷却过程中复合不充分。而Ni pv tabGAP和Stoller势则表现出更高的峰值损伤和相对紧凑的热峰但最终残留缺陷却较少说明其热峰冷却更快或者原子在热峰内更容易复合。排查技巧当比较不同势函数的级联结果时不要只看最终缺陷数。一定要分析级联的时空演化过程观察热峰的大小、形状、寿命。可以使用原子温度云图、径向分布函数随时间变化等工具。有时最终缺陷数量相近的两个势函数其级联的微观物理过程可能截然不同这会影响对后续现象如缺陷迁移、聚集的预测。5. 关键测试三重叠级联与长期缺陷演化单次级联像是“闪电战”而重叠级联模拟则是“持久战”它更接近材料在反应堆内经受长期辐照的真实情况。我们在一个模拟晶胞中依次注入数千个5 keV的PKA让新产生的级联与之前残留的缺陷结构发生相互作用。5.1 缺陷浓度的饱和现象随着注入级联数量对应辐照剂量的增加所有势函数模拟的体系其缺陷浓度都逐渐趋于饱和。这是一个典型的辐照损伤现象新缺陷的产生与旧缺陷的恢复复合、湮灭、聚集形成大缺陷达到了动态平衡。然而令人震惊的是达到饱和的缺陷浓度水平在不同势函数间存在巨大差异。Ni pv tabGAP势预测的饱和缺陷浓度显著低于其他三种势函数。即使我们将模拟剂量加倍从2000次级联增加到4000次它的缺陷浓度依然稳定在低水平没有向其他势函数的结果靠拢。而其他三种势函数Ni tabGAP Stoller Bonny的饱和浓度则彼此接近。这个结果非常反直觉因为在前面的单次级联测试中Ni pv tabGAP和Stoller势的行为颇为相似但在累积损伤的长期竞争中它们却分道扬镳。5.2 最终缺陷结构的微观分析观察2000次级联后的原子构型所有势函数都预测形成了由肖克利不全位错链构成的大尺寸间隙型位错环以及由空位聚集而成的层错四面体。这是面心立方金属中典型的辐照缺陷结构。定性上看缺陷类型相似。但定量上Ni pv tabGAP势体系中的缺陷密度明显更低大尺寸的位错环结构似乎也更“稀疏”。这说明决定长期辐照损伤演化的不仅仅是单次级联产生缺陷的“数量”更是这些缺陷的“命运”。缺陷的稳定性、迁移能力、以及它们之间相互作用的势垒这些都由势函数在近平衡区域即缺陷周围弛豫后的原子构型的能量面精确性所决定。Ni pv tabGAP势可能更准确地描述了间隙原子团簇的迁移势垒或空位团簇的稳定性使得缺陷更容易在后续级联的扰动下复合或发生改变从而抑制了缺陷的累积。重要提示重叠级联模拟是检验势函数预测长期辐照行为的“试金石”但计算成本极高。在进行此类模拟前务必用单次级联、阈值位移能、点缺陷性质等测试对势函数进行初步筛选。如果势函数在这些基础测试中表现怪异那么它在重叠级联中给出错误结果的可能性极大。不要盲目相信一个势函数仅仅因为它标榜“机器学习”或“基于DFT”。6. 关键试四与实验对比的RBS/C模拟模拟的终极目标是解释和预测实验。卢瑟福背散射/沟道谱技术是实验上表征近表面区域晶格损伤的常用手段。我们通过模拟计算了经过重叠级联损伤后的镍晶体在不同深度处的RBS/C谱并与实验数据直接对比。6.1 模拟与实验谱的匹配度我们将不同势函数产生的最终缺陷结构作为输入计算了相应的模拟RBS/C谱。对比发现在采用了更合理的核能量沉积剖面进行采样而非简单假设一个统一的TDE值后Ni pv tabGAP势模拟出的谱线与实验测量结果最为接近。其他势函数包括另一套tabGAP和两套EAM势模拟出的背散射产额均显著高于实验值。这意味着Ni pv tabGAP势所预测的缺陷类型、浓度及其空间分布综合起来更符合实验观测到的晶格无序度。这是一个强有力的证据表明在镍的案例中采用更精确的赝势Ni pv来生成训练数据所获得的机器学习势函数在预测真实辐照损伤累积方面可能具有更高的保真度。6.2 采样方法的影响我们对比了两种从模拟中提取结构用于RBS/C计算的方法一种是传统的基于NRT-dpa模型并假设固定TDE如40 eV的方法另一种是直接基于模拟中实际沉积的核能量剖面进行采样。结果表明后一种方法虽然增加了计算复杂度但能更真实地反映不同势函数在TDE上的差异对最终损伤剂量估算的影响从而使模拟与实验的对比更具物理意义。忽略这种差异可能会错误地归因或掩盖势函数之间的性能差别。7. 势函数优化尝试调整短程排斥相互作用既然发现了短程排斥相互作用是导致不同势函数行为差异的关键一个很自然的想法是能否对训练好的机器学习势函数进行“微调”专门优化其短程部分由于tabGAP势函数形式明确分离了对相互作用项和多体项这为我们进行此类操作提供了便利。我们尝试了三种修改方案1) 对“较软”的Ni tabGAP势用其原始外部排斥势重新拟合其两体项2) 将Ni tabGAP势的两体项在短程部分平滑过渡到Ni pv tabGAP势的“较硬”两体项3) 反向操作将Ni pv tabGAP势的两体项在短程部分过渡到Ni tabGAP的“较软”形式。7.1 修改后的效果验证将这些修改后的势函数重新投入重叠级联模拟结果颇具启发性。对于原版Ni tabGAP势缺陷产生多通过方法1和方法2使其短程部分“硬化”后其预测的饱和缺陷浓度显著下降变得与未修改的Ni pv tabGAP势结果非常接近。反之将Ni pv tabGAP势“软化”后其缺陷浓度有所上升但仍远低于原版Ni tabGAP势。这系列实验有力地证明了两点首先短程排斥相互作用确实是控制级联损伤演化特别是最终缺陷累积浓度的关键杠杆之一。其次对于表格化的机器学习势在训练后对其两体排斥部分进行有针对性的修正或替换在技术上是可行的并且能显著改变其在非平衡模拟中的行为。7.2 优化策略的局限与展望然而这种“外科手术式”的修改必须极其谨慎。我们只是简单地基于二聚体曲线的视觉匹配进行平滑过渡缺乏系统的优化。修改后的势函数虽然在重叠级联测试中表现改变但其在点缺陷形成能、表面能、弹性常数等平衡性质上是否依然保持精度需要重新全面验证。粗暴的修改可能会破坏势函数在训练数据域内的内在一致性。一个更稳健的思路是将短程排斥的准确性作为目标从一开始就纳入机器学习势的训练过程中。例如在训练数据集里不仅包含平衡构型、弹性变形、点缺陷、声子等信息还应 explicitly 加入高能二聚体、小间距晶格压缩、甚至是从第一性原理分子动力学级联模拟中提取的高能非平衡原子构型。让机器学习模型在训练阶段就“见识”并学习到这些极端情况下的相互作用从而生成一个在宽能量范围内都可靠的势函数。8. 总结与对从业者的建议通过这一系列从简到繁、从静态到动态的测试我们对机器学习势函数在镍辐射损伤模拟中的性能有了更立体和深刻的认识。核心结论是一个适用于辐射损伤模拟的优秀势函数必须同时在两个“战场”上表现优异——既要能精确描述平衡态和近平衡态的物理性质如晶格常数、弹性模量、缺陷形成能也必须能准确刻画短程高能排斥相互作用。对于从事相关模拟的研究者我个人的体会是切勿迷信“黑箱”不要因为一个势函数标榜“机器学习”或“第一性原理精度”就全盘接受。务必进行系统的基准测试特别是针对你的研究问题所涉及的关键物理过程。对于辐射损伤二聚体曲线、阈值位移能计算和单次级联模拟是最低限度的验证。理解训练数据的“基因”仔细阅读势函数原始文献弄清楚其训练数据集是如何构建的。使用了哪种赝势是否包含了高能数据训练集对你想研究的性质如沿特定方向的离位是否有足够的覆盖本次工作中Ni和Ni pv赝势带来的差异就是活生生的例子。长期行为是终极裁判对于研究辐照损伤累积、微结构演化等长期效应单次级联的结果可能具有误导性。只要计算资源允许尽量进行重叠级联模拟或者至少研究缺陷团簇的稳定性与迁移率。势函数在长期模拟中表现出的差异往往揭示了其在描述复杂缺陷相互作用能垒上的根本不同。实验对比是金标准尽可能地将模拟结果与实验观测进行定量对比如RBS/C谱、透射电镜观测的缺陷密度与类型、电阻率变化等。即使无法完全吻合趋势上的一致性能极大增强模拟结果的可信度。审慎进行势函数修改如果确实需要对现有势函数进行修改以改进其短程行为必须采用系统的方法并在修改后对势函数的全面性能进行重新评估。记录下所有的修改步骤和参数确保工作的可重复性。机器学习势函数无疑为材料辐照模拟带来了新的强大工具但它并非“银弹”。它放大了我们对高精度训练数据的需求也要求我们以更严谨、更系统的方式去验证和解读模拟结果。这条路没有捷径唯有通过细致扎实的基准测试和物理分析才能让这些强大的计算工具真正服务于对材料辐照行为更深刻、更准确的理解。