聚变堆金属材料中子辐照多尺度计算模拟

吕广宏

(北京航空航天大学)

中子辐照损伤是使聚变堆材料服役性能劣化的最主要环境因素之一，将严重影响聚变堆运行的安全性和稳定性[1]。但中子辐照下材料服役行为评价难度大、费用高、周期长，造成系统数据的匮乏，严重制约了其失效行为、规律评价与寿命预测技术的发展。在缺乏14.1 MeV高能强流中子装置进行辐照实验的情况下，聚变堆材料中子辐照行为研究主要聚焦于裂变中子实验模拟、离子辐照实验模拟和多尺度计算模拟[2]。由于嬗变速率、损伤速率、损伤空间分布情况以及He/dpa等条件与聚变条件不同，裂变中子实验无法完全准确预测聚变堆材料的中子辐照行为。因此，构建聚变堆金属材料中子辐照模拟平台是中国聚变工程实验堆(CFETR)工程设计的迫切需要。

美国和欧洲均大力支持核材料中子辐照模拟研究，如美国“Transformational Advances in Computational Multiscale Modeling of Radiation Effects in Materials”项目(2016—2021)、欧洲“Joint Programme on Nuclear Materials, the European Energy Research Alliance”项目(2016—2020)。国内也已开始裂变堆材料中子辐照的多尺度模拟研究，如中国原子能科学研究院开发的裂变中子辐照模拟软件已应用于模拟压力容器钢等结构材料中微结构的演化过程及辐照脆化趋势预测[3]。因此，国内外都已认识到材料中子辐照模拟平台的重要性，开始了这一方面的研究，这些材料中子辐照模拟平台均在建设或进一步完善中。同时，由于中子能谱不同、dpa损伤速率等差异，裂变堆材料中子辐照模拟平台被认为不适合直接应用于聚变堆材料中子辐照条件下结构与性能的预测。

目前国内外对于聚变堆材料的中子辐照模拟均处于起步阶段，尚缺乏系统的研究和完善的模拟平台，无法实现聚变堆材料中子辐照行为从纳米到毫米的跨空间尺度、从皮秒到月/年的全周期过程模拟。美国西北太平洋国家实验室和比利时核能研究中心开展的动力学蒙特卡罗(KMC)以及美国劳伦斯国家实验室开展的团簇动力学(CD)等多尺度模拟研究均依赖于法国里尔大学建立的缺陷信息数据库[4-9]，但该数据库仅包含空位、间隙原子、氦以及它们组成的小团簇，缺乏中等尺寸和大尺寸缺陷团簇以及位错环、晶界等复杂缺陷相关信息，因此无法进行长时间高剂量聚变中子辐照模拟研究。亟需开展第一性原理(FP)和分子动力学(MD)研究建立完善的缺陷性质数据库和级联碰撞数据库，开发准确高效的KMC和CD程序，实现不同模拟尺度间信息的有效传递。

在中子辐照损伤微观结构预测材料宏观力学性能方面，现有研究主要采用弥散障碍模型等经典理论。但由于聚变中子辐照产生的缺陷类型多、尺寸大、密度高，缺陷间相互作用复杂，现有的理论模型无法满足力学性能预测。位错动力学(DD)和晶体塑形有限元(CPFEM)是定量研究辐照缺陷对材料力学性能影响的两个重要方法，这两方面的研究均需要深入。并且国内外目前主要在这两个尺度上各自独立建模分析，两者之间的参数传递很少研究，使得CPFEM对辐照缺陷的硬化效应的预测能力不够，多为解释性模型。现亟需建立辐照微结构模拟程序MD/OKMC/CD与力学性能模拟程序DD/CPFEM及DD与CPFEM间的信息传递和有效衔接，实现从初级辐照损伤到力热性能响应的全链条模拟预测。

综上，在缺乏聚变中子源的情况下，聚变堆材料中子辐照行为尚无法准确预测，这严重制约了聚变堆材料服役行为评价与寿命评估。应用不同尺度模拟方法，针对聚变堆材料中子辐照行为预测的迫切需求，本项目以CFETR拟全面使用的钨材料作为研究对象，构建从微观(级联碰撞)→介观(缺陷结构演化)→宏观(力热性能)的中子辐照模拟平台，通过与典型裂变堆实验进行对比验证平台有效性，实现聚变中子辐照下钨材料的辐照微结构和力热性能的准确预测。

1 研究方案

通过多尺度计算模拟研究，建立钨材料中子辐照缺陷性质数据库和辐照初级损伤结构数据库，发展现有微结构演化和力热性能计算程序，构建金属材料聚变中子辐照计算模拟平台，模拟得到钨材料中子辐照条件下的微结构与力学性能，所得到的结果与现有钨材料裂变中子辐照实验进行对比验证，以确定模拟平台的准确性，籍此预测聚变中子辐照下钨材料的辐照微结构和力热性能。平台命名为NINUM3(Neutron Irradiation of Nuclear Fusion Metallic Materials: Multiscale Modeling & Simulation)。具体研究方案如图1所示。

首先建立关键参数数据库，包括辐照缺陷性质数据库和辐照初级损伤结构数据库。辐照缺陷性质数据库包含空位/自间隙原子/氢氦团簇结合能、扩散/转向/发射指前因子和势垒、作用半径等关键信息；辐照初级损伤结构数据库包含初级碰撞原子能量EPKA≤300 keV(相当于14.1 MeV聚变中子传递给钨原子的最大能量)与不同温度下的缺陷种类、数量与空间分布。这两个数据库是微结构模拟的重要输入。其次，应用OKMC方法开展短时间低dpa的缺陷演化计算，应用CD方法开展长时间高dpa的缺陷演化计算。通过多尺度模型与上述数据库的耦合，实现缺陷演化的跨尺度模拟，获得不同dpa下的缺陷种类、数量与空间分布。基于获得的缺陷微结构信息，应用DD与CPFEM方法计算不同dpa条件下的力学性能与热学性能。将计算得到的缺陷微结构和力学性能结果与裂变中子辐照实验结果进行比对[10-12](表1)，验证中子辐照模拟平台的有效性。在此基础上预测钨材料在聚变中子辐照条件下的微结构和力热性能变化。

图1 金属材料中子辐照平台建立与计算模拟总体方案Fig.1 Overall plan for establishing simulation platform for neutron irradiation in metals

表1 HFIR/JMTR裂变中子辐照条件下钨缺陷结构与力学性能Table 1 Microstructure defect and mechanical property of tungsten irradiated by neutron in HFIR/JMTR

2 研究进展

2.1 中子辐照损伤关键参数数据库建立

1) 辐照缺陷性质数据库

综合利用FP、分子静力学(MS)、MD方法研究了辐照缺陷热力学与动力学行为，获得了描述缺陷及缺陷间交互作用的特征参数。FP因具有精度高、物理意义明确但计算体系小的特点，可用来获得空位/自间隙原子/氢/氦等点缺陷及其构成的小型缺陷团簇的精确特征参数，构建相应缺陷模型。首先通过系统的FP研究发现纳米孔洞的稳定结构可通过最小化其表面积来准确确定，并揭示钨中纳米孔洞形成能与其表面积之间的线性关系(图2)。在此基础上，进一步开发了一种新的物理模型，可准确预测任意大小纳米孔洞的稳定结构和能量。该物理模型与DFT计算结果以及近期的纳米孔洞退火实验结果一致。

图2 钨中V1～V47空位团簇形成能随其Wigner-Seitz表面积的变化Fig.2 Formation energy of V1-V47 nanovoids as function of their Wigner-Seitz areas

MS和MD/AMD方法适用于较大型缺陷团簇、位错环、空洞等的热力学和动力学行为研究，但其结果准确性依赖于所用经验势，因此须根据DFT和基于量子力学的紧束缚势结果对现有钨和钨-氢/氦势函数的有效性进行评估，甄选并优化出适用于辐照损伤模拟的经验势。本项目对比了常用的6种势函数在点缺陷、空位团簇、间隙团簇、位错环和堆垛层错能等静态性质研究中的表现，并评估了这6种势函数在辐照级联下的表现，最终确定了适用于分子动力学辐照损伤模拟的钨经验势函数[13-14]。

聚变中子辐照会产生大量复杂缺陷(如位错环、空洞、氦泡等)，其运动涉及多个多原子协同跃迁过程，故采用加速分子动力学(AMD)和AKMC相结合的方法开展研究。AMD-AKMC结合dimer过渡态搜索方法可在无须事先了解体系跃迁机制下实时给出缺陷在原子尺度和微秒尺度下的演化信息。大型复杂缺陷热力学和动力学性质通常由多个物理参数共同决定，很难通过物理机制分析找出关键物理参量进行建模，因此尝试引入机器学习方法对大型复杂缺陷进行建模。该方法以大型复杂缺陷电子和结构组态信息为输入，热力学和动力学能量数据为输出，利用支持向量机进行多维回归建模，基于所得大量缺陷结构和能量数据进行训练，得到缺陷模型。应用AMD方法模拟研究了金属钨中的空位及空位团簇的动力学演化过程，获得了相应的动力学机制，确认了纳米尺度的位错环之间的旋转转化关系。

2) 辐照初级损伤结构数据库

采用MD方法研究了不同温度、不同入射方向PKA在各种钨体系中的级联损伤行为，获得了辐照初级损伤缺陷类型、数量、分布等关键信息，构建了辐照初级损伤结构数据库。MD模拟辐照级联损伤已有50多年的历史，方法较成熟。但在模拟高能PKA(最高达300 keV)级联损伤时需使用原子数在千万量级以上的超大体系，其中级联碰撞所产生的缺陷种类的判别和分析处理非常困难且效率低。采用“自动步长”法可提高计算效率，即模拟程序根据PKA能量、模型体系中原子受力和速度等自动调节步长。此外，目前级联损伤的模拟只考虑级联反应中具有较大动能原子导致的结构演化，不考虑辐照引起的电子激发效应的影响。本项目引入紧束缚模型来研究电子激发效应。

3) 数据库的构建程序和方法

所使用的软件及代码多在Linux操作系统下使用，为方便调用使用MySQL数据库管理系统。根据MySQL模式要求整合中子辐照的模拟计算和实验数据，建立中子辐照损伤跨尺度集成模拟数据库，并将数据库与模拟分析程序链接在一起，实现数据的有效衔接和调用。MySQL是一个关系型数据库管理系统，具有读取速度快、灵活性高、可管理性强且代码开源无版权限制等特点。

2.2 中子辐照计算模拟平台构建

1) 介观尺度缺陷演化行为的OKMC模拟

OKMC方法忽略晶格原子热振动，将辐照缺陷列为研究的主要对象，可有效扩展模拟的时间和空间尺度。中子辐照产生的初始损伤具有空间分布不均匀特征，且存在大量空间关联性缺陷(如弗伦克尔缺陷)，这些空间信息对微结构的演化至关重要。OKMC固有地保留了这些空间信息。本项目将采用OKMC方法模拟辐照初级损伤结构在介观尺度/微秒时间范围内的演化，消除或弱化空间相关性和不均匀性，为后续大尺度DD和CPFEM等力热性能模拟提供必要参数。计算效率低且难以并行是限制OKMC广泛应用的重要原因。本项目拟通过使用高效的二叉树和列表相结合的数据结构来提高程序运行效率，探索OKMC高效分区并行算法，并尝试GPU硬件加速。

图3 NINUM3-OKMC和KSOME模拟给出的75 keV PKA在300 K下产生的VAC和SIA缺陷存活率随时间的变化Fig.3 Surviving fraction of VAC and SIA as a function of time for a PKA of 75 keV at 300 K obtained by NINUM3-OKMC and KSOME

目前，已开发了具有自主知识产权的OKMC(NINUM3-OKMC)，并与欧美典型程序LAKIMOCA和KSOME进行了对比验证，研究结果吻合较好，证明该程序可用于钨中缺陷演化的准确模拟。图3示出了用不同OKMC程序模拟75 keV单级联中子辐照损伤在300 K下退火过程中空位(VAC)型和间隙(SIA)型缺陷存活率随时间的变化。结果显示，NINUM3-OKMC与KSOME[15]的模拟结果相一致。进一步，使用自主开发的NINUM3-OKMC程序模拟了低剂量中子辐照下的缺陷演化，并与美国橡树岭国家实验室HFIR实验结果进行了对比，发现缺陷数密度与实验结果在同一数量级。在模拟过程中，通过在体系中引入杂质来捕获快速迁移的自间隙原子，降低SIA迁移频率，提升OKMC模拟效率。

2) 缺陷长时间演化行为的CD模拟

辐照微结构演化是有效连接原子尺度模拟和力热性能计算的关键环节。CD方法是模拟微结构长时间演化的有效方法，在辐照损伤研究中得到了广泛应用。该方法考虑了缺陷之间的相互作用，直接求解一系列缺陷演化的ODE方程组。在高剂量聚变中子辐照下，必须考虑大尺寸缺陷(几百nm～μm)的演化行为以及嬗变氦的影响，导致模型的复杂度和方程数量(达百万量级)显著增大，提高了数值求解的难度。将确定性方法、随机方法以及缺陷分群方法有机结合，主要考虑氦-空位团簇、间隙性位错环等缺陷的形核、长大和粗化过程，尤其是对于大尺寸缺陷采用更为精确的分群方法来代替Fokker-Plank方程，实现微结构演化的精确模拟。

基于CD的基本理论，结合钨中辐照缺陷的特点，建立了间隙原子团簇、空位团簇、位错环等辐照缺陷演化的物理模型，尤其是考虑了钨中空位之间结合较弱甚至排斥的情况下如何处理的问题。同时考虑到后续不同He/dpa的计算，推导建立了He-空位团簇的group method公式。利用前期开发的确定性和随机性求解算法与程序(Radieff和MISA_SCD)，模拟计算了不同离子辐照剂量下晶粒尺寸钨中位错环的演化过程(图4)，与实验结果[16]基本一致，初步证明了方法的可行性。通过发展特征时间退火方法，实现了初级辐照损伤MC(IM3D)方法、短时间缺陷演化的OKMC(MMonCa)方法和长时间缺陷演化的CD(IRadMat)方法的耦合，建立了计及级联内缺陷空间关联的团簇动力学模型(CD-SC)，修正了传统CD对长时间缺陷动力学和尺寸分布行为的错误估计，实现了从离散到连续模型的有效过渡。

3) 辐照微结构力热性能计算

根据辐照缺陷性质数据库提供的位错与辐照缺陷交互作用结果，建立位错与辐照缺陷交互作用模型并编程整合到DD软件ParaDis中，研究单晶钨辐照微观结构演化，确定辐照缺陷对滑移系阻碍作用及位错密度演化规律，建立描述钨辐照后应力应变响应的晶体塑性单晶本构。基于有限元软件ABAQUS开发用户材料子程序建立多晶钨模型，根据上述研究所得不同辐照剂量下缺陷状态数据并考虑晶粒尺寸、晶体取向等信息，利用CPFEM模拟多晶钨辐照后的应力应变响应，建立样品尺度有限元预测不同辐照环境下辐照缺陷对材料硬化的影响。基于单/多晶中空洞/气泡分布形貌和傅里叶定律，编写热传导程序，计算不同辐照条件下材料微结构形貌的热导率。

图4 CD程序模拟不同离子辐照剂量下钨中位错环面积Fig.4 Loop area as a function of radiation damage dose (dpa) simulated by CD

本项目通过分析MD模拟的位错与空洞相互作用过程并结合弹性理论，建立了位错动力学框架下的空洞模型，模拟位错通过周期性排列的不同尺寸空洞，得到了与分子动力学模拟一致的通过过程和所需分切应力，验证了空洞模型的可靠性；模拟位错通过随机分布的不同尺寸、数密度空洞，统计其硬化作用，利用多种硬化模型进行了拟合，发现其中BKS模型能较好地描述空洞导致的硬化作用。建立了含声子-声子散射、声子-电子散射、声子-杂质散射、电子-电子散射和电子-杂质散射的热导率物理模型，研究了声子与位错间相互作用、晶粒尺寸对热导率的影响。给出了一定温度范围内未辐照金属钨的热导率，以及热导率随位错密度和晶粒尺寸的变化规律。计算结果与实验值吻合较好。

3 总结与展望

中子辐照是聚变堆材料开发面临的严峻挑战。本项目选择CFETR拟全面使用的钨材料作为研究对象，不需考虑合金化元素影响，因此金属钨中子辐照模拟研究兼具可操作性和工程意义。通过开发和升级不同尺度的计算方法和程序，实现从微观、介观到宏观的有效连接，搭建国内首个聚变堆金属材料中子辐照计算模拟平台。应用FP和MD等方法研究金属钨中本征缺陷、氢/氦缺陷和辐照初级损伤缺陷的热力学和动力学行为，获得描述缺陷性质、缺陷间及其与位错相互作用的关键特征参数，揭示缺陷行为机制及其所形成团簇的作用机制，构建完善的钨中辐照缺陷性质数据库。优化一套新的适用于高能中子辐照模拟的钨和钨-氢/氦的势函数，应用MD和AMD方法系统研究不同钨体系中辐照初始缺陷类型、数量和空间分布，阐明位错、晶界、氢/氦对辐照初始缺陷产生和分布的影响机制；发展高能PKA级联损伤模拟方法和超大体系中的缺陷分析方法，并开发相应的模拟和分析程序；构建完善的反映中子辐照信息的辐照初级损伤结构数据库。开发具有自主知识产权的高效OKMC程序。以MD级联碰撞相关缺陷信息为主要输入参数，搭建从微观到介观的桥梁。通过构建团簇动力学多元体系物理模型，完成对辐照缺陷时空尺度的进一步扩展，分别建立确定性方法、概率性方法以及分群方法求解主方程并发展其相应的高效求解算法，从介观尺度上模拟研究间隙原子、空位、位错环、空洞等辐照缺陷的演化过程。

本项目将通过发展自适应加速MD算法、自适应KMC算法、MD-KMC耦合、CD中引入原子结构信息等新模拟方法与引入机器学习缺陷建模和势函数优化等新模拟技术实现不同尺度方法间的有效连接；通过分析各程序对辐照缺陷产生和演化参数的敏感性，探明决定缺陷行为/力热性能的关键缺陷参数，建立辐照缺陷参数在不同尺度间的逐级提炼和传递模型以及交互反馈机制，实现参数的有效传递。这些关键技术不仅适用于钨基材料中子辐照损伤结构与性能研究，还可扩展应用于其他核材料的中子辐照研究。