李新梅,郑华庆,郝丽娟,宋 婧,胡丽琴,3,江 平
(1.合肥工业大学数学学院,安徽合肥230002;2.中国科学院核能安全技术研究所,中国科学院中子输运理论与辐射安全重点实验室,安徽合肥230031;3中国科学技术大学,安徽合肥230027)
网格权窗减方差技术及其在聚变堆屏蔽分析中应用研究
李新梅1,2,郑华庆2,郝丽娟2,宋 婧2,胡丽琴2,3,江 平1
(1.合肥工业大学数学学院,安徽合肥230002;2.中国科学院核能安全技术研究所,中国科学院中子输运理论与辐射安全重点实验室,安徽合肥230031;3中国科学技术大学,安徽合肥230027)
在聚变堆辐射屏蔽计算中,如何有效解决深穿透问题是近年来国际聚变辐射安全领域关注的焦点之一。针对该问题,本文研究了直角坐标系与圆柱坐标系下基于网格的权窗减方差技术。本文基于超级蒙卡核模拟软件系统SuperMC实现了该方法,并选取减方差技巧的基准例题进行测试与分析,初步得出“粗划真空或密度很小的区域、细分密度大的区域”的网格划分规律,能有效提高网格权窗计算效率。基于该规律对聚变屏蔽基准问题进行对比分析,新的网格划分与原始网格划分的计算效率相比,FOM因子提高了1.92倍。减方差技巧的基准例题和聚变屏蔽基准问题计算中,SuperMC通量计算结果与MCNP相比偏差均在0.5%以下,证明了本文中方法的正确性。
基于网格的权窗;聚变堆;屏蔽计算;减方差技巧;SuperMC
收敛速度慢是蒙特卡罗方法(简称蒙卡)在进行聚变反应堆[1-5]模拟的瓶颈性限制之一,尤其是对于聚变堆屏蔽计算中深穿透问题。减方差技巧[6]是解决此类问题的主要方法。目前,现有的减方差技巧包括几何分裂与赌、源偏倚、指数变换、权窗等,并在MCNP[6]、FLUKA[7]、Serpent[8]等蒙卡软件中实现。
多年的实际应用表明,目前权窗[9]是应用最广泛的减方差技巧。权窗分为基于栅元(cell-based)的权窗和基于网格(mesh-based)的权窗两种模式。基于栅元的权窗是根据栅元的重要性设置权窗上下限参数,引导重要粒子朝着感兴趣的区域输运。为了有效利用栅元权窗需将尺寸过大栅元划分成合适的小尺寸栅元。而基于网格的权窗是根据网格区域的重要性设置网格权窗上下限,无需划分真实模型,只需定义划分网格参数即可[10]。栅元权窗一般适用于简单几何问题,对于大型、复杂模型,其适用性较差,尤其针对带重复结构的复杂模型,而基于网格的权窗恰好可以弥补栅元权窗的这些缺点。
本文系统地研究了网格权窗的减方差技巧,并在中国科学院核能安全技术研究所·FDS团队自主研发的超级蒙卡核模拟软件系统SuperMC[11-17]的基础上进行实现,并选取减方差基准例题对程序进行测试与分析。
目前,国外实现的网格权窗技巧主要是基于直角坐标系和圆柱坐标系上建立的[18]。直角坐标系下的网格绘制比较方便、直观。圆柱坐标系下的网格可以旋转任意角度进行划分,比较复杂但使用广泛。相比于直角坐标系,圆柱坐标系的网格权窗中需将粒子当前位置坐标转换到直角坐标系上,其转换矩阵如下:
通过转换后的坐标或者实际坐标找到相应的网格权窗下限,继而进行网格权窗处理,其基本原理如下:
权重为wgt的粒子输运到某网格相应的能量区间时:
1) 当wgt 2) 当wgt>WU,进行分裂,粒子分裂成int(wgt/WS+ε)个,权重均为WS; 3) 当WL 网格权窗的使用能将粒子权重控制在合理的范围内,有效避免极高权重粒子对计数的扰动,最终提高计算效率。其流程见图1。 图1 网格权窗的流程图Fig.1 Flow Chart of Mesh-based Weight Window 由图1可知,在使用网格权窗时需考虑粒子当前能量和材料下的平均自由程,取抽样自由程、到下一边界距离和平均自由程的最小值作为输运长度。因此,使用网格权窗技巧修正粒子游动过程时需在下一边界处、碰撞点以及平均自由程处进行处理。 本文选取参考文献[19]中减方差技巧测试的基准例题,并使用MCNP权窗产生器[8]生成的网格权窗下限进行测试。 2.1 长真空圆柱基准测试例题 2.1.1 模型描述 长真空圆柱测试例题是减方差技巧基准例题[19-20],其几何由半径100cm、高2010cm的圆柱和半径210cm、高10cm的圆柱叠加而成(见图2)。其中,距底面180cm厚度的栅元是密度为2.03g/cm3的重混凝土;180~2000cm为真空栅元;2000~2010cm是密度为0.0203g/cm3的轻混凝土层。能量50% 2MeV与50% 14MeV、方向各向同性的中子点源位于几何体底层重混凝土中,距底面1.0×10-6cm,计数区域为图中最上层非真空栅元,类型为栅元体通量。由于模型所有栅元均是圆柱型的,因此采用圆柱坐标系下的网格划分。表1中给出了圆柱轴向上不同网格划分的参数值。 图2 长真空圆柱几何示意图Fig.2 Figure of Vacuum Cylinder Model 表1 网格划分方法 续表 2.1.2 测试结果与分析 表2中给出了不使用网格权窗和使用网格权窗1、2、3的计算结果,SuperMC与MCNP的相对偏差均小于0.5%。表2中在不使用网格权窗条件下无计算结果,这是因为模型中间有很长一段真空材料,属于典型的深穿透问题,如果不使用减方差技巧,无法得到稳定收敛的计算结果。 表2 长真空圆柱计算结果 表3中给出了使用网格权窗1的品质因子FOM为7.60,网格权窗2的FOM为7.54,网格权窗3的FOM为7.60。比较不使用与使用网格权窗的计算结果发现,有效利用网格权窗能够成倍提高计算效率。 表3 长真空圆柱计算效率比较 注:① FOM=1/(σ2T),σ为标准差,T为计算时间/min。 综合分析上述三组网格权窗划分对计算结果的影响,可以总结得到如下网格设置的规律:对于真空或者密度很小的材料区域划分尽量控制在一个粗网格内,无需再细分;对于材料密度大、密度大的材料区域划分则相反。 2.2 聚变屏蔽基准例题测试 2.2.1 模型描述 聚变屏蔽基准问题[19]的整个模型在899.16cm×690.85cm×678.18cm的长方体内,几何和材料如图3所示。中子点源位于中间真空栅元内,其坐标为(-356.87,232.02,157.40),能量为14MeV,方向各向同性。在源外部有55.88cm厚的含铁屏蔽层,计数区域为屏蔽层后的栅元,类型为栅元体通量。此模型采用两种方法划分:参考文献[19]中网格划分(网格1)、根据上述网格划分规律对整个模型进行划分(网格2)。 图3 聚变屏蔽模型Fig.3 Fusion Shielding Model 网格2(以y轴方向为例)划分方法如下:-29.21~0.00cm划分5份,0.00~178.85cm划分5份,178.50~196.10cm划分4份,196.15~223.79cm划分6份,223.79~259.39cm划分5份,259.39~279.75cm划分4份,279.75~407.33cm划分1份,407.33~661.64cm划分1份。 2.2.2 测试结果与分析 表4中给出不使用网格权窗和使用网格权窗1、2时的计算结果,SuperMC与MCNP计算结果的相对偏差均小于0.5%。可见,无论圆柱坐标系下还是直角坐标系下使用网格权窗的计算结果均与MCNP吻合程度良好。 表4 聚变热屏蔽计算结果 续表 表5中给出不使用网格权窗的品质因子FOM为7.15,使用网格权窗1后FOM为123.70,计算效率是不使用网格权窗时的17.30倍;使用网格权窗2后FOM为237.059,计算效率是不使用网格权窗时的33.16倍,是网格权窗1计算效率的1.92倍。可见,总结得出的网格划分规律在聚变堆屏蔽分析中是有效的。 表5 聚变热屏蔽计算效率比较 注:① FOM=1/(σ2T),σ为标准差,T为计算时间/min。计算效率为FOM的比值。 本文对蒙特卡罗粒子输运模拟中直角坐标系与圆柱坐标系下的网格权窗减方差技巧进行了研究。利用减方差技巧的基准例题进行测试与分析,初步得出“粗划真空或密度很小的区域、细分密度大的区域”的网格划分规律,能有效提高网格权窗的计算效率。基于该规律对聚变屏蔽基准问题进行对比分析,新的网格划分与原始网格划分的计算效率相比FOM因子提高了1.92倍。SuperMC通量计算结果与MCNP相比偏差均在0.5%以下,证明了本文中方法的正确性。本文网格划分方法是手工划分的,可进一步研究网格权窗的网格自动划分方法。 致谢 感谢FDS团队成员对本文工作的帮助与支持。 [1] Y. Wu, FDS Team. CAD-Based Interface Programs for Fusion Neutron Transport Simulation[J]. Fusion Engineering and Design, 2009, 84 (7-11), 1987-1992. [2] Y. Wu, J. Jiang, M. Wang, M. Jin, et al. A Fusion-Driven Subcritical System Concept Based on Viable Technologies[J].Nuclear Fusion, 2011, 51(10):103036. [3] S. Zheng, M. Chen, J. Li, et al. Neutronics Analysis for the Test Blanket Modules Proposed for EAST and ITER[J]. Nuclear Fusion, 47 (2007) 1053-1056. [4] Y.Chen, Y. Wu. Conceptual Study on High Performance Blanket in a Spherical Tokamak Fusion-DrivenTransmuter[J]. Fusion Engineering and Design, 2000, 49-50:507-512. [5] 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, Part B: 1305-1311. [6] X-5 Monte Carlo Team. MCNP-A General Monte Carlo N-Particle Transport Code, Version 5[R]. Los Alamos National Laboratory Report LA-UR-03-1987, April 24, 2003. [7] A. Ferrari, P. R. Sala, A. Fasso, et al. Fluka:A Multi-paricle Transport Code[R]. CERN-2005-010, Oct 12,2005. [8] J. L. Serpent -A Continuous Energy Monte Carlo Reactor Physics Burn up Calculation Code User’s Manual [M]. March 6, 2013: 131-133. [9] T. E. Booth. Genesis of the Weight Window and the Weight Window Generator in MCNP-A Personal History[J]. LA-UR-06-5807. 9, July 28, 2006. [10] 姜宏宇, 王汝赡. 网格权重窗方法及其在煤浆系统研究中的应用[J]. 核电子学与探测技术, 1998,18(3): 201-205. [11] 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.5.3. [12] 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, 82,161-168(2015). [13] Z.P. Chen, H.Q. Zheng, G.Y. Sun, et al. Preliminary Study on CAD-based Method of Characteristics for Neutron Transport Calculation[J]. Chinese Physics C. (2014) 38 (5), 058201-1 058201-7. [14] 吴宜灿, 李静惊, 李莹 等. 大型集成多功能中子学计算与分析系统VisualBUS的研究与发展[J]. 核科学与工程,2007, 27(4): 365-373. [15] 吴宜灿, 胡丽琴, 龙鹏程,等.先进核能软件发展与核信息学实践[M]. 中国科研信息化蓝皮书, 2013. 北京:科学出版社, 2013,12: 232-244. [16] Y. Wu, 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. [17] Y. Wu, X. Zhu, S. Zheng, et al. Neutronics Analysis of Dual-Cooled Waste Transmutation Blanket for the FDS[J]. Fusion Engineering and Design, 2002, 63-64: 133-138. [18] J. S. Hendricks. Superimposed Mesh Plotting in MCNP [R]. Los Alamos National Laboratory LA-UR-01-1033, September 2001. [19] C. V. Culbertson, J. S. Hendricks. An Assessment of the MCNP4C Weight Window[R]. Los Alamos National Laboratory LA-13668, December 1999. [20] T. E. Booth. A Sample Problem for Variance Reduction in MCNP[R]. Los Alamos National Laboratory, LA-10363-MS, October 1985. Mesh-basedWeightWindowVarianceReductionTechniquesanditsAppliedResearchonFusionReactorShieldingAnalysis LIXin-mei1,2,ZHENGHua-qing2,HAOLi-juan2,SONGJing2,HuLi-qin2,3,JIANGPing1 (1.Hefei University of Technology, Hefei, Anhui, 230002, China;2.Key Laboratory of Neutronics and Radiation Safety, Institute of Nuclear Energy Safety Technology, Chinese Academy of Sciences, Hefei, Anhui, 230031, China;3.University of Science and Technology of China, Hefei, Anhui, 230027, China) In the field of fusion radiation safety, it is a focus to solve deep penetration problem effectively for the radiation shielding calculation in recent years. In this paper, mesh-based weight window variance reduction techniques under Cartesian and cylindrical coordinate systems were studied and implemented in Super Monte Carlo Program for Nuclear and Radiation Simulation(SuperMC). During the verification with a variance reduction benchmark problem, a regulation of rough mesh division in vacuum or low density regions and fine division in high density regions for improving the computing efficiency was obtained. The fusion shielding benchmark problem was analyzed with a new mesh based on the regulation and the FOM increased by 1.92 times with the new mesh. The difference between the calculation results of SuperMC and MCNP was within 0.5% for the benchmark problem and fusion shielding problem. The accuracy was demonstrated. Mesh-based weight window;Fusion reactor;Shielding calculation;Variance reduction;SuperMC 2016-12-20 国家ITER 973计划(2014GB112001);中国科学院战略性先导科技专项(XDA03040000);国家自然科学基金(11605233、11305203) 李新梅(1988—),女,安徽人,硕士研究生,现主要从事蒙特卡罗粒子输运计算研究工作 胡丽琴:liqin.hu@fds.org.cn TL7 :A :0258-0918(2017)04-0577-062 程序测试与结果分析
3 总结