西安脉冲堆中子价值计算与分析

张 良，江新标，陈立新，杨 宁，张信一

（1.西北核技术研究所，陕西 西安 710024；2.强脉冲辐射环境模拟与效应国家重点实验室，陕西 西安 710024）

中子价值Φ＋（r，E，Ω）的定义［1-2］是：在临界反应堆中，在r处投入一个能量为E、运动方向为Ω 的中子引起的对稳定功率的贡献，其物理意义是中子对反应堆裂变功率的贡献大小。中子价值也称中子重要性或伴随中子通量，其概念在1945年由Wigner［3］在反应性方程中首次提出，它在反应堆微扰、敏感性分析和堆芯动态计算方面有着重要应用［4］。Bell和Glasstone［5］在1976年给出了伴随中子通量的玻耳兹曼输运方程，中子价值可直接通过求解该方程得到，但在一些几何较复杂的堆芯中，求解较为困难［6］，此时，在满足扩散近似条件的情况下，可通过求解中子价值的扩散方程得到。在一些蒙特卡罗方法计算动态参数的研究中，也提出了中子价值的多种近似或替代方法，例如Nauchi等［7-8］提出了用中子裂变产生的下一代中子数、Meulekamp等［9］提出了用中子下一代的诱发裂变数来近似代替中子价值，美国洛斯阿拉莫斯国家实验室的MCNP5-1.60 版本采用中子经若干代后产生的裂变中子总数来近似中子价值函数［10］，这些近似方法用于计算动态参数时均得到了较好的结果。

本文采用栅元计算程序WIMS 和扩散计算程序CITATION，计算西安脉冲堆的中子价值分布，研究中子价值随空间和能量的变化规律，并分析缓发中子有效份额与缓发中子份额、中子代时间与瞬发中子寿命存在差异的原因。

1 计算方法

首先采用WIMS程序计算西安脉冲堆的6群群常数，然后采用CITATION 程序进行全堆芯扩散计算，得到堆芯的中子价值分布。

1.1 西安脉冲堆堆芯模型

西安脉冲堆堆芯计算模型如图1所示。脉冲堆堆芯是正六边形模型，本文为了适应CITATION的几何描述，在左上和右下角加了两块黑体，相当于真空，中子进入黑体后不再返回堆芯。稳态控制棒共5根，最右边两根是安全棒，运行中安全棒均被提升至顶位，由于安全棒下部是燃料跟随体（与燃料棒材料相同），在计算中视为燃料棒。

图1 西安脉冲堆堆芯模型Fig.1 Core model of Xi'an Pulsed Reactor

1.2 WIMS程序计算6群群常数

WIMS程序应用中子输运理论来进行栅元计算。

中子输运方程［2］：

式中：φ（r，E，Ω）为中子在相空间点（r，E，Ω）处的中子通量密度；Σt（r，E）为中子在相空间点（r，E）处的裂变截面；Σs（r，E，Ω→E′，Ω′）为从（E，Ω）到（E′，Ω′）的散射截面；Sf（r，E，Ω）为裂变中子产生率；k为反应堆有效增殖因数。

WIMS程序采用SN 方法对各栅元求解式（1），得到φ（r，E），然后进行群常数计算，本文采用的6群分群结构及说明列于表1。

第g 群的群常数Σm，g可表示为：

式中，m 代表群常数的类型，包括总截面Σt、裂变截面Σf、群转移截面Σg-g′等。

表1 6群分群结构Table 1 Six-group structure

1.3 CITATION 程序扩散计算

CITATION 是采用有限差分扩散理论的堆芯分析程序，在堆芯计算中有着广泛的应用。该程序可处理一维、二维、三维扩散计算问题以及燃耗计算问题，计算精度高，且适用于多种几何形状的堆芯。

在由WIMS计算出群常数后，采用CITATION 求解多群共轭扩散方程［2］：

式中：Dg为群扩散系数；ν 为每次裂变产生的中子数。通过求解该方程可得到中子价值（r，E）的分布。本文的计算采用的是二维几何模型，如图1 所示，认为堆芯在高度Z 方向是均匀的。

2 CITATION 扩散计算的正确性验证

CITATION 程序可同时计算中子通量和中子价值，通过将其计算出的通量分布和连续能量截面的蒙特卡罗程序MCNP［11］计算出的通量分布进行比较，来验证CITATION 程序用于西安脉冲堆扩散计算的正确性。如图1所示，选取从左往右第14列，第7～21行的栅元（自上而下编号从1到15，1、14、15号栅元为石墨，2～4、6、10～13号栅元为燃料棒、5号栅元为控制棒，7～9号栅元为水），分别比较第2群（快群）和第6群（热群）在这15个栅元中的通量分布，结果如图2所示。两种程序的计算结果基本符合一致，说明CITATION 程序可较好地用于西安脉冲堆的扩散计算。

图2 MCNP和CITATION 计算得到的中子通量密度分布对比Fig.2 Comparison of neutron flux density calculated by MCNP and CITATION

3 西安脉冲堆中子价值的计算结果及分析

3.1 中子价值的空间分布

如图1所示，同样选取从左往右第14列，第7～21行的栅元，分别比较第2群（快群）和第6群（热群）的中子在这15个栅元中的价值分布，结果如图3所示。

图3 中子价值的空间分布Fig.3 Space distribution of neutron importance

对热群中子，中子价值在2、3号燃料栅元中先增长，然后在4号燃料栅元中有所下降，原因是4号栅元更靠近5号控制棒。然后经历在控制棒栅元中的急剧降低和在6号燃料栅元中的回升后，在7、8、9号水栅元中先降低再回升，原因是8号水栅元最远离燃料棒，其中的热中子最难引起裂变反应，因此中子价值较低。之后中子价值在11号燃料栅元中达到最大，然后依次降低，这是由于12、13号燃料棒逐渐远离堆芯。快中子的价值整体变化平缓，因为快中子的核反应截面较小，不会因剧烈吸收或引起裂变导致价值突然降低或升高。

3.2 中子价值随能量的分布

分别选择1、4、5、8号栅元，分别为石墨、燃料、控制棒和水栅元，研究中子价值在这4种栅元中的中子价值随能量的变化。由图4可知，石墨栅元的中子价值整体相对较低，原因是石墨栅元主要布置在堆芯外围。燃料栅元的中子价值随能量的降低而增大，这是因为燃料的裂变截面随中子能量的降低而升高。控制棒栅元的中子价值在能量降低时迅速下降，原因是控制棒材料对低能中子的吸收截面很大。水栅元的中子价值随能量的降低缓慢降低，原因是水对低能中子的吸收较强，而快中子能量更高，有更大的可能性经慢化后进入燃料栅元引起裂变。

图4 各种栅元在各能群中的中子价值Fig.4 Neutron importance via energy group in different cells

4 中子价值对动态参数计算结果的影响

缓发中子有效份额βeff和中子代时间Λ 是反应堆动态特性分析中常用到的重要动态参数，它们的定义式中均考虑了中子价值的影响，βeff和Λ 对应的不考虑中子价值影响的参数分别为缓发中子份额βo 和瞬发中子寿命lP，本文分析βeff与βo、Λ 与lP存在差异的原因。

4.1 βeff和βo 的比较

缓发中子有效份额βeff的定义式为：

4.2 Λ 和lP 的比较

中子代时间Λ 的定义式为：

中子代时间Λ 的物理意义是：由φ＋（r，E，Ω）作权重的中子数和由φ＋（r，E，Ω）作权重的裂变中子产生率之比。即由φ＋（r，E，Ω）作权重、每产生一代中子所需的时间。

瞬发中子寿命lP的定义式则为：

对于西安脉冲堆，lP≈93μs，Λ≈36μs。Λ和lP存在明显的差异，引起这种差异的最重要原因是燃料栅元中的中子价值较其他栅元中的大。如式（5）所示，Λ 定义式中分子中的中子价值包括了所有栅元的中子价值，而分母中的中子价值则只是燃料栅元中的中子价值。而lP的定义式中未考虑中子价值的影响，因此有Λ＜lP。

Λ 和lP的相对大小，主要是由燃料栅元的平均中子价值和全堆芯的平均中子价值的相对大小决定的，亦即主要是由中子价值的空间分布决定的。如对裸堆Godiva和Jezebel［12］，燃料的中子价值实际上就是全堆芯的中子价值，因此二者的Λ 和lP很接近，相对差别仅为5%左右［13］，而这个微小的差别是由于中子价值随能量的变化引起的，这也同时说明，只考虑中子价值随能量分布的因素，并不会引起Λ 和lP有很大的差异。

事实上，由于大多类型的反应堆中，燃料栅元的中子价值均大于或接近全堆芯的平均中子价值，因此，一般有Λ≤lP。

5 结论

本文采用栅元程序WIMS 计算了西安脉冲堆6群群常数，采用扩散程序CITATION 进行全堆芯扩散计算，得到了西安脉冲堆的中子价值在空间和能量上的分布，并获得以下分析研究结果：

1）得到了中子价值随空间的变化规律。在选取的15个栅元中，燃料栅元的中子价值较高，控制棒和水栅元中的中子价值相对较低。燃料栅元在远离堆芯或接近控制棒时，中子价值会降低。快中子价值在各栅元中的变化较热中子的价值变化平缓。

2）得到了各栅元的中子价值随能量的变化规律。燃料栅元的中子价值随能量的降低而增大，原因是燃料发生裂变的可能性随中子能量的降低而升高。控制棒栅元中的中子价值随能量的降低而降低，原因是控制棒材料对低能中子有很强的吸收作用。水栅元的中子价值同样随能量的降低有缓慢降低，原因是水对低能中子的吸收较强，高能中子则有更大的可能性经慢化后引发裂变。

3）分析了缓发中子有效份额βeff和缓发中子份额βo 存在差异的原因。βeff和βo 的相对大小，取决于中子价值随能量的分布。对于西安脉冲堆，缓发中子较瞬发中子能量低，其平均价值相对大，因此有βeff＞βo。

4）分析了中子代时间Λ 和瞬发中子寿命lP存在差异的原因。Λ 和lP的相对大小，主要是由中子价值随空间的分布决定的，对于包含西安脉冲堆在内的大多反应堆，燃料栅元中的中子价值较其他栅元中的高，因此有Λ≤lP。

［1］ 杜书华.输运问题的数值模拟［M］.湖南：湖南科技出版社，1988.

［2］ 谢仲生，邓力.中子输运理论数值计算方法［M］.西安：西北工业大学出版社，2005.

［3］ WIGNER E P.The collected works of Eugene Paul Wigner：Effect of small perturbations on pile period［M］.Berlin：Springer Verlag，1992.

［4］ HENRY A F.Nuclear reactor analysis［M］.Massachusetts：MIT Press，1975.

［5］ BELL G I，GLASSTONE S.Nuclear reactor theory［M］.New York：Van Nostrand Reinhold Company，1970.

［6］ FEGHHI S A H，SHAHRIARI M，AFARIDEH H.Calculation of neutron importance function in fissionable assemblies using Monte Carlo method［J］.Annals of Nuclear Energy，2007，34：514-520.

［7］ NAUCHI Y，KAMEYAMA T.Proposal of direct calculation of kinetic parametersβeffand Λ based on continuous energy Monte Carlo method［J］.Journal of Nuclear Science and Technology，2005，42（6）：503-514.

［8］ NAUCHI Y，KAMEYAMA T.Development of calculation technique for iterated fission probability and reactor kinetic parameters using continuous-energy Monte Carlo method［J］.Journal of Nuclear Science and Technology，2010，47（11）：977-990.

［9］ MEULEKAMP R K，van der MARCK S C.Calculating the effective delayed neutron fraction using Monte Carlo［J］.Nuclear Science and Engineering，2006，152：142-148.

［10］KIEDROWSKI B C，BOOTH T E，BROWN F B，et al.MCNP5-1.60feature enhancements &manual clarifications［R］.USA：Los Alamos National Laboratory，2010.

［11］BRIESMEISTER J F.Mont Carlo N-particle transport［R］.USA：Los Alamos National Laboratory，2000.

［12］OECD／NEA Nuclear Science Committee.International handbook of evaluated criticality safety benchmark experiments［M］.USA：OECD／NEA，2008.

［13］KIEDROWSKI B C，BROWN F.Adjoint wighted kinetics parameters with continuous energy Monte Carlo［J］.Tran Am Nucl Soc，2009，100：297-299.