Diamond Difference Research Articles

Efficient and Portable Vectorized Sweep Kernels for the Transport Equation on 3D Cartesian Grids Using the Kokkos Framework

This paper describes the implementation of efficient and portable vectorized sweep kernels as part of the resolution of the neutron transport equation on three-dimensional Cartesian grids using the discrete ordinates ( S n ) method for the angular variable and the diamond differencing (DD) scheme for the spatial discretization. Vectorization is set up along the directions within the same octant and is independent of the spatial discretization order; therefore, the extension of this technique to high-order DD or discontinuous Galerkin schemes is immediate. Our implementation is written in C++17 and relies on the Kokkos performance portability framework. This library allows one to express shared-memory parallelism (including vectorization) in a machine-independent way and supports many backends including CUDA and OpenMP. Our vectorization procedure relies on the portable single instruction multiple data types provided by Kokkos. The method has been implemented for DD schemes up to order 2 and yields promising results on CPUs supporting standard vector instructions.

Tensor networks for solving the time-independent Boltzmann neutron transport equation

Tensor network techniques, known for their low-rank approximation ability that breaks the curse of dimensionality, are emerging as a foundation of new mathematical methods for ultra-fast numerical solutions of high-dimensional Partial Differential Equations (PDEs). Here, we present a mixed Tensor Train (TT)/Quantized Tensor Train (QTT) approach for the numerical solution of time-independent Boltzmann Neutron Transport equations (BNTEs) in Cartesian geometry. Discretizing a realistic three-dimensional (3D) BNTE by (i) diamond differencing in space, (ii) multigroup-in-energy, and (iii) discrete ordinates collocation in angle leads to large generalized eigenvalue problems that generally require a matrix-free approach and large computer clusters. Starting from this discretization, we construct a TT representation of the PDE fields and discrete operators, followed by a QTT representation of the TT cores. We then solve the tensorized generalized eigenvalue problem using a fixed-point scheme with tensor network optimization techniques. We validate our approach by applying the method to two examples of 3D neutron transport problems, currently solved by the Los Alamos National Laboratory PARallel TIme-dependent SN (PARTISN) solver.1 We demonstrate that our TT/QTT method, executed on a standard desktop computer, leads to large compression. This allows for the storage of terrabyte-sized neutron angular flux eigenvectors in megabytes. Additionally, we create megabyte-sized full access TT representations of yottabyte-sized transport matrix operators. By leveraging the TT operators and solution methods, we obtain a 7500 times speedup when compared to the PARTISN solution time with an error of less than 10−5.

On Solution Regularity of the Two-Dimensional Radiation Transfer Equation and Its Implication on Numerical Convergence

In this paper, we deal with the differential properties of the scalar flux ϕ ( x ) defined over a two-dimensional bounded convex domain, as a solution to the integral radiation transfer equation. Estimates for the derivatives of ϕ ( x ) near the boundary of the domain are given based on Vainikko’s regularity theorem. The optimal pointwise error estimates in terms of the scalar flux are presented for the two classic finite difference methods: diamond difference (DD) and step difference (SD). Numerical results indicate the implication of the solution smoothness on the numerical convergence behavior.

Development of angle-dependent linear source approximation for three-dimensional method of characteristics transport analysis method in STREAM

This paper features the development of angle-dependent linear source (ADLS) approximation for the 2D/3D Method of Characteristics (MOC)/Diamond Differencing (DD) method and its performance analysis simulation. The ADLS MOC/DD applies the linear source approximation to two different types of sources, the sum of fission and scattering source and the angle-wise surface source appearing in the MOC/DD equation. It is demonstrated memory usage and simulation time can be reduced by increasing the source approximation order from the zeroth to the first order.

Spectral Nodal Methodology for Multigroup Slab-Geometry Discrete Ordinates Neutron Transport Problems with Linearly Anisotropic Scattering

In this paper, we propose a numerical methodology for the development of a method of the spectral nodal class that generates numerical solutions free from spatial truncation errors. This method, denominated Spectral Deterministic Method (SDM), is tested as a study of the solutions (spectral analysis) of neutron transport equations in the discrete ordinates (SN) formulation, in slab geometry, multigroup approximation, with linearly anisotropic scattering, considering a heterogeneous domain with fixed-source. The unknowns in the methodology are the cell-edge, and cell average angular fluxes, the numerical values calculated for these quantities coincide with the analytic solution of the equations. These numerical results are shown and compared with the traditional fine-mesh Diamond Difference (DD) method and the coarse-mesh spectral Green's function (SGF) method to illustrate the method's accuracy and stability. The solution algorithms problem is implemented in a computer simulator made in C++ language, the same that was used to generate the results of the reference work.

A NOVEL COARSE-MESH METHOD APPLIED TO NEUTRON SHIELDING PROBLEMS USING THE MULTIGROUP TRANSPORT THEORY IN DISCRETE ORDINATES FORMULATIONS

In this paper, we propose a new deterministic numerical methodology to solve the one-dimensional linearized Boltzmann equation applied to neutron shielding problems (fixed-source), using the transport equation in the discrete ordinates formulation (SN) considering the multigroup theory. This is a hybrid methodology, entitled Modified Spectral Deterministic Method (SDM-M), that derives from the Spectral Deterministic Method (SDM) and Diamond Difference (DD) methods. This modification in the SDM method aims to calculate neutron scalar fluxes with lower computational cost. Two model-problems are solved using the SDM-M, and the results are compared to the coarse-mesh methods SDM, Spectral Green's Function (SGF) and Response Matrix (RM), and the fine-mesh method DD. The numerical results were obtained in the programming language JAVA version 1.8.0_91.

COMPARING THE DISCONTINUOUS GALERKIN AND HIGH-ORDER DIAMOND DIFFERENCING METHODS FOR THE TRANSPORT EQUATION ON A LOZENGE-BASED HEXAGONAL GEOMETRY

This paper presents an implementation and a comparison of two spatial discretisation schemes over a hexagonal geometry for the two-dimensional discrete ordinates transport equation. The methods are a high-order Discontinuous Galerkin (DG) finite element scheme and a high-order Diamond Differencing (DD) scheme. The DG method has been, and is being, studied on the hexagonal geometry, also called a honeycomb mesh – but not the DD method. In this research effort, it was chosen to divide the hexagons into (at least) three lozenges. An affine transformation is then applied onto said lozenges to cast them into the reference quadrilaterals usually studied in finite elements. In practice, this effectively means that the equations used in Cartesian geometry have their terms and operators altered using the Jacobian matrix of the transformation. This was implemented in the discrete ordinates solver of the code DRAGON5. Two 2D benchmark problems were then used for the verification and validation, including one based on the Monju 3D reactor benchmark. It was found that the diamond-differencing scheme seemed better. It converged much faster towards the solution at comparable mesh refinements for first-order expansion of the flux. Even if this difference was not present for second-order, DG was slower, about two to four times slower.

Three-dimensional method of characteristics/diamond-difference transport analysis method in STREAM for whole-core neutron transport calculation

A new three-dimensional (3D) transport analysis method is developed and implemented in the light water reactor (LWR) whole-core analysis code STREAM. The method is named as 3D method of characteristics/Diamond-difference (MOC/DD). In this method, 3D variables such as flux and source are expressed as a combination of the two-dimensional (2D) radial component and the one-dimensional (1D) axial component. The 2D method of characteristics (MOC), which is widely used in 2D neutron transport analysis in reactor physics, is used to solve complicated radial geometries composed of fuel pin, burnable absorber, structural material, etc. The axial neutron source, which is discretized by the diamond-difference (DD) method, is merged into the 3D MOC source term. This paper presents a detailed derivation of the 3D MOC/DD method and a parallelization technique to apply the 3D MOC/DD method to a whole-core problem. The accuracy and performance are examined by solving the C5G7 benchmark and the BEAVRS whole-core benchmark problems

Least squares linear source approach in method of characteristics

Method of characteristics (MOC) has been widely used in reactor physics calculation due to its high computational precision and great ability to deal with complex geometry. The linear source approximation is proposed to reduce the number of the meshes to improve the efficiency of the MOC method. In this research, a new least squares linear source MOC method is proposed and implemented. The numerical results show that the new least squares linear source MOC code is cable to deal with tiny mesh, such as gap in the PWR cell and the comparison among diamond difference (DD) method and several linear source MOC methods is performed. Numerical results indicate that the least squares linear source MOC method achieves desire precision when using much less mesh and the new method is stable with the number of the mesh increasing.

An Extended Linear Discontinuous Method for One-group Fixed Source Discrete Ordinates Problems with Isotropic Scattering in Slab Geometry

Nowadays, the obtainment of an accurate numerical solution of fixed source discrete ordinates problems is relevant in many areas of engineering and science. In this work, we extend the hybrid Finite Element Spectral Green's Function method (FEM-SGF), originally developed to solve eigenvalue diffusion problems, for fixed source problems using as a mathematical model, the discrete ordinates formulation in one energy group with isotropic scattering in slab geometry. This new method, Extended Linear Discontinuous Discrete Ordinates (ELD-SN), is based on the use of neutron balance equations and the construction of a hybrid auxiliary equation. This auxiliary equation combines a linear discontinuous approximation and spectral parameters to approximate the neutron angular flux inside the cell. Numerical results for benchmark problems are presented to illustrate the accuracy and computational performance of our methodology. ELD-SN method is free from spatial truncation errors in S2 quadrature, and generate good results in the other quadrature sets. This method is more accurate than the conventional Diamond Difference (DD) and Linear Discontinuous (LD) methods, but surpassed by the Spectral Green's Function (SGF) method, for quadrature order greater than two.

Implementation and performance analysis of the massively parallel method of characteristics based on GPU

The method of characteristics (MOC) is one of the most common methods for solving the neutron transport equation in practical application. Researches have been focused on the acceleration techniques and the parallel algorithm for improving the efficiency of MOC. The Graphics Processing Unit (GPU) provides an alternative method of parallelizing the MOC neutron transport sweep. In this work, a GPU-paralleled 2D MOC code is implemented, which employs the diamond difference (DD) scheme and the step characteristics (SC) scheme. Different parallel schemes which are ray-level, energy-group-level, and polar-angle-level, are analyzed to choose the proper parallel scheme. The C5G7 2D benchmark is calculated to verify the accuracy and efficiency of the code in different schemes with single precision and double precision. The bottlenecks of the GPU code are identified and the code is classified into three categories, which are compute-bound, memory-bound, and latency-bound, according to the performance analysis model introduced in this paper. In addition, corresponding optimization strategies are applied to improve the performance according to the analysis. Moreover, the speed, power efficiency, and hardware cost are compared for CPU and GPU based on a fictitious quarter core PWR problem. Numerical results demonstrate that the energy group-level parallelization can obtain the optimal performance on GPU. Optimization strategies are effective to improve the efficiency of the calculation on GPU, which indicates that the performance analysis model is useful and effective to locate the limitation of the code. Moreover, the GPU-version code is about 30 times faster than the CPU-version code with double precision and about 100 times faster with single precision, while the desired accuracy is kept. And the GPU delivers superior performance in both speed, energy efficiency, and hardware cost.

A multi-group extended linear discontinuous method for fixed-source discrete ordinates problems in slab geometry

At present, neutron density calculation in non-multiplying media is relevant in many areas of engineering and science. In this paper, we propose the Extended Linear Discontinuous (ELD) method in multi-group discrete ordinates formulation, originally formulated for one-energy group fixed-source problems with isotropic scattering source in slab geometry. The proposed auxiliary equations are uncoupled on angular directions and combine the linear discontinuous approximation of the finite element method and the quasi-analytical general solution of the spectral nodal method. Thus, we can implement an efficient and simple algorithm using the conventional source iteration scheme for the sweeping equations. Numerical results for benchmark problems are presented to illustrate the accuracy and computational performance of the ELD method. The work shows that the main advantages of the proposed method are that the numerical scheme is stable for coarse-meshes, and its numerical results are more accurate than those generated by the Diamond Difference (DD) and Linear Discontinuous (LD) methods.

A multi-group extended linear discontinuous method for fixed-source discrete ordinates problems in slab geometry

At present, neutron density calculation in non-multiplying media is relevant in many areas of engineering and science. In this paper, we propose the Extended Linear Discontinuous (ELD) method in multi-group discrete ordinates formulation, originally formulated for one-energy group fixed-source problems with isotropic scattering source in slab geometry. The proposed auxiliary equations are uncoupled on angular directions and combine the linear discontinuous approximation of the finite element method and the quasi-analytical general solution of the spectral nodal method. Thus, we can implement an efficient and simple algorithm using the conventional source iteration scheme for the sweeping equations. Numerical results for benchmark problems are presented to illustrate the accuracy and computational performance of the ELD method. The work shows that the main advantages of the proposed method are that the numerical scheme is stable for coarse-meshes, and its numerical results are more accurate than those generated by the Diamond Difference (DD) and Linear Discontinuous (LD) methods.

Numerical analysis of spatial discretization schemes in ARES transport code for shielding calculation

The discrete ordinates method (SN) is one of the mainstream methods for radiation shielding calculation of nuclear systems. The accuracy of spatial discretization in the SN equations has a significant effect on the accuracy of the shielding calculations. This work numerically examines the error behavior of various spatial discretization schemes implemented in the ARES transport code for shielding calculations. The numerical diffusion effect, encountered when the uncollided flux component is significant, such as in fast neutron groups or multi-dimensional void regions, is studied by investigating the ray propagation problem and glancing void problem. Positivity is investigated in optically thick cells and multi-dimensional void regions using spatial discretization schemes that cannot preserve non-negativity. In addition, the local accuracy is examined through a series of one-cell problems with various exact flux smoothness conditions in the cell. The results indicate that local accuracy depends not only on the spatial discretization but also the smoothness of the exact flux. Globally, severe unphysical flux representations arise in the diamond difference (DD) scheme. Other zero-order schemes inhibit oscillation but present significant numerical diffusion, whereas finite element methods damp the oscillation quickly near the flux discontinuity and present better resilience against negativity than DD. All these results can provide guidelines on selection of spatial discretization in realistic applications.

Goal-Based Error Estimation, Functional Correction, h, p and hp Adaptivity of the 1-D Diamond Difference Discrete Ordinate Method

ABSTRACTThis paper uses local dual weighted residual (DWR) error indicators to flag cells for goal-based refinement in a 1-D diamond difference (DD) discretisation of the discrete ordinate (SN) neutron transport equations. Goal-orientated adaptive mesh refinement (GO-AMR) aims to produce a mesh that is optimal for a given goal or QoI (Quantity of Interest). h, p and hp refinement is implemented and applied to various test cases. A merit function is derived for the combined hp algorithm and is calculated for each refinement option within each flagged cell. The refinement option with the highest merit function is chosen for a given flagged cell. This paper also investigates the use of the DWR error estimation as a correction term for the originally calculated QoI. If error correction and GO-AMR are combined the DWR error indicators do not always give an optimal mesh for the corrected value of the QoI. Therefore, refinement indicators based on the error in the error correction term are used and tested in this work.

Goal-based h-adaptivity of the 1-D diamond difference discrete ordinate method

The quantity of interest (QoI) associated with a solution of a partial differential equation (PDE) is not, in general, the solution itself, but a functional of the solution. Dual weighted residual (DWR) error estimators are one way of providing an estimate of the error in the QoI resulting from the discretisation of the PDE.This paper aims to provide an estimate of the error in the QoI due to the spatial discretisation, where the discretisation scheme being used is the diamond difference (DD) method in space and discrete ordinate (SN) method in angle. The QoI are reaction rates in detectors and the value of the eigenvalue (Keff) for 1-D fixed source and eigenvalue (Keff criticality) neutron transport problems respectively. Local values of the DWR over individual cells are used as error indicators for goal-based mesh refinement, which aims to give an optimal mesh for a given QoI.

High-Order Linear Discontinuous and Diamond Differencing Schemes Along Cyclic Characteristics

We are investigating a new class of linear characteristics schemes along cyclic tracks for solving the transport equation for neutral particles with scattering anisotropy. These algorithms rely on linear discontinuous exact integration and diamond differencing, as implemented with the method of discrete ordinates. These schemes are based on linear discontinuous coefficients that are derived through the application of approximations describing the mesh-averaged spatial flux moments in terms of spatial source moments and of the beginning-of-segment and end-of-segment flux values. The linear discontinuous characteristics (LDC) and quadratic-order diamond differencing (DD1) schemes are inherently conservative. In this technical note, we intend to continue the development of the LDC and DD1 schemes by extending their application to cyclic trackings. This extension will make possible the representation of reflective or general albedo boundary conditions. We will present an improved and much shorter derivation of the LDC and DD1 schemes, compared to a previous presentation. Finally, we will implement the new schemes as Matlab scripts for solving a one-dimensional slab benchmark and in the DRAGON5 lattice code for solving a more representative two-dimensional eight-symmetry pressurized water reactor assembly mock-up.

Corner neutronic code

The neutronic code CORNER is based on the SN discrete ordinates method [1] and the PM scattering cross-section approximation. It is intended for high-precision deterministic neutronic calculations of fast-neutron reactors and can be used to solve two types of steady-state neutron transport and gamma quanta problems in a 3D hexagonal geometry: Keff problems (homogeneous) and source problems (inhomogeneous).The code is developed in the Fortran language and has a modular structure. Its key modules are a module for the preparation of neutron constants in the ANISN format; a geometric module containing a description of the core's loading map and fuel assembly types, including their axial meshing and material composition; a module for preparing angular quadrature sets; an input data module containing the parameters of the approximation used and the control parameters; a neutronic calculation module and a calculation data processing module.The Directional Theta-Weighted (DTW) difference scheme [2] has been built to approximate the spatial dependence. It has advantages over the DD (Diamond Difference) scheme broadly used in coarse-mesh problems.The energy dependence is represented by multigroup approximation. The angular variable is discretized by introduction of the angular quadrature set. Quadrature sets can be also defined by the user.An iterative solution process is used, including external iterations for the fission source and internal iterations for the scattering source. The paper presents the results of a cross-verification against the Monte Carlo MMK code [3] and on a model of the BN-800 reactor core.

High-Order Diamond Differencing along Finite Characteristics

We are investigating a new class of linear characteristics schemes along finite-length tracks for solving the transport equation for neutral particles with scattering anisotropy. These algorithms are based on diamond differencing, as implemented with the method of discrete ordinates. The quadratic-order diamond-differencing (DD1) scheme is based on linear discontinuous coefficients that are derived through the application of approximations describing the mesh-averaged spatial flux moments in terms of spatial source moments and of the beginning- and end-of-segment flux values. This DD1 linear characteristics scheme is inherently conservative. This approach is an improvement relative to other linear characteristics schemes because no information needs to be collected on internal surfaces. Consequently, the DD1 scheme is compatible with existing tracking files for the collision-probability method. The proposed scheme is verified in one-dimensional slab geometry where it is found to be equivalent to a discrete ordinates solution and on simple two-dimensional benchmarks made of regular squares or hexagons.

A three-dimensional high-order diamond differencing discretization with a consistent acceleration scheme

Abstract The extension of the diamond differencing scheme to high-order spatial approximation is described in the context of a regular, three-dimensional (3D) discrete ordinates method. This spatial discretization keeps the advantages of the well-known linear diamond differencing (DD) scheme in terms of rigorousness and simplicity of derivation, while extending it to high-order solutions. Along with the implementation of high-order diamond differencing comes the search for a consistent acceleration method for the S N source iteration. The proposed method relies on a Diffusion Synthetic Acceleration (DSA) scheme, combined with a Krylov Subspace Algorithm, GMRES . This two-level acceleration scheme has been proven efficient in case of realistic problems at any order. Numerical results are provided on 2D/3D legacy benchmarks and establish good properties of our solver.