基于MCNP5和ORIGEN2.1程序的停堆剂量程序开发与初步验证

马辉强,苑旭东,于 涛*,陈珍平,谢金森,吴菱艳

(1.南华大学 核科学技术学院,湖南 衡阳 421001;2.南华大学 湖南省数字化反应堆工程技术研究中心,湖南 衡阳 421001;3.中国核动力研究设计院 核反应堆系统设计技术重点实验室，四川 成都 610213)

0 引 言

反应堆结构件在反应堆运行过程中受中子照射而被活化成放射性核素，使得反应堆构件通常具有较高的放射性活度。核动力系统的运行维护和退役必须保证放射性物质对人员的辐照剂量限制以及对环境的污染最小化[1-2]，所以需要准确了解反应堆的放射性水平和各种核素的特性。

测量是源项活度信息以及剂量分布最直接的方法，但实际上反应堆的大部分构件由于所处位置难以进入或由于其辐射剂量率太大而使测量无法进行[3-4]，因此为了弥补测量方法的不足而较快地得到退役反应堆的源项以及剂量分布数据，对核动力系统中进行活化和剂量分布的计算评估非常有必要。本文将基于三维蒙特卡罗粒子输运程序MCNP5与点燃耗ORIGEN2.1，采用“严格两步法”(R2S)[5]开发出用于反应堆结构材料辐照停堆剂量程序MOCA，并进行初步验证。

1 方法基础

1.1 R2S方法

R2S方法基于以下计算步骤：1)中子输运计算：确定装置运行时的中子注量率在空间上的分布。2)材料活化计算：通过辐照历史条件及步骤1)的中子注量率分布，求出装置停机后材料的活化情况及衰变光子源在空间上的分布。3)光子输运计算：利用步骤2)求出的衰变光子分布，计算中子辐照停止以后，不同地点剂量率分布。R2S方法以通用的MCNP程序做输运计算，采用ORIGEN2.1程序做活化计算，通过自动接口处理程序自动完成各个计算过程的数据传递与处理，实现MCNP5与ORIGEN2.1的输运计算与活化计算以及活化后的空间剂量率计算的自动耦合。MCNP5和ORIGEN2.1程序构成R2S方法计算过程的关键部分。

1.2 中子输运计算

MCNP5[6]具有对复杂几何形状区域的适应性强等优点，对于辐照停堆剂量计算，MCNP5可以模拟中子、光子输运，计算中子通量分布、各种中子核反应率以及剂量率等物理参量。MCNP5及以下版本不能直接分析辐照停堆剂量的问题，需要耦合其他燃耗或者活化计算程序进行计算。

中子输运过程主要通过F4和Fm4计数卡，计算各个活化栅元中的平均中子通量密度以及重要核素的各种中子核反应率。F4输出的中子通量密度是归一化到每个中子源的通量密度，单位为：n/cm2，其计算方式为：

(1)

ORIGEN2.1要求替换的截面为微观单群截面，MCNP5中不能够直接输出核素的反应截面，而是计算每种核素在每裂变源中子下的平均核反应率，即每秒每裂变中子源每单位体积内的中子与介质原子核发生反应的总次数的统计平均值。MCNP5程序通过Fm4计数卡来获得各个活化栅元中的平均核反应率的，其计算方法为：

(2)

其中：i为所需求解的某一反应；j为某种需替换截面的核素；R为平均核反应率；σi,j(E)为核素裂变截面为指定能量下的反应截面；φ(r,E,t)为在位置r，能量E的中子通量。

1.3 活化计算

ORIGEN2.1[7]是美国橡树岭实验室开发的点燃耗计算程序，采用单群截面，所需计算时间短，能够模拟中子与核燃料中的核素发生各种中子核反应过程，计算核反应过程中放射性物质的积累、衰变等过程，给出核素的组成、放射性活度、中子和18群光子产额等源项参数，适合于不同类型反应堆的燃耗、活化和衰变过程的计算。ORIGEN2.1程序中的模型有几个假设，一是假设中子通量和反应截面不随核素成分变化；二是假设所有核素具有相同的空间点，经历同等中子辐照，故无法处理一些与系统几何外形相关的问题，如中子通量密度、剂量率的空间分布等，因此需要耦合输运程序实现剂量率的空间分布计算。

ORIGEN2.1程序中采用非齐次一阶常微分方程:

(3)

式中，Xi为核素i的核密度；N为核素的数量；lij为核素j的放射性衰变份额；N为核素的数量；λi为核素i的衰变常数；λj为核素j的衰变常数；φ为中子注量率；fik为核素k的中子吸收份额；δi为核素i的谱平均中子吸收截面；ri为核素i的连续移除率；Fi为核素i的连续供给率。

1.4 有效单群微观截面与绝对反应率

ORIGEN2.1中程序考虑了(n,γ)、(n,2n)、(n,3n)、(n,α)、(n,p)、(n,f)等中子核反应，具有一套与堆型、燃耗深度相对应的数据库，该数据库中中子反应截面采用单群微观截面(单位:cm2)。不同的反应堆内，中子能谱各有差异，同时反应堆内中子从堆芯输运到堆内构件、压力容器和主屏蔽墙，与金属、含硼水发生散射和各种核反应，在不同的空间位置中子能谱发生了显著变化。研究表明[8]：对距离堆芯较近处的构件活化计算，不做截面修正或许是可行的；但是距离堆芯较远处，中子能谱发生了显著变化，不进行截面修正可能造成计算结果偏差较大，不能简单地使用ORIGEN2.1自带的截面数据进行活化源项的计算，需要对活化栅元中重要的核素及其反应链上相关核素的反应截面进行修正。

(4)

ORIGEN2.1所需的截面号在MCNP5中分别为102(n,γ)、16(n,2n)、17(n,3n)/107(n,α)、-6(n,f)/103(n,p)，其中斜杠前代表锕系核素所需计算的截面，斜杠后代表裂变产物和活化产物所需计算的截面，而ORIGEN2.1的中子核反应数据库要求有六个截面，由于MCNP5输运计算过程中不考虑核素的激发态，对于中子核反应数据库后两个截面设为0或者不进行替换。

ORIGEN2.1计算时需要输入的是绝对中子通量密度，需要做如下方式转换处理：

φ绝对=F4×FMF

(5)

R绝对=Fm4×FMF×Ni

(6)

其中FMF为中子输运时的中子源强常数，由辐照功率归一化得出或者直接给定，Ni为核素或元素i的原子密度，单位为:1024/cm3。

1.5 自动耦合过程

根据配置文件中活化栅元的栅元号、中子辐照史以及需要进行核素截面替换的核素ID生成MCNP5中子输入文件进行中子输运计算，然后从MCNP5输出文件中获取活化栅元中的中子通量密度、重要核素的中子核反应率以及活化栅元中各个元素或核素的质量，并计算各个活化栅元中需要进行截面替换的核素的等效单群截面。然后再生成ORIGEN2.1的输入文件tape5.inp和tape3.inp，完成各个栅元的活化计算得出各活化栅元中的核子密度、核子活度及衰变光的能群子分布，最后根据衰变光子分布生成MCNP5输入文件进行光子输运计算，得出光子剂量率的空间分布。程序耦合计算的流程图如图1所示。

2 程序验证

基于MCNP5和ORIGEN2.1程序以及“严格两步法”，开发了反应堆结构材料辐照停堆剂量程序MOCA，为了验证MOCA程序的正确性，本文通过建立结构材料中子辐照停堆剂量模型，采用SuperMC[9]程序对MOCA程序计算的核子数、多群衰变光子的产额以及活化源附近剂量率的空间分布结果进行对比验证。SuperMC为中国科学院核能安全研究所自主研发的大型粒子输运软件，具备辐照停堆剂量计算功能。

2.1 压水堆燃料元件包壳验证算例

2.1.1 模型参数

压水堆燃料包壳是反应堆的重要结构材料之一，其材料的活化水平直接影响到燃料包壳的寿命与堆芯安全设计。选取M5材料包壳作为活化对象，建立燃料元件栅格模型，如图2所示，包壳的几何尺寸、核素成分及辐照方案分别如表1、2、3所示。中子源为裂变中子源，MOCA中ORIGEN2.1活化计算采用热中子库201、202、203，计算辐照后关键核素的核子数、光子能群产额，以及径向距离包壳5～50 cm、轴向-180～180 cm处的光子剂量率并与SuperMC计算的结果进行对比验证。

表1 包壳栅格尺寸

表2 包壳材料成分表

表3 中子辐照方案

2.1.2 结果与分析

1)重要核素产额

通过分析SuperMC计算的结果中各种核素对衰变光子释放率的贡献，如表4所示，发现M5包壳经过辐照后的衰变光子主要由59Ni、91Y、93Zr、94Nb产生，对比分析SuperMC与MOCA产生的这几种关键核素的核子数，结果表明MOCA计算的核子数与SuperMC的结果吻合较好，最大相对误差小于8%，对于辐照停堆剂量问题最大误差限制要求小于30%[2]，MOCA计算的放射性核素核子数与SuperMC最大相对误差小于8%是可以接受的。核素41094(94Nb)计算结果差异较大是由于MOCA在计算过程中，截面替换的核素种类考虑不够充分。

表4 放射性核素产额

2)衰变光子产额

MOCA中ORIGEN2.1使用自带的数据库并将衰变光子分为18个能群，SuperMC直接使用ENDF数据库并将衰变光子能群分为24群[10]。如图3所示，由于ENDF数据库中的衰变数据更为详细，包含了较多低能光子的衰变产额，导致在低能群中SuperMC计算得到的光子产额将大于MOCA计算得到的产额，但是在高能区的光子产额MOCA与SuperMC计算的结果符合较好。

3)空间剂量率分布

MOCA将活化计算得出的各个栅元中光子产额以及光子能谱分布生成MCNP5的输入文件进行光子输运计算，计算光子的空间剂量率分布。MOCA和SuperMC通过探测器计算距离包壳50 cm内，轴向-180 cm～180 cm范围内的光子剂量率。图4中(a)、(b)、(c)分别为轴向-180 cm、0 cm、180 cm处径向0～50 cm范围内的光子剂量率分布；图5中(a)、(b)、(c)分别为径向5 cm、25 cm、50 cm处轴向-180～180 cm范围内的光子剂量率分布。可以发现，MOCA与SuperMC计算的空间剂量率分布的趋势和数值吻合较好。

2.2 不锈钢辐照剂量验证算例

2.2.1 模型参数

铁合金作为主要的反应堆结构材料之一，如吊篮、压力容器等，经过长期的辐照之后，会产生大量的放射性核素，是反应堆厂房中主要的放射性源项，也是源项调查的主要对象之一。本例题通过对不锈钢的辐照活化、光子剂量率空间分布计算，然后与SuperMC程序计算的结果进行对比，验证MOCA程序的正确性。不锈钢材料模型结构示意图、几何参数以及成分表如图6、表5、表6所示，辐照方案如表7所示。

表5 不锈钢辐照模型尺寸表

表6 不锈钢材料成分表

表7 中子辐照方案

2.2.2 结果与分析

1)重要核素产额

根据SuperMC计算的结果中各种放射性核素对衰变光子释放率的贡献，发现不锈钢经过辐照后释放的衰变光子主要由14C、32P、55Fe、58Co、59Ni等核素产生，对比SuperMC与MOCA计算结果中这几种关键放射性核素的核子数量(如表8所示)，可以看出MOCA计算的放射性核素核子数与SuperMC的结果吻合较好。

2)衰变光子产额

MOCA与SuperMC计算不锈钢经过辐照和冷却后的衰变光子产额如图7所示，同样由于SuperMC使用的ENDF数据库中的衰变数据更为详细，导致SuperMC计算的结果中低能量衰变光子的数目多于MOCA计算结果中低能量衰变光子的数目。由于对衰变光子释放率的贡献主要集中在几种特定的核素，因此在特定的能量范围内MOCA计算的结果与SuperMC计算的结果的分布趋势基本一致，但由于SuperMC和MOCA中ORIGEN2.1对衰变光子能群边界的划分不同，所以在各能群范围内的光子数量具有一定的差异。

表8 放射性核素产额

3)空间剂量率分布

MOCA将活化计算得出的各个栅元中衰变光子产额以及光子能谱分布生成MCNP5输入文件进行光子输运计算，计算光子的空间剂量率分布。MOCA与SuperMC计算距离不锈钢表面40 cm范围内，轴向-60 cm～60 cm范围内的光子剂量率。图8中(a)、(b)、(c)分别为轴向-60 cm、0 cm、60 cm处径向5～40 cm范围内的光子剂量率分布；图9中(a)、(b)、(c)分别为径向5 cm、25 cm、50 cm处轴向-60～60 cm范围内的光子剂量率分布。可以发现，MOCA与SuperMC计算的空间剂量率分布吻合较好。

3 结 论

本工作基于通用蒙特卡罗输运程序MCNP5与点燃耗程序ORIGEN2.1，采用“严格两步法”开发了三维停堆剂量计算程序MOCA，通过建立辐照停堆剂量模型进行MOCA程序与SuperMC程序对比验证，得出以下结论:

1)由于ORIGEN2.1与SuperMC使用的数据库以及对光子能群的划分不同，在计算核子数、衰变光子数量及其能群分布存在细小的差异，但对衰变光子的空间剂量率影响不大。

2)经验证MOCA程序在计算核子数、衰变光子数量及其能群分布和光子剂量率的空间分布计算的结果与SuperMC计算的结果总体吻合较好，可以用于反应堆结构材料辐照后停堆剂量的计算。