Scalable angular adaptivity for Boltzmann transport

This paper describes an angular adaptivity algorithm for Boltzmann transport applications which for the first time shows evidence of $\mathcal{O}(n)$ scaling in both runtime and memory usage, where $n$ is the number of adapted angles. This adaptivity uses Haar wavelets, which perform structured $h$-adaptivity built on top of a hierarchical P$_0$ FEM discretisation of a 2D angular domain, allowing different anisotropic angular resolution to be applied across space/energy. Fixed angular refinement, along with regular and goal-based error metrics are shown in three example problems taken from neutronics/radiative transfer applications. We use a spatial discretisation designed to use less memory than competing alternatives in general applications and gives us the flexibility to use a matrix-free multgrid method as our iterative method. This relies on scalable matrix-vector products using Fast Wavelet Transforms and allows the use of traditional sweep algorithms if desired

Adaptive discontinuous Galerkin methods for the neutron transport equation

In this thesis we study the neutron transport (Boltzmann transport equation) which is used to model the movement of neutrons inside a nuclear reactor. More specifically we consider the mono-energetic, time independent neutron transport equation. The neutron transport equation has predominantly been solved numerically by employing low order discretisation methods, particularly in the case of the angular domain. We proceed by surveying the advantages and disadvantages of common numerical methods developed for the numerical solution of the neutron transport equation before explaining our choice of using a discontinuous Galerkin (DG) discretisation for both the spatial and angular domain. The bulk of the thesis describes an arbitrary order in both angle and space solver for the neutron transport equation. We discuss some implementation issues, including the use of an ordered solver to facilitate the solution of the linear systems resulting from the discretisation. The resulting solver is benchmarked using both source and critical eigenvalue computations. In the pseudo three--dimensional case we employ our solver for the computation of the critical eigenvalue for three industrial benchmark problems. We then employ the Dual Weighted Residual (DWR) approach to adaptivity to derive and implement error indicators for both two--dimensional and pseudo three--dimensional neutron transport source problems. Finally, we present some preliminary results on the use of a DWR indicator for the eigenvalue problem

Adaptive discontinuous Galerkin methods for the neutron transport equation

In this thesis we study the neutron transport (Boltzmann transport equation) which is used to model the movement of neutrons inside a nuclear reactor. More specifically we consider the mono-energetic, time independent neutron transport equation. The neutron transport equation has predominantly been solved numerically by employing low order discretisation methods, particularly in the case of the angular domain. We proceed by surveying the advantages and disadvantages of common numerical methods developed for the numerical solution of the neutron transport equation before explaining our choice of using a discontinuous Galerkin (DG) discretisation for both the spatial and angular domain. The bulk of the thesis describes an arbitrary order in both angle and space solver for the neutron transport equation. We discuss some implementation issues, including the use of an ordered solver to facilitate the solution of the linear systems resulting from the discretisation. The resulting solver is benchmarked using both source and critical eigenvalue computations. In the pseudo three--dimensional case we employ our solver for the computation of the critical eigenvalue for three industrial benchmark problems. We then employ the Dual Weighted Residual (DWR) approach to adaptivity to derive and implement error indicators for both two--dimensional and pseudo three--dimensional neutron transport source problems. Finally, we present some preliminary results on the use of a DWR indicator for the eigenvalue problem

Stabilized finite element methods for natural and forced convection-radiation heat transfer

Thermal radiation in forced and natural convection can be an important mode of heat transfer in high temperature chambers, such as industrial furnaces and boilers, even under non-soot conditions. Growing concern with high temperature processes has emphasized the need for an evaluation of the eect of radiative heat transfer. Nevertheless, the modelling of radiation is often neglected in combustion analysis, mainly because it involves tedious mathematics, which increase the computation time, and also because of the lack of detailed information on the optical properties of the participating media and surfaces. Ignoring radiative transfer may introduce signicant errors in the overall predictions. The most accurate procedures available for computing radiation transfer in furnaces are the Zonal and Monte Carlo methods. However, these methods are not widely applied in comprehensive combustion calculations due to their large computational time and storage requirements. Also, the equations of the radiation transfer are in non-dierential form, a signicant inconvenience when solved in conjunction with the dierential equations of ow and combustion. For this reason, numerous investigations are currently being carried out worldwide to assess computationally ecient methods. In addition ecient modelling of forced and natural convection-radiation would help to simulate and understand heat transfer appearing in various engineering applications, especially in the case of the heat treatment of high-alloy steel or glass by a continuously heating process inside industrial furnaces, ovens or even smaller applications like microwaves. This thesis deals with the design of such methods and shows that a class of simplied approximations provides advantages that should be utilized in treating radiative transfer problems with or without ow convection. Much of the current work on modelling energy transport in high-temperature gas furnaces or chemically reacting ows, uses computational uid dynamics (CFD) codes. Therefore, the models for solving the radiative transfer equations must be compatible with the numerical methods employed to solve the transport equations. The Zonal and Monte Carlo methods for solving the radiative transfer problem are incompatible with the mathematical formulations used in CFD codes, and require prohibitive computational times for spatial resolution desired. The main objectives of this thesis is then to understand and better model the heat treatment at the same time in the furnace/oven chamber and within the workpieces under specied furnace geometry, thermal schedule, parts loading design, initial operation conditions, and performance requirements. Nowadays, there is a strong need either for appropriate fast and accurate algorithms for the mixed and natural convection-radiation or for reduced models which still incorporate its main radiative transfer physics. During the last decade, a lot of research was focused on the derivation of approximate models allowing for an accurate description of the important physical phenomena at reasonable numerical costs. Hence, a whole hierarchy of approximative equations is available, ranging from half-space moment approximations over full-space moment systems to the diusion-type simplied PN approximations. The latter were developed and extensively tested for various radiative transfer problems, where they proved to be suciently accurate. Although they were derived in the asymptotic regime for a large optical thickness of the material, these approximations yield encouraging even results in the optically thin regime. The main advantage of considering simplied PN approximations is the fact that the integro-dierential radiative transfer equation is transformed into a set of elliptic equations independent of the angular direction which are easy to solve. The simplied PN models are proposed in this thesis for modelling radiative heat transfer for both forced and natural convection-radiation applications. There exists a variety of computational methods available in the literature for solving coupled convection-radiation problems. For instance, applied to convection-dominated ows, Eulerian methods incorporate some upstream weighting in their formulations to stabilize the numerical procedure. The most popular Eulerian methods, in nite element framework, are the streamline upwind Petrov-Galerkin, Galerkin/least-squares and Taylor-Galerkin methods. All these Eulerian methods are easy to formulate and implement. However, time truncation errors dominate their solutions and are subjected to Courant-Friedrichs-Lewy (CFL) stability conditions, which put a restriction on the size of time steps taken in numerical simulations. Galerkin-characteristic methods (also known by semi-Lagrangian methods in meteorological community) on the other hand, make use of the transport nature of the governing equations. The idea in these methods is to rewrite the governing equations in term of Lagrangian co-ordinates as dened by the particle trajectories (or characteristics) associated with the problem. Then, the Lagrangian total derivative is approximated, thanks to a divided dierence operator. The Lagrangian treatment in these methods greatly reduces the time truncation errors in the Eulerian methods. In addition, these methods are known to be unconditionally stable, independent of the diusion coecient, and optimally accurate at least when the inner products in the Galerkin procedure are calculated exactly. In Galerkin-characteristic methods, the time derivative and the advection term are combined as a directional derivative along the characteristics, leading to a characteristic time-stepping procedure. Consequently, the Galerkin-characteristic methods symmetrize and stabilize the governing equations, allow for large time steps in a simulation without loss of accuracy, and eliminate the excessive numerical dispersion and grid orientation eects present in many upwind methods. This class of numerical methods have been implemented in this thesis to solve the developed models for mixed and natural convection-radiation applications. Extensive validations for the numerical simulations have been carried out and full comparisons with other published numerical results (obtained using commercial softwares) and experimental results are illustrated for natural and forced radiative heat transfer. The obtained convectionradiation results have been studied under the eect of dierent heat transfer characteristics to improve the existing applications and to help in the furnace designs

Anisotropic Adaptivity and Subgrid Scale Modelling for the Solution of the Neutron Transport Equation with an Emphasis on Shielding Applications

This thesis demonstrates advanced new discretisation and adaptive meshing technologies that improve the accuracy and stability of using finite element discretisations applied to the Boltzmann transport equation (BTE). This equation describes the advective transport of neutral particles such as neutrons and photons within a domain. The BTE is difficult to solve, due to its large phase space (three dimensions of space, two of angle and one each of energy and time) and the presence of non-physical oscillations in many situations. This work explores the use of a finite element method that combines the advantages of the two schemes: the discontinuous and continuous Galerkin methods. The new discretisation uses multiscale (subgrid) finite elements that work locally within each element in the finite element mesh in addition to a global, continuous, formulation. The use of higher order functions that describe the variation of the angular flux over each element is also explored using these subgrid finite element schemes. In addition to the spatial discretisation, methods have also been developed to optimise the finite element mesh in order to reduce resulting errors in the solution over the domain, or locally in situations where there is a goal of specific interest (such as a dose in a detector region). The chapters of this thesis have been structured to be submitted individually for journal publication, and are arranged as follows. Chapter 1 introduces the reader to motivation behind the research contained within this thesis. Chapter 2 introduces the forms of the BTE that are used within this thesis. Chapter 3 provides the methods that are used, together with examples, of the validation and verification of the software that was developed as a result of this work, the transport code RADIANT. Chapter 4 introduces the inner element subgrid scale finite element discretisation of the BTE that forms the basis of the discretisations within RADIANT and explores its convergence and computational times on a set of benchmark problems. Chapter 5 develops the error metrics that are used to optimise the mesh in order to reduce the discretisation error within a finite element mesh using anisotropic adaptivity that can use elongated elements that accurately resolves computational demanding regions, such as in the presence of shocks. The work of this chapter is then extended in Chapter 6 that forms error metrics for goal based adaptivity to minimise the error in a detector response. Finally, conclusions from this thesis and suggestions for future work that may be explored are discussed in Chapter 7.Open Acces

Methods for Solving Discontinuous-Galerkin Finite Element Equations with Application to Neutron Transport

Cette thèse traite des méthodes d’éléments finis Galerkin discontinus d’ordre élevé pour la résolution d’équations aux dérivées partielles, avec un intérêt particulier pour l’équation de transport des neutrons. Nous nous intéressons tout d’abord à une méthode de pré-traitement de matrices creuses par blocs, qu’on retrouve dans les méthodes Galerkin discontinues, avant factorisation par un solveur multifrontal. Des expériences numériques conduites sur de grandes matrices bi- et tri-dimensionnelles montrent que cette méthode de pré-traitement permet une réduction significative du ’fill-in’, par rapport aux méthodes n’exploitant pas la structure par blocs. Ensuite, nous proposons une méthode d’éléments finis Galerkin discontinus, employant des éléments d’ordre élevé en espace comme en angle, pour résoudre l’équation de transport des neutrons. Nous considérons des solveurs parallèles basés sur les sous-espaces de Krylov à la fois pour des problèmes ’source’ et des problèmes aux valeur propre multiplicatif. Dans cet algorithme, l’erreur est décomposée par projection(s) afin d’équilibrer les contraintes numériques entre les parties spatiales et angulaires du domaine de calcul. Enfin, un algorithme HP-adaptatif est présenté ; les résultats obtenus démontrent une nette supériorité par rapport aux algorithmes h-adaptatifs, à la fois en terme de réduction de coût de calcul et d’amélioration de la précision. Les valeurs propres et effectivités sont présentées pour un panel de cas test industriels. Une estimation précise de l’erreur (avec effectivité de 1) est atteinte pour un ensemble de problèmes aux domaines inhomogènes et de formes irrégulières ainsi que des groupes d’énergie multiples. Nous montrons numériquement que l’algorithme HP-adaptatif atteint une convergence exponentielle par rapport au nombre de degrés de liberté de l’espace éléments finis. ABSTRACT : We consider high order discontinuous-Galerkin finite element methods for partial differential equations, with a focus on the neutron transport equation. We begin by examining a method for preprocessing block-sparse matrices, of the type that arise from discontinuous-Galerkin methods, prior to factorisation by a multifrontal solver. Numerical experiments on large two and three dimensional matrices show that this pre-processing method achieves a significant reduction in fill-in, when compared to methods that fail to exploit block structures. A discontinuous-Galerkin finite element method for the neutron transport equation is derived that employs high order finite elements in both space and angle. Parallel Krylov subspace based solvers are considered for both source problems and $k_{eff}$-eigenvalue problems. An a-posteriori error estimator is derived and implemented as part of an h-adaptive mesh refinement algorithm for neutron transport $k_{eff}$-eigenvalue problems. This algorithm employs a projection-based error splitting in order to balance the computational requirements between the spatial and angular parts of the computational domain. An hp-adaptive algorithm is presented and results are collected that demonstrate greatly improved efficiency compared to the h-adaptive algorithm, both in terms of reduced computational expense and enhanced accuracy. Computed eigenvalues and effectivities are presented for a variety of challenging industrial benchmarks. Accurate error estimation (with effectivities of 1) is demonstrated for a collection of problems with inhomogeneous, irregularly shaped spatial domains as well as multiple energy groups. Numerical results are presented showing that the hp-refinement algorithm can achieve exponential convergence with respect to the number of degrees of freedom in the finite element spac

Thermal modeling of short pulse collimated radiation in a participating medium”

In most traditional engineering applications, such as in the thermal analysis of boilers,furnaces, internal combustion engines, etc., as temporal variations in thermal quantities of interest are much slower than the time scale associated with the propagation of radiation,the transient term from the radiative transfer equation is usually neglected.The rapid progress in the field of short pulse laser in a participating medium has lead to some interesting applications such as materials processing, optical tomography, remote sensing,laser micro-surgery etc.The temporal radiative signals from a medium irradiated by short pulse lasers offer more useful information which reflects the internal structure and properties of the medium than that by the continuous light sources.The time scales of such applications are usually in the order of 10-12 to 10-15 seconds. Therefore, the consideration of the transient term in the radiation transport equation is necessary.To understand the physical phenomena involved in this complicated transport process, two cases are considered such as 1) one dimensional model and 2) a two dimensional model.The models assume a participating medium bounded by diffusely emitting and reflecting boundaries, one of the boundaries is irradiated with a short pulse laser beam.The gray wall assumption resembles more to the practical application as compared to the simplified assumption of black wall. Due to reflective characteristics of the gray boundary surface, the temporal spread changes significantly by the multiple reflections and partial transmissions at the surfaces, this draws attention to the present problem considered.The finite volume method is applied to solve the transient radiative transfer equation governing the above said physical phenomena.The fully implicit scheme is used to iii discretize the transient term.In the proposed approach, one does not have to split the intensity into diffused and collimated part unlike the case with other existing methods like DTM, DOM, and REM etc.Intensity of radiation can directly be evaluated by solving the governing transient radiative transfer equation.The proposed methodology is compared with other existing methods.The effects of optical thickness, scattering albedo, emissivity, and linear anisotropic scattering on the transmitted and reflected signals are studied.The performance of two different spatial schemes: STEP and CLAM also have been tested. It is seen that the CLAM scheme yields results to a greater accuracy in case of two dimensional problems and, hence, correctly predicts the speed of photon whereas STEP scheme overpredicts the same

Subgrid scale modelling of transport processes.

Consideration of stabilisation techniques is essential in the development of physical models if they are to faithfully represent processes over a wide range of scales. Careful application of these techniques can significantly increase flexibility of models, allowing the computational meshes used to discretise the underlying partial differential equations to become highly nonuniform and anisotropic, for example. This exibility enables a model to capture a wider range of phenomena and thus reduce the number of parameterisations required, bringing a physically more realistic solution. The next generation of fluid flow and radiation transport models employ unstructured meshes and anisotropic adaptive methods to gain a greater degree of flexibility. However these can introduce erroneous artefacts into the solution when, for example, a process becomes unresolvable due to an adaptive mesh change or advection into a coarser region of mesh in the domain. The suppression of these effects, caused by spatial and temporal variations in mesh size, is one of the key roles stabilisation can play. This thesis introduces new explicit and implicit stabilisation methods that have been developed for application in fluid and radiation transport modelling. With a focus on a consistent residual-free approach, two new frameworks for the development of implicit methods are presented. The first generates a family of higher-order Petrov-Galerkin methods, and the example developed is compared to standard schemes such as streamline upwind Petrov-Galerkin and Galerkin least squares in accurate modelling of tracer transport. The dissipation generated by this method forms the basis for a new explicit fourth-order subfilter scale eddy viscosity model for large eddy simulation. Dissipation focused more sharply on unresolved scales is shown to give improved results over standard turbulence models. The second, the inner element method, is derived from subgrid scale modelling concepts and, like the variational multiscale method and bubble enrichment techniques, explicitly aims to capture the important under-resolved fine scale information. It brings key advantages to the solution of the Navier-Stokes equations including the use of usually unstable velocity-pressure element pairs, a fully consistent mass matrix without the increase in degrees of freedom associated with discontinuous Galerkin methods and also avoids pressure filtering. All of which act to increase the flexibility and accuracy of a model. Supporting results are presented from an application of the methods to a wide range of problems, from simple one-dimensional examples to tracer and momentum transport in simulations such as the idealised Stommel gyre, the lid-driven cavity, lock-exchange, gravity current and backward-facing step. Significant accuracy improvements are demonstrated in challenging radiation transport benchmarks, such as advection across void regions, the scattering Maynard problem and demanding source-absorption cases. Evolution of a free surface is also investigated in the sloshing tank, transport of an equatorial Rossby soliton, wave propagation on an aquaplanet and tidal simulation of the Mediterranean Sea and global ocean. In combination with adaptive methods, stabilising techniques are key to the development of next generation models. In particular these ideas are critical in achieving the aim of extending models, such as the Imperial College Ocean Model, to the global scale

Synergies between Numerical Methods for Kinetic Equations and Neural Networks

The overarching theme of this work is the efficient computation of large-scale systems. Here we deal with two types of mathematical challenges, which are quite different at first glance but offer similar opportunities and challenges upon closer examination. Physical descriptions of phenomena and their mathematical modeling are performed on diverse scales, ranging from nano-scale interactions of single atoms to the macroscopic dynamics of the earth\u27s atmosphere. We consider such systems of interacting particles and explore methods to simulate them efficiently and accurately, with a focus on the kinetic and macroscopic description of interacting particle systems. Macroscopic governing equations describe the time evolution of a system in time and space, whereas the more fine-grained kinetic description additionally takes the particle velocity into account. The study of discretizing kinetic equations that depend on space, time, and velocity variables is a challenge due to the need to preserve physical solution bounds, e.g. positivity, avoiding spurious artifacts and computational efficiency. In the pursuit of overcoming the challenge of computability in both kinetic and multi-scale modeling, a wide variety of approximative methods have been established in the realm of reduced order and surrogate modeling, and model compression. For kinetic models, this may manifest in hybrid numerical solvers, that switch between macroscopic and mesoscopic simulation, asymptotic preserving schemes, that bridge the gap between both physical resolution levels, or surrogate models that operate on a kinetic level but replace computationally heavy operations of the simulation by fast approximations. Thus, for the simulation of kinetic and multi-scale systems with a high spatial resolution and long temporal horizon, the quote by Paul Dirac is as relevant as it was almost a century ago. The first goal of the dissertation is therefore the development of acceleration strategies for kinetic discretization methods, that preserve the structure of their governing equations. Particularly, we investigate the use of convex neural networks, to accelerate the minimal entropy closure method. Further, we develop a neural network-based hybrid solver for multi-scale systems, where kinetic and macroscopic methods are chosen based on local flow conditions. Furthermore, we deal with the compression and efficient computation of neural networks. In the meantime, neural networks are successfully used in different forms in countless scientific works and technical systems, with well-known applications in image recognition, and computer-aided language translation, but also as surrogate models for numerical mathematics. Although the first neural networks were already presented in the 1950s, the scientific discipline has enjoyed increasing popularity mainly during the last 15 years, since only now sufficient computing capacity is available. Remarkably, the increasing availability of computing resources is accompanied by a hunger for larger models, fueled by the common conception of machine learning practitioners and researchers that more trainable parameters equal higher performance and better generalization capabilities. The increase in model size exceeds the growth of available computing resources by orders of magnitude. Since $2012$, the computational resources used in the largest neural network models doubled every $3.4$ months\footnote{\url{https://openai.com/blog/ai-and-compute/}}, opposed to Moore\u27s Law that proposes a $2$-year doubling period in available computing power. To some extent, Dirac\u27s statement also applies to the recent computational challenges in the machine-learning community. The desire to evaluate and train on resource-limited devices sparked interest in model compression, where neural networks are sparsified or factorized, typically after training. The second goal of this dissertation is thus a low-rank method, originating from numerical methods for kinetic equations, to compress neural networks already during training by low-rank factorization. This dissertation thus considers synergies between kinetic models, neural networks, and numerical methods in both disciplines to develop time-, memory- and energy-efficient computational methods for both research areas