蒙特卡羅粒子輸運(yùn)計(jì)算網(wǎng)格計(jì)數(shù)方法研究

張 寬,宋 婧,陳珍平,孫光耀,郝麗娟,龍鵬程,胡麗琴?,FDS團(tuán)隊(duì)

(1.中國科學(xué)技術(shù)大學(xué),安徽 合肥230027;2.中國科學(xué)院核能安全技術(shù)研究所,中國科學(xué)院中子輸運(yùn)理論與輻射安全重點(diǎn)實(shí)驗(yàn)室,安徽 合肥230031)

蒙特卡羅粒子輸運(yùn)計(jì)算網(wǎng)格計(jì)數(shù)方法研究

張 寬1,2,宋 婧1,2,陳珍平1,2,孫光耀1,2,郝麗娟2,龍鵬程2,胡麗琴1,2?,FDS團(tuán)隊(duì)

(1.中國科學(xué)技術(shù)大學(xué),安徽 合肥230027;2.中國科學(xué)院核能安全技術(shù)研究所,中國科學(xué)院中子輸運(yùn)理論與輻射安全重點(diǎn)實(shí)驗(yàn)室,安徽 合肥230031)

在進(jìn)行反應(yīng)堆數(shù)值模擬時(shí),通過網(wǎng)格計(jì)數(shù)方法可以精細(xì)地統(tǒng)計(jì)整個(gè)堆芯中多種物理量的空間分布情況。本文基于超級蒙特卡羅核計(jì)算仿真軟件系統(tǒng)Super MC發(fā)展了一種反應(yīng)堆靜態(tài)物理參數(shù)的網(wǎng)格計(jì)數(shù)方法。通過日本原子能研究所臨界實(shí)驗(yàn)裝置TCA例題進(jìn)行了數(shù)值驗(yàn)證,計(jì)算結(jié)果與MCNP計(jì)算結(jié)果吻合較好,表明了本文方法的可行性與正確性。

蒙特卡羅;網(wǎng)格計(jì)數(shù);靜態(tài)參數(shù);Super MC

通量、反應(yīng)率(尤其是裂變率)是反應(yīng)堆物理分析中常用到的物理量[1],堆芯功率分布分析是核反應(yīng)堆安全運(yùn)行的基礎(chǔ)[2],是反應(yīng)堆設(shè)計(jì)、安全分析、故障診斷、實(shí)時(shí)控制等問題的基礎(chǔ)。在實(shí)際問題中,堆芯功率使用更方便,同時(shí)因堆芯中子通量分布與堆芯功率分布基本成比例,所以實(shí)際應(yīng)用中多分析堆芯功率分布[3]。分析反應(yīng)堆堆芯等重要區(qū)域的功率分布或通量、反應(yīng)率分布時(shí)需要進(jìn)行詳細(xì)的計(jì)數(shù)統(tǒng)計(jì)以獲得相應(yīng)物理量精細(xì)的空間分布。對于反應(yīng)堆用的蒙特卡羅程序而言,基于傳統(tǒng)計(jì)數(shù)方法的網(wǎng)格計(jì)數(shù)功能是一項(xiàng)重要功能,主要針對通量、反應(yīng)率、功率等靜態(tài)參數(shù)。例如:MCNP[4],Serpent[5]、MC21[6]等反應(yīng)堆蒙特卡羅計(jì)算程序也具有基本的網(wǎng)格計(jì)數(shù)功能。

本文在傳統(tǒng)的計(jì)數(shù)方法基礎(chǔ)上實(shí)現(xiàn)并改進(jìn)了通量、反應(yīng)率、功率等靜態(tài)參數(shù)的網(wǎng)格計(jì)數(shù)方法,在中國科學(xué)院核能安全技術(shù)研究所·FDS團(tuán)隊(duì)研發(fā)的超級蒙特卡羅計(jì)算仿真軟件系統(tǒng)Super MC[7]中進(jìn)行應(yīng)用實(shí)現(xiàn)并通過基準(zhǔn)例題校驗(yàn)。

1 基于徑跡長度法的靜態(tài)參數(shù)估計(jì)

其中:T是關(guān)心體積中所有粒子徑跡序列,li是i軌跡的長度,W是粒子的初始權(quán)重總和。

反應(yīng)率的計(jì)算一般是基于通量估計(jì)的基礎(chǔ)上實(shí)現(xiàn)的,反應(yīng)率的徑跡長度估計(jì)也可以通過公式(1)乘以微觀截面得到:

基于徑跡長度法的通量估計(jì)為:

堆芯任一點(diǎn)r處單位體積內(nèi)的功率,即r處的功率密度為:

2 網(wǎng)格計(jì)數(shù)方法

2.1 網(wǎng)格計(jì)數(shù)流程

本文設(shè)計(jì)實(shí)現(xiàn)的網(wǎng)格計(jì)數(shù)流程如圖1所示。首先根據(jù)輸入標(biāo)識判斷是直角坐標(biāo)系還是圓柱坐標(biāo)系,然后進(jìn)入相應(yīng)的計(jì)數(shù)模塊。對于圓柱坐標(biāo)系的計(jì)數(shù),首先要進(jìn)行坐標(biāo)轉(zhuǎn)換,后續(xù)步驟這兩種情況的計(jì)數(shù)方法一致:首先判斷該步的粒子類型與要計(jì)數(shù)的粒子類型是否一致,若一致則進(jìn)行計(jì)數(shù)運(yùn)算;然后判斷粒子當(dāng)前步是否在劃分網(wǎng)格的區(qū)域中,若在,則查找各個(gè)方向所在的小網(wǎng)格位置;下一步統(tǒng)計(jì)在此小網(wǎng)格中粒子徑跡長度并檢測是否需要計(jì)算反應(yīng)率或功率,若需要調(diào)用相應(yīng)計(jì)數(shù)模塊;最后把粒子位置移動(dòng)到此網(wǎng)格的邊界,進(jìn)行下一網(wǎng)格的統(tǒng)計(jì)。

圖1 網(wǎng)格計(jì)數(shù)流程圖Fig.1 Flow Chart of Mesh Tally

2.2 坐標(biāo)變換

圓柱坐標(biāo)系下,由于半徑方向不是直線邊界,所以圓柱坐標(biāo)系的計(jì)數(shù)與直角坐標(biāo)系有所不同,需要將粒子坐標(biāo)從世界坐標(biāo)系到圓柱坐標(biāo)系進(jìn)行轉(zhuǎn)換。此轉(zhuǎn)換需要一個(gè)輔助直角坐標(biāo)系P c s c＝x,其中s c＝ [s,t,u],寫成矩陣形式為:

式中:x＝[x o,y o,z o]是圓柱網(wǎng)格底面的中心;a＝[a1,a2,a3]是圓柱的軸向量;v＝ [v1,v2,v3]是θ＝0的參考方向;d＝ [d1,d2,d3]＝a×v。

圓柱坐標(biāo)系曲面方程可以表示為:

或者表示為:＝0

則直角坐標(biāo)系到圓柱坐標(biāo)系的轉(zhuǎn)換為:

得到轉(zhuǎn)換矩陣為[4]:

2.3 計(jì)數(shù)網(wǎng)格跨不同材料區(qū)域處理

蒙特卡羅粒子輸運(yùn)程序中網(wǎng)格功率計(jì)算一般基于公式(3)和上述網(wǎng)格計(jì)算方法進(jìn)行統(tǒng)計(jì),但該方法具有局限性:當(dāng)某一網(wǎng)格跨越不同材料的柵元且粒子的抽樣步長正好跨越不同材料區(qū)域時(shí),此步長功率計(jì)算并未準(zhǔn)確獲取每段徑跡的材料,而使用該步中起點(diǎn)所在柵元的材料作為整個(gè)步長計(jì)算的材料,與真實(shí)情況相比會(huì)造成一定的誤差。

為此本文進(jìn)行了功率計(jì)算功能的優(yōu)化:首先在每步輸運(yùn)中記錄該步長是否跨越不同材料的柵元,若不跨越則可用一般方法進(jìn)行處理;若跨越,則記錄穿越不同材料柵元的交點(diǎn)及相關(guān)材料信息,然后將該步長根據(jù)交點(diǎn)分成若干子步長,每個(gè)子步長取各自起始點(diǎn)所在柵元的材料進(jìn)行功率計(jì)算,然后進(jìn)行累加。具體流程見圖2。

圖2 功率計(jì)算跨網(wǎng)格處理流程圖Fig.2 Flow Chart of Crossing Mesh Processing for Power Calculation

3 數(shù)值驗(yàn)證

3.1 例題描述與計(jì)算過程

本文選用日本原子能研究所臨界實(shí)驗(yàn)裝置TCA(Tank-type Critical Assembly)[8]例題進(jìn)行基準(zhǔn)校驗(yàn)。TCA是日本原子能研究所的一個(gè)壓水堆臨界實(shí)驗(yàn)裝置,其燃料棒有兩種:MOX燃料棒和UO2燃料棒。MOX燃料棒中PuO2的富集度為4.91%,鈾為天然豐度,MOX燃料芯塊半徑為0.858 cm,燃料包殼內(nèi)半徑為:0.872 cm,外半徑為0.935 cm,活性區(qū)長度為90.93±0.5 cm;UO2燃料棒中235U的富集度為2.596%,芯塊直徑為1.25 cm,燃料包殼內(nèi)徑1.265 cm,外徑為1.341 cm,活性區(qū)長度144.15±0.3 cm,其燃料棒布置采用均勻擺放方式。

在Super MC中基于組件模型陣列自動(dòng)建模功能進(jìn)行TCA裝置建模[9-13],如圖3所示,并實(shí)現(xiàn)完整的材料、源項(xiàng)、計(jì)數(shù)建模,在Super MC中對模型全區(qū)域進(jìn)行網(wǎng)格劃分:徑向劃分10個(gè),軸向劃分10個(gè),角度方向劃分18個(gè)。本測試中均使用FDS團(tuán)隊(duì)自主研發(fā)的HENDL3.0[14-16]數(shù)據(jù)庫,使用網(wǎng)格計(jì)數(shù)方法計(jì)算通量、反應(yīng)率、功率的分布情況,計(jì)算結(jié)果與MCNP進(jìn)行對比。

圖3 TCA裝置模型圖Fig.3 TCA Criticality Facility Model(a)俯視圖;(b)剖面圖

3.2 測試結(jié)果及分析

對每個(gè)網(wǎng)格的通量、反應(yīng)率、功率計(jì)算結(jié)果進(jìn)行統(tǒng)計(jì)分析,結(jié)果如表1所示。由結(jié)果可知:所有網(wǎng)格計(jì)數(shù)結(jié)果與MCNP相對偏差都在1%以下,符合工程設(shè)計(jì)要求,同時(shí)至少有92.05%的網(wǎng)格結(jié)果在一個(gè)σ之內(nèi),符合統(tǒng)計(jì)學(xué)的要求,計(jì)算結(jié)果正確。使用Super MC的可視化模塊RVIS[15-20]對網(wǎng)格通量計(jì)算結(jié)果進(jìn)行可視化與對比,結(jié)果如圖4所示,由圖可見Super MC計(jì)算結(jié)果與MCNP結(jié)果變化趨勢一致。

表1 靜態(tài)參數(shù)測試結(jié)果Table 1 Static Parameter Test Results

圖4 網(wǎng)格通量結(jié)果可視化Fig.4 Visualization of Mesh Tally of Flux(a)MCNP俯視圖;(b)Super MC俯視圖(c)MCNP側(cè)視圖;(d)Super MC側(cè)視圖

4 結(jié)論

本文設(shè)計(jì)實(shí)現(xiàn)了蒙特卡羅粒子輸運(yùn)模擬中反應(yīng)堆靜態(tài)參數(shù)的網(wǎng)格計(jì)數(shù)方法,并考慮了計(jì)數(shù)網(wǎng)格跨不同材料區(qū)域的問題,通過日本原子能研究所臨界實(shí)驗(yàn)裝置TCA例題的測試結(jié)果表明,本文方法對網(wǎng)格通量、反應(yīng)率、功率的計(jì)數(shù)結(jié)果與MCNP吻合良好,表明了本文方法的可行性與正確性。

[1] 謝仲生,吳宏春等.核反應(yīng)堆物理分析[M]．西安:西安交通大學(xué)出版社.2004.

[2] 李樹,田東風(fēng),鄧力.中子時(shí)間常數(shù)的Monte Carlo計(jì)算方法[J]．北京:清華大學(xué)學(xué)報(bào)(自然科學(xué)版),2007:1057-1061.

[3] 李偉.基于堆外計(jì)數(shù)的堆芯功率分布重構(gòu)方法研究[D]．哈爾濱:哈爾濱工程大學(xué),2009:1-3.

[4] John S.Hendricks.Superimposed Mesh Plotting in MCNP[R]．LA-UR-01-1033,2001:3-5.

[5] Jaakko Lepp?nen.Serpent-a continuous energy Monte Carlo Reactor Physics Burnup Calculation Code User's Manual.March 6,2013:131-133.

[6] T.M.Sutton,T.J.Donovan et al The MC21 Monte Carlo Transport Code[J]．M&C＋SNA 2007:4.

[7] 曾勤,李瑩,盧磊 等.蒙特卡羅粒子輸運(yùn)計(jì)算自動(dòng)建模程序MCAM在ITER核分析建模中的應(yīng)用[J]．原子核物理評論,2006,23(2):138-141.

[8] Ma fi zur Rahman,T.Suzaki,et al.Validation study of the Monte Carlo code MVP for analysis of two-region TCA critical experiments with PWR-type MOX fuels[J]．Progress in Nuclear Energy,48(2006)703-726.

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

[10] 吳宜燦,李瑩,盧磊等.蒙特卡羅粒子輸運(yùn)計(jì)算自動(dòng)建模程序系統(tǒng)的研究與發(fā)展.核科學(xué)與工程,2006,26(1):20-27.

[11] 王國忠,黨同強(qiáng),熊健 等.MCAM4.8在ITER建筑大廳中子學(xué)建模中的應(yīng)用[J]．核科學(xué)與工程,2011,31(4):352-353.

[12] 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.

[13] Guozhong Wang,Jian Xiong,Pengcheng Long,et al.Progress and Applications of MCAM:Monte Carlo Automatic Modeling Program for Particle Transport Simulation[J]．Progress in Nuclear Science and Technology,2011,2:821-825.

[14] Dezheng Xu,Zhaozhong He,et al.Production and Testing of HENDL-2.1/CG Coarse-group Cross-section Library Based on ENDF/B-Ⅶ.0[J]．Fusion Engineering and Design,2010,85(10-12):2105-2110.

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

[16] WU YC,Xie ZS,et al.A discrete ordinates nodal method for one dimensional neutron transport calculation in curvilinear geometry[J]．Nuclear Science and Engineering,1999,133(3):350-357.

[17] 羅月童,龍鵬程,薛曄 等.面向中子學(xué)分析的集成可視化平臺SVIP的發(fā)展研究[J]．核科學(xué)與工程,2007,27(4):374-378.

[18] 吳宜燦,劉萍,胡麗琴等.大型集成概率安全分析軟件系統(tǒng)的研究與發(fā)展[J]．核科學(xué)與工程,2007,27(3):69-75。

[19] Pengcheng Long,Qin Zeng,Tao He,et al.Development of a Geometry-Coupled Visual Analysis System for MCNP[J]．Progress in Nuclear Science and Technology,2011,2:280-283.

[20] 龍鵬程,羅月童,鄒俊 等.基于可編程圖形處理器的可視化技術(shù)在中子學(xué)分析中的應(yīng)用研究[J]．核電子學(xué)與探測技術(shù),2010,30(8):1042-1045.

Mesh Tally Method Study for Monte Carlo Particles Transport Simulation

ZHANG Kuan1,2,SONG Jing1,2,CHEN Zhen-ping1,2,SUN Guang-yao1,2,HAO Li-juan2,LONG Peng-cheng2,HU Li-qin1,2?,FDS Team
(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,Hefei,Anhui,230031,China)

Mesh tally method can accurately count the distribution of various physical quantities in reactor numerical simulation.The mesh tally of static parameters of rectangular coordinate system and cylindrical coordinate system was studied,which based on Super Monte Carlo Simulation Program for Nuclear and Radiation Process Super MC.The TCA criticality facility of Japan Atomic Energy Research Institute was carried out for validation.The feasibility and the correctness of the method were demonstrated by the match of calculation results of Super MC and MCNP.

Monte Carlo;Mesh tally;Static parameters;Super MC

2016-02-20

國家ITER 973計(jì)劃(2011GB113006);中國科學(xué)院戰(zhàn)略性先導(dǎo)科技專項(xiàng)(XDA03040000);國家自然科學(xué)基金(91026004);安徽省自然科學(xué)基金(1308085QH138)

張 寬(1986—),男,河南人,碩士研究生,主要從事蒙特卡羅粒子輸運(yùn)計(jì)算研究工作

孫光耀:guangyao.sun＠fds．org．cn

TL329＋.2

A

0258-0918(2016)02-0200-05