Time Step Size Research Articles

A fast explicit time-splitting spectral scheme for the viscous Cahn–Hilliard equation with nonlocal diffusion operator

The viscous Cahn–Hilliard equation (VCH) is a parameter-dependent model encapsulated Cahn–Hilliard (CH) model and constrained Allen–Cahn (CAC) model. In this paper, we propose a fast explicit time-splitting spectral scheme for the VCH with nonlocal diffusion operator. The numerical scheme is constructed based on the second-order symmetric time-splitting method, which results in a splitting of the original problem into nonlocal linear part and nonlinear part, respectively. The nonlocal linear subproblem can be solved exactly by applying the Fourier spectral discretization in space and thus has no stability restriction on the time-step size. The nonlinear subproblem is solved via a second-order strong stability preserving Runge–Kutta method. The new scheme is time symmetric, fully explicit and conserves the discrete mass exactly. Theoretical analysis shows that the stability of the scheme depends on both the viscosity coefficient and the range of nonlocal interactions. In addition, a rigorously convergence analysis of the scheme is established, and the convergence rate is proved to be second-order in time and spectral-order in space, respectively. Four numerical experiments are performed to demonstrate the effectiveness of the proposed scheme.

A Semi-implicit Finite Volume Scheme for Incompressible Two-Phase Flows

AbstractThis paper presents a mass and momentum conservative semi-implicit finite volume (FV) scheme for complex non-hydrostatic free surface flows, interacting with moving solid obstacles. A simplified incompressible Baer-Nunziato type model is considered for two-phase flows containing a liquid phase, a solid phase, and the surrounding void. According to the so-called diffuse interface approach, the different phases and consequently the void are described by means of a scalar volume fraction function for each phase. In our numerical scheme, the dynamics of the liquid phase and the motion of the solid are decoupled. The solid is assumed to be a moving rigid body, whose motion is prescribed. Only after the advection of the solid volume fraction, the dynamics of the liquid phase is considered. As usual in semi-implicit schemes, we employ staggered Cartesian control volumes and treat the nonlinear convective terms explicitly, while the pressure terms are treated implicitly. The non-conservative products arising in the transport equation for the solid volume fraction are treated by a path-conservative approach. The resulting semi-implicit FV discretization of the mass and momentum equations leads to a mildly nonlinear system for the pressure which can be efficiently solved with a nested Newton-type technique. The time step size is only limited by the velocities of the two phases contained in the domain, and not by the gravity wave speed nor by the stiff algebraic relaxation source term, which requires an implicit discretization. The resulting semi-implicit algorithm is first validated on a set of classical incompressible Navier-Stokes test problems and later also adds a fixed and moving solid phase.

Stochastic modified equations for symplectic methods applied to rough Hamiltonian systems

Abstract We investigate stochastic modified equations to explain the mathematical mechanism of symplectic methods applied to rough Hamiltonian systems. The contribution of this paper is threefold. First, we construct a new type of stochastic modified equation. For symplectic methods applied to rough Hamiltonian systems, the associated stochastic modified equations are proved to have Hamiltonian formulations. Secondly, the pathwise convergence order of the truncated modified equation to the numerical method is obtained by techniques in rough path theory. Thirdly, if increments of noises are simulated by truncated random variables, we show that the error can be made exponentially small with respect to the time step size.

An analysis of high order FEM and IGA for explicit dynamics: Mass lumping and immersed boundaries

SummaryWe investigate the behavior of different shape functions for the discretization of hyperbolic problems. In particular, we consider classical Lagrange polynomials and B‐splines. The studies focus on the performance of the these functions as a spatial discretization approach combined with an explicit time marching scheme. In this regard, a major concern is the maximum eigenvalue that imposes restrictions on the critical time step size and suitable lumping techniques that yield a diagonal mass matrix. The accuracy of the discretization methods is assessed in an asymptotic manner in terms of the convergence of eigenvalues and eigenvectors. Further, the global accuracy is investigated in terms of the full spectrum. The results show that B‐spline discretization with a consistent mass matrix are more accurate than those based on Lagrange shape functions, which holds true in the boundary‐fitted as well as in the immersed setting. On the other hand, Lagrange shape functions are more robust with respect to standard lumping techniques, which cannot be directly applied for B‐splines without loss of accuracy. In general, we observe that none of the standard lumping schemes yields optimal results for B‐splines, even in the boundary‐fitted setting. For the immersed setting, also Lagrange shape functions show a drop in accuracy which depends on the position of the boundary that cuts the element. Several remedies are considered in order to overcome these issues, including interpolatory B‐spline bases as well as eigenvalue stabilization methods. While accuracy and stability can be improved using these remedies, we conclude from our study that it is still an open question, how to design a discretization method that achieves large critical time step sizes in combination with a diagonal mass matrix and high accuracy in the immersed setting. We note that these considerations primarily relate to linear structural dynamics applications, such as for example, structural acoustics. In nonlinear problems, such as automotive crash dynamics, other considerations predominate. An example of a one‐dimensional elastic‐plastic bar impacting a rigid wall is illustrative.

An implicit–explicit relaxation extrapolated Runge–Kutta and energy-preserving finite element method for Klein–Gordon–Schrödinger equations

An implicit–explicit (IMEX) relaxation extrapolated Runge–Kutta (RERK) and energy-preserving finite element method is designed for Klein–Gordon–Schrödinger (KGS) equations with periodic boundary conditions. First, the RERK method is employed to the discretization in time, which can keep the unconditionally stable and energy conservation after selecting the appropriate relaxation parameters γn. An IMEX time-marching discretization is constructed so that only a system of linear equations needs to be solved at any time step. Next, the finite element method is utilized to the discretization in space. By the temporal–spatial splitting technique, optimal error estimates in the L2-norm and H1-norm are obtained with the convergent order O(τs+hk+1) and O(τs+hk), without any time stepsize restrictions. At last, some numerical examples are presented to demonstrate the theoretical results.

A GPU based accelerated solver for simulation of heat transfer during metal casting process

The metal casting process, which is one of the key drivers of the manufacturing industry, involves several physical phenomena occurring simultaneously like fluid flow, phase change, and heat transfer which affect the casting yield and quality. Casting process modeling involves numerical modeling of these phenomena on a computer. In recent decades, this has become an inevitable tool for foundry engineers to make defect-free castings. To expedite computational time graphics processing units (GPUs) are being increasingly used in the numerical modeling of heat transfer and fluid flow. Initially, in this work a CPU based implicit solver code is developed for solving the 3D unsteady energy equation including phase change numerically using finite volume method which predicts the thermal profile during solidification in the metal casting process in a completely filled mold. To address the computational bottleneck, which is identified as the linear algebraic solver based on the bi-conjugate gradient stabilized method, a GPU-based code is developed using Compute Unified Device Architecture toolkit and was implemented on the GPU. The CPU and GPU based codes are then validated against a commercial casting simulation code FLOW-3D CAST® for a simple casting part and against in-house experimental results for gravity die casting of a simple geometry. Parallel performance is analyzed for grid sizes ranging from 10 × 10 × 10 to 210 × 210 × 210 and for three time-step sizes. The performance of the GPU code based on occupancy and throughput is also investigated. The GPU code exhibits a maximum speedup of 308× compared to the CPU code for a grid size of 210 × 210 × 210 and a time-step size of 2 s.

Nonsmooth Reduced Interface Models and Their Use in Co-Simulation of Mechanical Systems

Abstract In a co-simulation setup, the entire system is decomposed into a collection of individual subsystems that are interfaced together, with each subsystem being modeled and integrated separately according to its own requirements. To maintain the interconnectivity and consolidation of the primary system, these subsystems must communicate with each other through the interface and transfer certain information at the end points of a defined time interval termed macro time step. Inside the macro time step, the evolution of the interface variables has to be approximated as information about them will only be available again at the end of the step. In real-time simulations, the size of the macro time step and the accuracy of the approximated interface variables are critical factors; if the interface variables are approximated accurately, the size of the macro time step can be kept large enough to provide interactive rates without loss of accuracy and stability. This work focuses on systems where unilateral contact interactions are important and proposes reduced interface model concepts for such nonsmooth systems. The use of the proposed reduced interface model (RIM) is demonstrated in co-simulation to provide model-based approximation of the interface variables. The advantages of the proposed method are demonstrated through two representative case studies.

The Finite Element Method in Time for Multibody Dynamics

Abstract The generalized-α scheme has become the approach of choice for the time integration of the equations of motion of multibody systems. Despite its simplicity, this scheme presents drawbacks: the time-step size cannot be changed easily, making it difficult to implement time adaptivity, and the solution of periodic problems cannot be found easily. This paper explores an alternative approach based on the finite element method in time. The basic principles underpinning the approach are presented and both time-continuous and time-discontinuous approaches are investigated. Two types of Galerkin schemes will be presented here: the time-continuous and the time-discontinuous schemes. In the former, the displacement field is continuous across interelement boundaries, whereas discontinuities or “jumps” are allowed across interelement boundaries for the latter. Simple problems are treated to identify the best schemes. Families of schemes of various accuracies are presented. The first family, based on time-continuous elements, features schemes that do not present numerical dissipation. Asymptotic annihilation is achieved by the time-discontinuous elements that form the second family. The problem of kinematic constraints is treated within the framework of the finite element method in time. Special emphasis is devoted to the satisfaction of the kinematic constraints and their time derivative within a time element.

An eigenvalue stabilization technique for immersed boundary finite element methods in explicit dynamics

The application of immersed boundary methods in static analyses is often impeded by poorly cut elements (small cut elements problem), leading to ill-conditioned linear systems of equations and stability problems. While these concerns may not be paramount in explicit dynamics, a substantial reduction in the critical time step size based on the smallest volume fraction χ of a cut element is observed. This reduction can be so drastic that it renders explicit time integration schemes impractical. To tackle this challenge, we propose the use of a dedicated eigenvalue stabilization (EVS) technique.The EVS-technique serves a dual purpose. Beyond merely improving the condition number of system matrices, it plays a pivotal role in extending the critical time increment, effectively broadening the stability region in explicit dynamics. As a result, our approach enables robust and efficient analyses of high-frequency transient problems using immersed boundary methods. A key advantage of the stabilization method lies in the fact that only element-level operations are required.This is accomplished by computing all eigenvalues of the element matrices and subsequently introducing a stabilization term that mitigates the adverse effects of cutting. Notably, the stabilization of the mass matrix Mc of cut elements – especially for high polynomial orders p of the shape functions – leads to a significant raise in the critical time step size Δtcr.To demonstrate the efficiency of our technique, we present two specifically selected dynamic benchmark examples related to wave propagation analysis, where an explicit time integration scheme must be employed to leverage the increase in the critical time step size.

Linear complementarity formulations for contact analysis and wear prediction in gears

This work presents linear complementarity based formulations for contact analysis and wear prediction in gear trains. Assuming steady motion of gears transmitting constant input torque without transmission errors, three sets of linear complementarity formulations are proposed, first for rigid gear teeth, second for gear teeth with skin compliance, and third for gear teeth with full compliance, the last two being new in linear complementarity literature. The developed formulations are used for contact analysis of two examples of gear trains, a pair of spur gears in mesh, and a planetary gear train. The accuracy of the method is tested by considering various time step sizes and stability of the formulations demonstrated for wide range of time steps. Contact forces between meshing teeth are calculated at a pre-defined set of points on the path of contact, and then wear is calculated on gear tooth flank using Archard’s wear equation. These wear calculations are compared with existing results in literature which demonstrate the applicability of the developed LC formulations. Future work is projected to include more state-of-art friction models and extension to 3D contact analysis in helical gear trains.

Dynamics of the Fractional-Order Lorenz System Based on Adomian Decomposition Method and Its DSP Implementation

Dynamics and digital circuit implementation of the fractional-order Lorenz system are investigated by employing Adomian decomposition method (ADM). Dynamics of the fractional-order Lorenz system with derivative order and parameter varying is analyzed by means of Lyapunov exponents (LEs), bifurcation diagram, chaos diagram and phase diagram. Results show that the fractional-order Lorenz system has rich dynamical behavior and it is a potential model for application. It is also found that the minimum order is affected by numerical algorithm and time step size. Finally, the fractional-order system is implemented on DSP digital circuit. Phase diagrams generated by the DSP are consistent with that generated by simulation.

A high-order multi-time-step scheme for bond-based peridynamics

A high-order multi-time-step scheme (MTS) for the bond-based peridynamic (PD) model, an extension of classical continuous mechanics widely used for analyzing discontinuous problems like cracks, is proposed. The MTS scheme discretizes the spatial domain with a meshfree method and advances in time with a high-order Runge–Kutta method. To effectively handle discontinuities (cracks) that appear in a local subdomain in the solution, the scheme employs the Taylor expansion and Lagrange interpolation polynomials with a finer time step size, that is, coarse and fine time step sizes for smooth and discontinuous subdomains, respectively, to achieve accurate and efficient simulations. By eliminating unnecessary fine-scale resolution imposed on the entire domain, the MTS scheme outperforms the standard STS scheme for PD by significantly reducing computational costs, particularly for problems with discontinuous solutions, as demonstrated by comprehensive theoretical analysis and numerical experiments.

Determination of Maneuvering Force Coefficients for a Destroyer Model with OpenFOAM

_ Ship maneuvering performance can be predicted using various methods, including physical model tests, rapid simulations using coefficient-based forces, and high-fidelity computational fluid dynamics. This article considers the application of computational fluid dynamics for the prediction of hull force coefficients to be used as inputs to rapid maneuvering simulations. The open-source software OpenFOAM was used to simulate forces on the destroyer model DTMB 5415 for steady drift, oscillatory pure sway motion, and oscillatory yaw motion. Recommendations are provided regarding best practices in a number of areas, including domain size, mesh refinement, time step size, and turbulent modeling. Predicted forces for steady drift angles up to 12 degrees are typically within 10%of experimental values. For oscillatory sway and yaw motions, predicted forces are typically within 25%of experimental values. Keywords maneuvering; OpenFOAM; best modeling practices

Comprehensive investigation of the explicit, positivity preserving methods for the heat equation

In this paper-series, we investigate the performance of 12 explicit non-conventional algorithms. All of them have the convex combination property, thus they are unconditionally stable and preserve the positivity of the solution when they applied to the heat equation. In this part of the series, we construct several 2D systems to find how the errors depend on the time step size. Sweeps for other key parameters will be presented in the next part of the series.

Comprehensive investigation of the explicit, positivity preserving methods for the heat equation

In this paper-series, we investigate the performance of 12 explicit non-conventional algorithms in several 2D systems. All of them have the convex combination property, thus they are unconditionally stable and preserve the positivity of the solution when they are applied to the heat equation. In the first part of the series, we examined how the errors depend on the time step size and running times. Now we present additional numerical test results, where sweeps for parameters such as the stiffness and the wavelength of the initial function will be performed.

A highly parallel algorithm for simulating the elastodynamics of a patient-specific human heart with four chambers using a heterogeneous hyperelastic model

In this paper a highly parallel method is developed for simulating the elastodynamics of a four-chamber human heart with patient-specific geometry. The heterogeneous hyperelastic model is discretized by a finite element method in space and a fully implicit adaptive method in time, and the resulting nonlinear algebraic systems are solved by a scalable domain decomposition algorithm. The deformations of the cardiac muscles are quite complex due to the realistic geometry, the heterogeneous hyperelasticity of the cardiac tissue, and the myocardial fibers with active stresses. Moreover, the deformations in different chambers and at different phases of the cardiac cycle are very different. To simulate all the muscle movements including the atrial diastole, the atrial systole, the isovolumic contraction, the ventricular ejection, the isovolumic relaxation, and the ventricular filling, the temporal-spatial mesh needs to be sufficiently fine, but not too fine so that the overall computing time is manageable, we introduce a baseline mesh in space and a two-level time stepping strategy including a uniform baseline time step size to obtain the desired time accuracy and an adaptive time stepping method within a baseline time step to guarantee the convergence of the nonlinear solver. Through numerical experiments, we investigate the performance of the proposed method with respect to the material coefficients, the fiber orientations, as well as the mesh sizes and the time step sizes. For an unstructured tetrahedral mesh with more than 200 million degree of freedoms, the method scales well for up to 16,384 processor cores for all steps of an entire cardiac cycle.

One-Way Fluid–Structure Interaction Simulation of the Added Damping of T(0,1) Modes of Two Partially Overlapped Identical Plates in Still Water

Abstract Structures with a partially overlapped status in water can be seen in some engineering applications, and the fluid-structure coupling vibration behavior of two partially overlapped identical plates has been studied previously through the experimental method. In this study, the added damping of the T(0,1) modes of two submerged partially overlapped plates is numerically investigated through the one-way fluid-structure interaction (FSI) method. The relevant numerical settings, like the mesh size, calculated periods, time-step size, and vibration amplitude, were tested first. Then, the numerical results were compared with experimental results, and good agreements were found. Finally, numerical results were analyzed. The vibration status of two plates in the joint abutting area or the overlapped area has an important influence on the added mass variation. When the added mass is higher, the phase difference between modal force and vibration displacement is also greater, which is the main reason for the higher added damping. The relationship between the phase difference and the frequency in water can be approximately fitted to a straight line, which can probably be used to predict the added damping variations caused by fluid boundary changes of submerged structures.

Unconditional stability and error estimates of the modified characteristics FEMs for the Micropolar Navier–Stokes Equations

In this paper, the unconditional stability and optimal error estimate of the velocity, pressure and angular velocity for the modified characteristics FEMs of the unsteady Micropolar Naiver–Stokes Equations (MNSE) are presented. In this method, the nonlinear equation is linearized by the characteristic finite element method for dealing with the time derivative term and the convection term. Basing on the characteristic time-discrete system, the error function is split into a temporal error and a spatial error. With a rigorous analysis to the characteristic time-discrete system, we prove that the error between the numerical solution and the solution of the time-discrete system is τ-independent, where τ denotes the time stepsize. The stability results and optimal error estimates in L2 norm and H1 norm will be given. Finally, some numerical results will be provided to confirm our theoretical analysis.

Improved dynamic analysis for plane elasticity problems with a boundary meshfree method and distributed sources

In this paper, we present an improved dynamic analysis method to address plane elasticity problems. Distributed sources are employed in a boundary-free mesh-based approach known as Boundary Distributed Source (BDS). This method, which builds upon a modified Method of Fundamental Solutions (MFS), aims to eliminate singularities without the reliance on fictitious boundaries. Distributed sources are strategically positioned directly on the real boundary, centered over disks. The method's efficiency and simplicity are highlighted by its ability to operate without the need for precise data concerning adjacent points. The dual reciprocity method (DRM) is utilized to approximate the non-homogeneous inertial term. Furthermore, the present formulation investigates several factors that influence the convergence of the proposed method, including the number of collocation points, internal points, the radius of the circular disk, and the time step size. Based on the numerical results, it is evident that the method exhibits remarkable accuracy and efficiency when applied to plane elasticity problems including forced, harmonic, and free analyses. As compared to conventional numerical methods for dynamic analysis, the method outperforms in both accuracy and convergence rate.

Monolithic and local time-stepping decoupled algorithms for transport problems in fractured porous media

Abstract The objective of this paper is to develop efficient numerical algorithms for the linear advection-diffusion equation in fractured porous media. A reduced fracture model is considered where the fractures are treated as interfaces between subdomains and the interactions between the fractures and the surrounding porous medium are taken into account. The model is discretized by a backward Euler upwind-mixed hybrid finite element method in which the flux variable represents both the advective and diffusive fluxes. The existence, uniqueness, as well as optimal error estimates in both space and time for the fully discrete coupled problem are established. Moreover, to facilitate different time steps in the fracture-interface and the subdomains, global-in-time, nonoverlapping domain decomposition is utilized to derive two implicit iterative solvers for the discrete problem. The first method is based on the time-dependent Steklov–Poincaré operator, while the second one employs the optimized Schwarz waveform relaxation (OSWR) approach with Ventcel-Robin transmission conditions. A discrete space-time interface system is formulated for each method and is solved iteratively with possibly variable time step sizes. The convergence of the OSWR-based method with conforming time grids is also proved. Finally, numerical results in two dimensions are presented to verify the optimal order of convergence of the monolithic solver and to illustrate the performance of the two decoupled schemes with local time-stepping on problems of high Péclet numbers.