-
部分核素存在一段特殊的能量区间，在该能量区间内共振峰非常狭窄，由于使用的核测量仪器分辨率的限制以及物理现象本身的随机性，使得该能区内的共振参数数据无法直接测量得到，把这部分能量区间称为不可分辨共振能区[1]。在截面库的制作过程中，对不可分辨能区共振自屏效应的处理至关重要。为了给先进反应堆提供更精确的多群和连续点截面，打造拥有自主知识产权的截面制作程序，华北电力大学开发了先进快堆核截面处理程序AXSP[2]，旨在实现国产高精度截面制作程序的自动化和创新性。通过计算不可分辨能区共振自屏截面去考虑共振自屏效应，主要有两种计算方法，分别是概率表方法[3]和积分统计原理方法。概率表方法包括两种：第一种是“Ladder Sampling”方法，第二种是“Monte Carlo”方法[4]。“Ladder Sampling”方法最早由Levitt提出，在NJOY[5]中PURR模块、清华大学RXSP[6]中的PURC模块[7]以及西安交通大学的NECP-Atlas[8]中概率表模块[9]均采用了这种方法。“Monte Carlo”方法是由Dunn等提出，并在核数据处理程序AMPX[10]的概率表模块PURM[4]中实现。对于积分统计原理方法，这是一种数值积分求解方法，在NJOY[5]和快堆MC2[11]中使用，国内西安交通大学开发的核数据处理程序NECP-Atlas中也采用了这种方法并且对于近似项进行了精确计算，提高了计算精度。本工作主要针对不可分辨能区共振自屏效应的处理，分别使用“Ladder Sampling”方法和有效自屏截面近似方法在AXSP中开发两个模块，同时对概率表模块的计算效率进行改进，并通过基准题进行验证。
-
对ENDF/B-VII.1评价核数据库[12]中给出的不可分辨共振参数在选定的能量范围内的平均值及其满足的统计分布进行抽样，为特定的自旋序列构建“Ladder”。
先对共振间距抽样，共振间距符合Wigner公式：
$$ f\left(x\right) = \frac{\pi}{2}x{\mathrm{e}}^{-\tfrac{\pi}{4}{x}^{2}} \text{，} $$ (1) $$ x = \frac{{D}_{lj}}{ \left\langle {D}_{lj} \right\rangle } \text{，} $$ (2) 其中：Dlj为抽样共振间距；
${ \left\langle {D}_{lj} \right\rangle } $ 为平均共振间距；l表示轨道角动量；j表示总角动量。使用随机数发生器，其能自动产生符合(0,1)均匀分布的随机数R。并根据下列积分式(3)，将式(1)中f(x)的表达式代入积分式(3)中得到式(4)，对式(4)变形得到x的表达式，可求得x，求出x后代回式(2)可求出抽样共振间距Dlj：
$$ {\int }_{0}^{x}f\left(x\right)\mathrm{d}x = R \text{，} $$ (3) $$ {\int }_{0}^{x}\frac{\text{π}}{2}x{\mathrm{e}}^{-\frac{\text{π}}{4}{x}^{2}}\mathrm{d}x = R \text{，} $$ (4) $$ x = \sqrt{\frac{4}{\text{π}}\mathrm{ln}\frac{1}{1-R}} 。 $$ (5) 共振宽度的满足卡方分布函数：
$$ {p}_{v}\left(y\right) = \frac{v}{2\mathrm{\varGamma }\left(\dfrac{v}{2}\right)}{\left(\frac{vy}{2}\right)}^{\frac{v}{2}-1}{\mathrm{e}}^{-\frac{vy}{2}} , $$ (6) 其中：v为不同反应道对应的卡方分布的自由度，其为常数；
$ \varGamma $ 函数为数学中伽马函数；y可以表示为$$ y = \frac{{\varGamma }_{\lambda r}}{ \left\langle {\varGamma }_{\lambda r} \right\rangle } \text{，} $$ (7) ${\varGamma }_{\lambda r}$ 为反应道r(包括总截面、弹性散射、裂变、辐射俘获反应道)对应共振$ \lambda $ 的抽样宽度；$\left\langle \varGamma _{\lambda r} \right\rangle$ 为平均共振宽度。求解抽样共振宽度的方法与求解共振间距方法类似。在得到了共振参数的共振序列后，使用SLBW公式[5]计算点截面，采用了
$\psi\text{-}\chi$ 方法[5]计算多普勒展宽，给出了弹性散射截面的计算公式：$$ \begin{split} &{\sigma _{n}} = \sum\limits_{{l}} {\frac{{4{\pi}}}{{{k^2}}}} (2l + 1){\sin ^2}{\varphi _{\text{l}}}+ \\ &\sum\limits_{{l}} {\sum\limits_{{r}} {\frac{{4{\pi}}}{{{{{k}}^2}}}\frac{{2{{J}} + 1}}{{4{{I}} + 2}}\frac{{{\varGamma _{{nr}}}}}{{{\varGamma _{{tr}}}}}} } \Biggr\{ \left[\cos 2{\varphi _{\text{l}}} - \left(1 - \frac{{{\varGamma _{{nr}}}}}{{{\varGamma _{{tr}}}}}\right)\right]\times \\ &\psi (\theta ,X) + \sin 2{\varphi _{\text{l}}}\chi (\theta ,X)\Biggr\}\text{，} \end{split} $$ (8) 其中：l为轨道角通量；k表示中子波数；
$\varphi_{\mathrm{l}}$ 为相移因子；r为给定l的共振数；J为对应共振的自旋；I为总自旋；${\varGamma }_{nr}$ 和${\varGamma }_{tr}$ 为共振r的弹性散射宽度和总宽度；ψ和χ分别为复数误差函数的实部和虚部函数；其中$ \theta $ 和${X}$ 的表达式如下所示：$$ \theta = \sqrt {\frac{{{A}}}{{4{{kT}}{{{E}}_0}}}} {\varGamma _{{{tr}}}^{}} \text{，} $$ (9) $$ X = \frac{{2\left( {{{{E}}_0} - {{E}}_r'} \right)}}{{{\varGamma _{{{tr}}}^{}}}} \text{，} $$ (10) A为原子核质量与中子质量比值；k为波尔兹曼常数；T为温度；E0为入射中子能；
${E}'_{{r}}$ 为引发共振r的能量。通过上述步骤后得到了能量网格中一系列对应能量点下的各反应道截面，先对总截面进行排序，根据算法找出总截面最小值σmin和最大值σmax，以σmin对应的能量值作为第一个概率箱的能量，再根据式(11)：
$$ {\sigma _{{i}}} = {\sigma _{\min }} \boldsymbol\cdot {\left( {\frac{{{\sigma _{\max }}}}{{{\sigma _{\min }}}}} \right)^{\frac{{{{k}} - 2}}{{{{K}} - 2}}}} \text{，} $$ (11) 其中：k = 2, 3, ···, K得到σ2，σ3，σ4直到σn等截面，并以这些截面对应的能量值依次作概率箱的能量。同时，把其他截面也存入对应的概率箱中，通过每个概率箱中相应的不同反应道的截面值除以该概率箱的能量点数目可以得到不同反应道的截面平均值。并且每个概率箱的概率是通过该概率箱所有能量点数目与总的概率箱数目相除得到。
生成的概率表数据可以给蒙特卡罗程序直接使用，Bondarenko共振自屏截面可以作为多群截面使用，Bondarenko共振自屏截面可用概率表数据通过式(12)计算：
$$ {\sigma _x}\left( {{E}} \right) = \frac{{\displaystyle\sum\limits_{{i}} {\dfrac{{{{{P}}_{{i}}}\left( {{E}} \right){\sigma _{x{{i}}}}\left( {{E}} \right)}}{{{\sigma _0} + {\sigma _{{{ti}}}}\left( {{E}} \right)}}} }}{{\displaystyle\sum\limits_{{i}} {\dfrac{{{{{P}}_{{i}}}\left( {{E}} \right)}}{{{\sigma _0} + {\sigma _{{{ti}}}}\left( {{E}} \right)}}} }} {，} $$ (12) 其中：x表示不同的反应道(为t、n、f、γ时分别表示总截面、弹性散射、裂变和辐射俘获)；Pi(E)表示对应能量点E下概率表中第i个概率值(其中i为1到概率箱数目最大值)；σxi(E)表示对应能量点E下概率表中x反应道的第i个对应截面值；σ0为本底截面。针对式(12)中σx(E)的计算，分子将对应的概率表中第i个概率数值Pi(E)乘以对应x反应道的第i个截面值再与概率表中总截面值σti(E)加本底截面之和作商再依次求和，分母同理，得到分子分母的数值即可求出σx(E)。
-
对于某一能量点E*处的有效自屏截面，可以看成求解在该能量点附近能量区间的平均截面。各个反应道的0阶有效自屏截面可表示为
$$ {\overline \sigma _{0x}}\left( {{E^*}} \right) = \frac{{\displaystyle\int_{{{{E}}_1}}^{{{{E}}_2}} {{\sigma _x}\left( E \right){\phi _0}\left( E \right)} {\text{d}}E}}{{\displaystyle\int_{{{{E}}_1}}^{{{{E}}_2}} {{\phi _0}\left( E \right)} {\text{d}}E}} \text{，} $$ (13) 式中E*是能量范围[E1, E2]内的有效能量；
$ {\sigma }_{x}\left(E\right) $ 为x反应道截面值；${\phi }_{0}\left(E\right)$ 为0阶权重能谱。为了计算输运截面的有效值，还需计算中子通量的1阶分量${\phi }_{1}\left(E\right)$ 加权的1阶有效自屏总截面：$$ {\overline \sigma _{1t}^{}}\left( {{E^*}} \right) = \frac{{\displaystyle\int_{{{{E}}_1^{}}}^{{{{E}}_2^{}}} {{\sigma _x^{}}\left( E \right){\phi _1^{}}\left( E \right)} {\text{d}}E}}{{\displaystyle\int_{{{{E}}_1^{}}}^{{{{E}}_2^{}}} {{\phi _1^{}}\left( E \right)} {\text{d}}E}} \text{，} $$ (14) 式中
${\phi }_{1}^{}\left(E\right)$ 为1阶权重能谱。$$ {\overline \sigma _{0x}^{}}\left( {{E^ * }} \right) = {{{b}}_x^{}} + \frac{{\displaystyle\int_{{{{E}}_1^{}}}^{{{{E}}_2^{}}} {\frac{{{\sigma _{{{R}}x}^{}}\left( E \right)}}{{\overline \sigma + {\sigma _{{{Rt}}}^{}}\left( E \right)}}C\left( E \right)} {\text{d}}E}}{{\displaystyle\int_{{{{E}}_1^{}}}^{{{{E}}_2^{}}} {\frac{1}{{\overline \sigma + {\sigma _{{{Rt}}}^{}}\left( E \right)}}C\left( E \right){\text{d}}E} }} \text{，} $$ (15) $$ {\overline \sigma _{1t}^{}}\left( {{E^ * }} \right) = \frac{{\displaystyle\int_{{{{E}}_1^{}}}^{{{{E}}_2^{}}} {\frac{1}{{\overline \sigma + {\sigma _{{{Rt}}}^{}}\left( E \right)}}C\left( E \right){\text{d}}E} }}{{\displaystyle\int_{{{{E}}_1^{}}}^{{{{E}}_2^{}}} {\frac{1}{{{{\left[ {\overline \sigma + {\sigma _{{{Rt}}}^{}}\left( E \right)} \right]}^2}}}C\left( E \right){\text{d}}E} }} - {\sigma _0^{}} \text{，} $$ (16) 其中：bx是ENDF/B-VII.1评价核数据库[12]中给出的平滑本底截面；
$\overline{\sigma }$ 为平均截面；$ {\sigma }_{0} $ 为稀释截面；C(E)被假定为E的一个缓慢变化的函数，近似认为其为常数，即C(E)=1；${\sigma }_{{R}x}^{}\left(E\right)$ 的表达式为$$ {\sigma _{{{R}}x}}\left( E \right) = \sum\limits_{{s}} {\sum\limits_{{r}} {{\sigma _{xsr}}\left( {E - {{{E}}_{{{sr}}}}} \right)} } \text{，} $$ (17) 其中：
$ {\sigma }_{x\mathrm{s}\mathrm{r}} $ 表示x反应道下自旋s中共振r对应的截面；Esr是该共振的中心能量。$$ {\overline \sigma _{0x}} = {{{b}}_x} + \frac{{\bar \sigma {I_{0x}}}}{{1 - {I_{0t}}}} \text{，} $$ (18) $$ {\overline \sigma _{1t}} = {{{b}}_{t}} + \frac{{\bar \sigma {I_{{{1t}}}}}}{{1 - {I_{{{0t}}}} - {I_{{1t}}}}} \text{。} $$ (19) 0阶波动积分I0x和1阶波动积分I1t表达式如下：
$$ {I_{0x}} = \frac{1}{{{{{E}}_2} - {{{E}}_1}}}\int_{{{{E}}_{\text{1}}}}^{{{{E}}_{\text{2}}}} {\frac{{\displaystyle\sum\nolimits_{{{sr}}} {{\sigma _{x{{sr}}}}\left( {E - {{{E}}_{{{sr}}}}} \right)} }}{{\overline \sigma + \displaystyle\sum\nolimits_{{{sr}}} {{\sigma _{{{tsr}}}}\left( {E - {{{E}}_{{{sr}}}}} \right)} }}{\text{d}}E} \text{，} $$ (20) $$ {I_{{1t}}^{}} = \frac{1}{{{{{E}}_{\text{2}}^{}}-{{{E}}_{\text{1}}}}}\int_{{{{E}}_{\text{1}}}}^{{{{E}}_{\text{2}}^{}}} {\frac{{\overline \sigma \displaystyle\sum\nolimits_{{{sr}}}^{} {{\sigma _{x{{sr}}}^{}}\left( {E - {{{E}}_{{{sr}}}^{}}} \right)} }}{{{{[\overline \sigma + \displaystyle\sum\nolimits_{{{sr}}}^{} {{\sigma _{{{tsr}}}}\left( {E - {{{E}}_{{{sr}}}^{}}} \right)} ]}^2}}}{\text{d}}E} 。 $$ (21) 对波动积分做分解，得到三种不同的积分，其中孤立项Bxs和Dts的表达式如下：
$$ \begin{split} {B_{xs}^{}} = & \frac{1}{{{{{E}}_{\text{2}}^{}}-{{{E}}_{\text{1}}^{}}}}\int_{{{{E}}_{\text{1}}^{}}}^{{{{E}}_{\text{2}}^{}}} {\sum\limits_r {\frac{{{\sigma _{x{{sr}}}^{}}}}{{\overline \sigma + {\sigma _{{{tsr}}}^{}}}}} {\text{d}}E} \\ =&\frac{1}{{{{{D}}_{{s}}^{}}}}\int {P\left( \alpha \right)\int_{ - \infty }^\infty {\frac{{{\sigma _{x{\text{s}}\alpha }^{}}\left( \xi \right)}}{{\overline \sigma + {\sigma _{{{ts}}\alpha }^{}}\left( \xi \right)}}{\rm{d}}\xi {\rm{d}}\alpha } } \end{split} \text{，} $$ (22) $$ \begin{split} {D_{{ts}}} = & \frac{1}{{{{{E}}_{\text{2}}}-{{{E}}_{\text{1}}}}}\int_{{{{E}}_{\text{1}}}}^{{{{E}}_{\text{2}}}} {\sum\limits_{\text{r}} {\frac{{\overline \sigma {\sigma _{{{tsr}}}}}}{{{{\left[ {\overline \sigma + {\sigma _{{{tsr}}}}} \right]}^2}}}} {\text{d}}E} \\ = & \frac{1}{{{{{D}}_{\text{s}}}}}\int {P\left( \alpha \right)\int_{ - \infty }^\infty {\frac{{\overline \sigma {\sigma _{x{\text{s}}\alpha }}\left( \xi \right)}}{{{{\left[ {\overline \sigma + {\sigma _{{{ts}}\alpha }}\left( \xi \right)} \right]}^2}}}{\text{d}}\xi {\text{d}}\alpha } } \end{split} \text{，} $$ (23) 式中Ds为平均能级间距；
$ P\left(\alpha \right) $ 为共振参数$ \alpha $ 服从的概率分布；$\xi = E-{{{{E}}}}_{{sr}}$ ；第三种积分是序列内重叠积分：$$ {V_{0xs}} = \frac{1}{{{{{E}}_{\text{2}}}{-}{{{E}}_{\text{1}}}}} \int_{{{{E}}_{\text{1}}}}^{{{{E}}_{\text{2}}}} {\sum\limits_{{r}} {\sum\limits_{{{r}'} \ne {{r}}} {\frac{{{\sigma _{x{{sr}}}}}}{{\overline \sigma + {\sigma _{{{tsr}}}}}}\frac{{{\sigma _{x{{sr}}}}}}{{\overline{\sigma} + {\sigma _{{{tsr}}}} + {\sigma _{{{tsr}'}}}}}} } {\text{d}}E}。 $$ (24) 有效截面近似方法忽略了序列内重叠效应，因此波动积分(20)和(21)可以写成下式：
$$ {I_{0x}} = \sum\limits_{{s}} {B_{x{{s}}}}\left(1 - \sum\limits_{{{{s}}^{{'}}} \ne {{s}}} {{B_{{{t}}{{{s}}^{{'}}}}}} \right) \text{，} $$ (25) $$ {I_{1{{t}}}} = \sum\limits_{{s}} {D_{{{ts}}}}\left(1 - \sum\limits_{{{{s}}^{{'}}} \ne {{s}}} {B_{{{t}}{{{s}}^{{'}}}}}\right)^2 。 $$ (26) 由式(25)和(26)，0阶有效自屏截面与1阶有效自屏总截面近似方法的计算公式可以写成：
$$ {{\overline \sigma }_{0x}} = {{b}_x} + \frac{{\overline \sigma \displaystyle\sum\limits_{{s}} {{B_{x{{s}}}}\left(1 - \displaystyle\sum\limits_{{{{s}}'} \ne {{s}}} {B_{{{t}}{{{s}}'}}}\right) } }}{{1 - \displaystyle\sum\limits_{{s}} {{B_{{{ts}}}}\left(1 - \displaystyle\sum\limits_{{{{s}}'} \ne {{s}}} {B_{{{t}}{{{s}}'}}}\right) } }} \text{，} $$ (27) $$ {{\overline \sigma }_{1{{t}}}} = {{{b}}_{t}} + \frac{{\overline \sigma \displaystyle\sum\limits_{{s}} {{D_{{{ts}}}}\left(1 - \displaystyle\sum\limits_{{{{s}}^{{'}}} \ne {{s}}} {B_{{{t}}{{{s}}^{{'}}}}}\right)^2 } }}{{1 - \displaystyle\sum\limits_{{s}} {{B_{{{ts}}}}\left(1 - \displaystyle\sum\limits_{{{{s}}^{{'}}} \ne {{s}}} {B_{{{ts}'}}}\right) - \displaystyle\sum\limits_{{s}} {{D_{{{ts}}}}\left(1 - \displaystyle\sum\limits_{{{{s}}^{{'}}} \ne {\text{s}}} {B_{{{t}}{{{s}}^{{'}}}}}\right)^2 } } }} \text{，} $$ (28) 其中
${B}_{xs}$ ，${B}_{ts}$ 和${D}_{ts}$ 的计算公式为$$ \begin{gathered} {B_{xs}} = \frac{1}{{{{{D}}_{s}}}}\int {P\left( \alpha \right)\frac{{{\varGamma _x}}}{{2\cos 2{\varphi _{\text{l}}}}}} J\left( {\beta ,\theta ,\tan 2{\varphi _{\text{l}}},0} \right){\text{d}}\alpha , \\ \end{gathered} $$ (29) $$ \begin{gathered} {B_{ts}} = \frac{1}{{{{{D}}_{{s}}}}}\int {P\left( \alpha \right)\varGamma } J\left( {\beta ,\theta ,\tan 2{\varphi _{\text{l}}},\tan 2{\varphi _{\text{l}}}} \right){\text{d}}\alpha, \\ \end{gathered} $$ (30) $$ \begin{gathered} {D_{ts}} = \frac{1}{{{{{D}}_{{s}}}}}\int {P\left( \alpha \right)} \varGamma K\left( {\beta ,\theta ,\tan 2{\varphi _{\text{l}}},\tan 2{\varphi _{\text{l}}}} \right){\text{d}}\alpha, \\ \end{gathered} $$ (31) 式中
$ {\phi }_{\mathrm{l}} $ 为相移因子；$ {\mathrm{\varGamma }}_{x} $ 和$ \mathrm{\varGamma } $ 为x反应道对应共振宽度和总宽度；$ \theta $ 表达式见式(10)；其中，J，K积分[5]和$ \mathrm{\beta } $ 表示为$$ J\left( {\beta ,\theta ,a,b} \right) = \frac{1}{2}\int_{ - \infty }^\infty {\frac{{\psi \left( {\theta ,X} \right) + {b}\chi \left( {\theta ,X} \right)}}{{\beta + \psi \left( {\theta ,X} \right) + {a}\chi \left( {\theta ,X} \right)}}{\text{d}}X}, $$ (32) $$ K\left( {\beta ,\theta ,a,b} \right) = \frac{1}{2}\int_{ - \infty }^\infty {\frac{{\beta \big[ {\psi \left( {\theta ,X} \right) + {b}\chi \left( {\theta ,X} \right)} \big]}}{{{{\big[ {\beta + \psi \left( {\theta ,X} \right) + {a}\chi \left( {\theta ,X} \right)} \big]}^2}}}{\text{d}}X}, $$ (33) $$ \beta = \frac{{\overline \sigma }}{{{\sigma _{m}^{}}\cos 2{\varphi _{\text{l}}^{}}}} \text{，} $$ (34) $$ {\sigma _{m}} = \frac{{4{\text{π }}{g_{\text{J}}^{}}}}{{{{{k}}^2}}}\frac{{{\varGamma _{n}}}}{\varGamma } \text{，} $$ (35) $$ {g_{{J}}^{}} = \frac{{2{{J}} + 1}}{{2(2{{I}} + 1)}} \text{，} $$ (36) 式中ψ和χ分别为复数误差函数的实部和虚部函数；k表示中子波数；J为对应共振的自旋；I为总自旋。相移因子
$ {\varphi }_{l} $ 取0，对于式(29)、(30)和(31)中$ J(\beta ,\theta ,\mathrm{0,0}) $ ，$ K(\beta ,\theta ,\mathrm{0,0}) $ 积分的求解，采用Gauss-Jacobi求积方法，再对共振宽度的概率分布进行n次积分，对于具有1到4个自由度的χ2概率分布，采用了一组10个正交点和权重，将n次积分转化为n次求和。 -
基于上述理论方法，使用FORTRAN2003语言在AXSP中开发了PUnresXS模块和UnresXS模块。通过输入卡文件可以选择自己需要计算的模块，输入本底截面和温度的数量及数值，同时程序也实现多个核素一起计算的功能，只需填入核素名称即可根据评价库的路径依次进行读取和处理，方便用户使用。对于不同模块处理后的PENDF格式输出，存放在对应的文件夹下，易于查看，AXSP程序的计算流程如图1所示。
图 1 AXSP的计算流程
-
概率表模块PUnresXS开发完成后，发现其生成概率表的耗时较长，这使得整个程序的计算效率较低。在生成概率表的过程中，通过大量抽样共振参数计算得到了数量庞大的截面值，并且需要对截面值进行排序，其中使用的冒泡排序法效率较低。选取了五种排序方式与冒泡排序法生成概率表的时间做了对比，在同一台电脑上计算，处理器与占用内存均一样，测试了六个核素。
根据表1给出的六种排序算法生成概率表的时间对比可以看到，希尔排序与堆排序计算效率相当，计算效率高于其它排序算法，希尔排序平均时间复杂度为O(n^1.25)，堆排序时间复杂度为O(n*logn)，对于数据量非常大的情况，堆排序较希尔排序法更有优势。考虑到抽样产生了较大数据量的截面，所以在PUnresXS模块中采用了堆排序对算法做了改进。根据表2给出的计算效率对比可以看到，相较于之前的冒泡排序，堆排序使得概率表的计算效率提高了60%以上。
表 1 使用六种排序算法生成六个核素的概率表计算时间
单位：s 排序类型 235U 238U 241Pu 241Am 232Th 74As 冒泡排序 317.52 946.18 489.01 545.72 869.65 647.48 快速排序 187.96 384.67 229.23 288.69 361.28 277.63 直接插入排序 241.46 587.33 328.97 386.16 547.51 417.06 希尔排序 121.83 100.34 100.50 156.23 97.11 77.29 选择排序 313.36 899.60 475.70 536.03 821.42 612.96 堆排序 122.19 104.09 100.48 154.85 101.90 79.99 表 2 堆排序法较冒泡排序提高的计算效率
排序类型 235U 238U 241Pu 241Am 232Th 74As 冒泡排序/s 317.52 946.18 489.01 545.72 869.65 647.48 堆排序/s 122.19 104.09 100.48 154.85 101.90 79.99 提高效率/% 61.52 89.00 79.45 71.62 88.28 87.65 -
为了验证开发完成的不可分辨共振处理模块的正确性，基于ENDF/B-VII.1评价核数据库[12]，以NJOY2016[5]的结果为参考解，分别对PUnresXS模块与UnresXS模块进行了验证。
-
同时选取了高浓缩铀裸球GODIVA基准题[13]与钠冷快堆ZPR-6/7基准题[14]对其进行基准题验证。
GODIVA基准题模型由一个裸铀(富集度为93.71%的235U)球体组成，直径为6.883 in(1 in = 2.54 cm)，质量为52.42 kg，其中包括0.02 kg的236U和极少量的杂质(C，Si和Fe)。
ZPR-6/7基准题采用分体式方形组件设计，装载Pu-U-Zr燃料组件，反射层材料为贫化铀和不锈钢，外部为真空边界，临界测量温度为20 °C。ZPR-6/7按照区域组分和体积构建的均匀R-Z模型[15]，如图2所示。
图 2 ZPR-6/7堆芯RZ模型
基于ENDF/B-VII.1评价库，用NJOY制作了两套 ACE格式截面库。一套截面库用PUnresXS模块替换NJOY中概率表模块PURR制作，另一套直接全部使用NJOY中各模块制作，制作好的截面库提供给蒙特卡罗程序计算。
同时还计算了没有概率表的keff结果与带有概率表的keff计算结果并对它们进行了对比，分别给出了两道基准题的keff结果对比如表3和表4所列。可以看到，有概率表时两个程序计算的GODIVA基准题keff结果偏差为+62×10−5，ZPR-6/7基准题keff结果偏差为+30×10−5。没有概率表时计算的keff相较有概率表的情况偏小，且相较有概率表时两者keff绝对偏差大于400×10−5，这证明了概率表模块开发的必要性。有概率表时计算的keff结果表明在统计误差范围内，使用两个程序的计算结果是一致的。
表 3 GODIVA基准题keff计算结果
程序 有概率表 无概率表 偏差/10−5 AXSP 1.000 85+/-0.000 32 0.996 39+/-0.000 34 +446 NJOY 1.000 23+/-0.000 59 0.995 99+/-0.000 54 +424 表 4 ZPR-6/7基准题keff计算结果
程序 有概率表 无概率表 偏差/10−5 AXSP 0.992 71+/-0.000 32 0.988 44+/-0.000 32 +427 NJOY 0.992 41+/-0.000 33 0.988 28+/-0.000 35 +413 同时，给出了有概率表时计算的ZPR-6/7基准题内堆芯和外堆芯两个代表性区域的中子能谱的计算结果，如图3所示。
图 3 ZPR-6/7基准题内堆芯与外堆芯中子能谱对比(在线彩图)
可以看到，内堆芯和外堆芯中子通量密度的相对误差在绝大多数能量点小于2.5%；内外堆芯在中子能量大于5 MeV时其误差较大，是由于在此区间内中子通量极小，数量级为10−10，当投入粒子数量和计算代数不足时容易引入较大的误差，除能量较高的区域外，中子能谱都在三倍统计区域内，因此中子通量密度计算结果也是一致的。
-
基于ENDF/B-VII.1评价核数据库对有效自屏截面计算模块UnresXS做了核素截面对比与基准题验证，并对结果进行了分析。
设置相同的本底截面与温度，将UnresXS模块计算的有效自屏截面结果与NJOY2016中UNRESR模块计算结果进行了对比与误差分析。图4和图5给出了235U，238Pu，232Th，241Am四个锕系核素的有效截面对比与误差曲线。由图可见，UnresXS模块计算的四个锕系核素的截面结果与NJOY2016的结果符合良好，精度相当，且相对误差均在0.1%以内。
图 4 235U和238Pu有效截面对比与误差曲线(在线彩图)
图 5 232Th和241Am有效截面对比与误差曲线(在线彩图)
选取了ICSBEP手册[16]中的17道临界基准题对有效自屏截面计算模块UnresXS做进一步的验证，17道临界基准题包含主要的裂变核素(如232Th，233U，235U，238U，239Pu，241Pu等)，均为快中子谱基准题，材料温度均设置为300 K，对于235U还分为高富集铀(HEU)、中等富集铀(IEU)和低富集铀(LEU)三种情况。基于ENDF/B-VII.1核评价库，用NJOY制作MATXS格式库，能群用快堆199群，两套库的制作流程与PUnresXS模块类似。使用课题组开发的快堆多群截面程序MGGC[17]对两套库进行处理，产生问题相关的截面数据，再利用ANISN[18]程序去计算临界基准题。
临界计算结果如表5所列，由表5可见，AXSP和NJOY的计算结果绝对偏差均小于0.1×10−5，偏差很小。结果表明有效自屏截面计算模块UnresXS实现了与NJOY中有效自屏截面计算模块UNRESR相当的功能，并且计算结果准确可靠。
表 5 临界基准题计算结果对比
基准题 AXSP_keff NJOY_keff 偏差/10−5 Hmf1 1.000 100 1.000 100 0.00 Hmf2 1.002 140 1.002 140 0.00 Hmf41 1.004 480 1.004 480 0.00 Imf3 1.003 120 1.003 120 0.00 Imf4 1.007 390 1.007 390 0.00 Imf6 1.001 510 1.001 510 0.00 Pmf2 0.994 305 0.994 305 0.00 Pmf5 0.999 369 0.999 369 0.00 Pmf6 0.996 848 0.996 848 0.00 Pmf8 0.994 619 0.994 619 0.00 Pmf9 1.005 260 1.005 260 0.00 Pmf10 0.996 228 0.996 228 0.00 Pmf23 0.996 351 0.996 351 0.00 Pmf25 0.999 246 0.999 246 0.00 Umf1 0.999 957 0.999 957 0.00 Umf2 0.999 235 0.999 236 −0.10 Umf4 0.999 814 0.999 814 0.00 -
本工作通过对不可分辨能区共振自屏截面处理方法的研究，在先进核截面处理程序AXSP中开发了概率表模块PUnresXS和有效自屏截面计算模块UnresXS，通过改进的排序算法显著提高了概率表模块的计算效率。以NJOY2016的计算结果为参考解，概率表模块PUnresXS针对两道基准题的keff结果与参考解的绝对偏差均在100×10−5以内并且中子能谱符合良好，有效自屏截面计算模块235U，238Pu截面值与参考解对比显示误差均在0.1%以内，UnresXS针对17道临界基准题的keff结果与参考解对比绝对偏差在10×10−5以内。计算结果表明，两个模块的计算结果可靠，精度与NJOY2016相当，可以实现在AXSP程序所需承担的功能，能够为后续的模块计算提供可靠的截面结果。
在程序功能实现后，也将进一步进行完善和创新，包括对有效自屏截面近似计算方法的完善，考虑同一自旋序列中共振重叠的影响以及弹性散射截面和总截面中干涉的影响，让有效自屏截面计算模块的计算结果在低本底截面下给出更加精确合理的结果。
Development and Validation of Resonance Self-shielding Cross Section Calculation Module of the Unresolved Energy Region in AXSP Program
-
摘要: 在不可分辨能区，由于共振峰非常密集，受核测量仪器分辨率的限制，不能测得对应能量点下精确的共振参数数据。对该能区，评价核数据库只给出共振参数的平均值和分布函数。为了考虑不可分辨共振能区的共振自屏效应，需要根据这些具有概率分布性质的参数计算有效自屏截面。针对不可分辨共振能区的共振自屏效应的处理，采用“Ladder Sampling”和积分统计原理方法，在先进核截面处理程序AXSP中开发了概率表模块PUnresXS和有效自屏截面计算模块UnresXS，并对概率表计算模块中的算法进行了改进。以NJOY2016的计算结果作为参考解，分别对PUnresXS和UnresXS模块进行了验证，结果表明：概率表模块计算效率提高了60%以上，基准题keff与中子能谱结果与参考解符合良好，开发的共振自屏截面计算模块计算精度与NJOY2016相当，表明了这两个模块在AXSP程序中开发的正确性。
-
关键词:
- 不可分辨能区 /
- 共振自屏截面 /
- Ladder Sampling方法 /
- 积分统计原理近似法
Abstract: In the unresolved energy region, due to the dense distribution of resonance peaks, the precise resonance parameter data cannot be obtained due to the resolution limitation of the nuclear measuring instrument. For this energy region, the Evaluation Nuclear Data Files only provides the average values and distribution functions of the resonance parameters. In order to consider the resonance self-shielding effect in the unresolved resonance energy region, it is necessary to calculate the effective self-shielding cross-section based on these parameters with probability distribution properties. For the treatment of resonance self-shielding effect in the unresolved resonance energy region, a probability table module PUnresXS and an effective self-shielding cross-section calculation module UnresXS were developed in the advanced nuclear cross-section processing program AXSP based on the "Ladder Sampling" and integral statistics principle method. The algorithms in the probability table calculation module were also improved. By comparing with the reference solution obtained from NJOY2016 calculations, the PUnresXS and UnresXS modules were independently validated. The results showed that the efficiency of the probability table module was improved by more than 60%, the keff and neutron spectrum were in good agreement with the reference solution. The developed resonance self-shielding cross-section calculation module exhibited similar accuracy to NJOY2016, demonstrating the correctness of these two modules developed in the AXSP program. -
表 1 使用六种排序算法生成六个核素的概率表计算时间
单位：s 排序类型 235U 238U 241Pu 241Am 232Th 74As 冒泡排序 317.52 946.18 489.01 545.72 869.65 647.48 快速排序 187.96 384.67 229.23 288.69 361.28 277.63 直接插入排序 241.46 587.33 328.97 386.16 547.51 417.06 希尔排序 121.83 100.34 100.50 156.23 97.11 77.29 选择排序 313.36 899.60 475.70 536.03 821.42 612.96 堆排序 122.19 104.09 100.48 154.85 101.90 79.99 
下载: 导出CSV
表 2 堆排序法较冒泡排序提高的计算效率
排序类型 235U 238U 241Pu 241Am 232Th 74As 冒泡排序/s 317.52 946.18 489.01 545.72 869.65 647.48 堆排序/s 122.19 104.09 100.48 154.85 101.90 79.99 提高效率/% 61.52 89.00 79.45 71.62 88.28 87.65
下载: 导出CSV
表 3 GODIVA基准题keff计算结果
程序 有概率表 无概率表 偏差/10−5 AXSP 1.000 85+/-0.000 32 0.996 39+/-0.000 34 +446 NJOY 1.000 23+/-0.000 59 0.995 99+/-0.000 54 +424
下载: 导出CSV
表 4 ZPR-6/7基准题keff计算结果
程序 有概率表 无概率表 偏差/10−5 AXSP 0.992 71+/-0.000 32 0.988 44+/-0.000 32 +427 NJOY 0.992 41+/-0.000 33 0.988 28+/-0.000 35 +413
下载: 导出CSV
表 5 临界基准题计算结果对比
基准题 AXSP_keff NJOY_keff 偏差/10−5 Hmf1 1.000 100 1.000 100 0.00 Hmf2 1.002 140 1.002 140 0.00 Hmf41 1.004 480 1.004 480 0.00 Imf3 1.003 120 1.003 120 0.00 Imf4 1.007 390 1.007 390 0.00 Imf6 1.001 510 1.001 510 0.00 Pmf2 0.994 305 0.994 305 0.00 Pmf5 0.999 369 0.999 369 0.00 Pmf6 0.996 848 0.996 848 0.00 Pmf8 0.994 619 0.994 619 0.00 Pmf9 1.005 260 1.005 260 0.00 Pmf10 0.996 228 0.996 228 0.00 Pmf23 0.996 351 0.996 351 0.00 Pmf25 0.999 246 0.999 246 0.00 Umf1 0.999 957 0.999 957 0.00 Umf2 0.999 235 0.999 236 −0.10 Umf4 0.999 814 0.999 814 0.00
下载: 导出CSV
-
[1] 周峰. 不可分辨能区共振自屏截面处理方法研究[D]. 北京: 华北电力大学, 2022 ZHOU Feng. Research on Processing Method of Resonance Self-shield Cross Section in Unresolved Resonance Energy Region[D]. Beijing: North China Electric Power University, 2022. (in Chinese) [2] CULLEN D E. A Short History of ENDF/B Unresolved Resonance Parameters[R]. California: Lawrence Livermore National Laboratory, 2010. [3] LEVITT L B. Nuclear Science and Engineering, 1972, 49(4): 450. doi: 10.13182/NSE72-3 [4] DUNN M E, LEAL L C. Nuclear Science and Engineering, 2004, 148(1): 30. doi: 10.13182/NSE04-A2438 [5] MACFARLANE R, MUIR D W, BOICOURT R M, et al. The NJOY Nuclear Data Processing System, Version 2016[R]. New Mexico: Los Alamos National Laboratory, 2018. [6] 余健开. 堆用重要核数据评价处理及应用的关键方法与软件研发[D]. 北京: 清华大学, 2015. YU Jiankai. Research on the Key Methods of Evaluation, Processing and Application of Important Nuclear Data for Reactor and Code Development[D]. Beijing: Tsinghua University, 2015. (in Chinese) [7] YU J, HUO H, LI W, et al. Validation of a New PURC Module for Calculating Probability Table in Unresolved Resonance Region in RXSP Code[C]//22nd International Conference on Nuclear Engineering. Volume 5: Innovative Nuclear Power Plant Design and New Technology Application; Student Paper Competition. Prague: ASME. 2014: V005T17A034. [8] ZU Tiejun, XU Jialong, TANG Yongqiang, et al. Annals of Nuclear Energy, 2019, 123: 153. doi: 10.1016/j.anucene.2018.09.016 [9] 毕沪超, 祖铁军, 徐嘉隆, 等. 核动力工程, 2020, 41(3): 8. BI Huchao, ZU Tiejun, XU Jialong, et al. Nuclear Power Engineering, 2020, 41(3): 8. (in Chinese) [10] WIARDA D , DUNN M E, GREENE N M, et al. AMPX-6: A Modular Code System for Processing ENDF/B: ORNL/TM-2016/43[R]. Tennessee: Oak Ridge National Laboratory, 2016. [11] GREEN D M, PITTERLE T A. ETOE: A Program for ENDF/B to MC 2 Data Conversion[R]. Washington D C: United States Atomic Energy Commission, 1968. [12] CHADWICK M B, HERMAN M, OBLOŽINSKÝ P, et al. Nuclear Data Sheets, 2011, 112(12): 2887. doi: 10.1016/j.nds.2011.11.002 [13] LABAUVE R J. Bare, Highly Enriched Uranium Sphere (GODIVA), NEA/NSC/DOC (95)03, Volume II, HEU-MET-FAST-001[R]. New Mexico: Los Alamos National Laboratory, 2002. [14] SIMTH M A, LELL R M, Moneo P. ZPR-6 ASSENBLY 7: A Cylindrical Assembly with Mixed (Pu-U)-Oxide Fuel and Sodium with a Thick Depleted-Uranium Reflector: NEA/NSC/DOC(95)03/VI[R]. Chicago: Argonne National Laborary, 2003. [15] Methods and Issues for the Combined Use of Integral Experiments and Covariance Data: NEA/NSC/WPEC/DOC(2013)445[R]. Paris: OECD-Nuclear Energy Agency, 2013. [16] BESS J D. International Criticality Safety Benchmark Evaluation Project[R]. Paris: OECD-Nuclear Energy Agency, 2006. [17] 黄自锋. 共振峰重构与快堆超精细能群截面制作方法研究[D]. 北京: 华北电力大学, 2021. HUANG Zifeng. Research on Reconstruction and Linearization of Resonance Cross Section and Fabrication Method of Ultrafine Energy Group Cross Section of Fast Reactor[D]. Beijing: North China Electric Power University, 2021. (in Chinese) [18] ENGLE J R. A User’s Manual for ANISN, A One-Dimensional Discrete Ordinates Transport Code with Anisotropic Scattering[M]. Tennessee: Oak Ridge National Laboratory, 1967. -
点击查看大图
图(5) / 表 (5)
计量
- 文章访问数: 273
- HTML全文浏览量: 46
- PDF下载量: 7
- 被引次数: 0
甘公网安备 62010202000723号