锰浴法测量中子源发射率中各修正因子的蒙特卡罗计算

李春娟，刘毅娜，张伟华，王志强

(中国原子能科学研究院 计量测试部，北京 102413)

锰浴法测量中子源强度用于标定各种便携式中子源的发射率，是国际上普遍采用的绝对测量方法，有较高的准确度。它是将待测中子源置于体积很大的装有MnSO4溶液的容器中心，中子在水中充分慢化后被俘获，溶液中的55Mn吸收热中子后，生成放射性核素56Mn，通过56Mn→56Fe+β—+γ发生衰变，半衰期为(2.579±0.003) h。当中子源在锰浴中放置足够长时间后，锰池中的中子数将达平衡状态，通过测量此时锰液中56Mn发射出的γ射线，即可得出中子源发射率。

根据锰浴法测量中子源发射率的原理，需对由中子源发出的未被MnSO4溶液中55Mn俘获的中子组分进行修正，需考虑的修正项有：1) 锰池边界泄漏中子修正L；2) 溶液中氧和硫通过核反应引起的快中子损失修正O，简称氧硫吸收；3) 中子源及源容器对从溶液中反射来的热中子吸收修正S，简称源的自吸收；4) 被Mn俘获的热中子占被整个溶液俘获热中子的比例，即f。通常，只有不足50%的源中子会被溶液中的55Mn俘获，因此，修正因子准确与否直接影响测量结果的可靠性，而由于传统公式的局限性，不能详细考虑输运过程及各修正因子之间的相互影响，使得总修正因子的计算较困难，蒙特卡罗(MC)模拟方法在这些计算中发挥了至关重要的作用[1-4]，目前国际普遍采用MC模拟的方法[1]。

中国原子能科学研究院计量测试部(国防科技工业电离辐射一级计量站)建立了锰浴法测量中子源强度标准装置，并参与了由国际计量局(BIPM)组织的中子源强度国际比对[5]，比对中子源为241Am-Be(α,n)源，其为ISO标准源。本工作利用大型通用MC粒子输运程序MCNP[6]模拟计算测量过程中的各修正因子，对其不确定度进行详细评定，以提供一种用于MC模拟计算结果不确定度评定的方法。

1 计算模型的建立

1.1 几何模型的建立

锰浴法测量中子源强度标准装置如图1所示，系统主要由锰池、源容器、循环泵、流量计、测量室、两个NaI(Tl)探测器及相应的电子学系统组成。测量时将中子源置于锰浴中心36 h，待锰池中的锰液达饱和放射性后，保持硫酸锰溶液在锰池和屏蔽样品容器内连续循环，对56Mn溶液的饱和活化计数A进行测量。将已知活度(绝对测量得到)的56Mn溶液倒入锰池中，通过两路NaI(Tl)探测器测量锰液中56Mn的γ计数，确定锰浴探测效率ε[7]。

图1 中子源强度标准装置示意图

图2 模拟计算的几何模型

通过MC粒子输运程序MCNP对测量过程中的各修正因子进行模拟计算。对修正因子的计算有直接影响的结构主要包括不锈钢锰浴外壳、硫酸锰溶液、241Am-Be(α,n) ISO标准源(X3封装外壳)、Al桶及不锈钢容器，据此建立MC模拟计算几何模型(图2)。不锈钢锰浴外壳直径为100 cm、厚度为4 mm，通过烘干法测量硫酸锰溶液的浓度，计算得到溶液的密度为1.243 9 g/cm3，氢锰原子数比为58.08。241Am-Be(α,n)源为国际电离辐射咨询委员会组织的国际比对CCRI(Ⅲ)-K9.AmBe的传递标准源，源活性区为圆柱形，外径为17.6 mm，高度为16.3 mm，活性区内AmO2的质量分数为6.29%、Be的质量分数为93.71%，活性区密度为1.24 g/cm3。源采用双层不锈钢封装，封装外壳型号为X3，外观直径为22.4 mm，高度为31 mm。源的承托物包括铝桶及不锈钢容器，铝桶外径为33 mm，厚度为1.25 mm，不锈钢容器外径为38 mm，厚度为1 mm。

1.2 源的描述及中子反应截面数据库的选取

实验用241Am-Be(α,n)源的活性区内密度、发射α的氧化镅颗粒大小和封装外壳均满足ISO标准要求，故其封装外壳外的中子能谱与ISO 8529[8]中给出241Am-Be(α,n)源中子能谱一致。国际推荐将比对源置于球形容器中，并对容器抽真空，然后将其置于锰浴中心进行测量。此情况下计算时可将源近似描述为从球形容器内表面发射的各向同性球面源，源中子能谱可描述为ISO 8529中给出的241Am-Be(α,n)源中子能谱。

但本锰浴系统源的承托物为铝桶及不锈钢管容器，而源则置于承托物的底部，此结构会改变源中子入射的边界条件。为使计算更接近真实情况，确保计算结果的可靠性，计算中将源项描述为从源活性区发射的均匀分布的圆柱形体源。封装外壳内从源活性区发射的中子能谱是未知的，利用MCNP通过多次迭代的计算方法得到封装外壳内的裸体芯区中子能谱。具体计算过程如下。

首先利用ISO 8529给出的能谱作为裸体芯区体源的初始能谱，并采用相同的能群结构记录封装外壳外的能谱，计算出其与ISO 8529能谱的差别，将其叠加在初始能谱上作为裸体芯区体源源中子能谱进行重新计算，再次记录封装外壳外的能谱并与ISO 8529能谱进行比较，如此迭代计算，直到包壳外的能谱与ISO 8529能谱非常接近，此时的源能谱即可作为最终的源活性区体源中子能谱。修正过程中各次计算输入能谱的比较如图3所示，各次输出能谱与ISO 8529能谱的比较如图4所示。

图3 各次计算输入能谱的比较

图4 各次输出能谱与ISO 8529能谱的比较

另外，由于中子与源活性区材料及其不锈钢封装外壳的(n,2n)反应，计算得到源材料发射1个中子，封装外壳外的中子数约为1.016 6个，故将源中子的权重设为0.983 64，以使计算结果对源封装外壳外的1个源中子进行归一。

MCNP对中子输运过程的MC模拟基于中子与物质相互作用的反应截面，本计算中使用的中子反应截面数据均选自于ENDF/B-7[9]评价库，使用程序NJOY99[9]将其作为MCNP所需数据库的格式。

2 修正因子计算方法

为得到各修正因子并进行自洽性验证，模拟计算主要记录以下各量：锰池边界泄漏中子份额L，源活性区俘获中子份额S1，不锈钢封装外壳俘获中子份额S2，不锈钢源容器俘获中子份额S3，Al桶俘获中子份额S4，硫酸锰溶液俘获中子份额Nabs，硫酸锰溶液中55Mn通过56Mn(n,γ)反应所俘获的热中子份额f1，硫酸锰溶液中S通过S(n, γ)反应所俘获的热中子份额f2，硫酸锰溶液中S通过S(n, p)、S(n, α)反应所俘获的快中子份额O1、O2，硫酸锰溶液中O通过O(n,γ)反应所俘获的热中子份额f3，硫酸锰溶液中O通过O(n, α)反应所俘获的快中子份额O3，硫酸锰溶液中H通过H(n,γ)反应所俘获的热中子份额f4。其中，利用MCNP中的F1计算锰池边界泄漏中子，利用栅元注量计数器F4及乘子卡FM计算得到上述各核素通过核反应俘获的中子。硫酸锰溶液俘获中子份额计算方法如下。

设硫酸锰溶液中55Mn、32S、16O、1H的原子核数分别为n1、n2、n3、n4，它们与能量为E的中子的吸收截面分别为σ1abs(E)、σ2abs(E)、σ3abs(E)、σ4abs(E)，中子注量率谱为φ(E)，则硫酸锰溶液俘获中子总数Nabs即可通过式(1)得到。利用MCNP中的栅元注量计数器F4可记录硫酸锰溶液中的平均中子注量率谱，乘子卡FM可通过引用各核素的原子核数、中子与各核素相互作用的反应截面及由F4记录得到的中子注量率谱，得到硫酸锰溶液俘获中子总数。由于该计算对源封装外壳外的1个源中子进行归一，故模拟计算得到的硫酸锰溶液俘获中子总数等于硫酸锰溶液俘获中子份额。

n3σ3abs(E)+n4σ4abs(E))φ(E)dE

(1)

其他各量的计算方法与硫酸锰溶液俘获中子份额计算方法类似。通过以上各量即可得到各修正因子，其中锰池边界泄漏中子修正因子L可直接计算得到，源的自吸收修正因子S=S1+S2+S3+S4，55Mn俘获的热中子占被整个溶液俘获的热中子的比例f=f1/(f1+f2+f3+f4)，氧硫吸收修正因子O=O1+O2+O3。得到各修正因子后，通过式(2)，即可得到中子源发射率Q。

(2)

式中：A为测量得到的锰液饱和放射性活化计数；ε为锰浴系统的探测效率。

另外，若不考虑由中子源发出的未被MnSO4溶液中55Mn俘获的中子组分，而直接考虑由中子源发出的最终被硫酸锰溶液中55Mn通过56Mn(n, γ)反应俘获的中子份额F，也可直接由式(3)计算得到中子源发射率。由式(3)和式(1)计算得到的结果应该是一致的，F可通过MC模拟计算得到，即上述记录量中的f1。这种计算方法较简单，也可简化不确定度评价过程。但为对测量过程中各影响量有较详细的了解，并验证计算模型的可靠性，有必要对由中子源发出的未被MnSO4溶液中55Mn俘获的各中子组分进行计算。

Q=A/Fε

(3)

3 计算结果及实验验证

对10 000 000个粒子进行跟踪模拟，计算结果列于表1、2。从表1、2可看出：计算结果是自洽的，即直接计算得到的硫酸锰溶液俘获中子份额与分别计算得到的溶液中Mn、S、O、H所俘获的中子份额总和在其统计不确定度范围内一致，各材料的俘获中子份额及泄漏中子份额总和近似为1，F与f(1-L-O-S)的计算结果在0.6%内一致；锰池边界泄漏中子份额不足1.5%，可使源中子得到充分慢化并被锰浴吸收。

表1 各材料对中子吸收份额的MC模拟计算结果

表2 各修正因子的MC模拟计算结果

根据式(1)，将修正因子应用于溶液饱和放射性活度实验测量结果，得到待测241Am-Be(α,n)中子源的发射率，与参加比对的各机构所提供的结果进行比较，结果列于表3。比对结果表明，除法国LEN-LNHB外，本结果与各机构的实验室的结果在其不确定度范围内一致，验证了本计算模型的可靠性。

表3 各机构241Am-Be(α,n)中子源发射率对比[5]

4 模拟计算结果不确定度分析

不确定度评定是中子计量工作的重要组成部分，其表征了测量结果的分散性及可靠程度。根据式(1)及不确定度传播规律，中子源发射率Q测量结果的不确定度来源主要有56Mn溶液的饱和放射性活化计数A的测量不确定度μ(A)、锰浴系统探测效率的不确定度μ(ε)及各修正量的不确定度。这里主要对各修正量计算结果的不确定度及其对中子源发射率Q测量结果的不确定度的贡献进行分析。

MCNP程序计算各修正量的不确定度即锰池边界泄漏中子修正的不确定度μ(L)、溶液中氧和硫通过带电粒子反应引起的快中子损失修正的不确定度μ(O)、中子源及源容器对从溶液中反射回来的热中子吸收修正的不确定度μ(S)以及由锰所俘获的热中子份额f的不确定度μ(f)。根据不确定度的传播规律，各修正量的不确定度对中子源发射率不确定度的贡献可分别用式(4)～(7)计算得到。

(4)

(5)

μf(Q)=μ(f)

(6)

(7)

各修正量的不确定度的来源主要有统计不确定度、反应截面的不确定度、MnSO4溶液密度的不确定度、源和源的承托物材料密度的不确定度及锰浴半径的不确定度5个方面。首先评估得到各不确定来源引起的各修正量的标准偏差，认为其均满足矩型分布，得到各不确定度分量，然后通过不确定度的合成法则即可得到各修正量的不确定度，评估结果列于表4。各不确定来源引起的各修正量的标准偏差的评估方法分别如下。

表4 MCNP程序计算各修正量的不确定度

4.1 统计不确定度

统计抽样的好坏直接影响到计算结果的准确度，MCNP程序对其计算结果进行统计检查并给出相应统计标准偏差，结果在输出文件中给出。可通过增加计算粒子数和使用方差减小技术的方法减少统计标准偏差。这里取粒子数为10 000 000进行统计抽样。

4.2 反应截面的不确定度

MCNP程序计算是基于各反应截面库的，所用截面数据的不确定度直接影响计算结果。这里所使用的截面数据均选用美国的评价核数据库ENDF/B-7[9]的评价数据。通过改变不同截面数据库的办法评定反应截面的改变对各修正量计算结果的不确定度带来的贡献。计算发现，截面数据库的改变引起的L、O、S、f的最大改变分别为4.68%、21.8%、3.12%、0.04%。由于这里选用的截面数据是较新的数据，最终选定由截面带来的各修正量的标准偏差分别为4.0%、16%、3.0%、0.04%。

4.3 MnSO4溶液密度的不确定度

MnSO4溶液的浓度的不确定度为0.5%，可计算得到其密度的不确定度，通过MCNP程序的微扰计算功能(PERT卡)得到MnSO4溶液的密度在其不确定度范围内的改变对各修正量计算结果的影响，从而估计由MnSO4溶液密度的不确定度而引起的各修正量的标准偏差。

4.4 源和源的承托物材料密度的不确定度

源和源的承托物包括源活性区材料、不锈钢封装外壳、不锈钢源容器以及Al桶，估计各种材料的不确定度，通过MCNP程序计算得到源和源的承托物的密度在其不确定度范围内的改变对各修正量计算结果的影响，从而估计由源和源的承托物材料密度的不确定度引起的各修正量的标准偏差。

由于计算中源中子是从源封装外壳内发射出的体源，故改变源活性区材料及不锈钢封装外壳的密度均会改变封装外壳外的中子能谱，使其偏离ISO 8529所给出的能谱，故需在源活性区材料及不锈钢封装外壳的密度改变的情况下，重新对源中子能谱进行修正，而后得到各修正量，估计由源和源的承托物密度的不确定度而引起的各修正量的标准偏差。计算结果表明，源和源的承托物密度的改变主要对源的自吸收因子影响较大，对其他修正量的影响可忽略。

4.5 锰浴半径的不确定度

将锰浴的半径由50 cm改为50.5 cm，计算其对各修正量计算结果的影响，从而估计由锰浴半径的不确定度而引起的各修正量的标准偏差。

5 结论

利用MCNP程序对锰浴法测量中子源发射率中各修正因子进行计算是可行的，并能弥补传统计算方法的不足。本工作建立的计算方法和模型是合理的，通过多次迭代计算的方法得到合理的源活性区中子能谱有效地解决了计算中的一个难题，可将此计算方法推广应用于其他相关问题的计算。另外，本文也为MC模拟计算不确定度的评定提供了一较合理的方法。

参考文献：

[1] 许淑艳. 蒙特卡罗方法在实验核物理中的应用[M]. 北京：原子能出版社，1996.

[2] 武祯，李君利. 用于探测器校正因子计算的Monte Carlo方法[J]. 清华大学学报：自然科学版，2006，46(9)：1 585-1 588.

WU Zhen, LI Junli. Monte Carlo method for calculating particle radiation detector correction factors[J]. J Tsinghua Univ: Sci & Tech, 2006, 46(9): 1 585-1 588(in Chinese).

[3] 李春娟，陈军，吴海成，等.238U裂变电离室修正因子及探测响应的蒙卡模拟计算[J]. 原子能科学技术，2009，43(3)：208-214.

LI Chunjuan, CHEN Jun, WU Haicheng, et al. Monte Carlo simulation for the correction factors and detection response of238U fission chamber[J]. Atomic Energy Science and Technology, 2009, 43(3): 208-214(in Chinese).

[4] KHABAZ R. Appraisement of the correction factors for neutron reaction in the manganese bath using Monte Carlo calculation[J]. Radioanal Nucl Chem, 2012, 293: 455-462.

[5] ROBERTS N J, JONES L N, WANG Z, et al. International key comparison of measurements of neutron source emission rate (1999—2005): CCRI(Ⅲ)-K9.AmBe[J]. Metrologia, 2011, 48(Suppl.): 06018.

[6] BRIESMEISTER E. MCNP: A general Monte CarloN-particle transport code, Version 4C[R]. US: Los Alamos National Laboratory, 2000.

[7] 刘毅娜，王志强，李春娟，等. 锰浴法绝对测量中子源发射率[J]. 原子能科学技术，2013，47(6)：1 044-1 047.

LIU Yina, WANG Zhiqiang, LI Chunjuan, et al. Absolute measurement of neutron source emission rate with manganese bath method[J]. Atomic Energy Science and Technology, 2013, 47(6): 1 044-1 047(in Chinese).

[8] ISO 8529 Neutron reference radiations for calibrating neutron-measuring devices used for radiation protection purposes and for determining their response as a function of neutron energy[S]. Geneva, Switzerland: International Organization for Standardization, 1989.

[9] CHADWICKM B, OBLOZINSKY P, HERMAN M, et al. ENDF/B-Ⅶ.0: Next generation evaluated nuclear data library for nuclear science and technology[J]. Nuclear Data Sheets, 2006, 107(12): 2 931-3 060.

[10] MacFARLANE R E, MUIR D W, BOICOURT R M, et al. The NJOY nuclear data processing system, Volume Ⅰ: User’s manual, LA-9303-M[R]. US: Los Alamos National Laboratory, 1982.