Time Integration Research Articles

Cosmological Radiation Simulation Analysis: Taking Arepo-RT as an Example

With the fast-growing progress in computational physics, astrophysical simulations are made possible and continually upgraded with novel physical models. Contemporarily, many modern simulations feature the capability of cosmological radiation. In order to collate the recent improvements in this field, the numerical solver Arepo-RT is discussed as an important example. This study reorganizes and interprets the main physical models and numerical techniques used in Arepo-RT and similar solvers, with some details simplified to suit common readers. This paper describes simply and briefly the finite volume method, the radiative transfer model, the radiation source model, the Godunov-type schemes and the time integration, all of which are frequently used in modern radiation simulations. Apart from the realizations of the solvers, the results and limitations are also outlined in this paper, including the justification of reduced speed of light approximation, the results compared with JWST, and the issue of missing photons. These results shed light on guiding further improvements for Arepo-RT.

Characterization of Three-Mode Combination Internal Resonances in Electrostatically Actuated Flexible–Flexible Microbeams

Abstract Nonlinear intermodal coupling based on internal resonances in MEMS resonators has advanced significantly over the past two decades for various real-world applications. In this study, we demonstrate the existence of various three-mode combination internal resonances between the first five flexural modes of electrostatically actuated flexible–flexible beams and dynamic modal interaction between three modes via internal resonance. We first calculate the natural frequencies of the beam as a function of the stiffnesses of the transverse and rotational springs of the flexible supports, utilizing both analytical formulation and finite element analysis (FEA). Following this, we identify six combination internal resonances among the first five modes and use applied DC voltage to validate the exactness of one commensurable internal resonance condition (ω2=ω5−ω4). Subsequently, we studied a detailed forced vibration analysis corresponding to this resonance condition by solving the five-mode coupled governing equations through numerical time integration and the method of multiple scales. The results compellingly exhibit three-mode intermodal coupling among the second, fourth, and fifth modes as a function of excitation amplitude and frequency. Alongside this, intriguing nonlinear phenomena such as threshold behavior, saturation phenomena, and autoparametric instability are observed. Finally, this paper provides a systematic methodology for investigating three-mode combination internal resonances and related nonlinear dynamics, offering significant insights that could be used in observing phonon or mechanical lasing phenomena in MEMS resonators.

An elasto-visco-plastic constitutive model for snow: Theory and finite element implementation

Snow exhibits unique mechanical behaviour due to its evolving properties influenced by temperature, stress conditions, and viscous effects. This paper introduces a nonlinear constitutive model for snow, featuring new formulations for the yield function and strain rate potential, and incorporating viscosity, sintering, and degradation effects. The model is numerically integrated into the Abaqus/Standard FEM software using a fully implicit backward Euler method for time integration and Powell’s hybrid method for linearizing the nonlinear system of equations. The robustness and stability of the numerical scheme ensure accurate simulation of snow behaviour under various loading conditions. The model is finally validated against experimental data available in the literature, demonstrating its effectiveness and reliability in capturing the complex mechanical response of snow.

An explicit time integration method for Boussinesq approximation

AbstractThis paper presents a novel approach to the discretization of the Boussinesq equation using finite elements, coupled with an explicit pressure correction scheme. The entire numerical scheme is formulated exclusively in terms of matrix‐vector products, enabling a highly efficient numerical solution. The key advantage of this approach lies in its suitability for fast computations and seamless parallelization on graphics processing units (GPUs). We present the numerical scheme and discuss different computational examples.

Numerical study of bifurcation structures in reflected shock-wave/laminar-boundary-layer interaction within an end-wall tube

Abstract The primary aim of this study is to analyze the unsteady characteristics of the interaction between a reflected shock wave and a laminar boundary layer in an end-wall shock tube. Our direct numerical simulations at shock Mach numbers of Ms=1.9, 2.5, and 3.5 using a fifth-order WENO scheme and three-step Runge-Kutta time integration method revealed inhomogeneity and anisotropy in the shock bifurcation. Surprisingly, the upper and lower bifurcated structures maintain a notably asymmetric flow during the forward propagation of the reflected shock bifurcation. The inverse flow in the bifurcation resembles a crooked earthworm structure, exhibiting high-frequency oscillations indicative of instability. However, at higher shock intensities, the earthworm transforms into a stable strip-like configuration, facilitating the entrapment of inverse flow and leading to rapid bifurcation height growth and early convergence. Additionally, isolated islands with high density, temperature, and pressure emerge in the transitional region behind the bifurcated shocks, due to variations in wave propagation speed.

A thermodynamically‐based fractional model combined viscoelastic‐viscoplastic‐ductile damage with application to fiber‐reinforced polymer composites

AbstractA thermodynamically‐based fractional viscoelastic‐viscoplastic‐damage constitutive model combined with continuous damage mechanics (CDM) theory was established, in order to describe the rate‐dependent nonlinear behavior of fiber‐reinforced polymer composites (FRPCs). The fractional Helmholtz free energy consists of four contributions: viscoelastic (VE), viscoplastic (VP), hardening and damage, in which the VE and VP parts are constructed by fractional Zener and Scott‐Blair (SB) element forms respectively. The constitutive equation is obtained through Helmholtz free energy for the fractional Zener model, and plastic flow and hardening evolution law are all derived in the process. The ductile damage, coupled to both VE and VP free energy parts, is introduced through fractional damage energy release rates to model the degradation of material properties. The corresponding strain energy release rate and dissipation contributions are also derived. The fractional implicit time integration algorithms of proposed model are presented. The model is applied to validate tests of FRPCs under various loading conditions. The model validation and comparison are presented by simulating experimental data and existing models in the literature. And the corresponding evolution of dissipated energy is discussed to further valid the characterization ability of the model.Highlights A thermodynamical fractional constitutive model was developed for FRPCs. The Helmholtz free‐energy potential for fractional Zener model is adopted. The physical significance of fractional order parameters is explored. Fractional implicit integration algorithm of proposed model is implemented. The validation and comparison of the model are presented under various loads.

A local differential quadrature method for the generalized nonlinear Schrödinger (GNLS) equation

A local differential quadrature method based on Fourier series expansion numerically solves the generalized nonlinear Schrodinger equation. For time integration, a Runge-Kutta fourth-order method is used. Matrix stability analysis is used to examine the method's stability. Three test problems involving the motion of a single solitary wave, the interaction of two solitary waves, and a solution that blows up in finite time, respectively, demonstrate the accuracy and efficiency of the provided method. Finally, the numerical results obtained from the presented method are compared with the exact solution and those obtained in earlier works available in the literature.

Deep cryogenic silicon etching for 3D integrated capacitors: A numerical perspective

One promising approach to increase the capacity density of integral microcapacitors, microsupercapacitors, and microbatteries is three-dimensional structure design, where electrodes are exposed in three dimensions instead of conventional in-plane electrodes. Such structures include nanowires, nanotubes, nanopillars, nanoholes, nanosheets, and nanowalls. In this work, a cryogenic silicon etching process suitable for fabrication of structures with high electrode area is proposed. A numeric model of this process is experimentally calibrated and used for pillar array structure sidewall area optimization. The use of adaptive Runge–Kutta–Fehlberg time integrator allows to achieve almost linear overall computation complexity as a function of simulated etching time, despite the linear increase in conductance computation complexity with depth. A rule for choosing optimal geometric structure parameters under technological constraints is formulated. An optimized trefoil-like structure is proposed, resulting in a total 5.5% increase in sidewall area with respect to the hexagonal array of circular pillars, resulting in 20.33 sidewall area per unit chip area for 30 min long etch or 31.80 for 60 min long etch.

A plastic correction algorithm for full-field elasto-plastic finite element simulations: critical assessment of predictive capabilities and improvement by machine learning

AbstractThis paper introduces a new local plastic correction algorithm that is aimed at accelerating elasto-plastic finite element (FE) simulations for structural problems exhibiting localised plasticity (around e.g. notches, geometrical defects). The proposed method belongs to the category of generalised multi-axial Neuber-type methods, which process the results of an elastic prediction point-wise in order to calculate an approximation of the full elasto-plastic solution. The proposed algorithm relies on a rule of local proportionality, which, in the context of J2 plasticity, allows us to express the plastic correction problem in terms of the amplitude of the full mechanical tensors only. This lightweight correction problem can be solved for numerically using a fully implicit time integrator that shares similarities with the radial return algorithm. The numerical capabilities of the proposed algorithm are demonstrated for a notched structure and a specimen containing a distribution of spherical pores, subjected to monotonic and cyclic loading. As a second point of innovation, we show that the proposed local plastic correction algorithm can be further accelerated by employing a simple meta-modelling strategy, with virtually no added errors. At last, we develop and investigate the merits of a deep-learning-based corrective layer designed to reduce the approximation error of the plastic corrector. A convolutional architecture is used to analyse the neighbourhoods of material points and outputs a scalar correction to the point-wise Neuber-type predictions. This optional brick of the proposed plastic correction methodology relies on the availability of a set of full elasto-plastic finite element solutions to be used as a training data-set.

Domain decomposition with local time discretization for the nonlinear Stokes–Biot system

This work presents a domain decomposition method for the fluid-poroelastic structure interaction (FPSI) system, which utilizes local time integration for subproblems. To derive the domain decomposition scheme, we introduce a Lagrange multiplier and define time-dependent Steklov–Poincaré-type operators based on the interface conditions. These operators are employed to transform the coupled system into an evolutionary nonlinear interface problem, which is then solved using an iterative algorithm. This approach provides the flexibility to use different time discretization schemes and step sizes in subdomains, making it an efficient method for simulating multiphysics systems. We present numerical tests for both non-physical and physical problems to demonstrate the accuracy and efficiency of this method.

Spatially mixed implicit–explicit schemes in hydro-mechanically coupled soil dynamics

Mixed implicit–explicit (IMEX) schemes combine elements of both implicit and explicit time integration methods. This allows to efficiently handle boundary value problems with severe differences in stiffness, such as often encountered in soil dynamics. A spatially mixed IMEX scheme, in which different parts of the computational domain may be treated implicitly or explicitly depending on their behaviour, is discussed in this work. Domain decomposition by means of a mortar contact discretisation is used for this purpose and multiple systems of equations are solved. The merits of the scheme in solving problems of soil dynamics are demonstrated, for which the governing equations are extended to hydro-mechanical coupling.

Flow dynamics over cylinder in two‐dimensional benchmark problem: A finite element approximation

This research investigates 2D benchmark flow around a circular cylinder, utilizing the incompressible Navier–Stokes equations alongside the continuity and energy equations. The numerical solution is achieved through finite element discretization for space variable, combined with a second‐order Crank–Nicolson scheme for time integration. The computational results are derived using the FEATFLOW finite element based library package. Our study focuses on the dimensionless form of the flow equations and examines three key dimensionless parameters: drag, lift, and pressure drop. Upon applying finite element method (FEM) discretization, the system is converted into a set of linear or nonlinear ordinary differential equations or algebraic equations for steady‐state scenarios. We then apply the Newton–Raphson method as the outer nonlinear solver, and the multigrid method for efficiently resolving the linear subproblems. To ensure numerical accuracy, we evaluated the and errors, confirming that the experimental order of convergence matches the theoretical predictions. Flow profiles were both graphically represented and tabulated, offering a detailed understanding of the simulation results.

A novel iterative algorithm for nonlinear inelastic dynamic analysis of truss structures; MCB-Spline+Sp time integration method

The focus of this work is on the development of a novel algorithm for the nonlinear dynamic analysis of structures. In the sake of a general and easy-to-implement solution method, the tangential stiffness matrix is calculated by combining the Cauchy-Riemann equations derived from the differentiation of complex functions and using the nodal displacement perturbation vector. An implicit direct time integration method based on a cubic B-Spline is integrated with the incremental-iterative method using a quadrature rule based on the interpolation of spline functions, which is referred to as the MCB-Spline+Sp time integration method. A step-by-step algorithm for developing a computer code to implement this method is provided. The effectiveness of the proposed approach is evaluated through several nonlinear time history analyses, which demonstrate its accuracy and efficiency by observing less computational efforts and fast convergence rate. Moreover, a comparison is made between the developed method and well-established techniques such as the Newmark and Standard-Bathe time integration methods in combination with the iterative Newton-Raphson method.

Magneto-viscoelastic rod model for hard-magnetic soft rods under 3D large deformation: Theory and numerical implementation

The main purpose of this work is to develop a three-dimensional (3D) viscoelastic rod model for hard-magnetic soft (HMS) rods under large deformation which are widely used active structures in soft robotics. To do so, the Simo’s viscoelasticity theory has been rationally incorporated into the geometrically exact 3D curved rod model. The proposed model includes the deformation modes of axial tension, shear, bending, and torsion, which is applicable to the HMS rods with arbitrarily initial curved and twisted geometries under 3D large deformation. The viscoelastic constitutive equations of the HMS rod in the present formulation are formulated, which include the general relaxation functions. To obtain the expression for the magnetic load, the rotation-based magnetic free energy density is introduced, and the governing equations of the HMS rod with magnetic load and body force are presented. To obtain the numerical implementation, an implicit time integration algorithm that simply extends the generalized-α method for the rotation group, and the corresponding variational formulation and its linearization of the rod model are derived. To validate the model, five numerical examples, including 2D dynamic buckling, 3D static, and 3D dynamic problem are presented. The dynamic problems include the dynamic snap-through behavior of a bistable HMS arch and damped oscillation of a quarter arc cantilever under 3D deformation. The simulation results show good agreement with the results reported in the literature.

SUNDIALS time integrators for exascale applications with many independent systems of ordinary differential equations

SUNDIALS time integrators for exascale applications with many independent systems of ordinary differential equations

Computational prediction of impact sound insulation of a timber floor excited by an ISO standard tapping machine

During the last decades, researchers have actively developed procedures to predict impact sound insulation of timber slabs and floors using the finite element method (FEM). Since the simulations have mainly been performed in the low frequencies, a full understanding of their suitability for research and development (R&D) purposes has not yet been created. The object of the study reported in this paper was to determine the applicability of a modern FEM simulation procedure for predicting the normalised impact sound pressure level Ln of a full timber floor, when information provided by material manufacturers is used as input data. Force excitation caused by the ISO standard tapping machine was determined utilizing explicit time integration and FEM. The modelling of sound radiation itself was performed in a frequency domain FEM analysis. The models of a timber floor in its construction stages were constructed before the laboratory measurements of impact sound insulation so that the set-up corresponded to a simulated R&D task. Based on the prediction results, the applied procedure was considered suitable for R&D purposes and needs for further studies were identified.

An open-source time-domain wave-based room acoustic software in Python based on the nodal discontinuous Galerkin method

In this work, we introduce an open-source implementation of a time-domain wave-based room acoustic modeling software package, named DG_RoomAcoustics. In this software, the linear acoustic equations are spatially discretized by the nodal discontinuous Galerkin method, and are integrated in time by either the explicit Runge-Kutta or the arbitrary high-order derivatives (ADER) integration schemes. Following the principles of object-oriented programming paradigm, the software is structured to ensure generic applicability and to facilitate future extensions with additional functionalities (e.g., different time integration schemes, boundary conditions). A comprehensive presentation of the physical and numerical aspects of the problem is provided. A detailed exposition of the code structure and components is presented, all of which are released under an open-source license, fostering community feedback and collaborative contributions for ongoing improvements. A brief overview of the current capabilities of the software is introduced. Future work regarding possible functionality extensions and performance optimizing are discussed.

DG-IMEX method for a two-moment model for radiation transport in the [formula omitted] limit

We consider neutral particle systems described by moments of a phase-space density and propose a realizability-preserving numerical method to evolve a spectral two-moment model for particles interacting with a background fluid moving with nonrelativistic velocities. The system of nonlinear moment equations, with special relativistic corrections to O(v/c), expresses a balance between phase-space advection and collisions and includes velocity-dependent terms that account for spatial advection, Doppler shift, and angular aberration. The model is conservative for the correct O(v/c) Eulerian-frame number density and is consistent, to O(v/c), with Eulerian-frame energy and momentum conservation. This model is closely related to the one promoted by Lowrie et al. [1] and similar to models currently used to study transport phenomena in large-scale simulations of astrophysical environments. The proposed numerical method is designed to preserve moment realizability, which guarantees that the moments correspond to a nonnegative phase-space density. The realizability-preserving scheme consists of the following key components: (i) a strong stability-preserving implicit-explicit (IMEX) time-integration method; (ii) a discontinuous Galerkin (DG) phase-space discretization with carefully constructed numerical fluxes; (iii) a realizability-preserving implicit collision update; and (iv) a realizability-enforcing limiter. In time integration, nonlinearity of the moment model necessitates solution of nonlinear equations, which we formulate as fixed-point problems and solve with tailored iterative solvers that preserve moment realizability with guaranteed global convergence. We also analyze the simultaneous Eulerian-frame number and energy conservation properties of the semi-discrete DG scheme and propose a “spectral redistribution” scheme that promotes Eulerian-frame energy conservation. Through numerical experiments, we demonstrate the accuracy and robustness of this DG-IMEX method and investigate its Eulerian-frame energy conservation properties.

Implicit Peer Triplets in Gradient-Based Solution Algorithms for ODE Constrained Optimal Control

AbstractIt is common practice to apply gradient-based optimization algorithms to numerically solve large-scale ODE constrained optimal control problems. Gradients of the objective function are most efficiently computed by approximate adjoint variables. High accuracy with moderate computing time can be achieved by such time integration methods that satisfy a sufficiently large number of adjoint order conditions and supply gradients with higher orders of consistency. In this paper, we upgrade our former implicit two-step Peer triplets constructed in [Algorithms, 15:310, 2022] to meet those new requirements. Since Peer methods use several stages of the same high stage order, a decisive advantage is their lack of order reduction as for semi-discretized PDE problems with boundary control. Additional order conditions for the control and certain positivity requirements now intensify the demands on the Peer triplet. We discuss the construction of 4-stage methods with order pairs (3, 3) and (4, 3) in detail and provide three Peer triplets of practical interest. We prove convergence of order $$s-1$$ s - 1 , at least, for s-stage methods if state, adjoint and control satisfy the corresponding order conditions. Numerical tests show the expected order of convergence for the new Peer triplets.

Hydrodynamics of multicomponent vesicles: A phase-field approach

In this paper, we introduce a thermodynamically-consistent phase-field model to investigate the hydrodynamics of inextensible multicomponent vesicles in various fluid flows with inertial forces. Our model couples the fluid field, surface concentration field representing chemical species on the membrane, and vesicle dynamics, while enforcing global area and volume constraints through a Lagrange multiplier method. Specifically, we employ full Navier–Stokes equations for the fluid field, the Cahn–Hilliard equations for the species concentration on the membrane, a nonlinear advection–diffusion equation to describe the evolution of the vesicle membrane, and an additional equation to enforce local inextensibility. We utilize a residual-based variational multiscale method for the Navier–Stokes equations and a standard Galerkin finite element framework for the remaining equations. The PDEs are solved using an implicit, monolithic scheme based on the generalized-α time integration method. We extend previous models for homogeneous vesicles (Aland et al., 2014; Valizadeh and Rabczuk, 2022) to multicomponent vesicles, introducing Cahn–Hilliard equations while maintaining thermodynamic consistency. Additionally, we employ isogeometric analysis (IGA) for higher accuracy. We present a variety of two-dimensional numerical examples, including multicomponent vesicles in shear flow and Poiseuille flow with and without obstructions, using the resistive immersed surface method to handle obstructions. Furthermore, we provide three-dimensional simulations of multicomponent vesicles in Poiseuille flow.