張 寬,宋 婧,陳珍平,孫光耀,郝麗娟,龍鵬程,胡麗琴?,FDS團隊
(1.中國科學技術(shù)大學,安徽 合肥230027;2.中國科學院核能安全技術(shù)研究所,中國科學院中子輸運理論與輻射安全重點實驗室,安徽 合肥230031)
蒙特卡羅粒子輸運計算網(wǎng)格計數(shù)方法研究
張 寬1,2,宋 婧1,2,陳珍平1,2,孫光耀1,2,郝麗娟2,龍鵬程2,胡麗琴1,2?,FDS團隊
(1.中國科學技術(shù)大學,安徽 合肥230027;2.中國科學院核能安全技術(shù)研究所,中國科學院中子輸運理論與輻射安全重點實驗室,安徽 合肥230031)
在進行反應堆數(shù)值模擬時,通過網(wǎng)格計數(shù)方法可以精細地統(tǒng)計整個堆芯中多種物理量的空間分布情況。本文基于超級蒙特卡羅核計算仿真軟件系統(tǒng)Super MC發(fā)展了一種反應堆靜態(tài)物理參數(shù)的網(wǎng)格計數(shù)方法。通過日本原子能研究所臨界實驗裝置TCA例題進行了數(shù)值驗證,計算結(jié)果與MCNP計算結(jié)果吻合較好,表明了本文方法的可行性與正確性。
蒙特卡羅;網(wǎng)格計數(shù);靜態(tài)參數(shù);Super MC
通量、反應率(尤其是裂變率)是反應堆物理分析中常用到的物理量[1],堆芯功率分布分析是核反應堆安全運行的基礎[2],是反應堆設計、安全分析、故障診斷、實時控制等問題的基礎。在實際問題中,堆芯功率使用更方便,同時因堆芯中子通量分布與堆芯功率分布基本成比例,所以實際應用中多分析堆芯功率分布[3]。分析反應堆堆芯等重要區(qū)域的功率分布或通量、反應率分布時需要進行詳細的計數(shù)統(tǒng)計以獲得相應物理量精細的空間分布。對于反應堆用的蒙特卡羅程序而言,基于傳統(tǒng)計數(shù)方法的網(wǎng)格計數(shù)功能是一項重要功能,主要針對通量、反應率、功率等靜態(tài)參數(shù)。例如:MCNP[4],Serpent[5]、MC21[6]等反應堆蒙特卡羅計算程序也具有基本的網(wǎng)格計數(shù)功能。
本文在傳統(tǒng)的計數(shù)方法基礎上實現(xiàn)并改進了通量、反應率、功率等靜態(tài)參數(shù)的網(wǎng)格計數(shù)方法,在中國科學院核能安全技術(shù)研究所·FDS團隊研發(fā)的超級蒙特卡羅計算仿真軟件系統(tǒng)Super MC[7]中進行應用實現(xiàn)并通過基準例題校驗。
其中:T是關(guān)心體積中所有粒子徑跡序列,li是i軌跡的長度,W是粒子的初始權(quán)重總和。
反應率的計算一般是基于通量估計的基礎上實現(xiàn)的,反應率的徑跡長度估計也可以通過公式(1)乘以微觀截面得到:
基于徑跡長度法的通量估計為:
堆芯任一點r處單位體積內(nèi)的功率,即r處的功率密度為:
本文設計實現(xiàn)的網(wǎng)格計數(shù)流程如圖1所示。首先根據(jù)輸入標識判斷是直角坐標系還是圓柱坐標系,然后進入相應的計數(shù)模塊。對于圓柱坐標系的計數(shù),首先要進行坐標轉(zhuǎn)換,后續(xù)步驟這兩種情況的計數(shù)方法一致:首先判斷該步的粒子類型與要計數(shù)的粒子類型是否一致,若一致則進行計數(shù)運算;然后判斷粒子當前步是否在劃分網(wǎng)格的區(qū)域中,若在,則查找各個方向所在的小網(wǎng)格位置;下一步統(tǒng)計在此小網(wǎng)格中粒子徑跡長度并檢測是否需要計算反應率或功率,若需要調(diào)用相應計數(shù)模塊;最后把粒子位置移動到此網(wǎng)格的邊界,進行下一網(wǎng)格的統(tǒng)計。
圖1 網(wǎng)格計數(shù)流程圖Fig.1 Flow Chart of Mesh Tally
圓柱坐標系下,由于半徑方向不是直線邊界,所以圓柱坐標系的計數(shù)與直角坐標系有所不同,需要將粒子坐標從世界坐標系到圓柱坐標系進行轉(zhuǎn)換。此轉(zhuǎn)換需要一個輔助直角坐標系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。
圓柱坐標系曲面方程可以表示為:
或者表示為:=0
則直角坐標系到圓柱坐標系的轉(zhuǎn)換為:
得到轉(zhuǎn)換矩陣為[4]:
蒙特卡羅粒子輸運程序中網(wǎng)格功率計算一般基于公式(3)和上述網(wǎng)格計算方法進行統(tǒng)計,但該方法具有局限性:當某一網(wǎng)格跨越不同材料的柵元且粒子的抽樣步長正好跨越不同材料區(qū)域時,此步長功率計算并未準確獲取每段徑跡的材料,而使用該步中起點所在柵元的材料作為整個步長計算的材料,與真實情況相比會造成一定的誤差。
為此本文進行了功率計算功能的優(yōu)化:首先在每步輸運中記錄該步長是否跨越不同材料的柵元,若不跨越則可用一般方法進行處理;若跨越,則記錄穿越不同材料柵元的交點及相關(guān)材料信息,然后將該步長根據(jù)交點分成若干子步長,每個子步長取各自起始點所在柵元的材料進行功率計算,然后進行累加。具體流程見圖2。
圖2 功率計算跨網(wǎng)格處理流程圖Fig.2 Flow Chart of Crossing Mesh Processing for Power Calculation
本文選用日本原子能研究所臨界實驗裝置TCA(Tank-type Critical Assembly)[8]例題進行基準校驗。TCA是日本原子能研究所的一個壓水堆臨界實驗裝置,其燃料棒有兩種: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中基于組件模型陣列自動建模功能進行TCA裝置建模[9-13],如圖3所示,并實現(xiàn)完整的材料、源項、計數(shù)建模,在Super MC中對模型全區(qū)域進行網(wǎng)格劃分:徑向劃分10個,軸向劃分10個,角度方向劃分18個。本測試中均使用FDS團隊自主研發(fā)的HENDL3.0[14-16]數(shù)據(jù)庫,使用網(wǎng)格計數(shù)方法計算通量、反應率、功率的分布情況,計算結(jié)果與MCNP進行對比。
圖3 TCA裝置模型圖Fig.3 TCA Criticality Facility Model(a)俯視圖;(b)剖面圖
對每個網(wǎng)格的通量、反應率、功率計算結(jié)果進行統(tǒng)計分析,結(jié)果如表1所示。由結(jié)果可知:所有網(wǎng)格計數(shù)結(jié)果與MCNP相對偏差都在1%以下,符合工程設計要求,同時至少有92.05%的網(wǎng)格結(jié)果在一個σ之內(nèi),符合統(tǒng)計學的要求,計算結(jié)果正確。使用Super MC的可視化模塊RVIS[15-20]對網(wǎng)格通量計算結(jié)果進行可視化與對比,結(jié)果如圖4所示,由圖可見Super MC計算結(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è)視圖
本文設計實現(xiàn)了蒙特卡羅粒子輸運模擬中反應堆靜態(tài)參數(shù)的網(wǎng)格計數(shù)方法,并考慮了計數(shù)網(wǎng)格跨不同材料區(qū)域的問題,通過日本原子能研究所臨界實驗裝置TCA例題的測試結(jié)果表明,本文方法對網(wǎng)格通量、反應率、功率的計數(shù)結(jié)果與MCNP吻合良好,表明了本文方法的可行性與正確性。
[1] 謝仲生,吳宏春等.核反應堆物理分析[M].西安:西安交通大學出版社.2004.
[2] 李樹,田東風,鄧力.中子時間常數(shù)的Monte Carlo計算方法[J].北京:清華大學學報(自然科學版),2007:1057-1061.
[3] 李偉.基于堆外計數(shù)的堆芯功率分布重構(gòu)方法研究[D].哈爾濱:哈爾濱工程大學,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] 曾勤,李瑩,盧磊 等.蒙特卡羅粒子輸運計算自動建模程序MCAM在ITER核分析建模中的應用[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] 吳宜燦,李瑩,盧磊等.蒙特卡羅粒子輸運計算自動建模程序系統(tǒng)的研究與發(fā)展.核科學與工程,2006,26(1):20-27.
[11] 王國忠,黨同強,熊健 等.MCAM4.8在ITER建筑大廳中子學建模中的應用[J].核科學與工程,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] 吳宜燦,胡麗琴,龍鵬程等.先進核能系統(tǒng)設計分析軟件與數(shù)據(jù)庫研發(fā)進展[J].核科學與工程,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] 羅月童,龍鵬程,薛曄 等.面向中子學分析的集成可視化平臺SVIP的發(fā)展研究[J].核科學與工程,2007,27(4):374-378.
[18] 吳宜燦,劉萍,胡麗琴等.大型集成概率安全分析軟件系統(tǒng)的研究與發(fā)展[J].核科學與工程,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ù)在中子學分析中的應用研究[J].核電子學與探測技術(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計劃(2011GB113006);中國科學院戰(zhàn)略性先導科技專項(XDA03040000);國家自然科學基金(91026004);安徽省自然科學基金(1308085QH138)
張 寬(1986—),男,河南人,碩士研究生,主要從事蒙特卡羅粒子輸運計算研究工作
孫光耀:guangyao.sun@fds.org.cn
TL329+.2
A
0258-0918(2016)02-0200-05