基于分段高斯積分的在線多溫度核截面生成方法研究

陳 銳 郝麗娟 宋 婧 吳 斌 何 鵬

基于分段高斯積分的在線多溫度核截面生成方法研究

陳 銳1,2郝麗娟1宋 婧1吳 斌1何 鵬1

1（中國科學院核能安全技術(shù)研究所 中國科學院中子輸運理論與輻射安全重點實驗室 合肥 230031）
2（中國科學技術(shù)大學 合肥 230027）

反應堆運行過程中溫度不斷變化，在模擬中常采用在線多普勒展寬方法生成各種溫度下的中子核截面。已有的在線截面生成方法中，SIGMA1方法精度較高，但由于其使用了誤差函數(shù)及泰勒級數(shù)展開方法，截面的生成效率偏低。本文基于FDS團隊自主研發(fā)的超級蒙特卡羅核模擬軟件系統(tǒng)SuperMC，針對不同能段截面的特點，發(fā)展了基于分段高斯積分的在線多溫度核截面生成方法，在多普勒展寬共振峰較密集的區(qū)域使用高斯-厄米特積分方法，在低能區(qū)域使用高斯-勒讓德積分方法，在保證核截面精度的同時提高了截面生成效率。通過典型核素截面的對比以及臨界安全基準例題與多普勒反應性系數(shù)基準例題的測試，本文方法與SIGMA1方法相比平均計算效率提高5倍以上，且展寬溫度越高，效率提升越明顯。證明了該方法能夠快速并準確地生成各種溫度下的中子核截面，可用于反應堆多物理耦合計算。

中子截面，高斯積分，在線，多普勒展寬，SuperMC

靶核的熱運動引起的多普勒效應會導致反應截面隨溫度發(fā)生變化，直接影響中子輸運計算的結(jié)果。由于反應堆在實際運行過程中，堆內(nèi)溫度不斷變化且變化范圍廣，為了準確模擬反應堆內(nèi)中子與材料的相互作用，需要事先制作并存儲大量不同溫度點的核截面數(shù)據(jù)[1]，導致粒子輸運計算時內(nèi)存占用極大。為此，蒙特卡羅輸運程序常采用在線(On-the-fly)核截面生成方法獲得所需溫度的核截面數(shù)據(jù)，即在輸運模擬中根據(jù)核素溫度直接計算中子的反應截面，避免了直接存儲大量核截面數(shù)據(jù)，減少了內(nèi)存占用。由于在線核截面生成方法需實時計算中子各溫度的反應截面，對方法的效率要求較高。

目前常用的在線多溫度核截面生成方法主要有在線擬合方法[2]、靶核運動抽樣方法[3]、近似多極點方法[4]及SIGMA1方法[5]。其中SIGMA1方法生成的核截面精度較高，且占用內(nèi)存少[6-7]。該方法使用解析積分形式求解多普勒展寬方程，可將某一溫度截面精確展寬至目標溫度。但由于該方法中使用誤差函數(shù)及泰勒級數(shù)展開方法，導致核截面的在線生成效率偏低。Dean等[8-9]對SIGMA1方法進行了優(yōu)化，提出在共振峰較密集能區(qū)使用高斯-厄米特(Gauss-Hermite)積分，提高了多普勒展寬共振區(qū)域的核截面在線生成效率，但低能區(qū)域仍使用SIGMA1方法，綜合計算效率仍較低。

本文基于超級蒙特卡羅核模擬軟件系統(tǒng)SuperMC開展了基于分段高斯積分的在線多溫度核截面生成方法研究，在多普勒展寬共振區(qū)域使用高斯-厄米特積分的基礎上，提出在低能區(qū)域使用高斯-勒讓德(Gauss-Legendre)積分，提高在線展寬的綜合計算效率。SuperMC是一套通用、智能、精準的核設計與安全評價軟件[10-13]，已在聚變堆[14-16]、聚變裂變混合堆[17-18]以及鉛冷快堆[19-21]等先進核能系統(tǒng)中得到了廣泛的使用。

1 SIGMA1多溫度核截面生成方法

在自由氣體模型假設下，一般認為靶核熱運動速度滿足麥克斯韋-玻爾茲曼(Maxwell-Boltzmann)分布，當前溫度T2的有效核反應截面可由原始溫度T1的核截面遵循以下形式推導得出：

式中：T1、T2為熱力學溫度；v為入射中子速度；v′為中子與靶核的相對速率；σT(ν)代表溫度T時中子速度為v的截面；β=M /[2k(T2-T1)]，M為靶核質(zhì)量，k為玻爾茲曼常數(shù)。

引入無量綱變量x2=βν′2，y2=βν2，簡化式(1)得：

SIGMA1方法采用解析積分形式計算式(2)的積分值。為保證計算精度，需使用誤差函數(shù)及泰勒級數(shù)展開方法，導致此方法效率較低。

2 基于分段高斯積分的多溫度核截面生成方法

本文在保證精度要求的前提下（與SIGMA1方法生成截面的相對偏差小于0.1%），針對SIGMA1方法效率問題，根據(jù)其不同能量區(qū)間多溫度截面的特征，發(fā)展了基于分段高斯積分的在線多溫度核截面生成方法。在多普勒展寬共振峰較密集區(qū)域使用高斯-厄米特積分方法，在低能區(qū)域使用高斯-勒讓德積分方法。在保證核截面精度的同時提高了截面生成效率。

2.1 高斯-厄米特積分方法

對于指定積分項f( z)e-z2，可使用高斯-厄米特積分公式：

其中：zk為厄米特多項式Hn的節(jié)點，定義如下：

權(quán)重wk定義為：

令式(3)中z=x-y，則式(3)可表示為：

式(7)還可分解為：

由于z的范圍為-4≤z＜4[5]，故當y＞4時即可保證式(8)中第二項求積結(jié)果近似為0。此時式(4)中的高斯-厄米特積分公式可直接用于式(8)。同時，考慮到y(tǒng)＞4時式(2)中的e-（x+y）2項貢獻很小，可忽略不計。因此對于y＞4的多普勒展寬共振峰較密集區(qū)域（入射能量E＞16kT/A），式(2)可以表示為：

式中：f( z)=σ(z+y)3。

T1

高斯-厄米特積分中積分節(jié)點的選取會對計算結(jié)果的正確性與計算效率產(chǎn)生影響。圖1中以典型核素238U為例，給出了不同的積分節(jié)點的高斯-厄米特積分方法與SIGMA1方法生成的核截面的最大相對偏差，并給出了高斯-厄米特方法的截面生成時間。從圖1可以看出，節(jié)點越多，計算結(jié)果越精確，但效率越低。當節(jié)點N≥16時，即可保證相對偏差ε＜0.1%。由于在線截面生成效率隨節(jié)點數(shù)目增加而降低，為保證效率，多普勒展寬共振峰較密集區(qū)域(E＞16kT/A)選取16節(jié)點的高斯-厄米特積分方法。

圖1 不同高斯-厄米特積分節(jié)點生成截面的最大相對偏差及截面生成時間Fig.1 Max relative error and broadening time with different Gauss-Hermite quadrature orders.

2.2 高斯-勒讓德積分方法

德積分方法：

此時式(2)可表示為：

k德多項式Pn(x)的節(jié)點，定義如下：

wk定義為：

高斯-勒讓德積分中積分節(jié)點的選取也會對計算結(jié)果的正確性與計算效率產(chǎn)生影響。同樣以典型核素238U為例，分析不同積分節(jié)點對高斯-勒讓德積分方法結(jié)果的正確性與計算效率的影響。圖2給出了高斯-勒讓德積分方法與SIGMA1方法生成的核截面數(shù)據(jù)的最大相對偏差，以及高斯-勒讓德方法的截面生成時間。從圖2可以看出，當節(jié)點N≥2時，高斯-勒讓德積分方法最大相對偏差ε＜0.01%。由于在線截面生成效率隨節(jié)點數(shù)目增加而降低，為保證效率，低能區(qū)域(E≤16kT/A)選取兩節(jié)點的高斯-勒讓德積分方法。

圖2 不同高斯-勒讓德積分節(jié)點生成截面的最大相對偏差及展寬時間Fig.2 Max relative error and broadening time with different Gauss-Legendre quadrature orders.

3 方法驗證

本文通過典型核素多溫度截面在線生成測試及基準例題測試，從核截面數(shù)據(jù)本身和輸運計算結(jié)果兩方面驗證本文方法的正確性及效率。

3.1 典型核素多溫度截面在線生成測試

本文基于典型核素(16O、56Fe、90Zr、232Th、234U、235U、238U及239Pu) 300 K溫度時的總截面，在線生成了600 K、900 K、1200 K溫度時的總截面。

圖3 不同溫度下本文方法與SIGMA1方法生成238U總截面的相對偏差Fig.3 Relative error of 238U total cross sections between the proposed method and SIGMA1 method at different temperatures.

圖3 給出了使用本文方法與SIGMA1方法生成的238U總截面的相對偏差。從圖3可以看出，當能量升高至共振能區(qū)時，受計算精度影響，兩種方法的相對偏差隨著溫度升高逐漸增大，但在可分辨共振區(qū)范圍內(nèi)(E＜0.02 MeV)始終在0.1%以內(nèi)。由于不可分辨共振能區(qū)的反應截面無法通過準確計算公式得到，在連續(xù)能量核數(shù)據(jù)庫中以概率表方式存儲。與SIGMA1處理方法類似，本文直接復制該能區(qū)的概率表信息。由圖3可見，當E≥0.02 MeV時，本文方法與SIGMA1方法的誤差為0。

圖4給出了使用本文方法與SIGMA1方法在線生成各種溫度下典型核素總截面的計算效率比，文中計算效率比是指相同溫度下分別使用本文方法與SIGMA1方法展寬對指定核素在相同能量框架內(nèi)展寬效率之比，與展寬時間比呈倒數(shù)關系。

從圖4可以看出，本文方法的計算效率與SIGMA1方法相比，輕核的在線截面生成效率提高約1.3倍，而重核的在線截面生成效率提高10倍以上，平均效率提高5倍以上。這是由于輕核（圖4中16O）截面出現(xiàn)的共振峰較少，能量框架較為稀疏，此時SIGMA1方法無需過多使用泰勒級數(shù)展開便可得到展寬結(jié)果。因此本文方法與SIGMA1方法效率較為接近，約1.3倍。對于共振核素（如圖4中238U及239Pu），因其截面出現(xiàn)的共振峰較多，SIGMA1方法為保證精度需頻繁調(diào)用泰勒級數(shù)展開，導致其效率下降，而本文方法只需計算不同節(jié)點的求積函數(shù)即可完成截面處理，因此加速效果可高達20倍。從圖4中還可以看出，展寬所需目標溫度越高，加速效果越明顯。

圖4 本文方法與SIGMA1方法生成典型核素截面的展寬效率比Fig.4 Ratio of broadening efficiency between the proposed method and SIGMA1 method for typical nuclides’ cross section.

3.2 基準例題測試

為了進一步驗證本文的方法，基于本文方法生成的截面進行了大量的輸運計算基準例題的測試。本文給出Godiva及反應性多普勒系數(shù)基準例題的驗證結(jié)果。

3.2.1 Godiva基準例題測試

Godiva模型[22]為一個金屬鈾燃料球形裸堆，如圖5所示，燃料成分為235U、238U及234U。

圖5 Godiva模型1/8剖面圖Fig.5 1/8 section of Godiva model.

測試中使用原始溫度為300 K的截面，分別由本文方法及SIGMA1方法在線展寬至600 K，計算Godiva模型的有效增殖因子及中子能譜。其中每代計算粒子數(shù)為20000，循環(huán)500代并舍去前50個非活躍代。本文方法有效增殖因子計算結(jié)果為0.99996±0.00019，SIGMA1方法有效增殖因子計算結(jié)果為0.99993±0.00019，能譜計算結(jié)果見圖6。

圖6 Godiva能譜對比Fig.6 Energy spectra comparison of Godiva.

從圖6結(jié)果可以得出，本文方法與SIGMA1方法所產(chǎn)生的截面計算得到的有效增殖因子相對偏差僅為0.003%；能譜的最大偏差小于0.5%。

3.2.2 反應性多普勒系數(shù)基準模型測試

反應性多普勒系數(shù)基準模型旨在驗證反應性多普勒系數(shù)計算結(jié)果的正確性，要求計算結(jié)果與基準值的相對偏差不超過10%[23]。該基準模型為一個簡單的柵元結(jié)構(gòu)，燃料為UO2，慢化劑為輕水，包殼材料為天然鋯。柵元軸向無限長，慢化劑外表面使用反射邊界條件。慢化劑的溫度固定為600 K，燃料區(qū)的溫度包含600 K和900 K，分別代表熱態(tài)零功率(Hot Zero Power, HZP)和熱態(tài)滿功率(Hot Full Power, HFP)兩個運行狀態(tài)。該基準模型中通過分別計算燃料溫度為600 K和900 K時的有效增殖因子，得到溫度引起的反應性多普勒系數(shù)Dc，其表達形式如下：

式中：ΔρDop=(kHFP-kHZP)/( kHFP×kHZP)；ΔTFuel=300 K。

本文使用原始溫度為300 K的核截面，在線展寬得到600 K和900 K的核截面，計算Dc。每代計算粒子數(shù)為20000，循環(huán)500代并舍去前50個非活躍代。Dc的計算結(jié)果如圖7所示。

使用本文方法計算得到的不同富集度下Dc均在基準值誤差范圍內(nèi)，且與基準值的最大相對偏差為3.77%，低于基準模型中要求的10%，表明計算結(jié)果滿足基準模型的精度要求。

圖7 反應性多普勒系數(shù)計算結(jié)果Fig.7 Results of Doppler defect coefficients.

4 結(jié)語

本文發(fā)展了基于分段高斯積分的在線多溫度截面生成方法，在多普勒展寬共振峰較密集區(qū)域使用高斯-厄米特積分方法，在低能區(qū)域使用高斯-勒讓德積分方法。為了驗證本文發(fā)展的方法，進行了典型核素截面的對比以及Godiva與多普勒反應性系數(shù)基準例題的測試。計算結(jié)果表明，本文方法能夠快速并準確地生成各種溫度下的中子截面，在保證精度的前提下，典型核素的在線生成效率平均提高了5倍以上，某些核素如238U和239Pu提速高達 20倍，驗證了本文方法的正確性以及高效性，可用于反應堆多物理耦合計算。

致謝 本文開展研究工作中得到了FDS團隊其他成員的大力幫助和支持，在此深表感謝！

1 Trumbull T. Treatment of nuclear data for transport problems containing detailed temperature distribution[R]. Niskayuna, NY: Knolls Atomic Power Laboratory (KAPL), 2006.

2 Yesilyurt G, Martin W, Brown F. On-the- fl y Doppler broadening for Monte Carlo codes[J]. Nuclear Science and Engineering, 2012, 171(3): 239-257. DOI: 10.13182/ NSE11-67.

3 Viitanen T, Leppanen J. Explicit treatment of thermal motion in continuous-energy Monte Carlo tracking routines[J]. Nuclear Science and Engineering, 2012, 171(2): 165-173. DOI: 10.13182/NSE11-36.

4 Forget B, Xu S, Smith K. Direct Doppler broadening in Monte Carlo simulations using the multipole representation[J]. Annals of Nuclear Energy, 2014, 64: 78-85. DOI: 10.1016/j.anucene.2013.09.043.

5 Cullen D. Exact Doppler-broadening of tabulated cross section[J]. Nuclear Science and Engineering, 1976, 60(3): 199-229. DOI: 10.13182/NSE76-1.

6 Macfarlane R. The NJOY nuclear data processing system[R]. LA-UR-12-27079, Los Alamos National Laboratory, 2012.

7 Cullen D. PREPRO 2015-2015 ENDF/B pre-processing codes[R]. IAEA-NDS-39, International Atomic Energy Agency, 2015.

8 Dean C, Perry R, Neal R, et al. Validation of run-time Doppler broadening in MONK with JEFF3.1[C]. International Conference on Nuclear Data for Science and Technology, Jeju Island, South Korea, 2010.

9 Romano P, Trumbull T. Comparison of algorithm for Doppler broadening pointwise tabulated cross sections[J]. Annals of Nuclear Energy, 2015, 75: 358-364. DOI: 10.1016/j.anucene.2014.08.046.

10 Wu Y C, Song J, Zheng H Q, et al. CAD-based Monte Carlo program for integrated simulation of nuclear system SuperMC[J]. Annals of Nuclear Energy, 2015, 82: 161-168. DOI: 10.1016/j.anucene.2014.08.058.

11 吳宜燦, 胡麗琴, 羅月童, 等. 先進核能系統(tǒng)設計分析軟件與數(shù)據(jù)庫研發(fā)進展[J]. 核科學與工程, 2010, 30(1): 55-64.

WU Yican, HU Liqin, LUO Yuetong, et al. Development of design and analysis software for advanced nuclear systems[J]. Chinese Journal of Nuclear Science and Engineering, 2010, 30(1): 55-64.

12 吳宜燦, 宋婧, 胡麗琴, 等. 超級蒙特卡羅核計算仿真軟件系統(tǒng)SuperMC[J]. 核科學與工程, 2016, 36(1): 62-71.

WU Yican, SONG Jing, HU Liqin, et al. Super Monte Carlo simulation program for nuclear and radiation process SuperMC[J]. Chinese Journal of Nuclear Science and Engineering, 2016, 36(1): 62-71.

13 Wu Y C, Xie Z S, Fischer U. A discrete ordinates nodal method for one-dimensional neutron transport calculation in curvilinear geometries[J]. Nuclear Science and Engineering, 1999, 133(3): 350-357.

14 Wu Y C, FDS Team. Conceptual design and testing strategy of a dual functional lithium-lead test blanket module in ITER and EAST[J]. Nuclear Fusion, 2007, 47(11): 1533-1539. DOI: 10.1088/0029-5515/47/11/015.

15 Wu Y C, FDS Team. Design analysis of the China dual-functional lithium lead (DFLL) test blanket module in ITER[J]. Fusion Engineering and Design, 2007, 82: 1893-1903. DOI: 10.1016/j.fusengdes.2007.08.012.

16 Wu Y C, FDS Team. Conceptual design of the China fusion power plant FDS-II[J]. Fusion Engineering and Design, 2008, 83: 1683-1689. DOI: 10.1016/j.fusengdes. 2008.06.048.

17 Wu Y C, Jiang J Q, Wang M H, et al. A fusion-driven subcritical system concept based on viable technologies[J]. Nuclear Fusion, 2011, 51(10): 103036. DOI: 10.1088/ 0029-5515/51/10/103036.

18 Wu Y C, Zheng S, Zhu X, et al. Conceptual design of the fusion-driven subcritical system FDS-I[J]. Fusion Engineering and Design B, 2006, 81: 1305-1311. DOI: 10.1016/j.fusengdes.2005.10.015.

19 Wang M H, Huang H, Lian C, et al. Conceptual design of lead cooled reactor for hydrogen production[J]. International Journal of Hydrogen Energy, 2015 40: 15127-15131. DOI: 10.1016/j.ijhydene.2015.04.009.

20 Huang Q Y, Li J, Chen Y. Study of irradiation effects in China low activation martensitic steel CLAM[J]. Journal of Nuclear Materials, 2004, 329-333: 268-272. DOI: 10.1016/j.jnucmat.2004.04.056.

21 Huang Q Y, Gao S, Zhu Z Q, et al. Progress in compatibility experiments on lithium-lead with candidate structural materials for fusion in China[J]. Fusion Engineering and Design, 2009, 84: 242-246. DOI: 10.1016/j.fusengdes.2008.12.038.

22 Briggs J, Michael A, Yolanda R. International handbook of evaluated criticality safety benchmark experiments[R]. NEA/NSC/DOC (95)03, Nuclear Energy Agency, 2006.

23 Mosteller R. ENDF/B-V, ENDF/B-VI, ENDF/B-VII.0 results for the Doppler-defect benchmark[R]. LA-UR-07-0922, Los Alamos National Laboratory, 2007.

On-the-fly Doppler broadening method based on piecewise Gauss quadrature for generating neutron cross section

CHEN Rui1,2HAO Lijuan1SONG Jing1WU Bin1HE Peng1
1(Key Laboratory of Neutronics and Radiation Safety, Institute of Nuclear Energy Safety Technology,
Chinese Academy of Sciences, Hefei 230031, China)
2(University of Science and Technology of China, Hefei 230027, China)

Background: As temperatures are constantly changing in nuclear reactor operation, on-the-fly Doppler broadening methods are commonly adopted for generating nuclear cross section for various temperatures in transport simulation. Among the existing methods, the SIGMA1 approach is exact but inefficient due to involving of error function and Taylor series expansion. Purpose: In this paper, we proposed a piecewise Gauss quadrature method for on-the-fly Doppler broadening based on SuperMC to improve efficiency with given accuracy. Methods: According to the cross section features of different energy regions, Gauss-Legendre quadrature was used in low-energy region, while the Gauss-Hermite quadrature used in high-energy region. Results: The comparison of typical nuclides cross section and benchmarking with Godiva and Doppler coefficient of reactivity were presented in detail, which indicated

CHEN Rui, male, born in 1981, graduated from East China Institute of Technology in 2009, doctoral student, focusing on nuclear data processing

HE Peng, E-mail: peng.he@fds.org.cn

that the new method could generate neutron cross section rapidly and precisely at various temperatures. Compared with SIGMA1 method, the proposed method improved the computing efficiency by an average of 5 times, and higher temperature could promote the broadening efficiency. Conclusion: It showed the method could be applied in multi-physical coupling calculation of reactor.

Neutron cross section, Gauss quadrature, On-the-fly, Doppler broadening, SuperMC

TL329.2

10.11889/j.0253-3219.2017.hjs.40.040503

No.11305203、No.11405204)資助

陳銳，男，1981年出生，2009年畢業(yè)于東華理工大學，現(xiàn)為博士研究生，研究領域為核數(shù)據(jù)處理

何鵬，E-mail: peng.he@fds.org.cn

2016-12-08，

2016-12-26

Supported by National Natural Science Foundation of China (No.11305203, No.11405204)

Received date: 2016-12-08, accepted date: 2016-12-26