Time Step Size Research Articles

A sharp interface Lagrangian-Eulerian method for flexible-body fluid-structure interaction

This paper introduces a sharp-interface approach to simulating fluid-structure interaction (FSI) involving flexible bodies described by general nonlinear material models and across a broad range of mass density ratios. This new flexible-body immersed Lagrangian-Eulerian (ILE) scheme extends our prior work on integrating partitioned and immersed approaches to rigid-body FSI. Our numerical approach incorporates the geometrical and domain solution flexibility of the immersed boundary (IB) method with an accuracy comparable to body-fitted approaches that sharply resolve flows and stresses up to the fluid-structure interface. Unlike many IB methods, our ILE formulation uses distinct momentum equations for the fluid and solid subregions with a Dirichlet-Neumann coupling strategy that connects fluid and solid subproblems through simple interface conditions. As in earlier work, we use approximate Lagrange multiplier forces to treat the kinematic interface conditions along the fluid-structure interface. This penalty approach simplifies the linear solvers needed by our formulation by introducing two representations of the fluid-structure interface, one that moves with the fluid and another that moves with the structure, that are connected by stiff springs. This approach also enables the use of multi-rate time stepping, which allows us to use different time step sizes for the fluid and structure subproblems. Our fluid solver relies on an immersed interface method (IIM) for discrete surfaces to impose stress jump conditions along complex interfaces while enabling the use of fast structured-grid solvers for the incompressible Navier-Stokes equations. The dynamics of the volumetric structural mesh are determined using a standard finite element approach to large-deformation nonlinear elasticity via a nearly incompressible solid mechanics formulation. This formulation also readily accommodates compressible structures with a constant total volume, and it can handle fully compressible solid structures for cases in which at least part of the solid boundary does not contact the incompressible fluid. Selected grid convergence studies demonstrate second-order convergence in volume conservation and in the pointwise discrepancies between corresponding positions of the two interface representations as well as between first and second-order convergence in the structural displacements. The time stepping scheme is also demonstrated to yield second-order convergence. To assess and validate the robustness and accuracy of the new algorithm, comparisons are made with computational and experimental FSI benchmarks. Test cases include both smooth and sharp geometries in various flow conditions. We also demonstrate the capabilities of this methodology by applying it to model the transport and capture of a geometrically realistic, deformable blood clot in an inferior vena cava filter.

A linearly stabilized convolution quadrature method for the time-fractional Allen–Cahn equation

This paper proposes a linearly stabilized convolution quadrature method for solving the time-fractional Allen–Cahn equation. The stability condition is explicitly given such that the method is unconditionally stable for any time step size. The space is discretized by the central difference method. We prove that the fully discrete scheme preserves the discrete maximum principle, the discrete energy is bounded, and the modified discrete energy decays monotonically. Numerical simulations support the theoretical analysis.

Bending fluctuations in semiflexible, inextensible, slender filaments in Stokes flow: Toward a spectral discretization.

Semiflexible slender filaments are ubiquitous in nature and cell biology, including in the cytoskeleton, where reorganization of actin filaments allows the cell to move and divide. Most methods for simulating semiflexible inextensible fibers/polymers are based on discrete (bead-link or blob-link) models, which become prohibitively expensive in the slender limit when hydrodynamics is accounted for. In this paper, we develop a novel coarse-grained approach for simulating fluctuating slender filaments with hydrodynamic interactions. Our approach is tailored to relatively stiff fibers whose persistence length is comparable to or larger than their length and is based on three major contributions. First, we discretize the filament centerline using a coarse non-uniform Chebyshev grid, on which we formulate a discrete constrained Gibbs-Boltzmann (GB) equilibrium distribution and overdamped Langevin equation for the evolution of unit-length tangent vectors. Second, we define the hydrodynamic mobility at each point on the filament as an integral of the Rotne-Prager-Yamakawa kernel along the centerline and apply a spectrally accurate "slender-body" quadrature to accurately resolve the hydrodynamics. Third, we propose a novel midpoint temporal integrator, which can correctly capture the Ito drift terms that arise in the overdamped Langevin equation. For two separate examples, we verify that the equilibrium distribution for the Chebyshev grid is a good approximation of the blob-link one and that our temporal integrator for overdamped Langevin dynamics samples the equilibrium GB distribution for sufficiently small time step sizes. We also study the dynamics of relaxation of an initially straight filament and find that as few as 12 Chebyshev nodes provide a good approximation to the dynamics while allowing a time step size two orders of magnitude larger than a resolved blob-link simulation. We conclude by applying our approach to a suspension of cross-linked semiflexible fibers (neglecting hydrodynamic interactions between fibers), where we study how semiflexible fluctuations affect bundling dynamics. We find that semiflexible filaments bundle faster than rigid filaments even when the persistence length is large, but show that semiflexible bending fluctuations only further accelerate agglomeration when the persistence length and fiber length are of the same order.

Finite-Element Method for the Simulation of Lipid Vesicle/Fluid Interactions in a Quasi–Newtonian Fluid Flow

We present a computational framework for modeling an inextensible single vesicle driven by the Helfrich force in an incompressible, non-Newtonian extracellular Carreau fluid. The vesicle membrane is captured with a level set strategy. The local inextensibility constraint is relaxed by introducing a penalty which allows computational savings and facilitates implementation. A high-order Galerkin finite element approximation allows accurate calculations of the membrane force with high-order derivatives. The time discretization is based on the double composition of the one-step backward Euler scheme, while the time step size is flexibly controlled using a time integration error estimation. Numerical examples are presented with particular attention paid to the validation and assessment of the model’s relevance in terms of physiological significance. Optimal convergence rates of the time discretization are obtained.

Dimensional analysis of vapor bubble growth considering bubble–bubble interactions in flash boiling microdroplets of highly volatile liquid electrofuels

Electrofuels (e-fuels) produced from renewable electricity and carbon sources have gained significant attention in recent years as promising alternatives to fossil fuels for the transportation sector. However, the highly volatile e-fuels, such as short-chain oxymethylene ethers (OMEx) are prone to flash vaporization phenomena, which is associated with the formation and growth of vapor bubbles, followed by explosive bursting of the liquid jet. This phenomenon is important in many practical applications, for example, superheated liquid sprays in gasoline direct injection engines as well as cryogenic engines. The simulation of a flash boiling spray of such highly volatile liquid fuels is numerically challenging due to several reasons, including (1) the complexity of the bubble growth process in the presence of multiple vapor bubbles and (2) the need to use an extremely small time step size to accurately capture the underlying physics associated with the flash boiling process. In this paper, we first present a bubble growth model in flash boiling microdroplets considering bubble–bubble interactions along with the finite droplet size effects. A dimensional analysis of the newly derived Rayleigh–Plesset equation (RPE) with bubble–bubble interactions is then performed for Reynolds numbers of different orders of magnitude to estimate the relative importance of different forces acting on the bubble surface. Based on this analysis, a simplified nondimensional semi-analytical solution for bubble growth, which also includes the bubble–bubble interactions, is derived to estimate the bubble growth behavior with reasonable accuracy using the larger time step sizes for a wide range of operating conditions. The derived semi-analytical solution is shown to be a good approximation for describing the bubble growth rate over the whole lifetime of the bubble, thus making it useful for simulations of superheated sprays with large numbers of droplets and even more bubbles. The bubble–bubble interactions are found to have a significant impact on the bubble growth dynamics and result in delaying the onset of droplet bursting due to the slower growth of the vapor bubble compared to the bubble growth without bubble–bubble interactions. Furthermore, in a comparison with DNS results, the proposed bubble growth model is shown to reasonably capture the impact of bubble interactions leading to smaller volumetric droplet expansion.

Numerical investigation of the effects of symmetric and eccentric earthquake-induced pounding on accelerations and response spectra for two-storey structures

This paper examines the effect of symmetric and eccentric earthquake-induced pounding between adjacent two-storey structures on acceleration times histories and response spectra. The numerical investigation benefits from experimental data issued from an extensive shake table campaign involving two-storey steel–concrete composite structures and slab-to-slab contacts. It is shown that the high frequency branch of the floor response spectra is only related to the maximum impulse acceleration spike during time, characterized by its impulse value, contact duration and time-history shape. In order to numerically reproduce such characteristics, three-dimensional detailed models are set up by modeling the slab with shell elements or hexahedral elements. Previously developed explicit time integrator, the CD-Lagrange scheme, based on non-smooth contact dynamics approach, is employed for ensuring the contact conditions expressed in terms of velocity and impulse, without requiring any contact parameters such as contact stiffness, local damping and restitution coefficient. The eccentric pounding between slabs is simulated by introducing a small defect in the parallelism of the colliding slabs, as measured in experience. The results of eccentric and symmetric pounding, assuming elastic behavior for the steel and concrete materials, are compared with the experimental data in terms of acceleration time-histories and pseudo-acceleration response spectra. It is concluded that symmetric pounding predicts an acceleration spike with a more abrupt increase to the peak and a shorter contact duration in comparison to experimental data, leading to an overestimation of the high frequency branch of the response spectra. On the contrary, eccentric pounding achieves to reproduce the experimental acceleration time histories as well as response spectra up to 400 Hz. Finally, the robustness of the proposed numerical approach is underlined by changing the time step size of the CD-Lagrange scheme as well as by prospecting different relevant choices concerning the Rayleigh damping.

A new absorbing layer approach for solving the nonlinear Schrödinger equation

To simulate waves on unbounded domain, absorbing boundary conditions are usually needed to bound the computational domain and avoid boundary reflections as much as possible. In this paper, a new absorbing layer approach is presented to simulate the soliton propagation based on the cubic nonlinear Schrödinger (NLS) equation on unbounded domain, and a two-level finite difference scheme is constructed for solving this NLS problem. Both the analytical solution and the numerical solution are proved to be stable in L2-norm and l2-norm, respectively, and they decay exponentially in the absorbing (or called lossy) layers. Furthermore, the scheme is shown to be convergent with order O(τ2+h2), where τ and h are the time step and grid size, respectively. Numerical examples are given to illustrate the method and verify its effectiveness.

A uniformly accurate method for the Klein-Gordon-Dirac system in the nonrelativistic regime

A multiscale time integrator Fourier pseudospectral (MTI-FP) method is presented for discretizing the massive Klein-Gordon-Dirac (KGD) system which involves a small dimensionless parameter 0<ε≤1. In the nonrelativistic limit regime, the KGD system admits rapid oscillations in time as ε→0+. In addition, the nonlinear Yukawa interaction and the indefinite Dirac operator bring other significant difficulties. The main idea of the MTI-FP method is to construct a precise multiscale decomposition by the frequency (MDF) to the solution of the KGD system at each time step and then employ the Fourier pseudospectral discretization for the spatial derivatives followed with the exponential wave integrator (EWI) for the time marching. This approach is explicit, easy to implement and preforms significantly better than the classical methods in the literature. More specifically, we rigorously establish the uniform error bounds at O(τ+hm0−1) for all ε∈(0,1] and optimal quadratic temporal error bounds at O(τ2) in the ε=O(1) regime, where τ is the time step size, h is mesh size and m0 depends on the regularity of the solution. Extensive numerical results demonstrate that our error bounds are optimal and sharp. Finally, we apply the MTI-FP method to numerically study the nonrelativistic limit behaviors of the KGD system when ε→0+.

Heat exchanger with semi-circular tubes and mountings for effective heat transfer

Compact heat exchangers are widely being used in various devices owing to its minimized size and efficient operation. The traditionally fitted circular tubes in heat exchanger have been modified as semi-circular tubes with slit. This newly proposed design offers added tubes surface area of 63.69% for better fluid flow and heat exchange. To investigate, two dimensional, unsteady, incompressible and isothermal computational studies have been carried out. The different confinements (BR= 1.54d, 2d and 3d) and gap ratios (g = 0.07, 0.14 and 0.21)have been preferred in order to mimic confined fluid environment in compact heat exchangers over the Reynolds numbers range of 50-200. Ansys Workbench 18 has been employed in modeling, meshing and carrying out the simulations of the problem formulated. Pressure and velocity were coupled using the SIMPLE algorithm. The convective and diffusive terms in the governing equations were discretized using second order upwind and central difference scheme. The time step size was fixed as 0.005 and the simulation was run until flow maturity. Results claim that flow sheds well nearer to tubes with slit at BR 1.54, which is not happening for circular tubes at same BR. The arrangement with the minimum BR and maximum gap ratio offers best performance. Semi-circle combination separated by 0.21d offers best heat transfer with an improvement of 10.51% as compared to circular components at a blockage (BR) ratio of 1.54 and Re of 200.The concept of tube splits instead of circular tubes is more effective in dissipating heat whenever the flow shedding is deteriorating or completely eliminated. The overall performance parameter (ɛ) for tube split is highest when BR 3 and g 0.21.

A two-grid immersed finite element method with the Crank-Nicolson time scheme for semilinear parabolic interface problems

In this paper, we propose and analyze a new two-grid partially penalized immersed finite element method for solving the semilinear parabolic interface problems with meshes independent of the coefficient discontinuity. Based on the corresponding time-discrete system, we can unconditionally derive the optimal error estimates in both the L2 norm and semi-H1 norm, while previous works always require the coupling condition of time step and space size (e.g. condition τ=O(H)). Then, we design a two-grid algorithm based on Newton iteration to deal with nonlinear source term. It is shown, both theoretically and numerically, that the algorithm can achieve asymptotically optimal approximation in L2 norm (or semi-H1 norm) when the mesh size satisfies H=O(h3/2) (or H=O(h3)).

Optimal explicit Runge–Kutta time stepping for density-based finite-volume solvers

In this paper, we generate optimal Runge–Kutta stability polynomials for several finite-volume spatial discretizations. From these stability polynomials we generate Butcher tableaus, and compare the accuracy and efficiency of these new schemes to commonly used existing classical temporal schemes. Theoretical analysis shows that these optimal schemes are expected to be significantly more efficient, allowing larger time-step sizes to be used relative to the number of Runge–Kutta stages. These predictions are then verified in practice for several test cases including an isentropic vortex, Taylor–Green vortex, and a T106A turbine cascade. Significant performance improvements were observed for all of these cases, with negligible impact on solution accuracy. Importantly, these improvements are achieved with only minor modifications to the solver.

Hybrid asynchronous isogeometric Perfectly Matched Layer for transient elastodynamics

This paper presents an isogeometric formulation (IGA) of the Perfectly Matched Layer (PML) for two-dimensional time-domain elastodynamics in unbounded media. The displacement-based unsplit-field PML formulation is first discretized in space using B-spline basis functions of increasing polynomial order p, from 2 to 10, able to deal with the strong solution gradients in the thickness of the PML thanks to the high order continuity provided by the B-spline basis functions. In order to keep the high accuracy of IGA, implicit time integration is chosen in the PML subdomains, whereas the physical domain is integrated with an explicit time integration scheme. The domain decomposition is based on a dual formulation and Heterogeneous Asynchronous Time Integrator (HATI) approach. The coupling between non-conforming IGA-patches for the PML subdomains and the physical domain is carried out using a standard mortar method. It enables us to select for the IGA-PML subdomains both the appropriate space discretization and the appropriate time discretization with an implicit time integrator and a larger time step size. The hybrid multi-time step IGA-PML performance in terms of computation time, is assessed for two-dimensional elastodynamic problems in single and multi-layer media depending on the chosen discretization parameters in space and time.

Numerical Analysis of an SEIR Epidemic Model

This paper is concerned for the numerical methods of susceptible exposed infectious recovered (SEIR) epidemic model of Coronavirus disease factor. The model is handled numerically with three extraordinary numerical scheme, forward Euler, Runge-Kutta Plan (RK-4) and the proposed non-standard finite difference (NSFD) techniques. In epidemic model of infectious illnesses, energy is a demon omen property of the consistent framework since negative worth of a subpopulation is useless. The NSFD technique ends up being more important and trustable numerical system than forward Euler and RK-4 strategies. NSFD technique uses terrifically significant properties of tireless SEIR Coronavirus scourge model like energy and presence of equilibria while forward Euler and RK-4 systems can't hold these characteristics. Furthermore, the proposed NSFD is liberated from time step size while forward Euler and RK-4 depend upon the time step size. The numerical diversions with the aid of a numerical test is presented for the endorsements of the huge number of characteristics.

A general nonlinear order-reduction method based on the referenced nodal coordinate formulation for a flexible multibody system

An accurate and efficient formulation is favorable for dynamic analysis and control in the field of the flexible multibody system. This paper proposes a general nonlinear order-reduction method that can tackle overall motions and large deformations with a significant decrease in degrees of freedom, by incorporating the modal derivative techniques into the referenced nodal coordinate formulation (RNCF) developed earlier. To compute the modal derivatives straightforwardly, the closed-form expression of the tangent stiffness matrix's derivative is obtained. The elastic forces are expressed as cubic polynomials of modal coordinates, such that the geometric nonlinear deformations are explicitly expressed. The effectiveness of the proposed method is validated by several numerical examples. The presented geometrically nonlinear order-reduction method can achieve a great accuracy with a much fewer number of generalized coordinates, and it also inherits the large time step-sizes of the RNCF in dealing with large deformations and high-speed rotations.

Carleman linearization approach for chemical kinetics integration toward quantum computation

The Harrow, Hassidim, Lloyd (HHL) algorithm, known as the pioneering algorithm for solving linear equations in quantum computers, is expected to accelerate solving large-scale linear ordinary differential equations (ODEs). To efficiently combine classical and quantum computers for high-cost chemical problems, non-linear ODEs (e.g., chemical reactions) must be linearized to the highest possible accuracy. However, the linearization approach has not been fully established yet. In this study, Carleman linearization was examined to transform nonlinear first-order ODEs of chemical reactions into linear ODEs. Although this linearization theoretically requires the generation of an infinite matrix, the original nonlinear equations can be reconstructed. For the practical use, the linearized system should be truncated with finite size and the extent of the truncation determines analysis precision. Matrix should be sufficiently large so that the precision is satisfied because quantum computers can treat such huge matrix. Our method was applied to a one-variable nonlinear dot{y}=-{y}^{2} system to investigate the effect of truncation orders and time step sizes on the computational error. Subsequently, two zero-dimensional homogeneous ignition problems for H2–air and CH4–air gas mixtures were solved. The results revealed that the proposed method could accurately reproduce reference data. Furthermore, an increase in the truncation order improved accuracy with large time-step sizes. Thus, our approach can provide accurate numerical simulations rapidly for complex combustion systems.

Numerical study on the water entry of two-dimensional airfoils by BEM

Hydrodynamics of the traditional wedges, cylinders, and cones entering water have been heavily studied. However, the hydrodynamic characteristics of airfoils are practically unknown; the problem is inspired by the recently widely developed fixed-wing unmanned aerial-underwater vehicle, which can fly in the air, cruise in water and frequently enters and exits free water surface. In this study, the hydrodynamic characteristics and fully nonlinear free surface deformation for the water entry of airfoils have been investigated by the time domain boundary element method. First, convergence study based on the NACA0012 airfoil is considered. The time step and element size are adjusted, and the results show that the free surface elevation and pressure distribution agree well with those obtained through Reynolds-averaged Navier–Stokes equations. The results show that the higher the entry velocity, the greater the peak hydrodynamic force. Furthermore, the higher thickness ratio causes greater added mass. Gravity and velocity are found to significantly influence the hydrodynamic behavior, including the time-varying force, while thickness ratio has a great impact on added mass and pressure.

AI-aided multiscale modeling of physiologically-significant blood clots

We have developed an AI-aided multiple time stepping (AI-MTS) algorithm and multiscale modeling framework (AI-MSM) and implemented them on the AiMOS supercomputer. AI-MSM is the first of its kind to integrate multi-physics, including intra-platelet, inter-platelet, and fluid-platelet interactions, into one system. It has simulated a record-setting multiscale blood clotting model of 102 million particles, of which 70 flowing and 180 aggregating platelets, under dissipative particle dynamics to coarse-grained molecular dynamics. By adaptively adjusting timestep sizes to match the characteristic time scales of the underlying dynamics, AI-MTS optimally balances speeds and accuracies of the simulations.

Boosting the convergence of low-variance DSMC by GSIS

The low-variance direct simulation Monte Carlo (LVDSMC) is a powerful method to simulate low-speed rarefied gas flows. However, in the near-continuum flow regime, due to limitations on the time step and spatial cell size, it takes plenty of time to find the steady-state solution. Here we remove these deficiencies by coupling the LVDSMC with the general synthetic iterative scheme (GSIS) which permits the simulation at the hydrodynamic scale rather than the much smaller kinetic scale. As a proof of concept, we propose the stochastic-deterministic coupling method based on the Bhatnagar-Gross-Krook kinetic model. First, macroscopic synthetic equations are derived exactly from the kinetic equation, which not only contain the Navier-Stokes-Fourier constitutive relation, but also encompass the higher-order terms describing the rarefaction effects. Then, the high-order terms are extracted from LVDSMC and fed into synthetic equations to predict the macroscopic properties which are closer to the steady-state solution than LVDSMC. Finally, the state of simulation particles in LVDSMC is updated to reflect the change of macroscopic properties. As a result, the convergence to steady state is greatly accelerated, and the restrictions on cell size and the time step are removed. We conduct the Fourier stability analysis and simulate several canonical rarefied gas flows to demonstrate the advantages of LVDSMC-GSIS: when the Knudsen number is lower than 0.1, it can use the grid size about 10 times larger than that in traditional DSMC, and it can reduce the computational cost by two orders of magnitude in the flow regime.

Error Estimates of High-order Conservative Schemes for the Generalized Nonlinear Schrödinger Equation in One Dimension

In this paper, we use the high-order difference operator Dμ to compensate the central difference operator . Then, we discretize time by the time-midpoint method. Fully discrete schemes are obtained for the generalized nonlinear Schrödinger equation. For the high-order methods, conservation of mass and energy can be proved theoretically in the discrete sense. Then, to establish the optimal error estimates without any requirement on the grid ratio, we adopt the cut-off function technique. The error estimates of the schemes are proved to be of O(τ2 + h2(μ+1)) with time step size τ and mesh size h. Some numerical experiments are given to verify the accuracy order. In addition, the dynamics of the generalized nonlinear Schrödinger equation in one dimension are simulated.

Efficient simulation of waves in heterogeneous domains using the scaled boundary finite element method

AbstractSimulation of wave propagation phenomena is very challenging due to the restrictions placed on the spatial as well as on the temporal discretization. Image‐based domain discretization using hierarchical meshes (2D polygonal elements) in conjunction with the scaled boundary finite element method provides an efficient numerical simulation tool. In this setting, the computational effort is reduced as only a limited number of unique cell patterns are formed. This is beneficial while computing element matrices in the preprocessing stage. When solving the dynamic system, the advantages of similar elements can only be fully exploited for the explicit time‐stepping schemes in connection with mass lumping when the domain is discretized using low‐order elements (p ≤ 2). On the other hand, the time step size may be limited by the stability criterion. Moreover, high‐order elements require fewer degrees of freedom. Hence for a specific application, it is not straightforward to predict superiority of one time‐stepping scheme over the other. A comparative study of implicit and explicit time‐stepping schemes is presented using an artificial two‐material model.