基于CAD的点核剂量计算方法研究及初步应用

王 杰，石志勇，宋 婧，郝丽娟，胡丽琴,1，葛 鹏

(1. 中国科学技术大学，安徽合肥230027；2.中国科学院核能安全技术研究所，中国科学院中子输运理论与辐射安全重点实验室，安徽合肥230031；3. 火箭军装备研究院，北京市100094)

基于CAD的点核剂量计算方法研究及初步应用

王 杰1,2，石志勇3，宋 婧2，郝丽娟2，胡丽琴2,1，葛 鹏2

(1. 中国科学技术大学，安徽合肥230027；2.中国科学院核能安全技术研究所，中国科学院中子输运理论与辐射安全重点实验室，安徽合肥230031；3. 火箭军装备研究院，北京市100094)

点核积分方法是辐射剂量计算的基本方法之一，广泛应用于辐射防护领域。传统的点核剂量计算中采用文本方式描述计算模型，存在难以描述复杂几何、易出错且耗时的问题。针对该问题，本文基于FDS团队自主研发的超级蒙卡核模拟软件系统SuperMC，进行了基于CAD的点核剂量计算方法研究与程序开发，可基于实际问题的CAD模型直接进行支持多重源的光子点核剂量计算，提高了程序对复杂三维几何问题的处理能力，并包含较为完备的核数据库。使用ANSI/ANS6.6.1、ESIS和VisiPlan的基准例题对程序进行了测试验证，测试结果与VisiPlan4.0对比吻合良好。同时将该方法初步应用于ITER热室屏蔽的设计中，说明了本方法及程序处理复杂场景问题的能力。

剂量率；点核；CAD；ITER；SuperMC

点核积分法是一种格林函数积分法，适用于计算和处理γ射线在复杂几何空间的辐射屏蔽问题[1]。作为辐射剂量计算的重要方法之一，点核积分法不仅不受空间尺寸和屏蔽体厚度的限制，而且计算速度快，具有较高的计算效率。该方法采用线性衰减方法计算γ射线在物质中的衰减过程，考虑到散射通量的影响，引入积累因子来计算空间某一点散射光子所造成的影响的积累或增加[2]。

能量为E(MeV)，源强为1photon/sec的γ点源在空间R处产生的剂量率为:

(1)

式中:F(E)为通量剂量转换系数;B(E，Ld)为剂量积累因子;R为γ点源到剂量计算点处的直线距离;Ld(E)为光学距离：

(2)

式中：n是γ射线穿过的材料种类数；μi为第i种材料的线衰减系数；ti为射线在第i种材料中的穿透距离，光学距离的单位以平均自由程(mfp)表示。所谓的点核积分方法，就是在辐射屏蔽几何空间中计算以下积分：

D(r)=∫E∭VS·D(E,Ld,r)dEdV

(3)

式中：V为辐射源所在的几何区域；E为辐射源的能谱分布；S为辐射源的强度。由于几何区域和辐射源的复杂性，(3)式中的积分一般不能用解析方法求出，而是通过对辐射源进行空间和能谱离散使其成为点源，从而实现积分。

目前国内外广泛应用于屏蔽设计的点核积分程序有QAD[3]、PUTZ、VisiPlan[4]等。然而，随着先进反应堆[5-7]的发展，核设施的几何结构更加复杂，空间尺度也更大。QAD、PUTZ等一些采用文本形式描述计算模型的点核积分程序，难以方便地构建复杂的三维几何模型；而VisiPlan、Microshield等程序虽然增加了3D可视化建模功能，但仍需要对现实场景进行很大程度的简化才能计算，并且这些程序所使用的核数据库并不完备，缺少一些常用的核数据，不能很好地满足屏蔽设计的需求。为此，本文基于中国科学院核能安全技术研究所·FDS团队自主研发的超级蒙卡核模拟软件系统SuperMC[8-11]，进行了基于CAD点核剂量计算方法的研究和程序开发。

本文首先介绍了基于CAD点核剂量计算的关键方法，其次通过ANSI/ANS-6.6.1[12]、ESIS[13]和VisiPlan中的基准例题对程序进行了测试，最后将开发的程序应用于ITER热室屏蔽设计中。

1 基于CAD的点核剂量计算关键方法

1.1 CAD几何建模

由(2)式可以看出，光学距离Ld的计算需要预先计算出γ射线穿透的材料种类和距离。因此在使用点核积分法进行剂量计算之前，需要以构造实体几何法(CSG)形成相应的几何空间，采用射线追踪技术计算γ射线在几何空间中的穿透距离，然而用文本方式描述CSG模型存在着抽象易错并且耗时的问题。针对该问题，本文在SuperMC已有版本基础上，基于实际问题的CAD模型直接进行点核剂量计算。

基于边界表示法(BREP)描述的CAD几何模型可以在SuperMC中自动转换为基于构造实体几何法(CSG)表述的剂量计算模型[14]，算法的实质是把复杂模型分解为简单形体，然后基于基本的面(半空间)描述简单形体，用简单形体的布尔组合来完成对复杂形体的描述[15-16]。每个几何区域都被赋予了材料信息，并转换到CSG模型中用于点核剂量计算。SuperMC中采用层次树结构来描述几何和材料信息(组成成分、密度等)，方便了在射线追踪过程中对几何体的查找，层次树的根节点作为世界体，包含所有的几何区域，避免了出现某些区域未定义的情况。借助于CAD建模工具(CATIA、AutoCAD等)可以描述任意几何，从而能够方便的构建复杂的三维几何模型[17]。

1.2 辐射源

点核积分方法中对于辐射源的描述参数主要分有能谱、几何和源强。本文将辐射源的能谱离散为25个能群，能量范围为15keV～15MeV。每个放射源的能谱信息包含不同的能量值及该能量光子在所有出射光子中所占的比率。对于给定的能谱分布，(1-3)中对于能量的积分转换为：

(4)

式中：i为能群标号；Ei为第i个能群的能量值；FEi为Ei能量的光子在所有出射光子中所占的比率。

程序中支持的辐射源类型包括点源、线源(直线、圆圈)、面源(球面、圆柱面、圆盘、矩形面)和体源(球、圆柱、立方体等)。计算辐射源的体积分时，采用了基于统计误差的蒙卡积分方法：在抽样的过程中同时计算统计误差，并判断是否达到指定的值，从而自适应地确定蒙卡积分所需要的抽样点数，从而将线源、面源和体源离散为N个点源。对于均匀分布的源，在抽样点数为N的情况下，式(5)中的积分变为：

(5)

统计误差：

(6)

对于多重源的情况，则依次计算每个源产生的影响，所有辐射源产生的剂量总和作为该点的剂量率。

1.3 数据库

从点核积分公式(1)和公式(2)可以看出，点核剂量计算中主要需要三个参数：材料线衰减系数μ、剂量积累因子B和通量剂量转换系数F。线性衰减系数μ可表示为μm/ρ，μm为质量衰减系数。因此，进行点核剂量计算需要μm、B和F三个数据库。程序中内置了15keV～15MeV范围内25个能群的光子通量剂量转换系数F。

1.3.1 质量衰减系数

质量衰减系数表征了射线在材料中与物质的相互作用，对于γ射线，主要为光电效应、康普顿散射和电子对效应。本文在丹麦工业大学开发的WinXCom[18]程序和美国国家标准ANSI/ANS-6.4.3[19]的基础上设计制作了衰减系数数据库，包含了1-100号元素和水、空气、标准混凝土(NBS混凝土)的质量衰减系数。对于水、空气、NBS混凝土之外的混合物和化合物的质量衰减系数，本文采用下式计算：

(7)

式中：ωi为元素i的质量分数；μmi为元素i的质量衰减系数。程序将根据材料的组成成分，自动计算该材料的质量衰减系数。

1.3.2 积累因子

积累因子是进行γ射线屏蔽分析时所使用的一切点核公式的基础。常用的积累因子分为剂量积累因子和能量吸收积累因子，剂量积累因子又称为照射量积累因子，是剂量计算点实际剂量与未经碰撞粒子产生的剂量的比值[2]。常用的积累因子计算公式有G-P公式、泰勒公式、凯波多项式等。本文采用泰勒公式计算γ射线的剂量积累因子:

B(E,Ld)=A1e-α1Ld+(1-A1)e-α2Ld

(8)

式中参数A1、α1、α2是通过在矩方法解得的各向同性点源在无线均匀介质中输运结果的基础上拟合得到的，与材料和射线的能量相关。本文所使用的积累因子数据库包含了15种参考材料的泰勒公式参数，数据来源于ANSI/ANS-6.4.3[19]。对于不同的剂量计算点，程序支持选取不同的积累因子，一次计算中可选择多种积累因子。

2 基准验证

为了验证程序的正确性，本文采用ANSI[12]、ESIS[13]和VisiPlan的相关例题对程序进行了测试，测试结果与VisiPlan4.0[4]进行了对比。

2.1 ANSI/ANS例题测试

本文采用了ANSI/ANS 6.6.1[12]中的点源例题I.1、I.2和体源例题II.1、II.2对程序进行了测试，测试结果如图1所示。点源例题计算结果与VisiPlan4.0一致，体源例题所有计算点的计算结果与VisiPlan相比，偏差均在4.0%以内。

图1 ANSI/ANS例题测试结果Fig.1 The test results of ANSI/ANS problem

2.2 ESIS例题测试

本文采用ESIS发布的屏蔽例题对程序进行测试验证，其CAD模型如图2所示，详见文献[13]。测试结果如表1所示，和VisiPlan的相对偏差分别为0.87%和1.30%。

图2 ESIS例题的CAD模型Fig.2 The CAD model of ESIS problem表1 ESIS测试结果Table 1 The test results of ESIS

VisiPlan/(mSv/h)SuperMC/(mSv/h)deviationD149102×102486716×102-087%D283205×10-3842903×10-3130%

2.3 VisiPlan应用例题测试

本文采用了VisiPlan中的蒸汽发生器例题进行了测试，如图3。一回路和二回路中充满水，一回路传热管束简化为由铁和水的混合物组成的圆柱体，放射源在其中均匀分布。部分计算结果如图4所示，与VisiPlan相比，有两个点的相对偏差达到8%，其他均保持在5%以内。

图3 蒸汽发生器CAD模型Fig.3 The CAD model of steam generator

图4 蒸汽发生器例题计算结果Fig.4 The test results of steam generator

通过ANSI、ESIS和VisiPlan中基准例题的测试，本文计算结果与VisiPlan相比最大偏差为8%，出现偏差的主要原因是由于两个程序采用了不同的质量衰减系数数据库，并且对积分点的抽取方式不同而导致的。

3 初步应用：ITER热室屏蔽设计

为了说明程序对复杂场景问题的处理能力，本文对ITER热室屏蔽问题进行了初步分析，给出了热室中门7的设计方案。

ITER热室建筑位于托卡马克大厅[20]北边，如图5(a)所示。热室中存在7个混凝土(2.6g/cm3)门，位置如图5(b)所示，当4个被活化的第一壁组件被运输到热室中后，要求门2、3、4、7后的剂量应当低于10μSv/h，门5、6后的剂量应在1000μSv/h以下[21]。为此，需要计算出各个门的最佳厚度。

本文仅以门7为例，采用开发的程序对门7的厚度进行初步设计。为了简化计算过程，把每个第一壁组件等效为一个平板体源和圆柱体源的组合体。采用等效后的辐射源，可以计算出在门7前的γ剂量率为0.085Sv/h，与采用原始模型MCNP计算的结果(0.08Sv/h)符合很好，因此这样的简化是可行的。

图5 ITER HCB几何模型Fig.5 The geometry of ITER HCB (a) 热室CAD模型;(b) 热室中7个门的位置

在等效源的基础上，本文计算了在不同厚度下，门7后的剂量率。随着厚度的增加，γ剂量率的衰减如图6所示。

图6 不同厚度的门7对γ剂量率的衰减Fig.6 The attenuation of dose rate with different thickness

混凝土门的厚度X与衰减后的γ剂量率D可用下式表示：

(9)

b为常数，表征γ射线在混凝土中的衰减效果。根据图6，可以计算出b≈0.0585。

因此，为了使门7后的剂量率小于10μSv/h，门7的厚度至少为67cm。表1给出了在不同安全系数下，门7所需的厚度以及使用MCNP设计出的厚度。

表1 不同安全系数所要求的门7厚度Table 1 Minimum thickness of door 7

从表中看出，本文计算结果比MCNP结果大了5～7cm，这是因为点核积分法使用的是无限介质积累因子，包含了在边界处的反射贡献，导致计算结果偏高，因此本文得出的是一个偏保守的设计方案，能够确保实际的剂量率满足剂量限制。在计算时间上，使用MCNP进行一次计算需要数周的时间，而采用本文所提出的方法仅需要几个小时，大大缩减了计算时间。

4 结论

本文基于实际问题的CAD模型直接进行点核剂量计算，提高了程序对复杂三维几何问题的处理能力，可进行支持多重源的光子剂量计算，并包含较为完备的核数据库。采用ANSI/ANS-6.6.1、ESIS和VisiPlan中的基准例题对程序进行了测试验证，计算结果与VisiPlan4.0的计算结果吻合很好，证明了本文方法和程序的正确性。将开发的程序初步应用于ITER热室屏蔽问题中，给出了对于门7的初步设计方案，显示程序处理复杂场景问题的能力。在后续的研究中，将会结合实际反应堆应用场景开展应用。

致谢

本工作得到中科院核能安全技术研究所·FDS团队各位老师的帮助与指导，在此深表感谢。

[1] 万海霞. 辐射场可视化平台剂量分布快速计算程序开发[D]. 华北电力大学，2012.

[2] Schaeffer N M，Reactor Shielding for Nuclear Engineer[R]. Radiation Research Associates，Inc. Fort Worth，Tex. (USA),1973.

[3] Cain V R. QAD-CG: A Combinatorial Geometry Version of QAD-P5A Point-Kernel Code for Neutron and Gamma-ray Shielding Calculation[Z]. 1977，ORNL-CCC307.

[4] Vermeersch F，Van Bosstraeten C. Development of the VisiPlan ALARA planning tool[J]. Proceeding of the International Conference on Topical issues in Nuclear Radiation and Radioactive Waste Safety，Vienna Austria，August 31 to September 4，1998.

[5] Y. Wu，H. chen，S. Liu，et al. Fusion-Based Hydrogen Production Reactor and Its Material Selection[J]. Journal of Nuclear Materials，2009，386-388:122-126.

[6] Y. Wu，S. Zheng，X. Zhu，et al. Conceptual Design of the Fusion-driven Subcritical System FDS-I[J]. Fusion Engineering and Design，2006，81(8): 1305-1311.

[7] Y. Wu，FDS Team. Conceptual Design of the China Fusion Power Plant FDS-II[J]. Fusion Engineering and Design，2008，83: 1683-1689.

[8] J. Song，G. Sun，Z. Chen，et al. Benchmarking of CAD-based SuperMC with ITER Benchmark Model[J]. Fusion Engineering and Design，DOI: 10.1016/j. fusengdes.2014.05.003.

[9] 孙光耀，宋婧，郑华庆，等. 超级蒙特卡罗计算软件SuperMC2.0中子输运计算校验[J]. 原子能科学与技术，2013，第47卷增刊: 140-145.

[10] Y. Wu，J. Song and FDS Team. Development of Super Monte Carlo Calculation Program SuperMC 2.0. Proc. Int. Conf. ANS National Meeting-2013 ANS Winter Meeting and Technology Expo，American Nuclear Society，Washington D. C.，USA，November 10-14，2013.

[11] 吴宜灿，胡丽琴，龙鹏程，等. 先进核能系统设计分析软件与数据库研发进展[J]. 核科学与工程，2010，30(1): 55-64.

[12] Calculation and Measurement of Direct and Scattered Gamma Radiation from LWR Nuclear Power Plants，ANSI/ANS-6.6.1-1991，Amer. Nat. Standards Int.，New York，1991.

[13] Specification for Gamma Ray Shield Benchmark Applicable to A Nuclear Radioactive Waste Facility. ESIS Newsletter 5 (April 1973).

[14] Y. Wu，J. Song，H. Zheng，et al. CAD-Based Monte Carlo Program for Integrated Simulation of Nuclear System SuperMC[J]. Annals of Nuclear Energy，2015，82:162-168.

[15] 吴宜灿，李莹，卢磊 等. 蒙特卡罗粒子输运计算自动建模程序系统的研究与发展[J]. 核科学与工程，2006，26(1): 20-27.

[16] Y. Wu，FDS Team. CAD-based Interface Programs for Fusion Neutron Transport Simulation[J]. Fusion Engineering and Design，2009，84(7):1987-1992.

[17] 吴宜灿，李静惊，李莹 等. 大型集成多功能中子学计算与分析系统VisualBUS的研究与发展[J]. 核科学与工程，2007，27(4): 365-373.

[18] Gerward L，Guilber N，Jensen K B，et al. WinXCom-A Program for Calculation X-Ray Attenuation Cofficients[J]. Radiation Physics and Chemistry，2004，7(3): 653-654.

[19] Gamma-Ray Attenuation Cofficient and Buildup Factor for Engineering Materials，ANSI/ANS-6.4.3-1991，Amer. Nat. Standards Int.，New York，1991.

[20] Y. Li，L. Lu，A. Ding，et al. Benchmarking of MCAM 4.0 with the ITER 3D Model[J].Fusion Engineering and Design，2007，82(15): 2861-2866.

[21] Q. Yang. Dose Rates in Port Cells and Shielding of Port Cell Door，Nov. 2013，ITER_D_ LY6K73 v1.0.

Research and Preliminary Application of CAD-Based Point-KernelIntegral Method for Dose Calculation

WANG Jie1，2，SHI Zhi-yong3，SONG Jing2，HAO Li-juan2，HU Li-qin2，1，GE Peng2

(1. University of Science and Technology of China，Hefei，Anhui，230027，China;2. Key Laboratory of Neutronics and Radiation Safety，Institute of Nuclear Energy Safety Technology，Chinese Academy of Sciences，Hefei，Anhui，230031，China;3. EquipmentAcademy of the Rocket Force)

Point-kernel integral method is widely used in radiation shielding field and taken as a basic method for radiation dose calculation. However，for traditional point-kernel codes，describing the geometry of problems using text files is difficult to describe complex geometry and error-prone. Based on Super Monte Carlo Simulation Program for Nuclear and Radiation Process SuperMC，developed by FDS Team，a CAD-based point-kernel integral method was proposed in this study and implemented. The CAD model of problems can be directly used for calculation，which improves the ability of handing complex geometry. Gamma dose calculation can be performed with multi-source and complete database based on point-kernel method. Benchmark cases from ANSI/ANS-6.6.1，ESIS and VisiPlan were used to validate the code. The results were in consistent with those of VisiPlan4.0. Using SuperMC，the preliminary design of shielding design of ITER hot cell building was performed，which demonstrated the ability of this method for handling complex problems.

Dose Rate; Point-Kernel; CAD; ITER; SuperMC

2016-07-21

国家自然科学基金91026004；国家ITER 973计划2014GB112001

王 杰(1991—)，男，安徽人，硕士研究生，主要从事点核剂量计算方法研究工作

葛 鹏:peng.ge@fds.org.cn

329+.2

A

0258-0918(2016)05-0656-07