Large Time Steps Research Articles

Transient Friction Analysis of Pressure Waves Propagating in Power-Law Non-Newtonian Fluids

Modulated pressure waves propagating in the drilling fluids inside the drill string are a reliable real-time communication technology that transmit data from downhole to the surface during oil and gas drilling. In the analysis of pressure waves’ propagation characteristics, the modeling of transient friction in non-Newtonian fluids remains a great challenge. This paper establishes a numerical model for transient pipe flow of power-law non-Newtonian fluids by using the weighted residual collocation method. Then, the Newton–Raphson method is applied to solve the nonlinear equations. The numerical method is validated by using the theoretical solution of Newtonian fluids and is proven to converge reliably with larger time steps. Finally, the influencing factors of the wall shear stress are analyzed using this numerical method. For shear-thinning fluids, the friction loss of periodic flow decreases with the increase in flow rate, which is opposite to the variation law of friction with the flow rate for stable pipe flow. Keeping the amplitude of pressure pulsation unchanged, an increase in frequency leads to a decrease in velocity fluctuations; therefore, the friction loss decreases with the increase in frequency.

A Structure-Preserving Semi-implicit IMEX Finite Volume Scheme for Ideal Magnetohydrodynamics at all Mach and Alfvén Numbers

We present a divergence-free semi-implicit finite volume scheme for the simulation of the ideal magnetohydrodynamics (MHD) equations which is stable for large time steps controlled by the local transport speed at all Mach and Alfvén numbers. An operator splitting technique allows to treat the convective terms explicitly while the hydrodynamic pressure and the magnetic field contributions are integrated implicitly, yielding two decoupled linear implicit systems. The linearity of the implicit part is achieved by means of a semi-implicit time linearization. This structure is favorable as second-order accuracy in time can be achieved relying on the class of semi-implicit IMplicit–EXplicit Runge–Kutta (IMEX-RK) methods. In space, implicit cell-centered finite difference operators are designed to discretely preserve the divergence-free property of the magnetic field on three-dimensional Cartesian meshes. The new scheme is also particularly well suited for low Mach number flows and for the incompressible limit of the MHD equations, since no explicit numerical dissipation is added to the implicit contribution and the time step is scale independent. Likewise, highly magnetized flows can benefit from the implicit treatment of the magnetic fluxes, hence improving the computational efficiency of the novel method. The convective terms undergo a shock-capturing second order finite volume discretization to guarantee the effectiveness of the proposed method even for high Mach number flows. The new scheme is benchmarked against a series of test cases for the ideal MHD equations addressing different acoustic and Alfvén Mach number regimes where the performance and the stability of the new scheme is assessed.

The fewest switches surface hopping as an optimisation problem

Several ways of defining the probabilities of quantum state transitions in nonadiabatic trajectory surface hopping simulations have been proposed and successfully applied to date. Despite their success, there is a question about the uniqueness of the ways to define such probabilities – some formulations require mathematically-motivated but still ad hoc assumptions to resolve the otherwise under-determined problem. In this work, I present a new approach (termed FSSH-3) to define the hop proposal probabilities. Unlike other approaches, it does not require any ad hoc assumptions and is based on the least squares optimisation of certain functionals combined with particular choices of initial conditions for the target variables. A comprehensive comparative study of several established surface hopping methodologies is conducted using several multiple-state model Hamiltonians. The recently reported formulation of the FSSH-2 method is extended to the case of multiple states and validated. It is demonstrated that the performance of all studied methods, including the new FSSH-3 (in two variants) is robust with respect to integration time step used but only if local diabatisation approach is used. Under such conditions, all methods yield nearly identical results, in excellent agreement with fully quantum simulations, even with sufficiently large integration time steps.

Revisiting Lorenz’s Error Growth Models: Insights and Applications

This entry examines Lorenz’s error growth models with quadratic and cubic hypotheses, highlighting their mathematical connections to the non-dissipative Lorenz 1963 model. The quadratic error growth model is the logistic ordinary differential equation (ODE) with a quadratic nonlinear term, while the cubic model is derived by replacing the quadratic term with a cubic one. A variable transformation shows that the cubic model can be converted to the same form as the logistic ODE. The relationship between the continuous logistic ODE and its discrete version, the logistic map, illustrates chaotic behaviors, demonstrating computational chaos with large time steps. A variant of the logistic ODE is proposed to show how finite predictability horizons can be determined, emphasizing the continuous dependence on initial conditions (CDIC) related to stable and unstable asymptotic values. This review also presents the mathematical relationship between the logistic ODE and the non-dissipative Lorenz 1963 model.

Multi-time-step coupling of peridynamics and classical continuum mechanics for dynamic brittle fracture

Peridynamics (PD), as a nonlocal theory, is well-suited for solving problems with discontinuities, such as cracks. However, the nonlocal effect of peridynamics makes it computationally expensive for dynamic fracture problems in large-scale engineering applications. As an alternative, this study proposes a multi-time-step (MTS) coupling model of PD and classical continuum mechanics (CCM) based on the Arlequin framework. Peridynamics is applied to the fracture domain of the structure, while continuum mechanics is applied to the rest of the structure. The MTS method enables the peridynamic model to be solved at a small time step and the continuum mechanical model is solved at a larger time step. Consequently, higher computational efficiency is achieved for the fracture domain of the structure while ensuring computational accuracy, and this coupling method can be easily applied to large-scale engineering fracture problems.

Learned 1D Passive Scalar Advection to Accelerate Chemical Transport Modeling: A Case Study with GEOS-FP Horizontal Wind Fields

Abstract We developed and applied a machine-learned discretization for one-dimensional (1D) horizontal passive scalar advection, which is an operator component common to all chemical transport models (CTMs). Our learned advection scheme resembles a second-order accurate, three-stencil numerical solver but differs from a traditional solver in that coefficients for each equation term are output by a neural network rather than being theoretically derived constants. We subsampled higher-resolution simulation results—resulting in up to 16× larger grid size and 64× larger time step—and trained our neural-network-based scheme to match the subsampled integration data. In this way, we created an operator that has low resolution (in time or space) but can reproduce the behavior of a high-resolution traditional solver. Our model shows high fidelity in reproducing its training dataset (a single 10-day 1D simulation) and is similarly accurate in simulations with unseen initial conditions, wind fields, and grid spacing. In many cases, our learned solver is more accurate than a low-resolution version of the reference solver, but the low-resolution reference solver achieves greater computational speedup (500× acceleration) over the high-resolution simulation than the learned solver is able to (18× acceleration). Surprisingly, our learned 1D scheme—when combined with a splitting technique—can be used to predict 2D advection and is in some cases more stable and accurate than the low-resolution reference solver in 2D. Overall, our results suggest that learned advection operators may offer a higher-accuracy method for accelerating CTM simulations as compared to simply running a traditional integrator at low resolution. Significance Statement Chemical transport modeling (CTM) is an essential tool for studying air pollution. CTM simulations take a long computing time. Modeling pollutant transport (advection) is the second most computationally intensive part of the model. Decreasing the resolution not only reduces the advection computing time but also decreases accuracy. We employed machine learning to reduce the resolution of advection while keeping the accuracy. We verified the robustness of our solver with several generalization testing scenarios. In our 2D simulation, our solver showed up to 100 times faster simulation with fair accuracy. Integrating our approach to existing CTMs will allow broadened participation in the study of air pollution and related solutions.

A time-averaged method to analyze slender rods moving in tubes

The motion of slender rods in narrow tubes involves the coupling of axial movements and transverse vibrations, in which the transverse vibration exhibits high frequency and chaotic characteristics. These characteristics render the computational cost of the time integration to the motions of rods and tubes expensive. To solve this problem, this paper proposes a novel time-averaged method, which contains two innovations: 1) By time averaging the traditional dynamic equations, the time-averaged dynamic equations describing the rod-tube system are established based on the time-averaged concept used in turbulence problems; 2) Two time-averaged contact models are developed, by fitting the time-averaged sample set using the polishing function and the machine learning, respectively, to construct the relationship between the time-averaged contact gaps and the time-averaged contact forces. The time-averaged method can use a large time step, which is far larger than that required for the time integration of traditional methods, thus reducing the computational cost. Inclined and rolling experiments of the rod-tube system, as well as numerical experiments of the control rod drop in a nuclear reactor, are conducted to validate the accuracy of the time-averaged method. Comparative analyses are carried out between the time-averaged method and six continuous contact models, demonstrating the high computational efficiency of the proposed method. Moreover, the experimental results can provide benchmark data for algorithm verification for scholars studying dynamic contact algorithms.

Viscoelastic Fractional Model with a Non-Uniform Time Discretization for Laminated Glass: Experimental Validation

We discuss a novel approach, based on fractional calculus with a non-uniform time discretization, to numerically simulate interlayer viscoelastic behaviour and associated time-dependent deformation of laminated glass. Reference is made to the classic example of a simply supported laminated glass beam under long-duration loads. The fractional model is compared with some results obtained using the widely used finite element software ABAQUS 2021, which for the viscoelastic properties of the polymeric interlayer, utilizes the more traditional approach based on the Wiechert model and approximation via Prony series of the relaxation function and a uniform discretization of time for the numerical solution. The model is also validated through the comparison with experimental test. The novel approach based on fractional calculus presents two main advantages: 1) the definition of the model parameters from experimental data is simplified; and 2) the numerical implementation is easier and computationally more efficient. When a long observation time is considered, the use of a non-uniform time discretization presents the great advantage of not neglecting any part of the relaxation function. Use of traditional uniform time discretization requires the use of large time steps making it impossible to describe all the changes of the relaxation curve within the large time interval. Practical examples will be presented using viscoelastic models for Trosifol® Extra Stiff (PVB) and SentryGlas® interlayers. This methodology also shows potential to advance next generation standards for the design of structural laminated glass.

A comparative study of emerging material point method and FEM for forming simulation of textile reinforcements

For forming simulations of fabric composites, nonlinear Finite Element Method/FEM has been a long-standing tool to predict and mitigate defects such as wrinkling. However, small-time step requirements in explicit FEM codes, numerical instabilities, and large computational time are among challenges reported. This study presents an alternative fast forming simulation technique through an application of the so-called Material Point Method/MPM, which enables the use of much larger time steps along with fewer numerical instabilities. As a preliminary step towards assessment of this method, both standard 2D deformation modes and 3D hemispherical forming setups were employed, using a plain fabric weave at dry condition. The MPM results were compared to the conventional FEM simulations, as well as to the physical experiments. Notably, the MPM method showed a runtime 20 times faster than its FEM counterpart (under a comparable mesh size), yet with the same reliability in forming predictions as verified by experiments.

An Integral-like Numerical Approach for Solving Burgers’ Equation

The Burgers’ equation, commonly appeared in the study of turbulence, fluid dynamics, shock waves, and continuum mechanics, is a crucial part of the dynamical core of any numerical weather model, influencing simulated meteorological phenomena. While past studies have suggested several robust numerical approaches for solving the equation, many are too complicated for practical adaptation and too computationally expensive for operational deployment. This paper introduces an unconventional approach based on spline polynomial interpolations and the Hopf-Cole transformation. Using Taylor expansion to approximate the exponential term in the Hopf-Cole transformation, the analytical solution of the simplified equation is discretized to form our proposed scheme. The scheme is explicit and adaptable for parallel computing, although certain types of boundary conditions need to be employed implicitly. Three distinct test cases were utilized to evaluate its accuracy, parallel scalability, and numerical stability. In the aspect of accuracy, the schemes employed cubic and quintic spline interpolation perform equally well, managing to reduce the <i>ӏ</i><sub>1</sub>, <i>ӏ</i><sub>2</sub>, and <i>ӏ</i><sub>∞</sub> error norms down to the order of 10<sup>−4</sup>. Parallel scalability observed in the weak-scaling experiment depends on time step size but is generally as good as any explicit scheme. The stability condition is <i>ν</i>∆<i>t</i>/∆<i>x</i><sup>2</sup> > 0.02, including the viscosity coefficient <i>ν</i> due to the Hopf-Cole transformation step. From the stability condition, the schemes can run at a large time step size ∆<i>t</i> even when using a small grid spacing ∆<i>x</i>, emphasizing its suitability for practical applications such as numerical weather prediction.

Precise integration boundary element method for solving nonlinear transient heat conduction problems with temperature dependent thermal conductivity and heat capacity

The precise integration boundary element method (PIBEM) is used to solve nonlinear transient heat conduction problems with temperature dependent thermal conductivity and heat capacity in this article. At first, a boundary-domain integral equation can be obtained by using the fundamental solution for steady linear heat conduction problems. Secondly, the domain integrals are converted into boundary integrals by using the radial integration method to get a pure boundary integral equation, which can retain the advantage of dimension reduction of the boundary element method. Then, after discretization, a first-order nonlinear differential equation system in time is obtained, and the precise integration method with predictor-corrector technique is used to solve the system of nonlinear equations. In the end, several numerical examples are presented to verify the accuracy and stability of the presented method. The results obtained by the presented method are less affected by the time step size and satisfactory results can be obtained even for a large time step size.

An actuator sector model for wind power applications: a parametric study

Abstract. This paper investigates different actuator sector model implementation alternatives and how they compare to actuator line results. The velocity sampling method, tip/smearing correction, and time step are considered. A good agreement is seen between the line and sector model in the rotor plane and the wake flow. Using the sector model, it was possible to reduce the computational time by 75 % compared to the actuator line model as it is possible to run the simulations with a larger time step without compromising the accuracy considerably. The results suggest that the proposed velocity sampling method produces the closest results to the line model with different tip speed ratios. Moreover, the vortex-based smearing correction applied to the sector model results in the lowest error values, among the considered methods, to correct the radial load distributions. Also, it is seen that reducing the time step compared to the one used for the actuator disc/sector does not provide an advantage considering the increased computational time.

3‐D impulse‐based level‐set method for granular flow modeling

AbstractWe explore the viability of modeling dynamic problems with a new formulation of an impulse‐based Level‐Set DEM (LS‐DEM). The new formulation is stable, fast, and energy conservative. However, it can be numerically stiff when the assembly has substantial mass differences between particles. We also demonstrate the feasibility of modeling deformable structures in a rigid body framework and propose several enhancements to improve the convergence of collision resolution, including a hybrid time integration scheme to separately handle at rest contacts and dynamic collisions. Finally, we extend the impulse‐based LS‐DEM to include arbitrarily shaped topographic surfaces and exploit its algorithmic advantages to demonstrate the feasibility of modeling realistic behavior of granular flows. The new formulation significantly improves the performance of dynamic simulations by allowing larger time steps, which is advantageous for observing the full development of physical phenomena such as rock avalanches, which we present as an illustrative example.

Seismic modeling by combining the finite-difference scheme with the numerical dispersion suppression neural network

Seismic finite-difference (FD) modeling suffers from numerical dispersion including both the temporal and spatial dispersion, which can decrease the accuracy of the numerical modeling. To improve the accuracy and efficiency of the conventional numerical modeling, I develop a new seismic modeling method by combining the FD scheme with the numerical dispersion suppression neural network (NDSNN). This method involves the following steps. First, a training data set composed of a small number of wavefield snapshots is generated. The wavefield snapshots with the low-accuracy wavefield data and the high-accuracy wavefield data are paired, and the low-accuracy wavefield snapshots involve the obvious numerical dispersion including both the temporal and spatial dispersion. Second, the NDSNN is trained until the network converges to simultaneously suppress the temporal and spatial dispersion. Third, the entire set of low-accuracy wavefield data is computed quickly using FD modeling with the large time step and the coarse grid. Fourth, the NDSNN is applied to the entire set of low-accuracy wavefield data to suppress the numerical dispersion including the temporal and spatial dispersion. Numerical modeling examples verify the effectiveness of my proposed method in improving the computational accuracy and efficiency.

MODELING THE DYNAMICS OF DEFORMABLE OBJECTS BASED ON VOLUMETRIC PATCHES OF FREE FORMS

A method for solving nonlinear implicit numerical integration problems based on volumetric patches of free forms for modeling deformable objects with high-speed dynamics is proposed. This method improves the accuracy, consistency and controllability of deformation modeling and animation. The method is characterized by speed, accuracy, stability and uses optimization with region decomposition. The method is well suited for modeling deformable bodies with a large time step, in a wide range of deformation dynamics. We propose decomposed optimization, an optimization method based on free-form patches with region decomposition to minimize incremental potentials at each time step. The method uses quadratic matrix decomposition to combine non-overlapping subregions. The Hessian is evaluated once at the beginning of the time step. The advantages of the method are as follows. Geometric primitives and their mathematical models are proposed, which allow the reasonable application of these primitives to solve problems of volume-oriented modeling. Such requirements are met by volumetric patches of free forms based on analytical perturbation functions relative to the base triangles. The decomposed Hessian is constructed for each subregion and calculated using a set of vertices taken from a complete non-decomposed grid. The weights add the missing second-order Hessian data to the vertices of the subregions from the neighbors along the decomposition boundaries. Thanks to this, the descent is performed along the grid coordinates. There is no need to add gradients. During the descent, the gradient is determined. The Hessians under the region are calculated and factorized in parallel once per time step. They are used as an initializer at each iteration. Then the results are mixed together. This ensures stable and continuous high-quality modeling. An automated and reliable optimization method adapted for modeling nonlinear materials, high-speed dynamics and large deformations is proposed.

Quantized tensor networks for solving the Vlasov–Maxwell equations

The Vlasov–Maxwell equations provide an ab initio description of collisionless plasmas, but solving them is often impractical because of the wide range of spatial and temporal scales that must be resolved and the high dimensionality of the problem. In this work, we present a quantum-inspired semi-implicit Vlasov–Maxwell solver that uses the quantized tensor network (QTN) framework. With this QTN solver, the cost of grid-based numerical simulation of size $N$ is reduced from $O(N)$ to $O(\text {poly}(D))$ , where $D$ is the ‘rank’ or ‘bond dimension’ of the QTN and is typically set to be much smaller than $N$ . We find that for the five-dimensional test problems considered here, a modest $D=64$ appears to be sufficient for capturing the expected physics despite the simulations using a total of $N=2^{36}$ grid points, which would require $D=2^{18}$ for full-rank calculations. Additionally, we observe that a QTN time evolution scheme based on the Dirac–Frenkel variational principle allows one to use somewhat larger time steps than prescribed by the Courant–Friedrichs–Lewy constraint. As such, this work demonstrates that the QTN format is a promising means of approximately solving the Vlasov–Maxwell equations with significantly reduced cost.

Simulating charged particles in electromagnetic fields

Relativistic charged particles in electromagnetic fields follow paths that increase in complexity with an increasingly complex field. When traveling in these fields, the particles are acted upon by an electromagnetic force comprised of an electric component and a magnetic component. While solving for these paths in simple electromagnetic fields can be done analytically, the task becomes significantly more difficult when the fields get more complex. Thus, in these situations, numerical methods are required to find solutions. One of the most well-known such methods is the Boris method, which is explored in this paper. Before this method is applied to any complex situations, its accuracy must be ensured by testing on simple cases which are solvable by hand. These cases include a 1D electric field in the direction of the initial velocity of the particle, and a 1D magnetic field perpendicular to the particle's velocity. With the accuracy proven, the method was applied to the cases of a force-free field, a dipolar field, and a quadrupolar field. In the latter two cases the method produced very interesting results that could provide significant insight that would be very difficult to achieve analytically. In the case of the force-free field, however, the method shows some limitations, as a precise cancellation of the force produced by the electric and magnetic fields is required to produce a straight line and the Boris method has some difficulty achieving this, especially when using a large time step.

Efficiently linear and unconditionally energy-stable time-marching schemes with energy relaxation for the phase-field surfactant model

In this article, we consider consistent, efficient, and accurate numerical approximations for a recently developed Cahn–Hilliard type phase-field surfactant model enjoying the bounded-from-below condition of total energy. The phase-field surfactant model naturally satisfies the energy dissipation property in time. Based on this feature, we construct linearly decoupled and unconditionally energy-stable time-marching schemes based on an efficient variant of scalar auxiliary variable approach. To enhance the consistency of discrete modified energy and discrete original energy, a practical energy relaxation technique is adopted. We analytically proved the unique solvability and unconditional energy stability of the proposed schemes. In each time step, only two linear elliptic type equations with constants coefficients need to be solved and several algebraic operations need to be performed. Therefore, the proposed schemes are highly efficient and easy to implement. The space is discretized by the finite difference method and the linear multigrid algorithm is used to accelerate the convergence. The numerical results indicate that the proposed method not only has desired accuracy in time but also satisfies the energy dissipation law even if larger time steps are used. Furthermore, the surfactant laden spinodal decompositions and surfactant absorptions can be well simulated.

Extensions and investigations of space‐time generalized Riemann problems numerical schemes for linear systems of conservation laws with source terms

AbstractThe space‐time generalized Riemann problems method allows to obtain numerical schemes of arbitrary high order that can be used with very large time steps for systems of linear hyperbolic conservation laws with source term. They have been introduced in Berthon et al. (J. Sci. Comput. 55 (2013), 268–308.) in 1D and on 2D unstructured meshes made of triangles. The objective of this article is to complement them in order to answer some important questions arising when they are involved. The formulation is described in detail on quadrangle meshes, the choice of approximation basis is discussed and Legendre polynomials are used in practical cases. The addition of a limiter to preserve certain properties without compromising accuracy is also considered. Finally, the asymptotic behavior of the scheme in the diffusion regime is studied.

Asymptotic-preserving gyrokinetic implicit particle-orbit integrator for arbitrary electromagnetic fields

We extend the asymptotic preserving and energy conserving time integrator for charged-particle motion developed in Ricketson and Chacón (2020) [36] to include finite Larmor-radius (FLR) effects in the presence of electric-field length-scales comparable to the particle gyro-radius (the gyro-kinetic limit). We introduce two modifications to the earlier scheme. The first is the explicit gyro-averaging of the electric field at the half time-step, along with an analogous modification to the current deposition, which we show preserves total energy conservation in implicit PIC schemes. The number of gyrophase samples is chosen adaptively, ensuring proper averaging for large timesteps and the recovery of full-orbit dynamics in the small time-step limit. The second modification is an alternating large and small time-step strategy that ensures the particle trajectory samples gyrophases evenly. We show that this strategy relaxes the time-step restrictions on the scheme, allowing even larger speed-ups than previously achievable. We demonstrate the new method with several single-particle motion tests in a variety of electromagnetic field configurations featuring gyro-scale variation in the electric field. The results demonstrate the advertised ability to capture FLR effects accurately even when significantly stepping over the gyration time-scale.