Implicit Time Integration Research Articles

Smart modeling of soil-foundation interaction using coupled mechanisms: a numerical framework for liquefaction risk mitigation

This study investigates the impact of nearby structures on the cyclic settlement mechanisms of shallow foundations in liquefiable soils using a numerical model based on Biot’s porous media theory. The model predicts excess pore water pressure and settlement by coupling equilibrium and continuity equations, solved using an implicit time integration scheme. Soil nonlinearity under cyclic loading is represented using generalized plasticity, boundary surfaces, and non-associated models. Three scenarios are simulated to study the effect of spacing between light and heavy foundations and variation in acceleration intensity. Results show that as spacing between foundations increases, lateral displacement and settlement decrease. Excess pore water pressure generation also decreases with increased foundation spacing. Soil just below the foundation exhibits maximum settlement, decreasing with depth. When input acceleration increases from 0.1 g to 0.15 g and 0.2 g, settlement increases by 40%–55% and 90%–110% respectively for both light and heavy foundations, regardless of spacing. Excess pore water pressure also increases sharply with higher acceleration intensity. The findings highlight the importance of considering foundation-soil-foundation interaction effects in liquefaction-prone urban settings and provide insights for designing resilient shallow foundations. The advanced numerical modeling approach offers engineers a more informed way to mitigate liquefaction risk and build safer, more durable structures in earthquake-prone areas.

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 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.

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.

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.

A fully implicit, asymptotic-preserving, semi-Lagrangian algorithm for the time dependent anisotropic heat transport equation

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 χ∥/χ⊥⋙1 (with χ∥, χ⊥ the parallel and perpendicular heat diffusivities, respectively) in complicated magnetic field topologies in an accurate and efficient manner. 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 (χ∥/χ⊥=1010).

A dual-modified implicit time integration method for three-dimensional impact modelling within the framework of the SBFEM

This paper presents a dual-modified implicit time integration method (DMITIM) for three-dimensional impact modelling within the framework of the scaled boundary finite element method (SBFEM). Controlled by a single parameter, this method effectively dissipates high-frequency responses in impact problems without affecting low-frequency responses, ensuring time integration stability. To address the issue of false oscillations during impacts, additional Lagrange multipliers are introduced to supplement the velocity and acceleration corrections at each time integration step. For nonmatching meshes, by leveraging the mesh flexibility of the SBFEM, a rapid and straightforward remeshing method is utilized to establish a node-to-node contact scheme. The numerical results demonstrate that the DMITIM significantly reduces false oscillations in the velocity and contact force during impact events. It also showcases the flexibility of the SBFEM in handling nonmatching meshes. The combination of the SBFEM with the DMITIM provides enhanced simulation accuracy and stability for impact problems compared to the results from ABAQUS using the HHT-α time integration algorithm.

Timestep-dependent element interpolation functions in the method of matched sections on the example of heat conduction problem

The paper is devoted to further elaboration of the method of matched sections as a new technique within the finite element method. Like FEM it supposes that: a) the complex domain is represented as a mesh of nonintersecting simple elements; b) algebraic relations between the main parameters of an element are established from the governing differential equations; c) all relationships from all elements are assembled into one global matrix. On the other hand, it has two distinct features. The first one is that relations between kinematic and inner force parameters (called Connection equations) are derived from the approximate analytical solution of the governing equations rather than by the application of minimization techniques. The second one consists in that the conjugation between elements is provided between the adjacent sides (sections) rather than in the nodes of the elements. In application to the transient 2D heat conduction, it is assumed that for each small rectangular element, the 2D problem can be considered as the combination of two 1D problems – one is x-dependent, and another is y-dependent. Each problem is characterized by two functions – the temperature, T, and heat flux Q. In practical realization for rectangular finite elements, the method is reduced to the determination of eight unknowns for each element – two unknowns on each side, which are related by the connection equations, and the requirement of the temperature continuity at the center of the element. Another salient feature of the paper is an implementation of the original implicit time integration scheme, where the time step becomes the parameter of the element interpolation function within the element, i.e. it determines the behavior of the connection equations. This method was initially proposed by the first author for several 1D problems, and here for the first time, it is applied to 2D problems. The number of tests for rectangular plates exhibits the remarkable properties of the proposed time integration scheme concerning stability, accuracy, and absence of any restrictions as to increasing the time step.

Numerical model for time-dependent cracking and deflection of large-span concrete bridge structures

A new computational framework is proposed to analyse the simultaneous occurrence of multiple physical processes affecting the long-term behavior of large-scale concrete structures. Relying on the additivity of small strain, the framework comprehensively takes the nonlinearity of concrete creep, cyclic creep induced by traffic load, concrete shrinkage, concrete cracking, normal-strength steel rebars and prestressing strands relaxation into account, and other physical processes are further imagined. The algorithm is implemented with Abaqus/Standard with the aid of several user-supplied subroutines. The implicit time integration scheme is deemed appropriate for large-scale structures. The framework is systematically validated by several classical experimental results. Finally, a three-dimensional study of a three-span continuous rigid frame bridge with 270 m main span is performed. Comparison between the experimental and numerical results further reveals the contributions of various potential influencing factors. The proposed model shows promise for predicting the time-dependent performance of large-scale reinforcement concrete (RC) structures.

Fast implicit update schemes for Cahn–Hilliard-type gradient flow in the context of Fourier-spectral methods

This work discusses a way of allowing fast implicit update schemes for the temporal discretization of phase-field models for gradient flow problems that employ Fourier-spectral methods for their spatial discretization. Through the repeated application of the Sherman–Morrison formula we provide a rule for approximations of the inverted tangent matrix of the corresponding Newton–Raphson method with a selectable order. Since the representation of this inversion is exact for a sufficiently high approximation order, the proposed scheme is shown to provide a fixed-point-type iterative solver for gradient flow problems that require the solution of linear systems in the context of an implicit time-integration. While the proposed scheme is applicable to general gradient flow phase-field models, we discuss the scheme in the context of the Cahn–Hilliard equation, the Swift–Hohenberg equation, and the phase-field crystal equation for which we demonstrate the performance of the proposed method in comparison with classical solvers.

Krylov-based adaptive-rank implicit time integrators for stiff problems with application to nonlinear Fokker-Planck kinetic models

We propose a high order adaptive-rank implicit integrators for stiff time-dependent PDEs, leveraging extended Krylov subspaces to efficiently and adaptively populate low-rank solution bases. This allows for the accurate representation of solutions with significantly reduced computational costs. We further introduce an efficient mechanism for residual evaluation and an adaptive rank-seeking strategy that optimizes low-rank settings based on a comparison between the residual size and the local truncation errors of the time-stepping discretization. We demonstrate our approach with the challenging Lenard-Bernstein Fokker-Planck (LBFP) nonlinear equation, which describes collisional processes in a fully ionized plasma. The preservation of the equilibrium state is achieved through the Chang-Cooper discretization, and strict conservation of mass, momentum and energy via a Locally Macroscopic Conservative (LoMaC) procedure. The development of implicit adaptive-rank integrators, demonstrated here up to third-order temporal accuracy via diagonally implicit Runge-Kutta schemes, showcases superior performance in terms of accuracy, computational efficiency, equilibrium preservation, and conservation of macroscopic moments. This study offers a starting point for developing scalable, efficient, and accurate methods for high-dimensional time-dependent problems.

A monolithically implicit time-integration approach for a dislocation-density-based b.c.c. single crystal plasticity model

A monolithically implicit time-integration approach for a dislocation-density-based b.c.c. single crystal plasticity model

Nonlinear dynamic analysis of shear- and torsion-free rods using isogeometric discretization and outlier removal

AbstractIn this paper, we present a discrete formulation of nonlinear shear- and torsion-free rods introduced by Gebhardt and Romero (Acta Mechanica 232(10):3825–3847, 2021) that uses isogeometric discretization and robust time integration. Omitting the director as an independent variable field, we reduce the number of degrees of freedom and obtain discrete solutions in multiple copies of the Euclidean space $$\left( \mathbb {R}^3\right) $$ R 3 , which is larger than the corresponding multiple copies of the manifold $$\left( \mathbb {R}^3 \varvec{\times } S^2\right) $$ R 3 × S 2 obtained with standard Hermite finite elements. For implicit time integration, we choose the same integration scheme as Gebhardt and Romero in (2021) that is a hybrid form of the midpoint and the trapezoidal rules. In addition, we apply a recently introduced approach for outlier removal by Hiemstra et al. (Comput Methods Appl Mech Eng 387:114115, 2021) that reduces high-frequency content in the response without affecting the accuracy, ensuring robustness of our nonlinear discrete formulation. We illustrate the efficiency of our nonlinear discrete formulation for static and transient rods under different loading conditions, demonstrating good accuracy in space, time and the frequency domain. Our numerical example coincides with a relevant application case, the simulation of mooring lines.

Rocking block simulation based on numerical dissipation.

In this paper, a computational approach based on numerical dissipation is proposed to simulate rocking blocks. A rocking block is idealized as a solid body interacting with its foundation through a contact-based formulation. An implicit time integration scheme with numerical dissipation, set to optimally treat dissipation in contact problems, is employed. The numerical dissipation is ruled by the time step and the rocking dissipative phenomenon at impacts is accurately predicted without any damping model. A broad numerical campaign is conducted to define a regression law in analytic form for the setting of the time step, depending on the block size and aspect ratio, the contact stiffness, as well as the coefficient of restitution selected. The so-obtained regression law appears accurate and an a posteriori validation with cases not in the training dataset confirms the effectiveness of the approach. Finally, the comparison with available experimental tests highlights the approach efficacy for free rocking and harmonic loading cases (in a deterministic sense), and for earthquake-like loading cases (in a statistical sense). It is found that rocking blocks with sizes of interest for structural engineering (e.g., cultural heritage structures) can be simulated with time steps within 10-3 ÷ 10-1s, so allowing very fast computations.

A performance study of horizontally explicit vertically implicit (HEVI) time-integrators for non-hydrostatic atmospheric models

We conduct a thorough study of different forms of horizontally explicit and vertically implicit (HEVI) time-integration strategies for the compressible Euler equations on spherical domains typical of nonhydrostatic global atmospheric applications. We compare the computational time and complexity of two nonlinear variants (NHEVI-GMRES and NHEVI-LU) and a linear variant (LHEVI). We report on the performance of these three variants for a number of additive Runge-Kutta methods ranging in order of accuracy from second through fifth, and confirm the expected order of accuracy of the HEVI methods for each time-integrator. To gauge the maximum usable time-step of each HEVI method, we run simulations of a nonhydrostatic baroclinic instability for 100 days and then use this time-step to compare the time-to-solution of each method. The results show that NHEVI-LU is 3x faster than NHEVI-GMRES, and LHEVI is 8x faster than NHEVI-GMRES, for the idealized cases tested. The baroclinic instability and inertia-gravity wave simulations indicate that the optimal choice of HEVI method is LHEVI. For the fastest time-to-solution, the second and third order methods are best although the better accuracy of the high-order methods (particularly the fourth order method) should be considered.In the future, we will report on whether these results hold for more complex problems using, e.g., real atmospheric data and/or a higher model top typical of space weather applications.

New implicit time integration schemes combining high frequency damping with high second order accuracy

It is shown that weighted linear combinations of the generalised-α method and certain related higher order schemes allow for the formulation of unconditionally stable single step time integration methods with improved second order accuracy and more targeted high-frequency damping. It is also shown that, if the user controlled high frequency damping parameter ρ∞ is set to zero, the new schemes can be expressed as linear multistep backward difference formulae and, in a particular case, recover Park's method. The performance of the proposed methods is illustrated in terms of mathematical analysis and a number of linear and nonlinear numerical examples including finite element based solutions of the incompressible Navier-Stokes equations.

Line-Based High-Order Hybrid Scheme on Unstructured Grids for Compressible Flows

This study explores the use of line-based finite difference techniques on unstructured grids employing high-order stencils for practical aerodynamic applications. The primary focus is on finite-difference explicit and compact spatial discretization up to sixth order. Problems studied include compressible flows containing discontinuities, which necessitated the use of a hybrid finite difference–finite volume approach, where the former is shown to be effective and accurate for smooth regions, while the latter is robust in capturing discontinuities and shocks. To differentiate between discontinuous and smooth zones, two distinct shock indicators are used, one based on pressure gradients and the other on the divergence of velocity. To achieve faster convergence, spatial discretization is applied with line-based alternating direction implicit (ADI) time integration. A key focus of this work is to demonstrate the proposed methodology in flow scenarios relevant to aircraft applications that also test various aspects of the methodology. Consequently, the cases studied are isentropic vortex convection, standing normal shock, transonic NACA0012 airfoil, supersonic cylinder, unsteady flow past tandem circular cylinders, and supersonic flow past a forward-facing step. Comparisons were also done against an existing line-based finite-volume technique for unstructured grids that employed Monotonic Upstream-Centered Scheme for Conservation Laws (MUSCL) and fifth-order Weighted Essentially NonOscillatory (WENO5) spatial reconstruction schemes. It is observed that the hybrid finite difference–finite volume approach performed well with excellent agreement against available results, sometimes with far fewer degrees of freedom. Furthermore, the implicit ADI scheme is effective over unstructured grids and provides a three-time speed-up over the traditional Runge–Kutta fourth-order scheme.

A component-splitting implicit time integration for multicomponent reacting flows simulations

A component-splitting implicit time integration for multicomponent reacting flows simulations

Linearly implicit time integration scheme of Lagrangian systems via quadratization of a nonlinear kinetic energy. Application to a rotating flexible piano hammer shank.

This paper presents a general and practical approach for nonlinear energy quadratization based on the Euler-Lagrange formulation of the physical equations. A Scalar Auxiliary Variable -like method based on a phase formulation of the equations is applied. The proposed scheme is linearly implicit, reproduces a discrete equivalent of the power balance. It is applied to a rotating and flexible piano hammer shank. An efficient solving strategy leads to a quasi explicit algorithm which shows quadratic space/time convergence.