不同温度下α-Fe中级联碰撞分子动力学模拟研究

曹 晗，贺新福，王东杰，吴 石，贾丽霞，豆艳坤，杨 文

(中国原子能科学研究院 反应堆工程技术研究部，北京 102413)

BCC铁基材料(α-Fe)是反应堆中的重要结构材料，低温辐照脆化是这类材料服役面临的一个主要问题[1]。辐照导致材料中形成的基体损伤是引起辐照脆化的关键因素之一。在辐照材料中，基体损伤起源于入射粒子与晶格点阵原子相互作用引发的级联碰撞[2]，此过程会产生大量的缺陷。由于级联碰撞过程发生在原子尺度，并在约10 ps时间内完成，实验上无法直接观察，在此时间和空间尺度下，分子动力学方法是最有效、应用最广泛的一种模拟级联碰撞过程中缺陷产生与演化的方法。与此同时，级联碰撞的分子动力学模拟对于辐照机理认识，如缺陷退火机制、位错环形核长大机制，具有重要的指导意义，并能为动力学蒙特卡罗(kinetic Monte Carlo, KMC)和速率理论(rate theory, RT)等其他微观结构演化模拟方法提供如缺陷的类型、空间分布、存活率和成团率等必要的参数，其中，存活率为缺陷存活数量NFP占离位原子数量Nd的比例，成团率PC为间隙原子/空位团簇中含有的缺陷数量NC占缺陷存活数量的比例。

自20世纪90年代以来，利用分子动力学模拟方法研究α-Fe的级联碰撞一直是辐照损伤模拟的热点问题：王宁等[3]通过对α-Fe的级联碰撞模拟，得到了自间隙原子的取向与迁移规律；杨莉等[4]研究了He参与下α-Fe在低温时的级联碰撞过程，分析得到了He泡的形成机制；Peng等[5]研究了100 keV以上初级离位原子(PKA)引发的级联碰撞，并发现了位错环的生成，提出了〈100〉间隙型位错环形成的punch-out机制。然而，截至目前，在α-Fe的级联碰撞模拟研究中，依然存在着亟待解决的问题：业界一直尝试在级联碰撞产生的缺陷存活数量与PKA能量之间建立准确的数量关系，1955年，Kinchin和Pease建立Kinchin-Pease模型[6]，首次预测了离位原子数量与PKA能量的数量关系，经Norgett等[7]于1975年进行修正后，建立了NRT模型；对于缺陷存活数量，国内外多项研究显示[8-10]，级联碰撞产生的稳定缺陷数量与PKA能量并不是NRT模型的正比例关系，而是幂次关系；然而，有研究显示，由于PKA能量达20 keV以上时将会发生子级联[11]，导致不能用同一个拟合公式描述全部能量范围内存活缺陷数量与PKA能量的数值关系，不同PKA能量范围的模型存在区别的机理尚不明确，这种区别也会极大地影响模型的准确性和简洁性。

本文将利用分子动力学方法进行一系列不同温度、PKA能量条件的α-Fe级联碰撞模拟研究，研究级联碰撞后产生的缺陷构型和分布，拟为团簇、位错环等缺陷在微观尺度上演化的RT和KMC等模拟方法提供必要的参数，以确定温度与不同PKA能量范围下NFP-EPKA拟合曲线存在区别的相关性，为修正该模型奠定基础，进而完善反应堆用BCC铁基材料的多尺度模拟。

1 计算方法

本工作采用Marinca等在2007年建立的M07纯铁势函数[12-13]，该势函数与常用的M03和A04两种势函数相比，模拟结果能更准确地与第一性原理以及实验结果符合[12]。利用分子动力学软件Large-scale Atomic/Molecular Massively Parallel Simulator(LAMMPS)[14]研究了500、600和700 K 3个温度条件下BCC纯铁材料的级联碰撞过程，模拟盒子尺寸根据PKA能量大小分为80a0×80a0×80a0与150a0×150a0×150a0两种，其中a0表示晶格常数。计算过程中采用三维周期边界条件，在级联碰撞开始前将系统升温至模拟辐照温度并保温20～50 ps以确保系统达到热稳定状态。级联碰撞过程须在模拟初始状态下引入PKA，PKA的能量与模拟的盒子内原子数量和模拟时间的具体对应关系列于表1。为避免产生沟道效应，并尽量使入射方向随机分布，根据前人研究经验[11]，本工作针对每种PKA能量分别沿〈122〉、〈135〉和〈235〉 3个方向引入PKA。为使模拟结果达到统计效果，每种模拟条件计算3次。引入PKA后，每隔0.05 ps输出1次所有原子的位置信息，0.5 ps后，每隔0.25 ps输出1次，5.5 ps后，每隔0.5 ps输出1次，用以分析级联碰撞过程中缺陷的演化过程。为保证级联碰撞充分退火，根据PKA能量的不同，模拟总时长也有区别，PKA能量低于40 keV时，模拟时间为15 ps，PKA能量为40 keV和50 keV时，模拟时间分别为30 ps和40 ps。

表1 不同PKA能量的级联碰撞模拟盒子内原子数量Table 1 Number of atoms in displacement cascadessimulation boxes with different PKA energy

级联过程模拟完成后，采用Open Visualization Tool(OVITO)软件中的Wigner-Seitz晶胞方法判断间隙原子和空位[15]，缺陷类型根据晶胞内原子的数量判定。对级联碰撞模拟结果进行缺陷团簇分析，采用的间隙团簇截断半径为第二近邻距离，空位团簇截断半径为第三近邻距离，进而利用式(1)计算缺陷的成团率。

(1)

利用OVITO软件中的dislocation analysis(DXA)方法对模拟完成后产生的团簇进行位错环分析。

2 结果与讨论

2.1 级联碰撞模拟过程

一次完整的由中子引发的级联碰撞包括4个过程，分别为缺陷的产生、离位峰的出现、缺陷的湮灭以及退火形成稳定缺陷[2]。图1为600 K时由30 keV的PKA在含有1 024 000个原子的α-Fe中导致的级联碰撞缺陷演化过程可视化截图。从图中可看出，30 keV的PKA形成后，在0.1 ps的时间内便使PKA附近产生大量点缺陷；此算例在1 ps前后，Frenkel缺陷数量达到峰值；之后，经过14 ps的缺陷 复合，体系内最终产生了3个1/2〈111〉间隙型位错环、1个〈100〉空位型位错环、不同尺寸的间隙和空位团簇以及多个分布于基体中的点缺陷。在此过程中，由于空位的迁移能较间隙原子的迁移能高[16]，导致空位型缺陷不易迁移，所以空位型的位错环、团簇和单空位缺陷较间隙型缺陷更靠近PKA区域，从而呈现出间隙原子包围着空位的现象。在此算例中，产生的间隙原子快速扩散逃离了级联区域，而空位由于扩散能力较差，大部分仍留在级联区域，从而导致该区域形成空位型位错环[17]。

图1 600 K时30 keV PKA在α-Fe中导致的级联碰撞过程中产生的缺陷演化过程可视化截图Fig.1 Visual screenshot of defect evolution during displacement cascades induced by 30 keV PKA in BCC iron at 600 K

2.2 缺陷构型

对于级联碰撞后产生的稳定缺陷(图2)，主要分为3类。1) 弥散分布于基体中的点缺陷，包括间隙原子和空位两种，间隙原子以哑铃状的自间隙原子(dumbbell)存在，如图2b所示，对dumbbell中的两原子位置坐标进行分析，考虑到热振动会导致原子的位置发生微小偏移，可判断出dumbbell为〈110〉取向，点缺陷的迁移能相对较低，在长时间(s量级以上)的演化中，易相互结合或扩散至晶界、位错、表面等界面处而被湮灭，实际上，中子的不断入射和点缺陷的湮灭会使基体中的点缺陷浓度维持在较稳定的水平。2) 由少量间隙原子或空位组成的相应类型的小型团簇，间隙型小团簇中最近邻原子的排列方向与哑铃状自间隙原子一致，如图2c所示；而在空位型小团簇中两个最近邻空位的排列方向有〈100〉和〈111〉两种，如图2d所示。3) 由大量间隙原子或空位组成的位错环，在所有的算例中，主要形成了两种类型的位错环，分别为1/2〈111〉间隙型位错环和〈100〉空位型位错环，在间隙型位错环中，如图2e统计显示，两个最近邻原子的相对方向大部分为〈111〉，少量为〈112〉和其他方向，这是因为自间隙原子大量聚集时，〈111〉方向的dumbbell组成的团簇结合能更高[18]；在空位型位错环中，两个最近邻空位的相对方向则大多为〈100〉方向。

a——600 K时30 keV PKA在α-Fe中引发的级联碰撞产生的缺陷整体形貌；b——dumbbell中两原子取向矢量(由最近邻两原子位置坐标相减得到，余同)；c——间隙型小团簇中最近邻原子取向矢量；d——空位型小团簇中最近邻空位取向矢量；e——1/2〈111〉间隙型位错环中最近邻原子取向分布柱形图；f——空位型位错环构型图2 级联碰撞后产生的缺陷Fig.2 Defects produced by displacement cascades

2.3 缺陷存活数量

图3定量给出3种温度下不同PKA能量级联碰撞产生的缺陷数随时间的演化曲线。经分析得知，在误差允许范围内，随着PKA能量或温度的升高，出现缺陷数量峰值以及最终达到稳定状态所需的时间均会变长。图4为3种温度下不同PKA能量引起的级联碰撞导致的缺陷存活数量图。

图3 不同温度下级联碰撞缺陷数随时间的演化Fig.3 Evolution of number of displacement cascades’ defects with time at different temperatures

图4 不同温度下级联碰撞缺陷存活数量随PKA能量变化Fig.4 Number of survival defects at different temperatures vs energy of PKA

1975年，Norgett等对Kinchin-Pease模型进行了修正，建立了式(2)所示的NRT模型，定量描述了在级联碰撞中产生的离位原子数量与PKA能量的数量关系。

(2)

其中，Ed为靶材料的离位阈能。经多年的模拟研究，在NRT模型的基础上，业界普遍认为级联碰撞产生的稳定缺陷数量与PKA能量并不是NRT模型的正比例关系，而是幂次关系，即：

(3)

其中，a、b为系数，b取0.8左右，a和b均由级联碰撞模拟后得到的缺陷数量和PKA能量拟合公式得出。

图5 EPKA=30 keV前后拟合曲线中的bFig.5 b in fitting curve before and after EPKA=30 keV

2.4 团簇分析

图6为3种温度下不同PKA能量引起的级联碰撞导致的间隙和空位团簇尺寸、数量三维柱形图。从图可看出，PKA能量越高，出现大尺寸团簇的概率越大，团簇总数量越多；温度的升高则会使低能PKA导致的级联碰撞发生间隙团簇的形核长大行为变得容易，空位团簇则无明显区别。图7展示了3种温度下不同PKA能量引起的级联碰撞导致的间隙原子和空位的成团率。可看出，在误差允许的范围内，PKA能量的升高会使间隙原子和空位的成团率增大，温度的升高会使这种增大的趋势变缓，即在700 K时缺陷的成团率基本不随PKA能量变化，分别稳定于50%和40%左右，因此，在PKA能量变化时，高温时b的变化量与低温时的相比较小。

a——T=500 K，间隙团簇；b——T=600 K，间隙团簇；c——T=700 K，间隙团簇；d——T=500 K，空位团簇；e——T=600 K，空位团簇；f——T=700 K，空位团簇图6 级联碰撞间隙和空位团簇尺寸和数量三维柱形图Fig.6 3D cylindrical map of size and quantity of interstitial and vacancy clusters of displacement cascades

图7 不同温度下间隙原子和空位成团率Fig.7 Cluster fraction of interstitials and vacancies for different temperatures

3 结论

通过对α-Fe在500、600和700 K下的级联碰撞进行一系列的分子动力学模拟研究，得到以下结论。

1) 级联碰撞会使体系内产生大量辐照缺陷，包括点缺陷、间隙和空位团簇等，缺陷的尺寸决定了缺陷内最近邻原子的排列方向：间隙型小团簇中最近邻原子的排列方向与哑铃状自间隙原子一致，均为〈110〉方向；间隙型位错环中两个最近邻原子的相对方向大部分为〈111〉。

2) 高能PKA使得产生的位错环级别的大尺寸团簇的概率升高，从而有更多数量的缺陷存活，导致同一温度下不同PKA能量阶段的NFP-EPKA拟合曲线的b存在差别。

3) 在分段NFP-EPKA拟合曲线中，700 K时，EPKA=30 keV前后两阶段b相差最小，原因在于：温度的升高会使低能PKA导致的级联碰撞发生间隙团簇的形核长大行为变得容易，使得700 K时缺陷的成团率基本不随PKA能量变化，从而导致不同PKA能量引发的级联碰撞产生的缺陷数量增长趋势未发生明显变化。