Conserving energy and momentum in nonlinear dynamics: A simple implicit time integration scheme

A coupled implicit-explicit time integration method for compressible unsteady flows

Dynamical low-rank approximation allows for solving large-scale matrix differential equations (MDEs) with significantly fewer degrees of freedom and has been applied to a growing number of applications. However, most existing techniques rely on explicit time integration schemes. In this work, we introduce a cost-effective Newton’s method for the implicit time integration of stiff, nonlinear MDEs on low-rank matrix manifolds. Our methodology is focused on MDEs resulting from the discretization of random partial differential equations (PDEs), where the columns of the MDE can be solved independently. Cost-effectiveness is achieved by solving the MDE at the minimum number of entries required for a rank- r approximation. We present a novel cross low-rank approximation that requires solving the parametric PDE at r strategically selected parameters and O ( r ) grid points using Newton’s method. The selected random samples and grid points adaptively vary over time and are chosen using the discrete empirical interpolation method or similar techniques. The proposed methodology is developed for high-order implicit multi-step and Runge–Kutta schemes and incorporates rank adaptivity, allowing for dynamic rank adjustment over time to control error. Several analytical and PDE examples, including the stochastic Burgers’ and Gray–Scott equations, demonstrate the accuracy and efficiency of the presented methodology.

SummaryTo improve the computational efficiency of the numerical manifold method for discontinuous deformation simulations, a spatial‐domain coupled explicit‐implicit time integration algorithm is proposed. A subdomain partition algorithm based on a super manifold element is developed for the numerical manifold method to simulate dynamic motions of blocky rock mass. In different subdomains, explicit or implicit time integration method is employed respectively based on its contact and motion status. These subdomains interact through assembling the corresponding explicit or implicit time integration‐based matrices of different rock blocks. The computational efficiency of the discontinuity system under dynamic loading is improved by partially diagonalizing the global matrices. Two verification examples of a sliding block along an inclined plane under a horizontal acceleration excitation and a multiblock system acted on by dynamic forces are studied to examine the accuracy of the proposed numerical method, respectively. A highly fractured rock mass situated on an inclined slope subjected to seismic excitations is then studied to show the computational efficiency of the developed algorithm. The simulated results are in good agreement with those from the versions using purely implicit or explicit time integration algorithm for the numerical manifold method. The computational efficiency is shown to be higher using the proposed algorithm, which demonstrates its potential for application in dynamic analysis of highly fractured rock masses.

Methods for solving the stiff system of ordinary differential equations that constitute nuclear reaction networks are surveyed. Three semi-implicit time integration algorithms are examined; a traditional first-order-accurate Euler method, a fourth-order-accurate Kaps-Rentrop method, and a variable-order Bader-Deuflhard method. These three integration methods are coupled to eight different linear algebra packages. Four of the linear algebra packages operate on dense matrices (LAPACK, LUDCMP, LEQS, GIFT), three of them are designed for the direct solution of sparse matrices (MA28, UMFPACK, Y12M), and one uses an iterative method for sparse matrices (BiCG). The scaling properties and behavior of the 24 combinations (3 time integration methods times 8 linear algebra packages) are analyzed by running each combination on seven different nuclear reaction networks. These reaction networks range from a hardwired 13 isotope α-chain and heavy-ion reaction network, which is suitable for most multidimensional simulations of stellar phenomena, to a 489 isotope reaction network, which is suitable for determining the yields of isotopes lighter than technetium in spherically symmetric models of Type II supernovae. Each of the time integration methods and linear algebra packages are capable of generating accurate results, but the efficiency of the various methods—evaluated across several different machine architectures and compiler options—differ dramatically. If the execution speed of reaction networks that contain less than about 50 isotopes is an overriding concern, then the variable-order Bader-Deuflhard time integration method coupled with routines generated from the GIFT matrix package or LAPACK with vendor-optimized BLAS routines is a good choice. If the amount of storage needed for any reaction network is a concern, then any of the sparse matrix packages will reduce the storage costs by 70%-90%. When a balance between accuracy, overall efficiency, and ease of use is desirable, then the variable-order Bader-Deuflhard time integration method coupled with the MA28 sparse matrix package is a strong choice.

Numerical analysis of time integration errors for nonequilibrium radiation diffusion

A shape design sensitivity analysis (DSA) and optimization of structural transient dynamics are proposed for the finite deformation elastoplastic materials under impact with a rigid surface. A shape variation of the structure is considered using the material derivative approach in continuum mechanics. Hyperelasticitybased multiplicatively decomposed elastoplasticity is used for the constitutive model. The implicit Newmark time integration scheme is used for the structural dynamics. The design sensitivity equation is solved at each converged time step with the same tangent stiffness matrix as response analysis without iteration. The cost of sensitivity computation is more efficient than the cost of response analysis for the implicit time integration method. The efficiency and the accuracy of the proposed method are shown through the design optimization of a vehicle bumper. *Communicated by D. Tortorelli.

This article develops high‐order implicit time‐stepping methods combined with the fourth‐order central essentially‐non‐oscillatory (CENO) scheme for stiff three‐dimensional computational fluid dynamics problems having disparate characteristic time scales. Both aerodynamic and magnetohydrodynamic problems are considered on three‐dimensional multiblock body‐fitted grids with hexahedral cells. Several implicit time integration methods of third‐ and fourth‐order accuracy are considered, including the multistep backward differentiation formulas (BDF4), multistage explicitly singly diagonally implicit Runge‐Kutta (ESDIRK4), and Rosenbrock‐type methods (ROS34POW2). The resulting nonlinear algebraic system of equations is solved via a preconditioned Jacobian‐free inexact Newton–Krylov method with additive Schwarz preconditioning using block‐based incomplete LU decomposition. The performance of the high‐order implicit time‐stepping methods on smooth and stiff problems is compared with a standard fourth‐order explicit Runge‐Kutta (RK4) method. It is shown that the Rosenbrock methods, despite their advantage of only requiring the solution of linear systems, have significant drawbacks in terms of robustness issues for highly nonlinear compressible flows. The implicit BDF4 and ESDIRK4 methods are found to be much more efficient than the explicit fourth‐order RK4 method for a stiff resistive magnetohydrodynamic (MHD) problem discretized with the fourth‐order CENO method. When applied to the problem of vortex shedding governed by the Navier–Stokes equations, an A‐stable ESDIRK4 scheme proved to be the more robust and accurate implicit time‐marching scheme and was able to offer significant speedup compared with the RK4 method. Initial results are also shown for high‐order implicit time integration applied to two problems with discontinuities. The current study represents the first to achieve high‐order implicit time integration for MHD, enabling large time steps and substantial speedups for stiff MHD problems with high‐order accuracy, and it also represents the first to establish high‐order implicit time integration for high‐order CENO in space.

Implicit–explicit time integration of a high-order particle-in-cell method with hyperbolic divergence cleaning

Describes simulations of impact forging processes. Uses the explicit time integration finite element method, which is based on direct time integration of equation of motion, to compute the deformation of the workpiece and the dies. Uses the program developed to simulate the copper blow test performed on a 350,000J counter‐blow hammer. The calculated result reveals a good agreement in the final deformed configurations between the experiment and the explicit simulation. In order to compare this with the explicit method, the implicit time integration rigid‐plastic finite element program considering the inertia effect is also applied to the copper blow test simulation. As a result of the copper blow test simulation using the explicit program and the implicit program, finds that the calculated results have good agreements in available plastic deformation energy, forging load and equivalent plastic strain distribution. Finally, applies the developed program to simulations of multi‐blow forging processes. Presents three major findings from the multi‐blow forging simulations: (1) the continuous analysis technique used for the multi‐blow forging simulations works well; (2) the blow efficiency and the forging load generated by blow operations can be analysed efficiently and simulated results coincide with previous experimental and analytical ones; (3) the geometrical configuration of the workpiece is closely related to blow efficiency.

本文重点研究了基于指数变换的瞬态扩散方程求解方法，时间上采用指数变换减弱刚性问题的影响，并采用时间积分法提高计算精度。几何上采取细网有限差分以提高计算精度，并且采用平方共轭梯度法进行最后的方程组求解。本文基于上述方法开发了三维多群中子时空动力学计算程序。 Transient diffusion equation solution method based on exponential transform is studied mainly. Time integral aspect exponential transform weakens the influence of the stiff problem, and the time integration method is used to improve calculation accuracy. Fine mesh finite difference method is taken to improve the geometry calculation accuracy, and conjugate gradient squared method is used for solving final equations. Based on the above methods, multi-group neutron kinetics model calculation program is developed.

A new implicit-explicit local differential transform method (IELDTM) is derived here for time integration of the nonlinear (2 + 1)-dimensional advection–diffusion-reaction (ADR) equations. The IELDTM is adaptively constructed as a stability preserved and high order time integrator for spatially discretized ADR equations. For spatial discretization of the model equation, the Chebyshev spectral collocation method (ChCM) is utilized. A robust stability analysis and global error analysis of the IELDTM are presented with respect to the direction parameter \(\theta \). With the help of the global error analysis, adaptivity equations are derived to minimize the computational costs of the algorithms. The produced method is shown to eliminate the accuracy disadvantage of the classical \(\theta \)-method and the stability disadvantages of the existing differential transform-based methods. Two examples of the Burgers equation in one and two space dimensions and the Chapman oxygen-ozone ADR model are solved via the ChCM-IELDTM hybridization. The present time integrator is proven to provide more efficient numerical characteristics than the various multi-step and multi-stage time integration methods. The IELDTM is extensively compared with the widely used MATLAB solvers, ode45 and ode15s. The adaptive IELDTM has been proven to integrate the stiff Chapman ADR equations with optimum costs over relatively long-time intervals.

The role of time integration in energy conservation in Smoothed Particle Hydrodynamics fluid dynamics simulations

Two- and three-dimensional fluid-structure interaction computer programs for the simulation of nonlinear dynamics were developed and applied to a number of problems. The programs were created by coupling Arbitrary Lagrangian-Eulerian finite volume fluid dynamics programs with strictly Lagrangian finite element structural dynamics programs. The resulting coupled programs can use either fully explicit or implicit time integration. The implicit time integration is accomplished by iterations of the fluid dynamics pressure solver and the structural dynamics system solver. The coupled programs have been used to solve problems involving incompressible fluids, membrane and shell elements, compressible multiphase flows, explosions in both air and water, and large displacements. In this paper, we present the approach used for the coupling and describe test problems that verify the two-dimensional programs against an experiment and an analytical linear problem. The experiment involves an explosion underwater near an instrumented thin steel plate. The analytical linear problem is the vibration of an infinite cylinder surrounded by an incompressible fluid to a given radius.

Stability of semidiscrete formulations for parabolic problems at small time steps

We present a high-order difference method for problems in elastodynamics involving the interaction of waves with highly nonlinear frictional interfaces. We restrict our attention to two-dimensional antiplane problems involving deformation in only one direction. Jump conditions that relate tractions on the interface, or fault, to the relative sliding velocity across it are of a form closely related to those used in earthquake rupture models and other frictional sliding problems. By using summation-by-parts (SBP) finite difference operators and weak enforcement of boundary and interface conditions, a strictly stable method is developed. Furthermore, it is shown that unless the nonlinear interface conditions are formulated in terms of characteristic variables, as opposed to the physical variables in terms of which they are more naturally stated, the semi-discretized system of equations can become extremely stiff, preventing efficient solution using explicit time integrators. The use of SBP operators also provides a rigorously defined energy balance for the discretized problem that, as the mesh is refined, approaches the exact energy balance in the continuous problem. This enables one to investigate earthquake energetics, for example the efficiency with which elastic strain energy released during rupture is converted to radiated energy carried by seismic waves, rather than dissipated by frictional sliding of the fault. These theoretical results are confirmed by several numerical tests in both one and two dimensions demonstrating the computational efficiency, the high-order convergence rate of the method, the benefits of using strictly stable numerical methods for long time integration, and the accuracy of the energy balance.