分子动力学辐照模拟：级联碰撞损伤与缺陷演化分析

分子动力学辐照模拟在核材料研究中是核心手段，但级联碰撞初始能量设置争议大、缺陷识别方法选择影响统计结果、长时间尺度缺陷演化模拟计算量巨大，如何输出可靠损伤数据评估材料抗辐照性能，是核材料计算领域的关键技术挑战。

项目背景是一个钨(W)材料抗辐照性能评估任务。客户需要评估钨在14 MeV中子辐照下的初级离位损伤(PKA能量5 keV)特征。分子动力学辐照模拟在这种级联碰撞问题中几乎是唯一能从原子尺度直接观察缺陷形成过程的手段，但从PKA设置到缺陷识别到统计分析，技术细节密集。

体系构建是第一步。钨是BCC结构，晶格常数3.165Å。体系选取40×40×40个晶胞的立方盒子，含约128000个原子。PKA设为一个W原子，初始位置在体系中心，初始速度方向沿[135]方向（非密排方向，避免沟道效应使PKA穿透太远）。分子动力学辐照模拟中PKA能量不能直接设为5keV——如果给PKA原子速度v=sqrt(2E/m)，对应的动能5keV对应的速度约5×10⁵ m/s，需要极小的时间步长（~0.001 fs）才能稳定积分运动方程。用了变步长方案：PKA动能高于100eV时Δt=0.001 fs，50-100eV时Δt=0.1 fs，低于50eV时Δt=1 fs。整个级联过程约15 ps，总步数约20000步。

势函数选择是核心技术决策。钨的EFS势（Embedded-atom Method Finnis-Sinclair）是分子动力学辐照模拟中最常用的选择。但EFS势在描述高能碰撞阶段的N体相互作用不够精确——PKA能量>1keV时，近距离碰撞的散射截面偏差可达15%。因此用ZBL（Ziegler-Bistritzkhan-Littmark）普适屏蔽势在近距离（<1Å）替换EFS势，用样条函数平滑过渡到EFS描述范围。这个混合势方案在5keV级联模拟中验证过，离位原子数与OKMC模拟结果偏差在8%以内。分子动力学辐照模拟中势函数在近距离区域的处理精度直接决定级联碰撞的可靠性。

缺陷识别用Wigner-Seitz方法。每个原子位置与理想晶格的Wigner-Seitz胞比较：如果空位位置没有原子而间隙位置有多余原子，判定为Frenkel对。但识别阈值需要调优——原子偏离晶格位置超过0.3Å才算离位。5 keV PKA级联结束时统计结果：离位原子（DPA）峰值出现在0.3 ps时，达到168个；到5 ps时部分间隙原子复合回空位，剩余稳定Frenkel对84对（42个空位+42个间隙原子）。分子动力学辐照模拟的缺陷统计中必须区分峰值损伤和稳定损伤——两者可以差一倍以上，只有稳定损伤对宏观性能有影响。

缺陷类型分布进一步分析。84个稳定缺陷中：单空位52个、双空位6个、三空位2个，单间隙原子48个、双间隙6个、三间隙2个。间隙原子倾向于形成<111>方向的哑铃型构象——这是BCC金属间隙原子的经典构象。NRT标准（修正Kinchin-Pease模型）预测5keV级联的离位原子数N_NRT=0.8·E_PKA/(2·E_d)=0.8×5000/(2×90)=22.2，其中E_d=90eV是钨的离位阈能。实际模拟得到84个稳定Frenkel对，N_NRT预测只有22——因为NRT模型没有考虑级联内部的原子重排和高效复合过程。分子动力学辐照模拟的这个发现对理解钨的抗辐照性能至关重要：W的高效缺陷复合率使其稳定损伤远低于理论预测值。

温度效应分析是项目的延伸内容。在300K、500K、800K三个温度下分别做了5keV级联模拟（每个温度5个随机种子取平均）。300K稳定Frenkel对84±12，500K稳定67±9，800K稳定51±7。温度升高，稳定缺陷减少——因为热涨落促进间隙原子迁移和复合。分子动力学辐照模拟揭示的温度依赖性为客户评估高温服役条件下钨的辐照性能提供了定量依据。

计算效率方面，128000原子的体系在32核CPU上跑了约48小时/次级联。5个温度×5个种子=25次计算，总CPU时间1200小时。分子动力学辐照模拟的计算量主要来自极小时间步长阶段——前0.5ps的2000步占了总计算时间的60%。

https://www.keyanxueshu.com/