使用NATF进行基于网格的计算停机剂量率流程

说明

本篇文章简要介绍使用NATF来进行基于网格的计算停机剂量率的流程。本流程需要提前安装DAGMC, ALARA, PyNE, NATF。

• DAGMC的安装可以参考这里。
• PyNE的安装可以参考这里。
• ALARA的安装可以参考这里。
• NATF可以通过pip install -U natf安装最新版。

使用NATF计算基于网格的停机剂量率，其内核是PyNE R2S流程，只是通过ptrac代替DAGMC的Ray tracing获取网格栅元体积份额。这样可以不用DAGMC模型，而是使用已经存在的MCNP模型进行计算。主要步骤包括：

• 中子通量计算：计算指定区域的175群中子通量
• ptrac计算：在指定区域均匀抽样粒子，统计粒子的位置及所在栅元信息，用于获取网格栅元体积份额
• 活化计算：包括活化文件生成及ALARA活化计算
• 光子输运计算：光子输运获取停机剂量率

中子通量计算

使用NATF mesh R2S流程，我们需要中子输运计算得到的两个文件：

• outp : MCNP输出文件，里面包含栅元、材料信息
• meshtal : 网格计数输出文件，包含fmesh结果

注意:我们需要outp中的栅元、材料信息，因此中子输运输入文件中，需要使用PRINT关键字

Ptrac计算

Ptrac计算用于在指定区域均匀抽样，记录每个源粒子的位置及栅元信息。然后根据源粒子所在网格来统计每个网格单元的栅元体积份额。注意事项：

• 每个网格单元尽量覆盖1e3以上粒子
• ptrac输出文件格式选择bin

活化计算

进行停机剂量率计算，我们需要的主要信息包括：

1. 每个网格单元的材料成分。必须项。由ptrac及outp计算获取。
2. 每个网格单元的中子通量能谱。必须项。由meshtal文件提供。
3. 活化计算能力。必须项。由Analytic and Laplacian Adaptive Radioactivity Analysis软件:ALARA提供。
4. 辐照及冷却方案。必须项。需要在config.ini中指定。

在我们的工作文件夹下新建一个子文件夹用来存放R2S计算相关的文件。

mkdir r2s_run
cd r2s_run

初始化配置:setup

初始化配置r2s。执行：

natf_r2s.py setup

执行之后我们应该能够看到下面的信息：

# 一些warning
File "config.ini" has been written
File "alara_params.txt" has been written
Fill out the fields in these filse then run ">> r2s.py step1"

修改配置文件

setup这一步会生成两个文件：config.ini和alara_params.txt，我们需要修改这两个文件中的一些设置，包括辐照方案的设置和一些计算细节的设置。

修改config.ini设置

一般情况下，我们需要注意的几个设置是：

1. structured: 是否是结构化网格。结构化网格就是指MCNP中用的xyz或cyl网格。我们应该选True
2. sub_voxel: 是否采用sub_voxel模式，sub_voxel模式是我正在添加的代码，还没有完成，暂时不能选用，因此是False
3. reverse: 反转中子通量顺序。MCNP输出的中子通量是按照从低能到高能排序，而用ALARA进行活化计算时需要按照从高能到低能读取，因此需要把中子通量逆序排列。这里填True。
4. outp: MCNP outp文件
5. ptrac: MCNP ptrac文件
6. geom: DAGMC加载材料后的几何文件。如果使用MCNP几何文件，则此处空着。

[general]
# Specify whether this problem uses structured or unstructured mesh
structured: True
# Specify whether this problem uses sub-voxel r2s
sub_voxel: False

[step1]
# Path to MCNP MESHTAL file containing neutron fluxes or a DAG-MCNP5
# unstructured mesh tally .h5m file.
meshtal: meshtal
# Tally number within the meshtal file containing the fluxes for activation.
tally_num: 4
# The name of the tag used to store flux data on the mesh. For unstructured
# mesh this tag must already exist within the file specified in <meshtal>.
flux_tag: n_flux
# Path to the DAGMC material-laden geometry file (.h5m).
geom:
# If True the fluxes in the fluxin file will be printed in the reverse
# order of how they appear within the flux vector tag. Since MCNP and
# the Meshtal class order fluxes from low energy to high energy, this
# option should be true if the transmutation data being used is
# ordered from high-energy to low-energy.
reverse: True
# Number of rays to fire down each mesh row in each direction to calculate
# cell volume fractions.
num_rays: 10
# If true, rays will be fired down mesh rows in evenly spaced intervals.
# In this case <num_rays> must be a perfect square. If false, rays are fired
# down mesh rows in random intervals.
grid: False

[step2]
# List of decays times, seperated by commas. These strings much match exactly
# with their counterparts in the phtn_src file produced in step1. No spaces
# should appear in this line except the space between the time and the time unit
# for each entry.
decay_times:1 h
# The prefix of the .h5m files containing the source density distributations for
# each decay time.
output: source
# The name of the output files containing the total photon source intensities for
# each decay time
tot_phtn_src_intensities : total_photon_source_intensites.txt

修改alara_params.txt

material_lib alara_matlib
element_lib data/nuclib
data_library alaralib data/truncated_fendl2bin

cooling
    1 h
end

output zone
       integrate_energy
       # Energy group upper bounds. The lower bound is always zero.
       photon_source  data/truncated_fendl2bin  phtn_src 24 1.00E4 2.00E4 5.00E4 1.00E5
       2.00E5 3.00E5 4.00E5 6.00E5 8.00E5 1.00E6 1.22E6 1.44E6 1.66E6
       2.00E6 2.50E6 3.00E6 4.00E6 5.00E6 6.50E6 8.00E6 1.00E7 1.20E7
       1.40E7 2.00E7
end

#     flux name    fluxin file   norm   shift   unused
flux  my_flux     alara_fluxin  1e10     0      default

# Specify the irradiation schedule below.
# Syntax is found in the ALARA user manual
# This example is for a single 1 y pulse
schedule    my_schedule
    1 y my_flux my_pulse_history 0  s
end
pulsehistory  my_pulse_history
    1    0.0    s
end

#other parameters
truncation 1e-12
impurity 5e-6 1e-3
dump_file dump.file

r2s.py step1

r2s.py step1是用来根据r2s.py setup生成的配置文件生成ALARA输入文件，中子通量文件和ALARA材料库文件的。

当我们修改好上面两个文件后，就可以使用r2s.py生成ALARA计算需要的文件了：

r2s.py step1

材料数据库处理

ALARA输入文件中有三行是材料数据库相关的信息。

1. material_lib alara_matlib 这个alara_matlib是alara需要读取的材料成功定义的文件名。
2. element_lib data/nuclib 这个文件是定义元素的丰度用的。
3. data_library alaralib data/eaf2010bin

注意：生成的alara_matlib文件中，可能包含数据库中不识别的k:40核素，需要手动去除。

alara计算

当step1完成后，程序会给出运行alara的提示，按照提示运行：

alara alara_inp > output.txt

等待alara运行完成，会生成几个文件：

1. output.txt: ALARA输出文件，里面有活化后材料的核素组成以及光子发射率数据
2. phtn_src: ALARA输出文件，里面有详细的各个各个材料各个核素缠身的各个能群的光子数据

r2s step2

ALARA计算完成后，运行：

r2s.py step2

执行完成后，我们可以看到工作文件夹中多了几个文件：

1. phtn_src.h5: HDF5格式的光子源数据
2. source_i.h5m: 第i个冷却时间点的网格光子源数据。
3. total_photon_source_intensites.txt: 里面是各个冷却时间节点的总光子发射率。

光子输运

将MCNP输入文件修改为适合光子输运的文件：

• 去掉原来的SDEF卡
• 添加idum卡，定义光子抽样参数
• 添加停机剂量率计数卡

然后进行光子输运计算得到停机剂量率结果。

备注：在计算得到结果之后，还可以把几何图形及计算结果的网格数据进行可视化绘图。
绘图方法可以参考DAGMC数据绘图。