基于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