蒙特卡罗程序LoongSTARS MCX在快堆屏蔽问题中的计算

Monte Carto Program LoongSTARS MCX in Fast Reactor Shielding Problem Calculations

摘要: 快中子反应堆是我国核电发展三步走战略的关键一环，对于快堆屏蔽问题的计算是快堆设计与发展的重要工作，采用我国自主研发的液态金属冷却快堆设计分析龙码系统中的蒙特卡罗–确定论耦合粒子输运计算软件LoongSTARS MCX，通过模拟选定的SINBAD基准实验，验证LoongSTARS MCX程序在快堆屏蔽问题上的计算能力。整个验证以两种方式进行，首先通过与MCNP程序的计算结果对比，验证LoongSTARS MCX程序的功能性和计算精度。其次，将计算结果与SINBAD实验结果进行对比，验证LoongSTARS MCX程序的有效性和计算准度。结果表明，LoongSTARS MCX程序结果与MCNP程序结果高度吻合，相对偏差在±20%以内，并且与实验值符合较好，偏差在合理可靠范围内，LoongSTARS MCX程序可以用于快堆屏蔽问题上应用研究，具有良好的计算准度与精度。

Abstract: Fast neutron reactor is a key link in the three-step strategy of nuclear power development in our country, and the calculation of shielding problem of fast reactor is an important work for the design and development of a fast reactor, using the Monte Carto-deterministic method coupled particle transport calculation software LoongSTARS MCX in the LONGCODE system for liquid metal cooled fast reactor design analysis developed by our country, we simulated the selected SINBND benchmark experiment, verify the computing power of LoongSTARS MCX program on fast reactor shielding problem. The whole validation was carried out in two ways. First, the functionality and accuracy of the LoongSTARS MCX program were verified by comparing the results with those of the MCNP program. Secondly, compared with the Sinbad experimental results, the validity and accuracy of the LoongSTARS MCX program are verified. The results show that the results of LoongSTARS MCX are in good agreement with those of MCNP, the relative deviation is less than ±20%, and the deviation is in a reasonable and reliable range, the LoongSTARS MCX program can be used to study the shielding problem of fast reactor with good calculation accuracy and precision.

1. 引言

快堆作为第四代核电技术的重要组成部分，是未来核电进一步发展的方向之一，也是我国核能发展三步走战略[1]的关键一步。相较于传统热堆，快堆对于铀资源的利用率是其数十倍，并且能使核废料的产生得到最大程度降低[2]。经过数十年的努力，我国快堆实现了从无到有，从实验堆到示范堆的巨大跨越，但是距离形成我国完全自主化、成熟的第四代商业化快堆技术仍然有一段距离。

快堆技术的高质量发展与突破离不开反应堆数值计算软件的研发与应用。目前主要应用于反应堆模拟计算的方法有确定论方法和蒙特卡罗方法，确定论法是直接采用数值方法来解析粒子输运方程[3]从而获得近似解，比较知名的基于确定论法研发的程序有QAD、ANISN、DORT、TORT等；而相比于确定论方法，蒙特卡罗方法能够处理复杂三维几何问题，物理模型建立精准度高，但是蒙特卡罗法的缺点是收敛速度慢，存在统计误差，因此在解决计算规模较大的问题时，通常需要耗费巨大的计算成本才能得到可靠的结果，因此受限于计算机硬件条件，蒙特卡罗方法过去一直未能广泛应用于工程领域，但随着计算机技术的快速发展，反应堆物理以及安全分析中对于数值仿真计算精度的需求[4]，蒙特卡罗方法因其独特优势在核领域受到了越来越多的关注与研究。

国内外多家研究机构都相继开展了基于蒙特卡罗方法的软件研发工作。国际上有美国LosAlamos国家实验室研发的MCNP [5]程序，能够很好的应用于三维复杂几何结构中的中子、光子、电子或者三者耦合输运问题的模拟计算中，MCNP程序具有灵活性高、通用性强、功能性全面的特点，可进行辐射屏蔽计算、辐射剂量测定、临界系统计算(包括超临界以及次临界)、探测器响应与设计分析计算、反应堆本体物理计算等；法国原子能委员会CEA开发的TRIPOLI-4 [6]程序，可模拟中子、光子、正负电子以及耦合情况下的三维输运问题，TRIPOLI-4程序可用于辐射屏蔽问题计算、临界问题计算、堆芯物理问题计算等，同时支持固定源次临界模式解决辐射屏蔽问题，并且可以实现CAD模型与TRIPOLI相互转换；日本BNFL实验室同AEA Technology实验室合作开发的MCBEND [7]程序，可用于屏蔽以及剂量计算，支持输入输出可视化。国内主要有由中科院FDS核能研究团队基于CAD开发的超级蒙卡核模拟软件系统SuperMC [8]，应用于反应堆物理及辐射屏蔽安全分析等，可进行燃耗计算、活化计算、屏蔽计算、堆芯物理计算等多种问题，同时支持多种格式数据库输入，可进行中子、光子、电子以及耦合中子和光子输运计算，并且支持自动几何建模、可视化、高性能并行计算，内置丰富的减方差方法，使得程序高效化、智能化；清华大学核能所开发的堆用蒙卡分析程序RMC [9]，具备多物理多尺度计算的能力，可进行临界问题本征值、燃耗、瞬态过程等计算，程序中具备几何处理、新燃耗算法、源收敛加速、并行计算以及温度相关截面处理等技术，可提高计算效率；国家核电技术有限公司北京软件技术中心与清华REA团队联合开发的自主化蒙特卡罗粒子输运程序cosRMC [10]，可应用裂变反应堆临界计算，燃耗计算和屏蔽计算，支持几何处理加速、计数器优化、核截面处理优化、并行算法计算等多种技巧，同时包含多种减方差方法提高计算效率，支持建模到结果全过程可视化；北京应用物理与计算数学研究所开发的三维中子–光子耦合输运蒙特卡罗模拟软件JMCT [11]，该软件具备CAD可视化输入输出界面，支持连续和多群能量模式，考虑了包括热化在内的各种核反应，可精细计算反应堆全堆芯pin功率及时空分布，能够模拟固定源、临界本征值及伴随输运问题。

快堆屏蔽计算是快堆设计、安全分析以及辐射防护中一项非常重要的工作，但是基于快堆进行全堆试验从而获得满足实际测量需求具有很强的局限性，因此通常采用计算机模拟计算的方式来获得较为满意的测量需求。本文基于我国自主研发的液态金属冷却快堆设计分析龙码系统中蒙卡源项屏蔽分析软件LoongSTARS MCX，选取了三例OECD国际基准例题TUD SiC、FNG BLKT以及Winfrith Water，使用LoongSTARS MCX进行模拟计算，将计算结果与MCNP计算结果以及实验测量值进行对比分析，以验证LoongSTARS MCX对于快堆屏蔽问题的计算能力，为龙码软件系统应用于快堆屏蔽计算提供验证支持。

2. 蒙卡程序LoongSTARS MCX

LoongSTARS MCX是西安交通大学自主研发的蒙特卡罗粒子输运模拟软件，该软件目前具备中子–光子耦合输运计算、输运–燃耗耦合计算、光子点核积分计算、几何建模可视化、计算结果可视化、生成NECP-Hydra输入、大规模并行等功能。对于屏蔽深穿透问题，可以采用蒙特卡罗–确定论耦合方法，提高计算效率。软件可应用于裂变堆堆芯临界计算，各类辐射屏蔽计算及聚变堆包层屏蔽和增殖计算。与MCNP等蒙特卡罗粒子输运程序类似，LoongSTARS MCX软件具有输入编写方式简单，功能更为全面等优点，同时LoongSTARS MCX支持支持多种输入方式，其中包括：1) 直接在XML格式的输入卡片中定义求解问题的几何、材料等信息，作为LoongSTARS MCX的输入；2) 在XML格式的输入卡片中，指定已有的MCNP输入卡片的路径和名称，LoongSTARS MCX可直接读取已有的MCNP输入卡片中的几何、材料等信息，进行问题的求解；3) 将CAD模型导入SALOME可视化建模平台，基于该平台可自动生成MCNP输入卡片，之后通过方式2即可实现可视化三维模型的自动建模；4) 用户可通过网页客户端设置求解问题的参数，并可以在远程后台上自动生成方式1的XML格式的输入卡片，可以方便用户多选择使用。

3. SINBAD屏蔽基准题

SINBAD屏蔽基准题库[12]建立的主要目的是构建和维护一套国际通用基准实验标准库，以供屏蔽设计用户使用，包括进行数据验证和计算机程序验证，该数据库目前包括102个屏蔽基准，涵盖的材料包括：空气、氧、水、铝、铍、铜、石墨、混凝土、铁、铅、锂、镍、铌、碳化硅、钠、不锈钢、钨、钒及其混合物，基准实验可分为三类：裂变反应堆屏蔽(48个基准)，聚变包层中子(31个)和加速器屏蔽(23个)。该数据库在许多国家和国际项目范围内被广泛使用，如压水堆压力容器监测、核聚变计划(ITER反应堆研究)、核数据验证、原子能机构核数据项目等，同时SINBAD基准题被广泛用于反应堆数值模拟计算程序验证，本文选取三例快中子源基准题进行程序验证。

3.1. TUD SiC

TUD基准实验[13]由德累斯顿工业大学(TUD)团队构建，涉及的屏蔽材料包括铁、碳化硅和钨。本文选取了碳化硅基准实验进行了程序模拟验证，SiC屏蔽实验模型如图1 [14]所示，实验利用能量为125 KeV的氘核撞击加速器中的T-Ti靶得到14 MeV D-T中子源，测量伴随α粒子得到D-T中子源的强度，源的时间分布与exp[−(t/1.4 ns)2]成正比，D-T中子源位于组件前5.3 cm处，SiC组件尺寸为45.7 × 45.7 × 71.1 cm³，探测器采用高度和直径为3.8 cm圆柱形NE213闪烁光谱仪，该材料的质量密度为0.874 g/cm³，元素组成为54.8% H和45.2% C。探测器位于P1、P2、P3、P4 (a = 12.70, b = 27.94, c = 43.18, d = 58.42 cm)四个位置，分别使用LoongSTARS MCX以及MCNP测量了四处位置E > 1 MeV的中子能谱并与实验值进行了对比。

Figure 1. TUD SIC calculation model

图1. TUD SiC计算模型

3.2. FNG BLKT

该基准实验为了验证ITER屏蔽系统设计的正确性，模拟了ITER内部屏蔽结构并进行了一系列的中子学实验[15]，模型的结构和材料组成与实际结构和材料组成相同，计算几何结构如图2所示，它由三部分组成：第一部分(屏蔽块)是由有机水玻璃板、AISI-316不锈钢和1 cm厚的铜层组成，横截面长宽均为100 cm，源与组件的第一层铜距离为5.3 cm，第二部分是由厚度均为2.2 cm的AISI-316和Cu交替构成，此部厚度为30.8 cm，截面为47 cm × 47 cm。第三部分是由聚乙烯制成的屏蔽物，以防止背景辐射，该屏蔽层部分地围绕第一部分并完全包围第二部分，它的外形尺寸为127 cm × 127 cm × 82 cm。利用活化箔技术测量了中子积分通量随穿透深度的函数，所有使用的箔片直径为18 mm，厚度从1到3 mm不等，具体取决于位置深度。

3.3. Winfrith Water

基准实验示意结构如图3 [16]，由一个轻型结构支撑一个垂直铝管组成，支撑结构位于一个横截面为228 cm × 177 cm，高度为172 cm的充水铝缸内，中子源悬挂在支撑结构的臂上，一个或多个(最多为

Figure 2. FNG BLKT computational mode

图2. FNG BLKT计算模型

8个)中子源等距对称分布在测量管周围，实验中子源由²⁵²Cf瞬发裂变源提供。实验测量了³²S (n, p)³²P探测器的反应率，测量位置位于铝管源平面以及在源平面上下15 cm和30 cm处，硫探测器直径和高均为28 mm，实验同时需要使用NE213液体闪烁体探测器测量1 MeV以上快中子能谱，测量位于源平面，中子源与闪烁体探测器距离分别为10.16、15.24、20.32、25.40、30.48、35.56 cm。

Figure 3. The Winfrith Water Schematic of the benchmark experiment

图3. Winfrith Water基准实验示意图

4. 计算结果与分析

图4为TUD SiC基准题实验与计算中子能谱，LoongSTARS MCX计算结果与MCNP程序计算结果相比，P1处两种程序相对偏差在±1.9%以内，P2处相对偏差在±2.7%以内，P3处相对偏差在±8.2%以内，P4处相对偏差在±13%以内，整体而言四处位置计算结果平均相对偏差普遍在±4%以内，最大相对偏差绝对值不超过13%，表明LoongSTARS MCX计算精度与主流、成熟蒙卡程序MCNP相当；而对于计算结果与实验结果，当探测器处于P1位置时，计算结果远大于实验结果，这种偏差较大的现象在参考文献[17] MCS程序模拟计算中同样出现，但与基准题库给出的程序结果相同，另外三个位置的计算结果与实验结果基本一致。

(a) P1位置 (b) P2位置

(c) P3位置 (d) P4位置

Figure 4. TUD SiC benchmark experiment and calculation of neutron energy spectrum

图4. TUD SiC基准题实验与计算中子能谱

表1~5为五种类型活化箔探测器测得的反应率计算值与实验值，由于计算模型属于深穿透问题，从源到探测器位置中子需要经历数十个自由程，中子通量密度梯度大，因此采用粒子重要性减方差方法将中子输运到深穿透部位。在计算中，活化箔尺寸按照实际大小给出，所有剂量学反应截面数据均取自IRDFF-v1.05文件。对于⁵⁶Fe(n, p)⁵⁶Mn，27^Al(n, a)²⁴Na，⁹³Nb(n, 2n)^92mNb以及⁵⁸Ni(n, 2n)⁵⁷Ni探测器LoongSTARS MCX计算结果，在穿透深度较浅时，计算值与实验值相对偏差普遍在±10%以内，在深穿透位置，计算值与实验值相对偏差明显增大，最大相对偏差绝对值接近40%，这种现象同样出现在参考文献中[18]，对于¹¹⁵In(n, n’)^115mIn探测器，计算值与实验值相对偏差均在±10%以上，在深穿透处，最大偏差可达54%，参考文献中也证实了这一现象，但并未给出原因，整体来说，5种探测器计算结果与实验结果普遍吻合较好，在深穿透处偏差较大，由参考文献使用MCNP以及FENDL-1和EFF-3数据库模拟结果来看，这种偏差在合理范围内，整体来看五种探测器反应率计算值与实验值较为接近，并且LoongSTARS MCX与MCNP计算结果也符合的较好，计算结果平均相对偏差绝对值在3%以内，偏差较大处可能与MCNP采用的核数据库为ENDE/B-VII.1而LoongSTARS MCX采用的是Loong ATLAS以ENDF/B-VIII.0库制作的数据库有关。

图5为WinfrithWater基准题中子能谱图，根据文献中所给出的不同的源项配置方案，在相应探测器测量点位置进行中子通量计数。可以看出MCNP和LoongSTARS MCX的计算结果与实验测量值在探测器距离源较近时符合较好，在距离较远时，在低能区计算结果与实验值偏差较大，整体来看，在1.5至

Table 1. ⁵⁶Fe(n, p)⁵⁶Mn reaction rates of experimental and calculated values

表1. ⁵⁶Fe(n, p)⁵⁶Mn实验值与计算值反应率

深度 Depth/cm	实验值 Experiment	MCNP	C/E	LoongSTARS MCX	C/E
3.43	8.47E−05	8.70E−05	1.03	8.74E−05	1.03
10.32	1.41E−05	1.44E−05	1.02	1.42E−05	1.01
17.15	3.51E−06	3.52E−06	1.00	3.44E−06	0.98
23.95	1.03E−06	1.02E−06	0.99	9.96E−07	0.97
30.80	3.31E−07	3.17E−07	0.96	3.14E−07	0.95
41.85	6.61E−08	6.37E−08	0.96	6.54E−08	0.99
46.85	2.80E−08	2.67E−08	0.95	2.71E−08	0.97
53.80	1.00E−08	9.06E−09	0.91	8.54E−09	0.85
60.55	3.67E−09	3.26E−09	0.89	3.26E−09	0.89
67.40	1.30E−09	1.17E−09	0.90	1.55E−09	1.19
74.40	4.50E−10	4.17E−10	0.93	4.11E−10	0.91

Table 2. ²⁷Al(n, a)²⁴Na reaction rates of experimental and calculated values

表2. ²⁷Al(n, a)²⁴Na实验值与计算值反应率表

深度 Depth/cm	实验值 Experiment	MCNP	C/E	LoongSTARS MCX	C/E
3.43	8.50E−05	8.88E−05	1.04	8.86E−05	1.04
10.32	1.47E−05	1.46E−05	0.99	1.45E−05	0.99
17.15	3.60E−06	3.58E−06	0.99	3.53E−06	0.98
23.95	1.07E−06	1.04E−06	0.97	1.04E−06	0.97
30.80	3.44E−07	3.26E−07	0.95	3.18E−07	0.92
41.85	7.06E−08	6.76E−08	0.96	6.60E−08	0.93
46.85	2.94E−08	2.83E−08	0.96	2.72E−08	0.93
53.80	1.09E−08	1.01E−08	0.93	9.22E−09	0.85
60.55	3.71E−09	3.49E−09	0.94	3.43E−09	0.92
67.40	4.72E−10	5.07E−10	1.07	4.75E−10	1.01
74.40	8.50E−05	8.88E−05	1.04	8.86E−05	1.04

Table 3. ⁹³Nb(n, 2n)^92mNb reaction rates of experimental and calculated values

表3. ⁹³Nb(n, 2n)^92mNb实验值与计算值反应率表

深度 Depth/cm	实验值 Experiment	MCNP	C/E	LoongSTARS MCX	C/E
3.43	3.33E−04	3.68E−04	1.11	3.58E−04	1.08
10.32	5.48E−05	5.80E−05	1.06	5.60E−05	1.02
17.15	1.34E−05	1.37E−05	1.02	1.32E−05	0.99
23.95	3.80E−06	3.84E−06	1.01	3.79E−06	1.00
30.80	1.21E−06	1.17E−06	0.97	1.14E−06	0.94
41.85	2.69E−07	2.25E−07	0.84	2.28E−07	0.85
46.85	1.04E−07	9.45E−08	0.91	9.35E−08	0.90
53.80	3.64E−08	3.17E−08	0.87	3.15E−08	0.87
60.55	1.21E−08	1.13E−08	0.93	1.15E−08	0.95
67.40	4.51E−09	3.95E−09	0.88	4.13E−09	0.92
74.40	1.44E−09	1.42E−09	0.99	1.61E−09	1.12
81.10	5.10E−10	5.06E−10	0.99	5.25E−10	1.03
87.75	2.27E−10	1.90E−10	0.84	2.46E−10	1.08
92.15	1.44E−10	8.81E−11	0.61	1.04E−10	0.72

Table 4. ¹¹⁵In(n, n^’)^115mIn reaction rates of experimental and calculated values

表4. ¹¹⁵In(n, n^’)^115mIn实验值与计算值反应率表

深度 Depth/cm	实验值 Experiment	MCNP	C/E	LoongSTARS MCX	C/E
3.43	2.15E−04	1.91E−04	0.89	1.91E−04	0.89
10.32	5.97E−05	5.25E−05	0.88	5.24E−05	0.88
17.15	1.89E−05	1.66E−05	0.88	1.65E−05	0.87
23.95	6.59E−06	5.56E−06	0.84	5.55E−06	0.84
30.80	2.30E−06	1.90E−06	0.83	1.86E−06	0.81
41.85	6.06E−07	4.02E−07	0.66	3.96E−07	0.65
46.85	2.30E−07	1.99E−07	0.87	1.91E−07	0.83
53.80	8.87E−08	6.88E−08	0.78	6.66E−08	0.75
60.55	3.34E−08	2.50E−08	0.75	2.36E−08	0.71
67.40	1.61E−08	9.35E−09	0.58	8.76E−09	0.54
74.40	7.00E−09	3.49E−09	0.50	3.28E−09	0.47
81.10	2.51E−09	1.40E−09	0.56	1.23E−09	0.49

Table 5. ⁵⁸Ni(n, 2n)⁵⁷Ni reaction rates of experimental and calculated values

表5. ⁵⁸Ni(n, 2n)⁵⁷Ni实验值与计算值反应率表

深度 Depth/cm	实验值 Experiment	MCNP	C/E	LoongSTARS MCX	C/E
3.43	2.87E−05	2.90E−05	1.01	2.89E−05	1.01
10.32	4.15E−06	4.26E−06	1.03	4.22E−06	1.02
17.15	9.73E−07	9.57E−07	0.98	9.44E−07	0.97
23.95	2.57E−06	2.58E−06	1.00	2.55E−06	0.99
30.80	8.18E−06	7.63E−06	0.93	7.37E−06	0.90

8.8 MeV的能量范围内，计算值与实验值符合较好，在其余区间由参考文献[19]给出的计算结果知，这种偏差同样出现在程序模拟计算结果上，并且本论文模拟结果相比于参考文献整体更接近实验值。同时将LoongSTARS MCX与MCNP计算结果进行对比，在探测器与源距离25.4 cm之前，不同程序计算结果相对偏差在±10%以内，在之后最大相对偏差绝对值接近20%。

图6是不同源项配置方案下硫探测器的反应率测量结果，横坐标为铝量管内探测器与源所在水平面的距离。总体而言模拟计算结果与实验值符合较好，计算值与实验值比值接近于1，而不同程序计算结果平均相对偏差绝对值为6.5%，最大相对偏差绝对值为16.4%。

(a) R = 10.16 cm (b) R = 15.24 cm

(c) R = 20.32 cm (d) R = 25.40 cm

(e) R = 30.48 cm (f) R = 35.56 cm

Figure 5. Winfrith Water Benchmark experiment and calculation of neutron energy spectrum

图5. Winfrith Water基准题实验与计算中子能谱

(a) R = 10.16 cm

(b) R = 15.24 cm (c) R = 25.40 cm

(d) R = 30.48 cm (e) R = 35.56 cm

Figure 6. Winfrith Water reference ³²S(n, p) ³²P experiment and calculation of reaction rate

图6. Winfrith Water基准题³²S (n, p) ³²P实验与计算反应率

5. 结论

本文使用LoongSTARS MCX程序对TUD SiC、FNG BLKT以及Winfrith Water三例具备快中子屏蔽问题特点的基准题进行了模拟计算，将计算结果与参考程序MCNP计算结果进行了对比，TUD SiC基准题四处位置中子能谱相对偏差普遍在±10%以内，最大偏差绝对值不超过18%；FNG BLKT基准题五种类型探测器反应率计算结果与MCNP计算结果平均相对偏差绝对值在3%以内，最大偏差绝对值为13.4%；对于Winfrith Water基准题，6种源分布方式中子能谱计算结果与MCNP结果相比，相对偏差在±20%以内，硫探测器计算结果平均相对偏差绝对值为6.5%，最大偏差绝对值为16.4%。以上程序结果对比相对偏差均满足在±20%允许范围以内，初步验证了LoongSTARS MCX程序可以用于快堆屏蔽问题计算，计算精度与主流、成熟的蒙卡程序MCNP相当。

将LoongSTARS MCX计算结果与实验结果进行比较，确实存在一定的偏差。但经过与参考文献给出的计算结果以及基准题报告给出的测试结果对比，LoongSTARS MCX计算结果在合理、可靠的范围内。结果表明，在多个实验模型中，LoongSTARS MCX程序的计算结果是可靠的，能够满足快堆屏蔽问题中计算有效性与准度的要求。

基金项目

本研究得到了国家重点研发项目(2022YFB1902700)、国家教育部装备预研联合基金(8091B042203)、国家自然科学基金(11875129)、国家强脉冲辐射模拟与效应重点实验室基金(SKLIPR1810)、辐射应用创新中心基金(KFZC2020020402)、北京大学核物理与技术国家重点实验室基金(NPT2023KFY06)、中国铀业有限责任公司与华东理工大学核资源与环境国家重点实验室联合创新基金(2022NRE-LH-02)、中央高校基础研究基金(2023JG001)等支持。

NOTES

^*通讯作者。