Fully Implicit Time Integration: Fast Solvers and Implicit-Explicit Schemes [Slides]
Abstract not provided.
SEARCH · Search NASA
Search indexed NASA NTRS and DOE OSTI research on propulsion, heat transfer, battery materials and energy systems. Follow report and document links to the original sources.
Quote a phrase for an exact phrase match. Source license links do not imply unrestricted reuse.
Abstract not provided.
Splitting-based time integration approaches such as fractional step, alternating direction implicit, operator splitting, and locally one dimensional methods partition the system of interest into components, and solve individual components implicitly in a cost-effective way. Here this work proposes a unified formulation of splitting time integration schemes in the framework of general-structure additive Runge–Kutta (GARK) methods. Specifically, we develop implicit-implicit (IMIM) GARK schemes, provide the order conditions for this class, and explain their application to partitioned systems of ordinary differential equations. We show that classical splitting methods belong to the IMIM GARK family, and therefore can be studied in this unified framework. New IMIM-GARK splitting methods are developed and tested using parabolic systems.
Here, we consider the adaptive-rank integration of multi-dimensional time-dependent advection-diffusion partial differential equations (PDEs) with variable coefficients. We employ a standard finite-difference method for spatial discretization coupled with high-order diagonally implicit Runge-Kutta temporal schemes. The discrete equation is a generalized Sylvester equation (GSE), which we solve with a projection-based adaptive-rank algorithm structured around two key strategies: (i) constructing dimension-wise subspaces using a novel atypical extended Krylov strategy, and (ii) efficiently solving the basis coefficient matrix with a preconditioned GMRES solver. The low-rank decomposition is performed in 2D using SVD and with high-order SVD (HOSVD) in 3D to represent the tensor in a compressed Tucker format. For d-dimensional problems (here, d = 2 or 3), the computational complexity and memory storage of the approach are found numerically to scale as and $\mathscr{O}(Nr^2) + \mathscr{O} (r^{d+1})$ and $\mathscr{O}(Nr) + \mathscr{O} (r^{d})$, respectively, with the one-dimensional resolution and the maximal rank during the Krylov iteration (which we find to be largely independent of on our numerical examples). We present numerical examples that illustrate the advertised properties of the algorithm.
Not Available
Here, this work describes a crystal plasticity formulation combining several mathematical, numerical, and implementation choices to produce a highly efficient model. Specifically, the key choices in the implementation are (1) representing orientations with modified Rodrigues parameters, (2) implementing a fully coupled implicit time integration for the elastic stretch, the crystal orientations, and the model internal variables, (3) implementing the model in the NEML2 constitutive modeling framework, based on PyTorch, to vectorize the calculations and port the computation to GPUs and other hardware accelerators, and (4) an exact implementation of the consistent tangent matrix, even for arbitrary coupling to other field variables beyond the displacements, like temperature, neutron fluence, etc. The first two features of the model are, to our knowledge, novel. The paper considers each of these choices individually as well as the final model as a whole. This includes a full description of modified Rodrigues parameters, their advantages over other representations of orientations, the mathematical formulae and tools required to implement a model with modified Rodrigues parameters, and a detailed description of the geometry of the space of modified Rodrigues parameters (in an appendix). It also includes a description of a fully implicit time integration scheme for the orientations and the advantages in representing orientations with modified Rodrigues parameters in implementing such a model. The work then assess, via numerical examples, the advantages of fully coupled implicit time integration versus more common decoupled and explicit time integration schemes. These studies demonstrate the computational advantages of fully coupled integration versus other time integration algorithms, though the performance of the competing models depends on the complexity of the underlying single crystal model. The study concludes by demonstrating that the choice of time integration method affects the sharpness of the predicted texture, with explicit methods for integrating the orientations overestimating texture sharpness and implicit methods underestimating texture sharpness.
The standard particle-in-cell (PIC) method employs explicit finite-difference (FD) methods (e.g. the leap-frog scheme) for both spatial and temporal integrations. Here, we employ a pseudospectral method for solving the Poisson equation and a fully implicit time integration to achieve exact energy conservation. The advantage of a pseudospectral field solver is its spectral accuracy in solving field solutions. Earlier studies of implicit time integration of PIC FD equations can enforce exact energy exchange between field and particles, resulting in exact energy-conserving schemes. Here, we prove that the exact energy conservation property can be carried over to the pseudospectral scheme. Simultaneously, we provide a solution to ensure a pseudospectral charge continuity equation. We demonstrate the new scheme in a 2D electrostatic PIC code. In conclusion, theoretical results are confirmed via numerical examples.
Abstract not provided.
In this paper, we extend the operator-split asymptotic-preserving, semi-Lagrangian algorithm for time dependent anisotropic heat transport equation proposed in Chacón et al. (2014) [18] to use a fully implicit time integration with backward differentiation formulas. The proposed implicit method can deal with arbitrary heat-transport anisotropy ratios $\mathcal{X}$∥ /$ \mathcal{X}$⟂ $\ggg$ 1 (with $\mathcal{X}$∥, $ \mathcal{X}$⟂ the parallel and perpendicular heat diffusivities, respectively) in complicated magnetic field topologies in an accurate and efficient manner. Further, the implicit algorithm is second-order accurate temporally and demonstrates an accurate treatment at boundary layers (e.g., island separatrices), which was not ensured by the operator-split implementation. The condition number of the resulting algebraic system is independent of the anisotropy ratio, and is inverted with preconditioned GMRES. We propose a simple preconditioner that renders the finite-dimensional linear operator compact, resulting in mesh-independent convergence rates for topologically simple magnetic fields, and convergence rates scaling as ~ (NΔt) 1/4 (with N the total mesh size and Δt the timestep) in topologically complex magnetic-field configurations. We demonstrate the accuracy and performance of the approach with test problems of varying complexity, including an analytically tractable boundary-layer problem in a straight magnetic field, and a topologically complex magnetic field featuring magnetic islands with extreme anisotropy ratios $\mathcal{X}$∥ /$ \mathcal{X}$⟂ = 10 10 ) .
Multigrid (MG) is widely recognized as a highly effective solver for the model problem, the Laplacian, but textbook MG fails on most problems of interest. MG methods have been applied to complex, real-world applications with careful consideration of the physical model and discretization. In this work we develop the first step in applying MG methods to science and engineering relevant magnetohydrodynamics (MHD) tokamak models in the M3D-C1 (https://m3dc1.pppl.gov) fusion energy science code. The semi-implicit time integrator in M3D-C1 is composed of many linear solves. The implicit advance of the momentum equation is the most challenging and is the focus of this work. The current production solver in M3D-C1 is a block Jacobi (BJ) preconditioner within a Krylov solver, where blocks group degrees of freedom on planes of constant toroidal coordinate. BJ convergence degrades as the number of planes increases due to the spectral properties of the matrix preconditioned with BJ. The partially magnetic field-aligned, regular toroidal grid structure in M3D-C1 is amenable to semi-coarsening geometric MG in the toroidal direction. This paper develops such a solver and demonstrates competitive performance on a runaway electron model of a SPARC (https://cfs.energy/technology/sparc) disruption, and superior robustness on a stellarator model on which the BJ solver fails to converge.
The development of implicit time integration capabilities for axisymmetric full-F continuum simulations of ion neoclassical transport is reported. The approach involves the implicit treatment of the gyrokinetic Vlasov equation coupled to the nonlinear Fokker–Planck collision model in the long-wavelength limit approximation. To facilitate implicit simulations, advanced preconditioning of individual physics operators is developed, and a global multi-physics preconditioner is constructed by adopting an operator splitting methodology. The algorithm is implemented in the finite-volume code COGENT and is applied to study neoclassical transport properties for both the main ion species and the lithium impurity species in the closed-field-line region of the LTX- β tokamak. The implicit COGENT simulations elucidate the role of non-local transport effects, while demonstrating substantial speedup over the corresponding explicit approach.
BISON is a modern finite-element based nuclear fuel performance code that has been under development at the Idaho National Laboratory (USA) since 2009 [1]. The code is applicable to both steady and transient fuel behavior and can be used to analyze 1D (spherically symmetric), 2D (axisymmetric and generalized plane strain) or 3D geometries. BISON is the fuel performance code used within CASL for LWR fuel under both normal operating and accident conditions. BISON is built using the INL Multiphysics ObjectOriented Simulation Environment, or MOOSE [2, 3]. MOOSE is a massively parallel, finite element-based framework to solve systems of coupled non-linear partial differential equations using the Jacobian-Free Newton Krylov (JFNK) method [4]. This enables investigation of computationally large problems, for example a full stack of discrete pellets in a LWR fuel rod, or every rod in a full reactor core. MOOSE supports the use of complex two and three-dimensional meshes and uses implicit time integration, important for the widely varied time scale in nuclear fuel simulation. An object-oriented architecture is employed which greatly minimizes the programming effort required to add new material and behavioral models. The flexibility of the implicit and fully coupled multiphysics approach comes with a need for constructing suitable approximations for the Jacobian matrix of the coupled system used for either preconditioning a Krylov solve or in a direct Newton solve. Preconditioning options for Bison problems need to be revisited with new preconditioning methods becoming available.
This work presents a new multiscale method for coupling the 3D Maxwell's equations to the 1D telegrapher's equations. While Maxwell's equations are appropriate for modeling complex electromagnetics in arbitrary-geometry domains, simulation cost for many applications (e.g. pulsed power) can be dramatically reduced by representing less complex transmission line regions of the domain with a 1D model. By assuming a transverse electromagnetic (TEM) ansatz for the solution in a transmission line region, we reduce the Maxwell's equations to the telegrapher's equations. Here, we propose a self-consistent finite element formulation of the fully coupled system that uses boundary integrals to couple between the 3D and 1D domains and supports arbitrary unstructured 3D meshes. Additionally, by using a Lagrange multiplier to enforce continuity at the coupling interface, we allow for an absorbing boundary condition to also be applied to non-TEM modes on this boundary. We demonstrate that this feature reduces non-physical reflection and ringing of non-TEM modes off of the coupling boundary. By employing implicit time integration, we ensure a stable coupling, and we introduce an efficient method for solving the resulting linear systems. We demonstrate the accuracy of the new method on two verification problems, a transient O-wave in a rectilinear prism and a steady-state problem in a coaxial geometry, and show the efficiency and weak scalability of our implementation on a cold test of the Z-machine MITL and post-hole convolute.
Abstract Peridigm is a meshfree peridynamics code written in C++ for use on large-scale parallel computers. It was originally developed at Sandia National Laboratories and is currently managed as an open-source, community driven software project. Its primary features include bond-based, state-based, and non-ordinary state-based constitutive models, bond failure laws, contact, and support for explicit and implicit time integration. To date, Peridigm has been used primarily by methods developers focused on solid mechanics and material failure. Peridigm utilizes foundational software components from Sandia’s Trilinos project and was designed for extensibility. This paper provides an overview of the solution methods implemented in Peridigm , a discussion of its software infrastructure, and demonstrates the use of Peridigm for the solution of several example problems.
This work presents the design of nonlinear stabilization techniques for the finite element discretization of Euler equations in both steady and transient form. Implicit time integration is used in the case of the transient form. A differentiable local bounds preserving method has been developed, which combines a Rusanov artificial diffusion operator and a differentiable shock detector. Nonlinear stabilization schemes are usually stiff and highly nonlinear. This issue is mitigated by the differentiability properties of the proposed method. Moreover, in order to further improve the nonlinear convergence, we also propose a continuation method for a subset of the stabilization parameters. The resulting method has been successfully applied to steady and transient problems with complex shock patterns. Numerical experiments show that it is able to provide sharp and well resolved shocks. Furthermore, the importance of the differentiability is assessed by comparing the new scheme with its non-differentiable counterpart. Numerical experiments suggest that, for up to moderate nonlinear tolerances, the method exhibits improved robustness and nonlinear convergence behavior for steady problems. Additionally, in the case of transient problem, we also observe a reduction in the computational cost.
In this work, an implicit and conservative numerical scheme is proposed for the isotropic quantum Fokker-Planck equation describing the evolution of degenerate electrons subject to elastic collisions with other electrons and ions. The electron-ion and electron-electron collision operators are discretized using a discontinuous Galerkin method, and the electron energy distribution is updated by an implicit time integration method. The numerical scheme is designed to satisfy all conservation laws exactly. Numerical tests and comparisons with other modeling approaches are shown to demonstrate the accuracy and conservation properties of the proposed method.
The assembly of nematic colloids relies on long-range elastic interactions that can be manipulated through external stimuli. Confinement and the presence of a hydrodynamic field alter the defect structures and the energetic interactions between the particles. Here, the assembly landscape of nanoparticles embedded in a nematic liquid crystal confined in a nanochannel under a pressure-driven flow is determined. The dynamics of the liquid crystal tensor alignment field is determined through a Poisson-Bracket framework, namely the Stark–Lubensky equations, coupled with the zero-Reynolds momentum equations and the liquid crystal Landau-de Gennes free energy functional. A second order semi-implicit time integration and a three-dimensional Galerkin finite element method are used to resolve flow and nematic fields under several conditions. In general, the zero Reynolds flow displaces the defects around the particles in the upstream direction and renders the surface anchoring ineffective when the flow strength dominates over the nematic elasticity. More importantly, the potential of mean force for particle assembly is non-monotonic independent of surface anchoring. Our results show that the confinement length scale determines the repulsion/attraction transition between colloids, while the flow strength modifies the static defect structure surrounding the particles and determines the magnitude of the energetic barrier for successful assembly. In the attractive regime, the particles move at different rates through the nematic until one particle eventually catches up with the other. This process occurs against or along the direction of flow depending on the flow strength. Ultimately, these results provide a template for engineering and controlling the transport and assembly of nanoparticles under far-from equilibrium conditions in anisotropic media.
We report this paper aims to develop numerical approximations of the Keller–Segel equations that mimic at the discrete level the lower bounds and the energy law of the continuous problem. We solve these equations for two unknowns: the organism (or cell) density, which is a positive variable, and the chemoattractant density, which is a non-negative variable. We propose two algorithms, which combine a stabilized finite element method and a semi-implicit time integration. The stabilization consists of a nonlinear artificial diffusion that employs a graph-Laplacian operator and a shock detector that localizes local extrema. As a result, both algorithms turn out to be nonlinear and can generate cell and chemoattractant numerical densities fulfilling lower bounds. However, the first algorithm requires a suitable constraint between the space and time discrete parameters, whereas the second one does not. We design the latter to attain a discrete energy law on acute meshes. We report some numerical experiments to validate the theoretical results on blowup and nonblowup phenomena. In the blowup setting, we identify a locking phenomenon that relates the L ∞ (Ω)-norm to the L 1 (Ω)-norm limiting the growth of the singularity when supported on a macroelement.