Coarse Mesh Finite Difference Research Articles

Methods and performance of a GPU-based pinwise two-step nodal code VANGARD

The methods and performance of a GPU-based pinwise two-step nodal core calculation code VANGARD are presented which is targeted for use in the commercial nuclear design analyses. In order to realize high-speed core cycle depletion with a single consumer-grade GPU mounted on ordinary personal computers, various elaborated computational methods are introduced as follows. The computationally intense pinwise nodal calculation is optimized by adaptive employment of low-order expansion. The two-level surface current correction method is used for updating pinwise incoming partial currents within the assembly-level coarse mesh finite difference framework. A cross section data compression technique is introduced to fit the enormous pinwise group constant data into the limited GPU memory. A massively parallelized depletion scheme based on the Chebyshev Rational Approximation Method is implemented. The neighbor-informed burnup correction method is developed to ensure the accuracy of the pinwise depletion for gadolinia-bearing pins. Performance examinations are performed for the realistic APR1400 and AP1000 core problems to demonstrate that three-dimensional pinwise cycle depletion can be completed within 3 minutes with high accuracy and resolution.

High-fidelity neutronics simulation of channel-type molten salt reactors

As one of the generation IV reactors, molten salt reactors (MSRs) have attracted much attention in recent years. Generally, molten salt acting as both fuel and coolant circulates through graphite channels and external primary loop in a liquid-fueled MSR with a channel-type core. To accurately simulate the neutronics of channel-type MSRs, a GPU paralleled high-fidelity neutron transport code ThorMOC has been improved by implementing a 3D rasterized coarse mesh finite difference (RCMFD) method and a quasi-2D delayed neutron precursor (DNP) drift model. A Molten Salt Reactor Experiment (MSRE) model and a small Thorium-based MSR (sTMSR) model with stationary and flowing fuel are used for the verification of the improved ThorMOC. The reference results of the two stationary reactor models are obtained from Monte Carlo code OpenMC, which indicates that ThorMOC performs a high accuracy for stationary MSR simulations. The reactivity loss caused by fuel flowing in the steady state MSRE is calculated by ThorMOC and agrees well with the measured value from the MSRE. Moreover, the influences of geometry simplifications of fuel channels and downcomer on the reactivity loss are evaluated by comparing their results with those of the high-fidelity model. The reactivity loss in the sTMSR core is computed, and the effects of fuel velocity in core channels and residence time outside the core on the reactivity loss are also analyzed. The whole-core transport calculations of the 3D MSRE and sTMSR core models with 8 energy groups cost several hundred seconds.

Implementation and optimization of hybrid parallel strategy in HNET

For the self-developed three-dimensional whole-core High-fidelity NEutron Transport calculation program HNET, although the numerical acceleration algorithms can improve the computation performance of high-fidelity neutron transport in terms of algorithms and models, it still faces some critical issues such as long computing time and enormous memory requirement. The rapid development of high-performance clusters provides a foundation for the application of massively parallel computing. Most current MOC programs are based on a single-type variable to achieve efficient parallelism, and spatial domain decomposition methods are the most common parallel schemes. However, its parallelism is limited and cannot fully utilize the current state-of-the-art computer resources. To solve this problem, hybrid parallel strategies are implemented in HNET to further expand the parallel degree, improve the speed of computation, and reduce the memory requirement. A hybrid MPI/OpenMP method based on domain decomposition and characteristic rays is proposed for the method of characteristics (MOC). For domain decomposition, the simulation domain is divided into spatial subdomains, with each subdomain handled by different processes. On this basis, the characteristic ray parallelism is implemented taking advantage of the inherent parallelism of the characteristic rays. Meanwhile, the optimization of the hybrid parallel strategy further improves the computation speed by eliminating atomic operations and using private pointer arrays and other techniques. In addition, in the framework of generalized equivalence theory (GET) based two-level coarse mesh finite difference method (CMFD), there is also a certain time consumption in solving CMFD linear system. Hence, for CMFD, a hybrid MPI/OpenMP method based on domain decomposition and secondary domain decomposition can be used. Using the secondary domain decomposition method, each subdomain is divided into sub-subdomains, with each sub-subdomain handled by threads, which allows CMFD to utilize the resources of characteristic ray parallelism in MOC, and also increases the speed of CMFD computation. Numerical results show that for both steady-state and transient calculations, the hybrid MPI/OpenMP in HNET can further expand the parallelism and accelerate the computation. It can take full use of parallel resources and achieve large-scale parallelism.

The NILO-CMFD Method for Iteratively Solving Coupled Neutron Transport–Thermal Hydraulics Problems

The coarse-mesh finite difference (CMFD) method is commonly used to accelerate the iterative convergence of single-physics neutron transport problems. For multiphysics problems, the neutron cross sections depend on the temperature and density, both of which depend on the fission heat source; the resulting nonlinear feedback can significantly degrade the performance of CMFD and even cause instability. In this paper, we propose, for a class of one-dimensional (1-D) model multiphysics problems, a new nonlinearly implicit low-order (NILO) CMFD (NILO-CMFD) acceleration method to improve the performance of CMFD-based methods for solving loosely coupled multiphysics problems. Our numerical testing and Fourier analysis show that for the 1-D model problems, the new NILO-CMFD method achieves the same rapid convergence rate that CMFD achieves for single-physics problems.

The progress on the 3D whole-core neutronics calculation of PANGU code

PANGU code is a modern computer code for the physics analysis of pebble-bed high temperature gas-cooled reactors (HTGRs). This paper presents the recent progress on the development of PANGU’s capability in three-dimensional (3D) whole-core neutronics calculations. Firstly, diffusion and SP3 solvers based on nodal Green function method (NGFM) have been developed for the 3D whole-core calculation. Secondly, the 3D discrete-ordinate (Sn) neutron transport solver with coarse mesh finite difference (CMFD) acceleration has been developed to provide reference solutions, as well as being used for pre-homogenization calculations. Thirdly, two methods have been proposed for treating the strong absorbers in the framework of 3D whole-core diffusion calculations. Finally, the performances of the newly developed methods and solvers are tested by typical numerical cases based on HTR-10.

Development and verification of the hexagonal core simulation capability of the NECP-X code

The whole-core high-fidelity calculation capability for hexagonal cores is developed in the NECP-X code for VVER and other hexagonal reactors, including the 2D/1D based transport and Global-Local resonance calculation methods. Furthermore, the spatial domain decomposition-based parallel strategy and the unstructured coarse-mesh finite difference (CMFD) acceleration method is developed to improve the efficiency, which can realize rectangular, hexagonal, and other polygon grid-based pin-by-pin acceleration. The C5G7 3D hexagonal benchmark problem is calculated to test the performance of the transport solver. The results agree well with the Monte Carlo solutions and demonstrate good acceleration performance. The code is then applied to the “X-2” VVER-1000 benchmark, and the results show that the bias of eigenvalue is less than 30pcm, and the maximum bias of assembly power is 0.61 %.

Jacobian-free Newton Krylov coarse mesh finite difference algorithm based on high-order nodal expansion method for three-dimensional nuclear reactor pin-by-pin multiphysics coupled models

The nuclear reactor core is a large-scale, complicated and tight-coupled nonlinear multiphysics coupled system including neutron transport, heat conduction of nuclear fuel pins, coolant flow, heat transfer and so on. To meet the demands of the advanced nuclear reactor analysis and design, there is an urgent need to reduce the computational costs and improve the convergence rate for the three-dimensional (3D) nuclear reactor core multiphysics coupled models, especially for a high-fidelity pin-by-pin (at the level of fuel pin by fuel pin) full core simulation. In this work, we develop CNJ, an efficient Jacobian free Newton Krylov (JFNK) coarse mesh finite difference (CMFD) algorithm based on the nodal expansion method (NEM) to simultaneously solve the complicated pin-by-pin neutronics-thermal hydraulic (N-TH) coupled models. The high-order NEM can improve the numerical accuracy of the CMFD method on a coarse mesh and the JFNK method can ensure the rapid convergence and high efficiency on large-scale and nonlinear coupled discrete systems. By combining the CMFD, NEM and JFNK methods and making full use of their respective advantages, the CNJ algorithm can obtain high accuracy and efficiency even on a coarse mesh to greatly reduce the number of solution variables and the computational cost of the 3D coupled models. Finally, two representative and complicated 3D pin-by-pin N-TH coupled models are analyzed to evaluate the numerical accuracy and efficiency advantage of CNJ in detailed.

The Legendre Polynomial Axial Expansion Method

This work presents a new formulation of the axial expansion transport method explicitly using Legendre polynomials for arbitrarily high-order expansions. This new formulation also features an alternative method of axial leakage calculation to allow for nonextruded flat source region meshes. This alternative axial leakage is introduced alongside a balance equation requirement to ensure that neutron balance is preserved in the coarse mesh for a given axial leakage formulation, which allows for effective coarse mesh finite difference acceleration. A matrix exponential table method is derived to allow for fast computations of arbitrarily high-order matrix exponentials for this work and precludes the need for further research into matrix exponential calculations for this method. Numerical results are presented that demonstrate the stability of the axial expansion method in systems with voidlike regions, showcase the speedup from matrix exponential tables, and investigate the axial convergence of the method in terms of both expansion order and mesh size.

ARCHER - a new Three-Dimensional method of characteristics neutron transport code for Pebble-bed HTR with coarse mesh finite difference acceleration

A new deterministic code ARCHER using Three-Dimensional (3-D) Method of Characteristics (MOC) has been developed by Institute of Nuclear and new Energy Technology (INET), Tsinghua University, aiming to obtain high-fidelity transport solution of pebble-bed High Temperature gas-cooled Reactor (HTR) with explicit pebble-bed geometry. However, 3-D MOC usually suffers from the slow convergence rate without efficient acceleration methods. In this work, the Coarse Mesh Finite Difference (CMFD) based on regular mesh is implemented to accelerate 3-D MOC calculation, which is available for the problems with cylinder and cuboid geometry. The complex geometry of reactor core with randomly stacked spherical fuel elements and irregular helium flow area, as well as the topological relationships between fine Flat Source Region (FSR) mesh for MOC and coarse CMFD mesh, are the key but challenging tasks for high-fidelity transport solution. Constructive Solid Geometry (CSG) and Boolean operation are constructed to describe problem geometry, which has the ability to deal with complex topological relationships between meshes. In addition, tree grid structure is utilized to reduce the time complexity of ray tracing and obtain the grid mapping relationship between coarse mesh and fine mesh. Spatial domain decomposition parallel based on MPI is utilized to alleviate the memory pressure and pursue high computational performance. Ray parallel is further achieved through OpenMP with shared memory. Several problems have been used to verify the geometric modeling capability and computational performance of this new 3-D MOC code for pebble-bed HTR.

Frequency Transform Method for Transient Analysis of Nuclear Reactors in Monte Carlo

The Frequency Transform method is used for the first time to efficiently model a multiple-second transient problem with Monte Carlo (MC). This is achieved by coupling MC with a time-dependent coarse mesh finite difference (TD-CMFD) diffusion solver. TD-CMFD presents a large advantage over commonly used point kinetics equations since it preserves spatial resolution during the transient and provides equivalence with the high-order method through nonlinear diffusion coefficients. As TD-CMFD computes time-dependent and spatially dependent neutronics information, it also computes frequencies that describe the rate of change of neutron and delayed precursor concentrations. These frequencies are used in MC shape function calculations as an approximation for the time derivatives. As the simulation proceeds, MC calculations update the multigroup cross sections, currents, and diffusion coefficients that are needed in TD-CMFD, and in turn, TD-CMFD updates the frequencies. Our results show the success of the Frequency Transform method in prescribed transient problems on the C5G7 geometry and on a fuel pin geometry. The Frequency Transform method showed significant improvement compared to the Adiabatic approximation, which does not use any frequency information in the MC calculation. The improvements in spatial resolution are shown to be a direct result of frequencies.

A Stable Condition and Adaptive Diffusion Coefficients for the Coarse-Mesh Finite Difference Method

Coarse-mesh finite difference (CMFD) method is a widely used numerical acceleration method. However, the stability of CMFD method is not good for the problems with optically thick regions. In this paper, a stability rule named the “sign preservation rule” in the field of numerical heat transfer is extended to the scheme of CMFD. It is required that the disturbance of neutron current is positively correlated with that of the negative value of flux gradient. A necessary condition for stability of the CMFD method is derived, an adaptive diffusion coefficient equation is proposed to improve the stability of CMFD method, and the corresponding revised CMFD method is called the rCMFD method. With a few modifications of the code, the rCMFD method was implemented in the hexagonal-Z nodal SN (discrete-ordinates) solver in the NECP-SARAX code system. The rCMFD method and other similar acceleration methods were tested by three fast reactor problems which were obtained by modifying the hexagonal pitches of a benchmark problem. The numerical results indicated that the rCMFD method showed better stability than the traditional CMFD method and the artificially diffusive CMFD (adCMFD) method and a better convergence rate than the adCMFD method and the optimally diffusive CMFD (odCMFD) method for these fast reactor problems.

Development and verification of a GPU-based MOC transient solver

The rapidly developed computational capability of the Graphic Processing Unit (GPU) provides an alternative platform for transport transient calculation. In this work, a GPU-based 2D method of characteristic (MOC) transient fixed source problem (TFSP) solver is implemented on lattice physics code ALPHA (Advanced Lattice Physics code based on Heterogeneous Architecture), which employs the full implicit method (FIM) and pin-based coarse mesh finite difference (CMFD) acceleration. C5G7-TD benchmark is calculated to verify the transient solver built-in ALPHA. To exploit the potential capability of GPU, sensitivity studies are performed on convergence criteria of MOC, which demonstrate single float-point precision (FP32) can be applied in the transient simulation. The numerical results of the ALPHA code show a good agreement with the results of MPACT and the ALPHA code with FP32 can provide enough accuracy and good performance on executing time.

Development of a GPU-based three-dimensional neutron transport code

A GPU-based three-dimensional (3D) neutron transport code, named ThorMOC, is developed in this work. A GPU-parallelized transport sweeping algorithm based on a planar method of characteristics (MOC) is implemented in ThorMOC and a two-level coarse mesh finite difference (CMFD) method is used to accelerate the convergence of the MOC transport solution. To reduce the memory burden of 3D calculation, the transport sweeping algorithm is improved by azimuthal angle decomposition. This developed ThorMOC code is verified by a Takeda benchmark and a C5G7 3D extension benchmark, and the performance of the azimuthal angle decomposition algorithm is evaluated by the whole core calculation of a C5G7 hexagonal benchmark. The numerical results obtained from ThorMOC show a good agreement with those from the Monte Carlo code OpenMC. Besides, the azimuthal angle decomposition algorithm can reduce the GPU storage requirement by more than 40% for the C5G7 hexagonal benchmark.

Investigation of Nodal Numerical Adjoints From CMFD-Based Acceleration Methods

In this work, the attribute of mathematical adjoints acquired from various CMFD (coarse-mesh finite difference) accelerated nodal methods based on the nodal expansion method (NEM) is presented. The direct transposition is implemented to the NEM-CMFD system matrix that includes correction factors to acquire the adjoint flux. Three different acceleration schemes are considered in this paper, which are, namely, CMFD, pCMFD, and one-node CMFD methods, and the self-adjoint attribute of migration operator for each acceleration scheme is studied. With regard to the one-node CMFD acceleration, a mathematical description for encountering negative adjoint flux values is given alongside an adjusted one-node CMFD scheme that stifles such an occurrence. The overall features of the aforementioned acceleration methods are recognized through analyses of 2D reactor problems including the KWU PWR 2D benchmark problem. Further systematic assessment is conducted based on the first-order perturbation theory, where the obtained adjoint fluxes are applied as weighting functions. It is clearly shown that the adjusted one-node CMFD scheme results in an improved reactivity estimation by excluding the presence of negative adjoint flux values.

Triangular Polynomial Expansion Nodal Method for VVER Core Analysis

Due to the low computational cost, nodal diffusion methods are still commonly used to simulate full-core reactor problems. This work represents the developmental effort to build an accurate nodal kernel to treat hexagonal geometry in the core simulator code PARCS. An innovative method called TriPEN-9 has been developed by splitting a hexagonal assembly into six triangular nodes and solved using cubic polynomial expansion for the scalar flux with nine-term expansion coefficients. The nodal diffusion calculation is further accelerated with the multilevel coarse-mesh finite difference method. The verification of the TriPEN-9 method on the VVER full-core problem is provided with the model based on the NURESIM (Nuclear Reactor Simulator)-SP1 V1000-2D-C1-tr benchmark problem. The Serpent Monte Carlo code is used as a reference solution for verification and to generate homogenized group-constants data for PARCS. Exact discontinuity factors were generated in GenPMAXS, a cross-section processing code, using a similar expansion method as the TriPEN-9 core solver method with the utilization of heterogeneous solutions from Serpent. Implementing the TriPEN-9 method in PARCS, this approach can exactly reproduce the solutions from the high-fidelity Serpent calculations.

A hybrid finite-difference/nodal method for benchmarking the HEXNEM3 coarse-mesh nodal method for time-dependent neutron diffusion applications

The solution accuracy of the recently developed modal analytical coarse-mesh finite-difference (ACMFD) formulation of the HEXNEM3 nodal flux expansion method for solving the time-dependent two-group neutron diffusion equation in hexagonal-z geometry is verified through a series of benchmark problems in both WWER-440 and WWER-1000 geometries. Reference solutions to these mathematical test problems are obtained by a hybrid finite-difference/nodal method developed specifically for the purpose of this study. The accuracy of the hybrid scheme for solving the neutron diffusion equation is also benchmarked against a set of test problems with published reference solutions. The coded implementations of the newly formulated HEXNEM3 method and the simpler but computationally expensive hybrid method demonstrated a very good mutual agreement in all test cases presented in this article.

Jacobian-free Newton Krylov two-node coarse mesh finite difference based on nodal expansion method

A Jacobian-Free Newton Krylov Two-Nodal Coarse Mesh Finite Difference algorithm based on Nodal Expansion Method (NEM_TNCMFD_JFNK) is successfully developed and proposed to solve the three-dimensional (3D) and multi-group reactor physics models. In the NEM_TNCMFD_JFNK method, the efficient JFNK method with the Modified Incomplete LU (MILU) preconditioner is integrated and applied into the discrete systems of the NEM-based two-node CMFD method by constructing the residual functions of only the nodal average fluxes and the eigenvalue. All the nonlinear corrective nodal coupling coefficients are updated on the basis of two-nodal NEM formulation including the discontinuity factor in every few newton steps. All the expansion coefficients and interface currents of the two-node NEM need not be chosen as the solution variables to evaluate the residual functions of the NEM_TNCMFD_JFNK method, therefore, the NEM_TNCMFD_JFNK method can greatly reduce the number of solution variables and the computational cost compared with the JFNK based on the conventional NEM. Finally the NEM_TNCMFD_JFNK code is developed and then analyzed by simulating the representative PWR MOX/UO2 core benchmark, the popular NEACRP 3D core benchmark and the complicated full-core pin-by-pin homogenous core model. Numerical solutions show that the proposed NEM_TNCMFD_JFNK method with the MILU preconditioner has the good numerical accuracy and can obtain higher computational efficiency than the NEM-based two-node CMFD algorithm with the power method in the outer iteration and the Krylov method using the MILU preconditioner in the inner iteration, which indicates the NEM_TNCMFD_JFNK method can serve as a potential and efficient numerical tool for reactor neutron diffusion analysis module in the JFNK-based multiphysics coupling application.

Convergence analysis for the CMFD accelerated 2D/1D neutron transport calculation method based on Fourier analysis

The two-dimension/one-dimension (2D/1D) transport method and the coarse mesh finite difference (CMFD) method are widely used for high-fidelity transport calculation. The CMFD method could dramatically reduce the calculating time while preserving accuracy. But the CMFD method accelerating the 2D/1D transport calculations may lead to divergence problems in some high scattering rate problems. The results of Fourier analysis indicate that the convergence is improved by increasing transport inner iterations. But more inner iterations lead to time burden in the method of characteristic (MOC) calculation. Fourier analysis is performed on the MOC based 2D/1D method to analyze the convergence feature, including increasing 1D transport inner iterations, iteration sequence of the 2D and 1D calculations, and so on. The results show that the effect of the 1D transport inner iteration on the convergence is related to the coarse-mesh height, and more 1D iteration could improve the convergence rate without introducing much more computation burden. Compared to the 1D2D computation sequence, the 2D1D sequence calculation has better convergence performance.

Two different methodologies for neutron flux calculation for the isotopic depletion problem in subcritical systems driven by an external neutron source in a coarse-mesh framework

The isotopic depletion calculation is typically solved by using the quasi-static approximation. In this approximation, the neutron flux is assumed as constant along the time steps in which the burn-up cycle is divided. This work proposes two different calculation methodologies to obtain the neutron flux, used in the isotopic depletion calculation in subcritical systems driven by an external neutron source, in a coarse mesh framework. The so-called depletion calculation methodology with Coarse Mesh Finite Difference (CMFD), and depletion calculation methodology with Pseudo-Harmonic Expansion (PHE). In these methodologies, the neutron flux at the beginning of each burn-up step is calculated with CMFD and the PHE methods, respectively. In order to verify the applicability of both methodologies to an ADS-type reactor, it is simulated a 45 MW small-reactor with burn-up cycle of 90 days, and an 800 MW large-reactor with burn-up cycle of 600 days. The subcriticality adjustment mechanism is based on the effective multiplication factor, solving the neutron diffusion equation with the Nodal Expansion Method (NEM). In addition, the averaged power density and the intensity of the neutron source are evaluated throughout the burn-up cycle. The results obtained by both methodologies are in good agreement with each other.

Neutron transport analysis of C5G7-TD benchmark with PANDAS-MOC

Here we present a newly developed neutronics code, PANDAS-MOC (Purdue Advanced Neutronics Design and Analysis System with Methods of Characteristics), to find the high fidelity solutions for the steady state and transient neutron transport analysis. Its essential method is the 2D/1D coupling method, in which the 2D planar solution is determined by the method of characteristics (MOC) along groups of characteristic tracks, the 1D axial solution is solved by the nodal expansion method (NEM), and they are coupled by the transverse leakage. The calculation is further accelerated by multi-level coarse mesh finite difference (CMFD) method. Moreover, its transport capability is verified by the 2D and 3D exercises of the C5G7-TD benchmark. The estimated eigenvalue of the steady-state, and the fractional core fission rate and reactivity over the transient activities are well-agreed with the published results of other codes, such as MPACT and nTRACER.