基于自动网格划分与权窗平滑的自适应减方差方法

郑华庆　李 斌　宋 婧　郝丽娟　龙鹏程　胡丽琴　吴宜灿

（中国科学院核能安全技术研究所 中子输运理论与辐射安全重点实验室　合肥 230031）

基于自动网格划分与权窗平滑的自适应减方差方法

郑华庆李 斌宋 婧郝丽娟龙鹏程胡丽琴吴宜灿

（中国科学院核能安全技术研究所 中子输运理论与辐射安全重点实验室合肥 230031）

如何更高效发挥减方差的作用是蒙特卡罗方法在先进核能系统屏蔽分析应用中的研究热点之一。本文发展了一种基于自动网格划分与权窗平滑的自适应减方差方法，在蒙特卡罗和确定论耦合的一致性伴随驱动的重要性抽样方法（Consistent Adjoint Driven Importance Sampling, CADIS）基础上，利用计算机辅助设计（Computer Aided Design, CAD）自动转换和自由程自动划分SN网格，通过确定论方法伴随预计算，实现基于伴随通量的区域权窗参数自动配置，并对伴随通量变化剧烈区域进行权窗平滑处理，保证了粒子在不同区域的有效偏倚，进一步提高计算的效率，从而解决大空间蒙特卡罗计算难以收敛的问题。该方法已初步应用于中国铅基反应堆（China Lead-based Reactor, CLEAR）堆顶盖的屏蔽计算分析，该案例具有结构复杂、屏蔽材料厚重的特点，测试结果表明本方法将计算效率提高近10倍。

蒙特卡罗，耦合，减方差，辐射屏蔽

屏蔽分析是核装置物理设计与分析的关键环节。蒙特卡罗方法具备处理复杂几何和模拟精细的物理过程的能力，因而在屏蔽分析方面得到广泛应用。减方差方法能够加快蒙特卡罗方法的收敛速度，可以有效解决应用蒙特卡罗方法进行屏蔽计算时收敛慢的问题。利用传统减方差方法对复杂核系统计算时，需要使用者具备丰富的经验并进行反复的试算，才能得到比较好的效果，对使用者要求较高且无法高效地发挥减方差方法的优势。如何降低减方差方法使用对经验的依赖，更高效地发挥减方差方法的作用，是蒙特卡罗方法在先进核能系统屏蔽分析应用中的研究热点之一。

美国橡树岭国家实验室Wagner等［1］先后发展了基于一致性伴随驱动的重要性系列抽样方法（Consistent Adjoint Driven Importance Sampling, CADIS；Forward-Weighted CADIS；Multi Steps CADIS），其原理是：先通过SN方法对模型进行粗网格预计算，利用其计算结果为后续的蒙特卡罗计算提供相应的减方差参数，从而使蒙特卡罗方法计算深穿透问题时的精度和效率大大提升。CADIS方法经过近20年的发展，已在辐射屏蔽分析中得到广泛的应用。但是其无法自动设置SN模型的网格参数，同时确定程序产生的权窗参数可能存在过度分裂的情况，因此需要对这两个问题进行更加深入的研究。

本文提出了基于自动网格划分与权窗平滑的自适应减方差方法，并集成在FDS团队［2-8］自主研发的超级蒙特卡罗核计算仿真软件系统SuperMC中［9-12］。该方法已在加速驱动次临界系统屏蔽分析中成功应用，可以更高效、方便地发挥减方差方法的作用，有效解决深穿透问题。

1　方法描述

针对传统的CADIS方法无法自动设置SN模型网格参数的问题，本文发展了一种基于计算机辅助设计（Computer Aided Design, CAD）自动转换和自由程预处理的SN网格参数自动生成方法。

该方法首先根据辐射源的分布和位置以及计算模型的几何复杂性进行蒙特卡罗的预计算，统计得到随着粒子类型、能量、材料变化的自由程，然后根据式（1）计算得到平均自由程。

式中：λ为平均自由程；Σ为宏观截面。

然后根据平均自由程自动有效地划分整个模型在空间中的网格边界参数值，将其作为SN自动建模的输入参数，将CAD模型自动生成SN计算模型［13-14］。该网格划分与蒙特卡罗计算的权窗网格对应，以保证粒子在每个网格中发生适当次数的分裂，从而降低估计量方差，提高计算效率。

根据自适应减方差方法，网格的权窗下限WL为：式中：R为响应量；S为正向计算源强；φ*为伴随通量；Cu为权窗上下限比值，参考MCNP说明书［15］，一般为5。

同时对伴随通量变化剧烈的区域引入了平滑因子ρ （0＜ρ≤1），按照式（3）：

式中：WL为平滑处理后的权窗下限；WL0为原始计算得到的权窗下限；ρ为平滑因子。

对权窗下限进行平滑处理，可有效缓解权重急剧变化的趋势，减少粒子经过这些区域的过度分裂，避免计算时间的额外增加，提高计算效率。

2　在ADS中的应用

加速器驱动次临界系统（Accelerator Driven subcritical System, ADS）目前被认为是一种较为理想的核废料嬗变处理装置，被国际原子能机构（International Atomic Energy Agency, IAEA）称为“新出现的核废料嬗变及能量产生的核能系统”。中国科学院核能安全技术研究所FDS团队在中国科学院战略性先导科技专项“未来先进核裂变能—ADS嬗变系统”［16］支持下，针对中国铅基反应堆（China Lead-based Reactor, CLEAR）开展全面设计与研发工作［17-19］。

本文选用中国铅基反应堆CLEAR堆顶盖屏蔽分析作为基于蒙特卡罗和确定论耦合的自适应减方差的应用例题，其几何结构如图1所示。

图1　CLEAR的CAD模型Fig.1　CAD model of CLEAR.

在次临界运行状态下，质子束流以能量250MeV、流强为10mA的状态轰击重金属散裂靶，产生散裂中子驱动反应堆运行。通过蒙特卡罗程序模拟质子轰击散裂靶，得到散裂中子产额及能谱［20-22］。如图2所示，散裂中子能量区间宽在1.0×10-5-180MeV，平均能量为7.97MeV，其中大于1MeV的所占份额为29.7%。

图2　质子轰击靶窗产生的散裂中子谱［20］Fig.2　Spectrum of spallation neutron［20］.

在分析堆顶盖（包括堆顶固定屏蔽、大小旋塞）对堆内中子、光子的衰减效果时，由于顶盖中含有较厚的重混凝土生物屏蔽层，使得在计算顶盖外侧的中子/光子通量（剂量）存在严重的“深穿透”效应，直接蒙特卡罗模拟计算统计误差难以满足收敛标准（参考MCNP说明书［15］，一般为＜5%）。

减方差方法研究中通常使用FOM （Figure of Merit）衡量蒙特卡罗软件的计算效率［15］，定义为FOM=1/σ2T，其中：σ是计数平均值的相对误差估计；T是计算时间。因σ2∝1/N，T∝N，故其理想情况下FOM应为常数。

由自适应减方差方法得到的中心轴网格权窗下限如图3所示，横坐标为中心轴网格的坐标值，纵坐标为权窗下限值，由于中子在放射源发射达到堆顶部的过程中，权窗下限值无明显突变，故不需要进行平滑处理。

图3　中心轴网格权窗下限Fig.3　Lower limit of mesh weight window at central axis.

SuperMC计算结果如表1所示，在没有添加权窗的情况下，计算品质因子FOM仅为0.0729，添加自适应减方差产生的源偏倚和权窗的情况FOM为0.7272，计算效率提升到原来的9.975倍。

表1　堆顶盖计算结果Table 1　Calculation results of reactor cover.

3　结语

本文实现了基于一致性伴随驱动的重要性抽样方法，并在此基础上进行了改进，发展了一种基于自动网格划分与权窗平滑的自适应减方差方法：利用CAD自动转换和自由程自动划分SN网格，通过确定论方法伴随预计算，实现基于伴随通量的区域权窗参数自动配置，并对伴随通量变化剧烈区域进行权窗平滑处理，保证了粒子在不同区域的有效偏倚，进一步提高计算的效率，从而解决大空间蒙特卡罗计算难以收敛的问题。通过CLEAR屏蔽计算分析的堆顶盖通量计算的应用，结果表明该自适应减方差方法将计算效率提高近10倍。

致谢感谢FDS团队其他成员对本文工作的支持。

1 Wagner J C, Peplow D E, Mosher S W, et al. Review of hybrid （deterministic/Monte Carlo） radiation transport methods, codes, and applications at Oak Ridge National Laboratory［J］. Progress in Nuclear Science and Technology, 2011, 2: 808-814

2 Wu Y, FDS Team. Conceptual design activities of FDS series fusion power plants in China［J］. Fusion Engineering and Design, 2006, 81（23-24）: 2713-2718. DOI: 10.1016/j.fusengdes.2006.07.068

3 Qiu L, Wu Y, Xiao B, et al. A low aspect ratio tokamak transmutation system［J］. Nuclear Fusion, 2000, 40: 629-633. DOI: 10.1016/S0920-3796（98）00312-3

4 Wu Y, Qian J, Yu J. The fusion-driven hybrid system and its material selection［J］. Journal of Nuclear Materials, 2002, 307-311: 1629-1636

5 Wu Y, Jiang J, Wang M, et al. A fusion-driven subcritical system concept based on viable technologies［J］. Nuclear Fusion, 2011, 51（10）: 103036

6 Wu Y, FDS Team. Conceptual design of the China fusion power plant FDS-II［J］. Fusion Engineering and Design, 2008, 83（10-12）: 1683-1689. DOI: 10.1016/j.fusengdes. 2008.06.048

7 Wu Y, FDS Team. Fusion-based hydrogen production reactor and its material selection［J］. Journal of Nuclear Materials, 2009, 386-388: 122-126. DOI: 10.1016/ j.jnucmat.2008.12.075

8 Wu Y, 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

9 Wu Y, FDS Team. CAD-based interface programs for fusion neutron transport simulation［J］. Fusion Engineering and Design, 2009, 84: 1987-1992. DOI: 10.1016/j.fusengdes.2008.12.041

10 Wu Y. Design status and development strategy of China liquid lithium-lead blankets and related material technology［J］. Journal of Nuclear Materials, 2007, 367-370: 1410-1415. DOI: 10.1016/j.jnucmat.2007. 04-031

11 Wu Y, 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

12 Wu Y, Song J, Zheng H, et al. CAD-based Monte Carlo program for integrated simulation of nuclear system SuperMC［J］. Annual of Nuclear Energy, 2015, 82: 161-168. DOI: 10.1016/j.anucene.2014.08.058

13 Song J, Sun G, Chen Z, et al. Benchmarking of CAD-based SuperMC with ITER Benchmark model［J］. Fusion Engineering and Design, 2014, 89: 2499-2503. DOI: 10.1016/j.fusengdes.2014.05.003

14 孙光耀, 宋婧, 郑华庆, 等. 超级蒙特卡罗软件SuperMC2.0中子输运计算校验［J］. 原子能科学技术, 2013, 47（增2）: 140-145

SUN Guangyao, SONG Jing, ZHENG Huaqing, et al. Benchmark of neutron transport simulation capability of super Monte Carlo calculation program SuperMC2.0［J］. Atomic Energy Science and Technology, 2013, 47（Suppl 2）: 140-145

15 Hu H, Wu Y, Chen M, et al. Benchmarking of SNAM with the ITER 3D model［J］. Fusion Engineering and Design, 2007, 82: 2867-2871. DOI: 10.1016/ j.fusengdes.2007. 06.015

16 Long P, Zou J, Huang S, et al. Development and application of SN auto-modeling tool SNAM 2.1［J］. Fusion Engineering and Design, 2010, 85: 1113-1116. DOI: 10.1016/j.fusengdes.2010.02.016

17 X-5 Monte Carlo Team. MCNP - a general Monte Carlo N-particle transport Code［R］. Version 5. Report LA-UR-03-1987, Los Alamos National Laboratory, April 24, 2003

18 詹文龙, 徐瑚珊. 未来先进裂变能—ADS嬗变系统［J］ .中国科学院院刊, 2012, 27（3）: 375-381

ZHAN Wenlong, XU Hushan. Advanced fission energy program - ADS transmutation system［J］. Bulletin ofChinese Academy of Sciences, 2012, 27（3）: 375-381

19 Wu Y, Bai Y, Song Y, et al. Development strategy and conceptual design of China lead-based research reactor［J］. Annals Nuclear Energy, 2016, 87: 511-516. DOI: 10.1016/j.anucene.2015.08.015

20 吴宜灿, 柏云清, 宋勇, 等. 中国铅基研究反应堆概念设计研究［J］. 核科学与工程, 2014, 34（2）: 56-63

WU Yican, BAI Yunqing, SONG Yong, et al. Conceptual design of China lead-based research reactor CLEAR-I［J］. Nuclear Science and Engineering, 2014, 34（2）: 56-63

21 Wu Y, Bai Y, Wang W, et al. Overview of China lead-alloy cooled rector development and ADS program in China［C］. POWER - the 9th International Topical Meeting on Nuclear Thermal-hydraulics, Operation and Safety （NUTHOS-9）,Kaohsiung, Taiwan, September 9-13, 2012

22 Yang Q, Li B, Chen C, et al. Primary shielding design and analysis for CLEAR-I［C］. 21st International Conference on Nuclear Engineering （ICONE-21）, Chengdu, Sichuan, China, July 29-August 2, 2013

Adaptive variance reduction method based on automatic mesh generating and weight window smoothing

ZHENG HuaqingLI BinSONG JingHAO LijuanLONG PengchengHU LiqinWU Yican
（Key Laboratory of Neutronics and Radiation Safety, Institute of Nuclear Energy Safety Technology, Chinese Academy of Sciences, Hefei 230031, China）

Background: In the radiation shielding analysis of advanced nuclear reactor with Monte Carlo method, it is one of the researches focuses on using variance reduction methods more effectively. Purpose: An adaptive variance reduction method based on automatic mesh generating and weight window smoothing is developed in this paper. Methods: With hybrid calculation of Monte Carlo method and deterministic theory, consistent adjoint driven importance sampling （CADIS） method was developed by automatic mesh generating with computer aided design（CAD） and pre-calculation of free path. The adjoint flux calculation results with deterministic method are treated as core values to generate source biasing parameters and weight window parameters in the form of mesh. And smooth treatment is introduced for alleviating the lower bound of weight window violent change to make sure the efficiency of particle bias in different areas and improve calculation efficiency. This method can solve the convergence problem in Monte Carlo simulation with large space geometry. Results: The adaptive method has been applied in the shielding

Monte Carlo, Hybrid calculation, Variance reduction, Radiation shielding

ZHENG Huaqing, male, born in 1984, graduated from University of Chinese Academy of Sciences with a doctor’s degree in 2015, engaged in

WU Yican, E-mail: yican.wu@fds.org.cn

TL7

10.11889/j.0253-3219.2016.hjs.39.090502

——

中国科学院战略性先导科技专项（No.XDA03040000）、国家磁约束核聚变能发展研究专项（No.2014GB112001）、国家自然科学基金（No.11305205、

No.11305203、No.11405204）资助

郑华庆，男，1984年出生，2015年于中国科学院大学获博士学位，现从事蒙特卡罗粒子输运计算研究工作

吴宜灿，E-mail: yican.wu@fds.org.cn

Supported by Strategic Priority Research Program of Chinese Academy of Sciences （No.XDA03040000）, the National Special Program for ITER

（No.2014GB12001）, National Natural Science Foundation of China （No.11305205, No.11305203, No.11405204）

the research of Monte Carlo particle transport simulation

2016-02-11，

2016-03-22

calculation of China lead-based reactor （CLEAR） cover with complex geometry and massive shielding body. The result showed that the computing efficiency is improved by nearly 10 times. Conclusion: The adaptive variance reduction method can be applied to the shielding calculation of advanced nuclear fission reactor.