Application of global variance reduction method to calculate a high-resolution fast neutron flux distribution for TMSR-SF1

Application of global variance reduction method to calculate a high-resolution fast neutron flux distribution for TMSR-SF1 1 2 first-author 1 1 1 1 2 1 1 2 corresp Corresponding author, mayuwen@sinap.ac.cn Corresponding author, mayuwen@sinap.ac.cn mayuwen@sinap.ac.cn Shanghai Institute of Applied Physics, Chinese Academy of Sciences, Shanghai 201800, China Shanghai Institute of Applied Physics, Chinese Academy of Sciences, Shanghai 201800, China University of Chinese Academy of Sciences, Beijing 100049, China University of Chinese Academy of Sciences, Beijing 100049, China The solid-fuel thorium molten salt reactor (TMSR-SF1) is a 10-MWth fluoride-cooled pebble bed reactor. As a new reactor concept, one of the major limiting factors to reactor lifetime is radiation-induced material damage. The fast neutron flux ( E > 0.1 MeV) can be used to assess possible radiation damage. Hence, a method for calculating high-resolution fast neutron flux distribution of the full-scale TMSR-SF1 reactor is required. In this study, a two-step subsection approach based on MCNP5 involving a global variance reduction (GVR) method, referred to as forward-weighted consistent adjoint driven importance sampling (FW-CADIS), was implemented to provide fast neutron flux distribution throughout the TMSR-SF1 facility. In addition, instead of using the general source specification cards, the user-provided SOURCE subroutine in MCNP5 source code was employed to implement a source biasing technique specialized for TMSR-SF1. In contrast to the one-step analog approach, the two-step subsection approach eliminates zero-scored mesh tally cells and obtains tally results with extremely uniform and low relative uncertainties. Furthermore, the maximum fast neutron fluxes of the main components in TMSR-SF1 are provided, which can be used for radiation damage assessment of the structural materials. TMSR-SF1 Fast neutron flux Global variance reduction MCNP 1 Introduction 1 Fluoride-salt-cooled high-temperature reactors (FHRs) [ 1 1 , 2 2 ] are a Gen-Ⅳ nuclear energy system that combines the graphite-matrix coated particle fuel (TRISO) adopted in high-temperature gas-cooled reactors and a molten fluoride coolant salt derived from molten salt reactors. This advanced reactor design has several advantages such as intrinsic safety features and potential economic advantages. In 2011, the Chinese Academy of Sciences launched the "Thorium Molten Salt Reactor Nuclear Energy System" (TMSR) project [ 3 3 , 4 4 ] , and a 10-MWth solid-fuel thorium molten salt reactor (TMSR-SF1) was developed by the Shanghai Institute of Applied Physics, which will be deployed in the next 5 to 10 years. As a new reactor concept, the radiation damage assessment of the structural materials is crucial to TMSR-SF1 design. Therefore, developing calculation methods to determine the high-resolution neutron flux (especially the fast neutron flux) distribution of the full-scale TMSR-SF1 reactor is necessary. Currently, Monte Carlo (MC) methods are the main type of computational method applied to model TMSR-SF1 [ 5 5 , 6 6 ] . MC methods are based on a statistical method which simulates individual particles and tracks their histories. These events and information are recorded to obtain their average statistical behaviors as tallies, whereas, the possibility of particle transport to different regions varies widely. Therefore, it is very difficult for MC methods to solve global problems that require obtaining tallies over the entire problem space with uniformly low relative uncertainty. However, the fast neutron ( E > 0.1 MeV) flux from the core is rapidly attenuated since TMSR-SF1 contains a large amount of neutron moderating material. Therefore, the calculation of the fast neutron flux distribution of TMSR-SF1 is a typical deep-penetration global optimization problem which is extremely computationally expensive for a MC simulation. To overcome the inherent statistical limitations in the global problem of MC simulation, numerous global variance reduction (GVR) methods have been proposed to simultaneously increase the computational efficiency and simulation performance for large models [ 7 7 - 10 10 7-10 ]. Among these GVR methods, forward-consistent adjoint driven importance sampling (FW-CADIS) is one of the most efficient and widely used methods. In 2007, FW-CADIS was proposed by Wagner to improve the performance of MC simulations for global targets [ 11 11 ] . The principle of the FW-CADIS method is to develop an importance function to generate optimal GVR parameters (i.e., source biasing and particle weight window parameters) for MC simulation. The importance function has a negative relevant correlation with the forward flux, which can be beneficial for biasing the MC particles which are evenly distributed over the entire tally region. The FW-CADIS method employs both the forward and the adjoint solutions from deterministic transport calculations to produce the importance function and, thus, to realize variance reduction (VR) for global response problems. In this study, a two-step subsection approach involving the FW-CADIS technique (Sub-GVR method) was applied to calculate the fast neutron flux distribution of the TMSR-SF1 facility. Firstly, a criticality calculation was performed to generate the space- and energy-dependent fission neutron source. Secondly, a fixed-source GVR calculation employing the FW-CADIS method was carried out to produce the detailed fast neutron flux distribution. A three-dimensional discrete ordinate (SN) code DENOVO [ 12 12 ] was employed to produce the importance function, and a utility code called SN2MCNP was developed to convert the importance function into source biasing and weight window parameters. Compared with other applications of the FW-CADIS method, the fission neutron source region of TMSR-SF1 has a complex geometry which contains a circular truncated cone. For this reason, MCNP5 was employed to implement our application since it provides a user-defined SOURCE subroutine which can be used to flexibly define a source by writing code in this subroutine. A dedicated sampling methodology for TMSR-SF1 was employed in the SOURCE subroutine to implement the source biasing technique. The remainder of this paper is organized as follows. The methodology including the fission source probability distribution function (PDF), adjoint theory, and FW-CADIS are introduced in Sect.2. The procedure used to implement the Sub-GVR method is described in Sect.3. The calculation models are presented in Sect.4. The simulation results and a detailed discussion are presented in Sect.5, and conclusions are drawn in Sect. 6. 2 Methodology 1 2.1 Source probability distribution function 2 The spatial distribution of the fission neutron source in the core can be expressed as <math display="block" id="M1"><mrow><mi>S</mi><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo><mo>=</mo><msub><mover accent="true"><mi>R</mi><mo>¯</mo></mover><mtext>f</mtext></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo><mstyle displaystyle="true"><munder><mo>∑</mo><mi>m</mi></munder><mrow><msub><mi>ν</mi><mi>m</mi></msub><mfrac><mrow><msub><mi>σ</mi><mrow><mi>m</mi><mtext>f</mtext></mrow></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo><msub><mi>N</mi><mi>m</mi></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo></mrow><mrow><mstyle displaystyle="true"><munder><mo>∑</mo><mi>m</mi></munder><mrow><msub><mi>σ</mi><mrow><mi>m</mi><mtext>f</mtext></mrow></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo><msub><mi>N</mi><mi>m</mi></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo></mrow></mstyle></mrow></mfrac></mrow></mstyle><mo>,</mo></mrow></math> where <math id="M2"><mi>m</mi></math> is the m th fissionable nuclide, <math id="M3"><mi>i</mi></math> , j , and <math id="M4"><mi>k</mi></math> are the indices used to specify a fuel element, <math id="M5"><mrow><msub><mover accent="true"><mi>R</mi><mo>¯</mo></mover><mtext>f</mtext></msub></mrow></math> stands for the total fission rate in a fuel element, <math id="M6"><mrow><msub><mi>σ</mi><mrow><mi>m</mi><mtext>f</mtext></mrow></msub></mrow></math> denotes the effective microscopic fission cross-section for nuclide <math id="M7"><mi>m</mi></math> , <math id="M8"><mrow><msub><mi>N</mi><mi>m</mi></msub></mrow></math> refers to the nuclide atomic density of nuclide <math id="M9"><mi>m</mi></math> within the fuel, and <math id="M10"><mrow><msub><mi>ν</mi><mi>m</mi></msub></mrow></math> represents the average number of neutrons per fission for nuclide <math id="M11"><mi>m</mi></math> . If one lets <math id="M12"><mrow><mover accent="true"><mi>ν</mi><mo>¯</mo></mover><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo><mo>=</mo><mstyle displaystyle="true"><munder><mo>∑</mo><mi>m</mi></munder><mrow><msub><mi>ν</mi><mi>m</mi></msub><mfrac><mrow><msub><mi>σ</mi><mrow><mi>m</mi><mtext>f</mtext></mrow></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo><msub><mi>N</mi><mi>m</mi></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo></mrow><mrow><mstyle displaystyle="true"><munder><mo>∑</mo><mi>m</mi></munder><mrow><msub><mi>σ</mi><mrow><mi>m</mi><mtext>f</mtext></mrow></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo><msub><mi>N</mi><mi>m</mi></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo></mrow></mstyle></mrow></mfrac></mrow></mstyle></mrow></math> , then the space PDF of fission neutrons can be given by <math display="block" id="M13"><mrow><msup><mi>S</mi><mo>′</mo></msup><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo><mo>=</mo><mfrac><mrow><msub><mover accent="true"><mi>R</mi><mo>¯</mo></mover><mtext>f</mtext></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo></mrow><mrow><mstyle displaystyle="true"><munder><mo>∑</mo><mrow><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi></mrow></munder><mrow><msub><mover accent="true"><mi>R</mi><mo>¯</mo></mover><mtext>f</mtext></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo><mstyle displaystyle="true"><munder><mo>∑</mo><mi>m</mi></munder><mrow><mover accent="true"><mi>ν</mi><mo>¯</mo></mover><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo></mrow></mstyle></mrow></mstyle></mrow></mfrac><mo>.</mo></mrow></math> In addition, the energy PDF of fission neutrons in a fuel element with indices <math id="M14"><mi>i</mi></math> , j , and <math id="M15"><mi>k</mi></math> is <math display="block" id="M16"><mrow><mi>χ</mi><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo><mo>=</mo><mfrac><mrow><mstyle displaystyle="true"><munder><mo>∑</mo><mi>m</mi></munder><mrow><msub><mi>χ</mi><mi>m</mi></msub><mo stretchy="false">(</mo><mi>E</mi><mo stretchy="false">)</mo><mfrac><mrow><msub><mi>ν</mi><mi>m</mi></msub><msub><mi>σ</mi><mrow><mi>m</mi><mtext>f</mtext></mrow></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo><msub><mi>N</mi><mi>m</mi></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo></mrow><mrow><mstyle displaystyle="true"><munder><mo>∑</mo><mi>m</mi></munder><mrow><msub><mi>ν</mi><mi>m</mi></msub><msub><mi>σ</mi><mrow><mi>m</mi><mtext>f</mtext></mrow></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo><msub><mi>N</mi><mi>m</mi></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo></mrow></mstyle></mrow></mfrac></mrow></mstyle></mrow><mrow><mstyle displaystyle="true"><mrow><mo>∫</mo><mrow><mstyle displaystyle="true"><munder><mo>∑</mo><mi>m</mi></munder><mrow><msub><mi>χ</mi><mi>m</mi></msub><mo stretchy="false">(</mo><mi>E</mi><mo stretchy="false">)</mo><mfrac><mrow><msub><mi>ν</mi><mi>m</mi></msub><msub><mi>σ</mi><mrow><mi>m</mi><mtext>f</mtext></mrow></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo><msub><mi>N</mi><mi>m</mi></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo></mrow><mrow><mstyle displaystyle="true"><munder><mo>∑</mo><mi>m</mi></munder><mrow><msub><mi>ν</mi><mi>m</mi></msub><msub><mi>σ</mi><mrow><mi>m</mi><mtext>f</mtext></mrow></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo><msub><mi>N</mi><mi>m</mi></msub><mo stretchy="false">(</mo><mi>i</mi><mo>,</mo><mi>j</mi><mo>,</mo><mi>k</mi><mo stretchy="false">)</mo></mrow></mstyle></mrow></mfrac></mrow></mstyle></mrow></mrow></mstyle><mtext>d</mtext><mi>E</mi></mrow></mfrac><mo>,</mo></mrow></math> where <math id="M17"><mrow><msub><mi>χ</mi><mi>m</mi></msub><mo stretchy="false">(</mo><mi>E</mi><mo stretchy="false">)</mo></mrow></math> is the fission spectrum corresponding to the m th fissionable nuclide. The fission spectrum of a fissionable nuclide can be characterized by a Watt fission spectrum which can be described as [ 13 13 ] <math display="block" id="M18"><mrow><msub><mi>χ</mi><mi>m</mi></msub><mo stretchy="false">(</mo><mi>E</mi><mo stretchy="false">)</mo><mo>=</mo><mi>C</mi><mi>exp</mi><mo stretchy="false">(</mo><mo>−</mo><mi>E</mi><mo>/</mo><mi>a</mi><mo stretchy="false">)</mo><mi>sinh</mi><msup><mrow><mo stretchy="false">(</mo><mi>b</mi><mi>E</mi><mo stretchy="false">)</mo></mrow><mrow><mn>1</mn><mo>/</mo><mn>2</mn></mrow></msup><mo>,</mo></mrow></math> where <math id="M19"><mi>a</mi></math> and <math id="M20"><mi>b</mi></math> are parameters associated with the fissionable nuclide and incident neutron energy, respectively. 2.2 GVR method 2 2.2.1 Background of GVR 3 To distinguish between different objectives of MC simulation problems, the terms ''source-detector’’, "regional" and "global" are used in this paper to represent problem objectives that are defined as follows. 1) Source-detector: The goal is to calculate the flux or response at a single location. 2) Regional: The goal is to calculate the flux or response over a portion of the problem space. 3) Global: The goal is to calculate the flux or response throughout the entire problem space. For the source-detector and regional problems, one’s interest is in determining a quantity with low statistical uncertainty at a localized region. Numerous VR techniques [ 14 14 ] have been applied successfully, such as the DXTRAN sphere, forced collision techniques, among others. For the global problem, the solution of the whole problem space with uniform low statistical uncertainty is required. Hence, we want to achieve uniform MC particle density over the entire phase space. Many GVR methods have been proposed to solve the global problems, and one of the most successful methods is to employ an adjoint function to speed-up convergence. 2.2.2 Adjoint theory 3 The adjoint function has physical significance as it denotes the importance of a particle to some objective function (e.g., flux, energy deposition, dose, among others) and the main task of most MC simulations is to calculate the objective function at some desired tally regions. Therefore, the adjoint function is well suited for the biasing MC simulations, which leads to the utilization of the adjoint function to generate VR parameters for MC simulations. The specific response at some location can be written as <math display="block" id="M21"><mrow><mi>R</mi><mo>=</mo><mstyle displaystyle="true"><mrow><munderover><mo>∫</mo><mrow><mn>4</mn><mi>π</mi></mrow><mrow/></munderover><mrow><mstyle displaystyle="true"><mrow><munderover><mo>∫</mo><mi>E</mi><mrow/></munderover><mrow><mstyle displaystyle="true"><mrow><munderover><mo>∫</mo><mi>V</mi><mrow/></munderover><mrow><mi>ϕ</mi><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mover accent="true"><mi>Ω</mi><mo>^</mo></mover><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo></mrow></mrow></mstyle></mrow></mrow></mstyle></mrow></mrow></mstyle><msub><mi>σ</mi><mtext>d</mtext></msub><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mover accent="true"><mi>Ω</mi><mo>^</mo></mover><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo><mtext>d</mtext><mi>r</mi><mtext>d</mtext><mover accent="true"><mi>Ω</mi><mo>^</mo></mover><mtext>d</mtext><mi>E</mi><mo>,</mo></mrow></math> where <math id="M22"><mrow><msub><mi>σ</mi><mtext>d</mtext></msub></mrow></math> denotes some response function and <math id="M23"><mi>ϕ</mi></math> is the particle flux. Based on the identity of the adjoint [ 15 15 ] and letting <math id="M24"><mrow><msub><mi>σ</mi><mtext>d</mtext></msub><mo>=</mo><msup><mi>s</mi><mo>†</mo></msup></mrow></math> , the following equivalent expression can be derived [ 11 11 ] , <math display="block" id="M25"><mrow><mi>R</mi><mo>=</mo><mstyle displaystyle="true"><mrow><munderover><mo>∫</mo><mrow><mn>4</mn><mi>π</mi></mrow><mrow/></munderover><mrow><mstyle displaystyle="true"><mrow><munderover><mo>∫</mo><mi>E</mi><mrow/></munderover><mrow><mstyle displaystyle="true"><mrow><munderover><mo>∫</mo><mi>V</mi><mrow/></munderover><mrow><mi>ϕ</mi><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mover accent="true"><mi>Ω</mi><mo>^</mo></mover><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo></mrow></mrow></mstyle></mrow></mrow></mstyle></mrow></mrow></mstyle><msup><mi>s</mi><mo>†</mo></msup><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mover accent="true"><mi>Ω</mi><mo>^</mo></mover><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo><mtext>d</mtext><mi>r</mi><mtext>d</mtext><mover accent="true"><mi>Ω</mi><mo>^</mo></mover><mtext>d</mtext><mi>E</mi><mo>=</mo><mstyle displaystyle="true"><mrow><munderover><mo>∫</mo><mrow><mn>4</mn><mi>π</mi></mrow><mrow/></munderover><mrow><mstyle displaystyle="true"><mrow><munderover><mo>∫</mo><mi>E</mi><mrow/></munderover><mrow><mstyle displaystyle="true"><mrow><munderover><mo>∫</mo><mi>V</mi><mrow/></munderover><mrow><msup><mi>ϕ</mi><mo>†</mo></msup><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mover accent="true"><mi>Ω</mi><mo>^</mo></mover><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo></mrow></mrow></mstyle></mrow></mrow></mstyle></mrow></mrow></mstyle><mi>s</mi><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mover accent="true"><mi>Ω</mi><mo>^</mo></mover><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo><mtext>d</mtext><mi>r</mi><mtext>d</mtext><mover accent="true"><mi>Ω</mi><mo>^</mo></mover><mi mathvariant="normal">d</mi><mi>E</mi><mo>,</mo></mrow></math> where <math id="M26"><mi>s</mi></math> stands for a fixed source and the superscript <math id="M27"><mrow><msup><mrow/><mo>†</mo></msup></mrow></math> indicates adjoint quantities. From Eq. (6) one can conclude that the adjoint function represents the expected contribution (i.e., importance) to the response R from a source particle in some phase space. Therefore, the adjoint function can be used to develop the variance reduction parameters. However, the space-, energy-, and angular-dependent adjoint function may generate a huge amount of variance reduction parameters (tens or hundreds of GB) which would be costly to use in the MC simulation. Consequently, the scalar adjoint function and fixed source are used in the GVR methodology as follows: <math display="block" id="M28"><mrow><mrow><mo>{</mo><mtable columnalign="left"><mtr><mtd><msup><mi>φ</mi><mo>†</mo></msup><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo><mo>=</mo><mstyle displaystyle="true"><mrow><munder><mo>∫</mo><mrow><mn>4</mn><mi>π</mi></mrow></munder><mrow><msup><mi>ϕ</mi><mo>†</mo></msup><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mover accent="true"><mi>Ω</mi><mo>^</mo></mover><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo></mrow></mrow></mstyle><mtext>d</mtext><mover accent="true"><mi>Ω</mi><mo>^</mo></mover></mtd></mtr><mtr><mtd><mi>S</mi><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo><mo>=</mo><mstyle displaystyle="true"><mrow><munder><mo>∫</mo><mrow><mn>4</mn><mi>π</mi></mrow></munder><mrow><msup><mi>s</mi><mo>†</mo></msup><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mover accent="true"><mi>Ω</mi><mo>^</mo></mover><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo></mrow></mrow></mstyle><mtext>d</mtext><mover accent="true"><mi>Ω</mi><mo>^</mo></mover></mtd></mtr></mtable></mrow><mo>.</mo></mrow></math> 2.2.3 Theoretical model of FW-CADIS 3 The FW-CADIS method is derived from adjoint theory. In the FW-CADIS methodology, the adjoint source is some response function <math id="M29"><mrow><msub><mi>σ</mi><mtext>d</mtext></msub></mrow></math> weighted with the inverse of the forward flux: <math display="block" id="M30"><mrow><msup><mi>S</mi><mo>†</mo></msup><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo><mo>=</mo><mfrac><mrow><msub><mi>σ</mi><mtext>d</mtext></msub><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo></mrow><mrow><mstyle displaystyle="true"><mrow><munder><mo>∫</mo><mi>E</mi></munder><mrow><mi>φ</mi><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo><msub><mi>σ</mi><mtext>d</mtext></msub><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo></mrow></mrow></mstyle><mtext>d</mtext><mi>E</mi></mrow></mfrac><mo>.</mo></mrow></math> Intuitively, this is useful for increasing the adjoint importance where the forward flux is low and decreasing the adjoint importance where the forward flux is high. Using the adjoint source mentioned above, the adjoint function <math id="M31"><mrow><msup><mi>φ</mi><mo>†</mo></msup><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo></mrow></math> can be derived [ 11 11 ] . Subsequently, the following equations are used to create GVR parameters: <math display="block" id="M32"><mrow><mover accent="true"><mi>S</mi><mo>^</mo></mover><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo><mo>=</mo><mfrac><mrow><msup><mi>φ</mi><mo>†</mo></msup><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo><mi>S</mi><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo></mrow><mrow><mstyle displaystyle="true"><mrow><munderover><mo>∫</mo><mi>E</mi><mrow/></munderover><mrow><mstyle displaystyle="true"><mrow><munderover><mo>∫</mo><mi>V</mi><mrow/></munderover><mrow><msup><mi>φ</mi><mo>†</mo></msup><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo><mi>S</mi><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo></mrow></mrow></mstyle></mrow></mrow></mstyle><mtext>d</mtext><mi>r</mi><mtext>d</mtext><mi>E</mi></mrow></mfrac><mo>=</mo><mfrac><mrow><msup><mi>φ</mi><mo>†</mo></msup><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo><mi>S</mi><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo></mrow><mi>R</mi></mfrac><mo>,</mo></mrow></math> <math display="block" id="M33"><mrow><msub><mi>w</mi><mn>0</mn></msub><mo>=</mo><mfrac><mi>S</mi><mover accent="true"><mi>S</mi><mo>^</mo></mover></mfrac><mo>=</mo><mfrac><mi>R</mi><mrow><msup><mi>φ</mi><mo>†</mo></msup><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo></mrow></mfrac><mo>,</mo></mrow></math> <math display="block" id="M34"><mrow><mi>w</mi><mo>=</mo><mfrac><mi>R</mi><mrow><msup><mi>φ</mi><mo>†</mo></msup><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo></mrow></mfrac><mo>.</mo></mrow></math> The biased source distribution is given by Eq. (9) which is derived from Eq. (6) and the concept of importance sampling [ 16 16 ] . The starting source weight of particles is expressed by Eq. (10). For the transport process, the objective statistical weight of the particles is given by Eq. (11). Notice that Eq. (10) and Eq. (11) have the same form, which means the source particles are created with a weight that is consistent with the objective weight for the position they are born into. This is an important feature of the FW-CADIS method because it eliminates the inconsistency between source biasing and transport biasing, and avoids decreases in the computing efficiency due to the waste of computational resources on unnecessary splitting/rouletting of the starting particles. 3 Implementation of Sub-GVR method 1 A flowchart of the Sub-GVR method is presented in Fig. 1 Fig. 1 . In the first step, a criticality calculation is executed with the eigenvalue mode (KCODE) of MCNP5 and an FMESH tally is used to record the effective microscopic fission cross-section for each fissionable nuclide and the total fission rate in every mesh cell. After this, a utility python code SN2MCNP is employed to covert the mesh tally (MESHTAL file) to the space- and energy-dependent fission neutron PDF according to Eq. (2), Eq. (3), and Eq. (4). The second step is divided into four sub-steps to implement GVR. The first sub-step employs a forward-mode DENOVO calculation using the fission neutron PDF derived from the previous step as a neutron source to generate the forward flux. It should be noted that the forward-mode DENOVO calculation step is intended to acquire an adjoint source to be used for S N adjoint-mode calculation, rather than to provide a solution to the problem. Therefore, a relatively coarse mesh is used in the DENOVO model. For the second sub-step, an adjoint source is produced according to Eq. (8) and subsequently an adjoint-mode DENOVO calculation is executed, and the scalar adjoint flux is produced. In sub-step 3, SN2MCNP is used again to read the scalar adjoint function from the DENOVO flux file, and then combine these values with the fixed neutron source used in the forward-mode DENOVO calculation to calculate the response according to Eq. (6). SN2MCNP is also applied to convert the scalar adjoint function into the source biasing and weight window parameters. For the source biasing parameters, SN2MCNP outputs a mesh-based biased source distribution file (based on Eq. (9)). For the weight window parameters, MCNP5 requires lower weight window bounds (WWBs) instead of the center of the weight windows (intervals) as given in Eq. (11). In MCNP5, the width of the window is controlled by the weight window parameter card (WWP card); the parameter <math id="M35"><mi>c</mi></math> as the ratio of upper and lower WWBs ( <math id="M36"><mrow><mi>c</mi><mo>=</mo><msub><mi>w</mi><mtext>u</mtext></msub><mo>/</mo><msub><mi>w</mi><mtext>l</mtext></msub></mrow></math> ) can be input in the WWP card. Therefore, the space- and energy-dependent lower WWBs for MCNP5 are given by <math display="block" id="M37"><mrow><msub><mi>w</mi><mtext>l</mtext></msub><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo><mo>=</mo><mfrac><mi>R</mi><mrow><msup><mi>φ</mi><mo>†</mo></msup><mo stretchy="false">(</mo><mi>r</mi><mo>,</mo><mi>E</mi><mo stretchy="false">)</mo><mo stretchy="false">(</mo><mfrac><mrow><mi>c</mi><mo>+</mo><mn>1</mn></mrow><mn>2</mn></mfrac><mo stretchy="false">)</mo></mrow></mfrac><mo>.</mo></mrow></math> As a consequence, SN2MCNP uses Eq.(12) to generate the standard MCNP5 weight windows file (WWINP file). In sub-step 4, the user-provided SOURCE subroutine in the MCNP5 source code is employed to apply the source biasing technique to the massive mesh-based source, which cannot be dealt with by the conventional method (using general source specification cards). Besides, the fission neutron source region of TMSR-SF1 has a complex geometry which is inconsistent with the geometry described by a mesh-based source. Therefore, this inconsistency must also be dealt with in the user-provided SOURCE subroutine. The algorithms used to implement the source sampling methodology in the SOURCE subroutine have four steps, as follows. First, a mesh cell is sampled directly from the mesh-based biased source distribution file. Second, the starting positions of particles are sampled uniformly within the mesh cell. The particle source region is described by a regular grid (cylindrical or rectangular grid) and therefore the area covered by the grid may exceed the actual active core. However, we can resample the starting positions of particles until they are strictly located in the active core or the resampling time exceeds the limit number N (defined by IDUM card). Third, starting energy is sampled according to the biased energy distribution of the chosen mesh cell. In addition, the flight directions of the starting particles are isotropic by default in the SOURCE subroutine. Finally, the statistical weight of the starting particles is modified accordingly based on Eq. (10). The source sampling algorithm is implemented in the Source subroutine. The flowchart of the Source subroutine is shown in Fig. 2 Fig. 2 . 4 TMSR-SF1 facility model description 1 The 10-MWth TMSR-SF1 is a graphite-moderated, salt-cooled thermal reactor. In this design, 6-cm-diameter pebbles fill a cylindrical cavity with upper and lower circular truncated cones. The diameter and height of the cylindrical cavity are 135.0 cm and 180.0 cm, respectively. The minimum diameter of the upper and lower circular truncated cones is 30.0 cm. There are two types of pebbles loaded into the core cavity, namely fuel pebbles containing thousands of TRISO which are located in the top area, and graphite pebbles placed in the bottom area. The pebbles are cooled by the 2LiF-BeF 2 mixed salt (FLiBe) that is kept at an operating temperature between 625 °C and 650 °C. As shown in Fig. 3 Fig. 3 , the core cavity is surrounded by thick graphite reflectors radially and axially. The other major components of TMSR-SF1 include a Hastelloy pressure vessel, aluminum silicate fiber thermal insulator, concrete biological shield, FLiBe plenum, boronated carbon bricks, stainless steel support plate, and stainless steel thermal shield. Table 1 Table 1 gives the geometric dimensions of the facility model. 1 1 Geometry 1 1 Size (cm) 1 1 Core radius 1 1 67.5 1 1 Core height 1 1 240.6 1 1  Reflector  radius 1 1 142.5 1 1 Reflector height 1 1 306.4 1 1 Pressure vessel radius (inner) 1 1 149.3 1 1 Pressure vessel radius (outer) 1 1 152.8 1 1 Pressure vessel height 1 1 681.4 1 1 Thermal insulator radius (inner) 1 1 270.0 1 1 Thermal insulator radius (outer) 1 1 310.0 1 1 Thermal insulator height 1 1 754.7 1 1 Concrete shield radius (inner) 1 1 310.0 1 1 Concrete shield radius (outer) 1 1 410.0 1 1 Concrete shield height 1 1 860.0 1 1 Stainless steel thermal shield radius 1 1 310.0 1 1 Stainless steel thermal shield thickness 1 1 35.0 1 1 Boronated carbon bricks shield radius 1 1 310.0 1 1 Boronated carbon bricks shield thickness 1 1 50.0 1 1 Steel support plate radius 1 1 149.3 1 1 Steel support plate thickness 1 1 21.5 1 1 Upper FLiBe plenum radius 1 1 67.5 1 1 Upper FLiBe plenum thickness 1 1 40.0 1 1 Upper Argon plenum radius 1 1 67.5 1 1 Upper Argon plenum thickness 1 1 60.0 1 1 Lower FiLiBe plenum radius 1 1 146.3 1 1 Lower FiLiBe plenum thickness 1 1 14.0 5 Results and discussion 1 To calculate the detailed spatial distribution of the fast neutron flux, a fine rectangular mesh tally was used in this research. All single mesh cells had a uniform size of 4.0 cm × 4.0 cm × 4.0 cm, and the X-Y-Z dimensions of the MCNP5 model of TMSR-SF1 were 820.0 cm × 820.0 cm × 860.0 cm. Therefore, the total number of mesh cells is 9,035,375. In addition, all our calculations were performed on a tower server with an Intel Xeon Processor E7-4830 (12-core, 2.1 GHz). A one-step eigenvalue mode analog MCNP5 calculation was implemented initially so as to make a comparison. This calculation was performed with 6950 active k eff -cycles each with 2×10 5 particles per cycle resulting in a total of 1.39×10 9 histories. The resulting fast neutron flux distribution and associated relative uncertainties are plotted in Fig. 4 Fig. 4 . As shown in Fig. 4a Fig. 4a , most regions have no tally scores, while the nonzero scored mesh cells are mainly concentrated within the pressure vessel because of the existing of large amount neutron-moderator material. Low relative uncertainties of the tallies were obtained near the core (source region). The mesh cells far from the core have very few scores or no scores at all. The first step of Sub-GVR method employed a mesh tally including 10,830 orthometric mesh cells and the size of each cell was 7.224 cm × 7.224 cm × 7.224 cm. In this work, only 235 U and 238 U needed to be considered in the fission neutron source definition, since the calculation was performed at the beginning-of-cycle condition. TMSR-SF1 is a low-power experimental reactor with a relatively low burn-up level which makes it impossible to use thorium on a large scale. Currently, only trace-level thorium is used as an irradiation material for radiation experimental research and the fuel in the fuel pebbles is UO 2 (17.0 wt% 235 U enrichment). Therefore, only 235 U and 238 U need to be accounted for in the neutron source definition, since all other fissionable nuclides make only negligible contributions to the fission in TMSR-SF1. The data for the fissionable nuclides 235 U and 238 U, including values for the parameters <math id="M38"><mi>a</mi></math> and <math id="M39"><mi>b</mi></math> in the Watt fission spectrum function and the number of neutrons per fission, are summarized in Table 2 Table 2 . The parameters <math id="M40"><mi>a</mi></math> and <math id="M41"><mi>b</mi></math> for 235 U correspond to fission reactions caused by thermal neutrons [ 17 17 ] , whereas the fission reactions of 238 U are almost all caused by fast neutrons. Therefore, the parameters corresponding to fast incident neutrons are used [ 18 18 ] . Figure 5 Figure 5 illustrates the spectral difference for the Watt fission spectra from the two nuclides. 1 1 Nuclide 1 1 a (MeV) 1 1 b (MeV -1 ) 1 1 ν 1 1 235 U 1 1 0.988 1 1 2.249 1 1 2.430 ± 0.008 1 1 238 U 1 1 0.920 1 1 3.121 1 1 2.810 ± 0.020 For our application of the Sub-GVR method, the DENOVO S N grid had approximately 1,113,228 cells over the TMSR-SF1 facility model, that is, 102 × 102 × 107 cells in the X-Y-Z directions. A relatively low-order S4 angular quadrature and P3 Legendre scattering cross-section expansion approximation (upscattering was deactivated) were employed to reduce the computational resource requirements. The forward and adjoint DENOVO calculations took 28.17 and 20.34 min of wall-clock time, respectively. The details of the calculated space PDF of the fission neutron source and the corresponding space PDF of the biased source are shown in Fig. 6 Fig. 6 . Compared with the PDF of the fission source, the source particles in the external region of the core will be sampled more frequently when using the biased source, which encourages the source particles to have a larger escape probability from the core. In addition, it should be noted that some of the areas in mesh cells are outside the core, which can be ignored by the sampling algorithm in the SOURCE subroutine. Fig. 7 Fig. 7 illustrates the energy PDF of the fission neutron source and the biased source in an arbitrarily selected mesh cell in the core. The energy PDF of the fission neutron source derived from Eq. (3) is approximately equal to the Watt fission spectrum of 235 U. For the biased source, the sampling probabilities in the high-energy region are increased. Furthermore, the sampling probability for the energy region below 0.1 MeV is set to zero since the source at this region makes no contribution to the fast neutron flux. The fixed-source MCNP5 run in the second step utilized 1×10 9 particle histories and required 8114.63 min of calculation time. Plots of the calculated fast neutron flux mesh tally and associated relative uncertainties are shown in Fig. 8 Fig. 8 . A dramatic improvement in the mesh tally result was obtained by the Sub-GVR method. Extremely uniform and low relative uncertainties were obtained throughout the entire TMSR-SF1 facility model. Fig. 9 Fig. 9 shows the detailed distributions of relative uncertainties for the one-step analog and Sub-GVR approaches. Fig. 9a Fig. 9a shows the PDF of relative uncertainty, giving the fraction of mesh cells that have relative uncertainties at certain values, and Fig. 9b Fig. 9b presents the cumulative distribution functions (CDF) of relative uncertainty, giving the fraction of mesh cells that have relative uncertainties below certain values. As shown in Fig. 9a Fig. 9a , the one-step analog calculation yielded relative uncertainties varying uniformly in the range of 0.0–0.98 and rising to 0.98–1.0 radically. For the Sub-GVR scheme, the relative uncertainties were mainly distributed in the range 0.0–0.1, and the fraction of relative uncertainties larger than 10% was almost zero. The CDF graph, presented in Fig. 9b Fig. 9b , is particularly useful in illustrating the performance of the two methods. With the one-step analog method, only 2.5, 3.3, and 4.2% of the total mesh cells had relative uncertainty values less than 5, 10, and 20%, respectively, while with the Sub-GVR method, 93.4, 99.7, and 99.9% of the mesh cells had relative uncertainties below these levels. In conclusion, the Sub-GVR method produces a fast flux mesh tally with far more uniform and lower relative uncertainties than the one-step analog method, thus dramatically improving the computational efficiency for this problem. To quantitatively compare the computational efficiencies of the two approaches, the mesh tally figure-of-merit (FOM) is introduced [ 19 19 ] , which can be defined as <math display="block" id="M42"><mrow><mi>F</mi><mi>O</mi><mi>M</mi><mo>=</mo><mfrac><mn>1</mn><mrow><msup><mover accent="true"><mi>r</mi><mo>¯</mo></mover><mn>2</mn></msup><mi>T</mi></mrow></mfrac><mo>,</mo></mrow></math> where <math id="M43"><mover accent="true"><mi>r</mi><mo>¯</mo></mover></math> is the mean relative uncertainty of the mesh tally cells and <math id="M44"><mi>T</mi></math> is the total execution time. In the Sub-GVR method, <math id="M45"><mi>T</mi></math> is the sum of total DENOVO execution time and MCNP5 execution time. The ratio of the FOM from Sub-GVR and the one-step analog calculation can be defined as speed-up [ 20 20 ] and is used to quantify the efficiency improvements obtained by using the Sub-GVR method. Table 3 Table 3 gives detailed information of the two calculation cases. The FOM obtained from applying the Sub-GVR method is 1.11e-1 while that of the analog approach is 1.37e-4, which provides a speed-up of 810.2. 1 1 Method 3 1 Execution Time (min) 1 1 Nonzero fraction (c) (%) 1 1 Mean Rel. Unc. (%) 1 1 FOM 1 1 Speed-up 1 1 1 1 Step 1 1 1 Step 2 5 1 Total 1 1 One-step Analog 1 1 10416.36 1 1 N/A (a) 1 1 10416.36 1 1 38.10 1 1 83.56 1 1 1.37E-4 1 1 1.0 1 1 Sub-GVR 1 1 5832.25 1 1 8163.14 (b) 1 1 13995.39 1 1 100.00 1 1 2.54 1 1 1.11E-1 1 1 810.2 (a) Not applicable. (b) Sum of the DENOVO execution time and the MCNP5 execution time. (c) The fraction of mesh cells with a nonzero score. The maximum fast flux results of the main components in TMSR-SF1 were also calculated by the Sub-GVR method. Fig. 10a Fig. 10a gives the maximum fast flux of the main components. From the maximum-flux region, which includes the center of the core, to the minimum-flux region, which includes the stainless steel thermal shield, the fast flux is attenuated by approximately nine orders of magnitude. Therefore, this is a typical deep-penetration problem that is difficult to calculate by analog MC simulation, as demonstrated above. However, by applying the Sub-GVR method, results with uniform and low relative uncertainties were obtained, as shown in Fig. 10b Fig. 10b . The relative uncertainties of all maximum fast flux results were less than 5%, which is considered as generally reliable based on years of experience using MCNP5 on a wide variety of problems [ 21 21 ] . 6 Conclusion 1 A two-step subsection approach involving a GVR technique, Sub-GVR, has been implemented and applied in the calculation of a detailed fast neutron distribution throughout a full-scale model of TMSR-SF1. The results obtained using the Sub-GVR method and an analog method were compared, and the former method was demonstrated to be extremely effective in producing the fast neutron flux mesh tally with relatively uniform and low relative uncertainties. The FOM value of the Sub-GVR method was 810.2 times greater than that of the analog method. Furthermore, the maximum fast flux results of the main components in TMSR-SF1 were given, which can be used for the radiation damage assessment of the structural materials. In conclusion, the utilization of the Sub-GVR method has efficiently improved the statistical certainty and considerably reduced the global variance of the tally results through the entire geometry, and enabled the calculation of high-fidelity distributions of the fast neutron flux throughout the full-scale TMSR-SF1 facility. References C.W. Forsberg , P.F. Peterson , P.S. Pickard . Molten-Salt-Cooled Advanced High-Temperature Reactor for Production of Hydrogen and Electricity . Nucl. Tech . 144 ( 3 ), 289 - 302 ( 2003 ). doi: 10.13182/nt03-1 http://doi.org/10.13182/nt03-1 Molten-Salt-Cooled Advanced High-Temperature Reactor for Production of Hydrogen and Electricity M. Fratoni and E. Greenspan . Neutronic Feasibility Assessment of Liquid Salt-Cooled Pebble Bed Reactors . Nucl. Sci. Eng . 168 ( 1 ), 1 - 22 ( 2011 ). doi: 10.13182/nse10-38 http://doi.org/10.13182/nse10-38 Neutronic Feasibility Assessment of Liquid Salt-Cooled Pebble Bed Reactors M.H. Jiang , H.J. Xu , Z.M. Dai , Advanced fission energy program—TMSR nuclear energy system . Bull. Chin. Acad. Sci . 27 , 366 - 374 ( 2012 ). doi: 10.3969/j.issn.1000-3045.2012.03.016 http://doi.org/10.3969/j.issn.1000-3045.2012.03.016 Advanced fission energy program—TMSR nuclear energy system TMSR group , TMSR-SF1 conceptual design overview, Chapter 2 , Shanghai Institute of Applied Physics , Shanghai, China ( 2015 ). TMSR-SF1 conceptual design overview, Chapter 2 R.M. Ji , M.H. Li , Y. Zou , et al . Impact of photoneutrons on reactivity measurements for TMSR-SF1 . Nucl. Sci. Tech . 28 , 76 . ( 2017 ). doi: 10.1007/s41365-017-0234-7 http://doi.org/10.1007/s41365-017-0234-7 Impact of photoneutrons on reactivity measurements for TMSR-SF1 R.M. Ji , Y. Dai , G.F. Zhu , et al . Evaluation of the fraction of delayed photoneutrons for TMSR-SF1 . Nucl. Sci. Tech . 28 , 135 . ( 2017 ). doi: 10.1007/s41365-017-0285-9 http://doi.org/10.1007/s41365-017-0285-9 Evaluation of the fraction of delayed photoneutrons for TMSR-SF1 M. A. COOPER and E. W. LARSEN , Automated Weight Windows for Global Monte Carlo Particle Transport Calculations . Nucl. Sci. Eng . 137 , 1 ( 2001 ). doi: 10.13182/nse00-34 http://doi.org/10.13182/nse00-34 Automated Weight Windows for Global Monte Carlo Particle Transport Calculations A. Davis , A. Turner . Comparison of global variance reduction techniques for Monte Carlo radiation transport simulations of ITER . Fusion Eng. Des . 86 , 2698 - 2700 . ( 2011 ). doi: 10.1016/j.fusengdes.2011.01.059 http://doi.org/10.1016/j.fusengdes.2011.01.059 Comparison of global variance reduction techniques for Monte Carlo radiation transport simulations of ITER X. C. Nie , J. Li , S.L. Liu , et al . Global variance reduction method for global Monte Carlo particle transport simulations of CFETR . Nucl. Sci. Tech . 28 , 115 ( 2017 ). doi: 10.1007/s41365-017-0270-3 http://doi.org/10.1007/s41365-017-0270-3 Global variance reduction method for global Monte Carlo particle transport simulations of CFETR T. Shi , J. Ma , H. Huang , et al . A new global variance reduction technique based on pseudo flux method . Nucl. Eng. Des . 324 : 18 - 26 . ( 2017 ). doi: 10.1016/j.nucengdes.2017.08.001 http://doi.org/10.1016/j.nucengdes.2017.08.001 A new global variance reduction technique based on pseudo flux method J.C. Wagner , E.D. Blakeman , D.E. Peplow . Forward-Weighted CADIS Method for Global Variance Reduction . Trans. Am. Nucl. Soc . 97 , 630 - 633 ( 2007 ). Forward-Weighted CADIS Method for Global Variance Reduction T.M. Evans , A.S. Stafford , K.T. Clarno . Denovo—A New Three-Dimensional Parallel Discrete Ordinates Code in SCALE . Nucl. Tech . 171 , 171 - 200 ( 2010 ). doi: 10.13182/nt10-05 http://doi.org/10.13182/nt10-05 Denovo—A New Three-Dimensional Parallel Discrete Ordinates Code in SCALE J.F. Briesmeister , (Ed.), 2000 . MCNPTM– A General Monte Carlo N-Particle Transport Code . Version 4C , LA-13709-M , Manual . Los Alamos National Laboratory . MCNPTM– A General Monte Carlo N-Particle Transport Code T.E. Booth , R.A. Forster , R.L. Martz . MCNP variance reduction developments in the 21st century . Nucl. Tech . 180 ( 3 ), 355 - 371 ( 2012 ). doi: 10.13182/nt12-a15349 http://doi.org/10.13182/nt12-a15349 MCNP variance reduction developments in the 21st century G.I. Bell and S. Glasstone , Nuclear Reactor Theory , ( Van Nostrand and Reinhold , New York 1967 ) Nuclear Reactor Theory M.H. Kalos and P.A. Whitlock , Monte Carlo Methods—Volume 1: Basics , John Wiley & Sons , New York ( 1986 ) Monte Carlo Methods—Volume 1: Basics L. Cranberg , G. Frye , N. Nereson , L. Rosen . Fission Neutron Spectrum of 235U , Phys. Rev . 103 , 662 - 670 ( 1956 ). doi: 10.1103/PhysRev.103.662 http://doi.org/10.1103/PhysRev.103.662 Fission Neutron Spectrum of 235U A. Vasiliev et al . Development of a CASMO-4/SIMULATE-3/MCNPX Calculation Scheme for PWR Fast Neutron Fluence Analysis and Validation against RPV Scraping test data . Annals of Nuclear Energy , 34 : 615 - 627 ( 2007 ). doi: 10.1016/j.anucene.2007.02.020 http://doi.org/10.1016/j.anucene.2007.02.020 Development of a CASMO-4/SIMULATE-3/MCNPX Calculation Scheme for PWR Fast Neutron Fluence Analysis and Validation against RPV Scraping test data J.C. Wagner , D.E. Peplow , S.W. Mosher . FW-CADIS Method for Global and Regional Variance Reduction of Monte Carlo Radiation Transport Calculations . Nucl. Sci. Eng . 176 ( 1 ), 37 - 57 ( 2014 ). doi: 10.13182/nse12-33 http://doi.org/10.13182/nse12-33 FW-CADIS Method for Global and Regional Variance Reduction of Monte Carlo Radiation Transport Calculations A. Haghighat , H. Hiruta , B. Petrovic , J. C. Wagner . Performance of the Automated Adjoint Accelerated MCNP (A3 MCNP) for Simulation of a BWR Core Shroud Problem . Proc. International Conference on Mathematics and Computation, Reactor Physics, and Environmental Analysis in Nuclear Applications , Madrid, Spain , Sep. 27–30, (1999 ). Performance of the Automated Adjoint Accelerated MCNP (A3 MCNP) for Simulation of a BWR Core Shroud Problem X-5 Team Monte Carlo . MCNP - A General Monte Carlo N-Particle Transport Code , Version 5 Volume I : Overview and Theory ( 2003 ). Los Alamos, NM : Los Alamos National Laboratory . Overview and Theory 10.1007/s41365-019-0650-y This work was supported by the Chinese TMSR Strategic Pioneer Science and Technology Project (No.XDA02010000) and the Frontier Science Key Program of Chinese Academy of Sciences (No.QYZDY-SSW-JSC016). 2019-01-10 2019-03-25 2019-04-28 2019-07-12 2019-08-01 2021-04-11 © China Science Publishing & Media Ltd. (Science Press), Shanghai Institute of Applied Physics, the Chinese Academy of Sciences, ChineseNuclear Society and Springer Nature Singapore Pte Ltd. 2019 This is a post-peer-review, pre-copyedit version of an article published in Nuclear Science and Techniques. The final authenticated version is available online at: http://dx.doi.org/10.1007/s41365-019-0650-y http://dx.doi.org/10.1007/s41365-019-0650-y . 2019 Nuclear Science and Techniques 08 30