Implicit-explicit Runge-Kutta Research Articles

Spatiotemporal pattern formations in stiff reaction-diffusion systems by new time marching methods

Reaction-diffusion systems are extensively used in the modeling of developmental biology and in chemical systems to explain the Turing instability are generally highly stiff in both reaction and diffusion terms. For numerical simulations of stiff reaction-diffusion systems, this paper introduces a novel family of implicit-explicit Runge-Kutta (IMEX RK) type methods. Despite being implicit in nature, the developed methods require no numerical or analytical inversion of the coefficient matrix, that is, computationally explicit. Fourier analysis is used to assess the stability results for the developed methods with the model two-dimensional reaction diffusion equation. The efficiency and robustness of the developed methods are validated by numerical simulations of spatiotemporal patterns for reaction-diffusion systems governing phase-separation, the Schnakenberg model, and the electrodeposition process. The results of computed solutions are also compared to those reported in the literature.

Adaptive local discontinuous Galerkin methods with semi-implicit time discretizations for the Navier-Stokes equations

In this paper, we present a mesh adaptation algorithm for the unsteady compressible Navier-Stokes equations under the framework of local discontinuous Galerkin methods coupled with implicit-explicit Runge-Kutta or spectral deferred correction time discretization methods. In both of the two high order semi-implicit time integration methods, the convective flux is treated explicitly and the viscous and heat fluxes are treated implicitly. The remarkable benefits of such semi-implicit temporal discretizations are that they can not only overcome the stringent time step restriction compared with time explicit methods, but also avoid the construction of the large Jacobian matrix as is done for fully implicit methods, thus are relatively easy to implement. To save computing time as well as capture the flow structures of interest accurately, a local mesh refinement (h-adaptive) technique, in which we present detailed criteria for selecting candidate elements and complete strategies to refine and coarsen them, is also applied for the Navier-Stokes equations. Numerical experiments are provided to illustrate the high order accuracy, efficiency and capabilities of the semi-implicit schemes in combination with adaptive local discontinuous Galerkin methods for the Navier-Stokes equations.

High order well-balanced asymptotic preserving finite difference WENO schemes for the shallow water equations in all Froude numbers

In this paper, high order semi-implicit well-balanced and asymptotic preserving finite difference WENO schemes are proposed for the shallow water equations with a non-flat bottom topography. We consider the Froude number ranging from O(1) to 0, which in the zero Froude limit becomes the “lake equations” for balanced flow without gravity waves. We apply a well-balanced finite difference WENO reconstruction, coupled with a stiffly accurate implicit-explicit (IMEX) Runge-Kutta time discretization. The resulting semi-implicit scheme can be shown to be well-balanced, asymptotic preserving (AP) and asymptotically accurate (AA) at the same time. Both one- and two-dimensional numerical results are provided to demonstrate the high order accuracy, AP property and good performance of the proposed methods in capturing small perturbations of steady state solutions.

On the Construction of Conservative Semi-Lagrangian IMEX Advection Schemes for Multiscale Time Dependent PDEs

This article is devoted to the construction of a new class of semi-Lagrangian (SL) schemes with implicit-explicit (IMEX) Runge-Kutta (RK) time stepping for PDEs involving multiple space-time scales. The semi-Lagrangian (SL) approach fully couples the space and time discretization, thus making the use of RK strategies particularly difficult to be combined with. First, a simple scalar advection-diffusion equation is considered as a prototype PDE for the development of a high order formulation of the semi-Lagrangian IMEX algorithms. The advection part of the PDE is discretized explicitly at the aid of a SL technique, while an implicit discretization is employed for the diffusion terms. In this way, an unconditionally stable numerical scheme is obtained, that does not suffer any CFL-type stability restriction on the maximum admissible time step. Second, the SL-IMEX approach is extended to deal with hyperbolic systems with multiple scales, including balance laws, that involve shock waves and other discontinuities. A conservative scheme is then crucial to properly capture the wave propagation speed and thus to locate the discontinuity and the plateau exhibited by the solution. A novel SL technique is proposed, which is based on the integration of the governing equations over the space-time control volume which arises from the motion of each grid point. High order of accuracy is ensured by the usage of IMEX RK schemes combined with a Cauchy–Kowalevskaya procedure that provides a predictor solution within each space-time element. The one-dimensional shallow water equations (SWE) are chosen to validate the new conservative SL-IMEX schemes, where convection and pressure fluxes are treated explicitly and implicitly, respectively. The asymptotic-preserving (AP) property of the novel schemes is also studied considering a relaxation PDE system for the SWE. A large suite of convergence studies for both the non-conservative and the conservative version of the novel class of methods demonstrates that the formal order of accuracy is achieved and numerical evidences about the conservation property are shown. The AP property for the corresponding relaxation system is also investigated.

Numerical approximation to semi-linear stiff neutral equations via implicit–explicit general linear methods

Among the initial value problems of semi-linear neutral equations, there are a class of so-called stiff problems, where the classical explicit methods do not work as the methods have only bounded stability regions, which confine the computational stepsize to be excessively small and thus leads to an unsuccessful calculation. For resolving this difficult issue, ones turned to develop the implicit methods with unbounded stability regions to solve the stiff problems. Nevertheless, it is well-known that the implementation of an implicit method needs a large computational cost. In order to improve the computational efficiency, in Refs. Tan and Zhang (2018, 2020), the authors adopted the implicit–explicit (IMEX) splitting technique to derive the extended IMEX one-leg methods and IMEX Runge–Kutta methods, respectively. Unfortunately, these two methods have the serious order barrier. So far, for stiff neutral equations (SNEs), no IMEX method with order more than two has been found. To improve the computational accuracy and efficiency of IMEX methods, in the present paper, we construct a class of extended implicit–explicit general linear (EIEGL) methods for solving semi-linear SNEs. Under some suitable conditions, an EIEGL method is proved to be stable and convergent of order p whenever the underlying implicit–explicit general linear (IEGL) method has order p and stage order p. With applications to several concrete problems of SNEs, the computational accuracy of EIEGL methods are further illustrated, where we also verify that the convergence order of EIEGL methods can exceed two, namely, third- and fourth-order EIEGL methods can be obtained. Moreover, based on a numerical comparison with the extended implicit general linear (EIGL) methods, the advantage of EIEGL methods in computational efficiency is shown.

A coupled high-order implicit-explicit flux reconstruction lattice Boltzmann method for nearly incompressible thermal flows

In this paper, based on the double distribution function (DDF) model and the Boussinesq approximation, a high-order implicit-explicit flux reconstruction thermal lattice Boltzmann method (FRTLBM) is proposed for simulating nearly incompressible thermal flows. Two discrete velocity Boltzmann equations (DVBEs), one governing the flow field and the other governing the temperature field, are solved by using a high-order flux reconstruction scheme for spatial discretization and a second-order implicit-explicit Runge-Kutta scheme (IMEX RK) for time discretization in the generalized curvilinear coordinates with uniform and non-uniform grids. Numerical validations of the proposed method are implemented by simulating (a) porous plate problem, (b) natural convection in a square cavity, (c) Rayleigh-Bénard convection, (d) natural convection in a concentric annulus and (e) mixed heat transfer from a heated circular cylinder. Numerical results demonstrate that the present method can achieve the high-order accuracy and it is stable and accurate even at the relatively high Rayleigh number (Ra = 107 and 108). The present method is efficient due to the small amount of mesh, the simple algorithm and the time step unlimited by the relaxation time. Since the rescontruction is conducted in a single cell, the present method is compact for parallel computing. Moreover, the present method has a good adaptability of complicated geometries.

An assessment of implicit-explicit time integrators for the pseudo-spectral approximation of Boussinesq thermal convection in an annulus

We analyze the behavior of an ensemble of time integrators applied to the semi-discrete problem resulting from the spectral discretization of the equations describing Boussinesq thermal convection in a cylindrical annulus. The equations are cast in their vorticity-streamfunction formulation that yields a differential algebraic equation (DAE). The ensemble comprises 28 members: 4 implicit-explicit multistep schemes, 22 implicit-explicit Runge-Kutta (IMEX-RK) schemes, and 2 fully explicit schemes used for reference. The schemes whose theoretical order varies from 2 to 5 are assessed for 11 different physical setups that cover laminar and turbulent regimes. Multistep and order 2 IMEX-RK methods exhibit their expected order of convergence under all circumstances. IMEX-RK methods of higher-order show occasional order reduction that impacts both algebraic and differential field variables. We ascribe the order reduction to the stiffness of the problem at hand and, to a larger extent, the presence of the DAE. Using the popular Crank-Nicolson Adams-Bashforth of order 2 (CNAB2) integrator as reference, performance is defined by the ratio of maximum admissible time step to the cost of performing one iteration; the maximum admissible time step is determined by inspection of the time series of viscous dissipation within the system, which guarantees a physically acceptable solution. Relative performance is bounded between 0.5 and 1.5 across all studied configurations. Considering accuracy jointly with performance, we find that 6 schemes consistently outperform CNAB2, meaning that in addition to allowing for a more efficient calculation, the accuracy that they achieve at their operational, dissipation-based limit of stability yields a lower error. In our most turbulent setup, where the behavior of the methods is almost entirely dictated by their explicit component, 13 IMEX-RK integrators outperform CNAB2 in terms of accuracy and efficiency.

A high-order implicit-explicit flux reconstruction lattice Boltzmann method for viscous incompressible flows

In this paper, a high-order implicit-explicit flux reconstruction lattice Boltzmann method (FRLBM) is proposed for simulating viscous incompressible flows. The discrete velocity Boltzmann equation (DVBE) is solved by using a high-order flux reconstruction scheme (FR) for spatial discretization and a high-order implicit-explicit Runge-Kutta scheme (IMEX RK) for time discretization in the generalized curvilinear coordinates with uniform and non-uniform grids to obtain high simulation efficiency. Numerical validations of the proposed method are implemented by simulating (a) Taylor-Green vortex problem, (b) doubly periodic shear layer flow, (c) lid-driven cavity flow, (d) cylindrical Couette flow and (e) flow over a circular cylinder. Numerical results demonstrate that the present method is stable and accurate even at high Reynolds number. Compared with the existing high-order off-lattice Boltzmann methods, the present method is more efficient. In addition, the present method holds the compactness of the standard LBM for parallel computing, but greatly improves the adaptability of complicated geometries.

A unified asymptotic preserving and well-balanced scheme for the Euler system with multiscale relaxation

The design and analysis of a unified asymptotic preserving (AP) and well-balanced scheme for the Euler Equations with gravitational and frictional source terms is presented in this paper. The asymptotic behaviour of the Euler system in the limit of zero Mach and Froude numbers, and large friction is characterised by an additional scaling parameter. Depending on the values of this parameter, the Euler system relaxes towards a hyperbolic or a parabolic limit equation. Standard Implicit–Explicit Runge–Kutta schemes are incapable of switching between these asymptotic regimes. We propose a time semi-discretisation to obtain a unified scheme which is AP for the two different limits. A further reformulation of the semi-implicit scheme can be recast as a fully-explicit method in which the mass update contains both hyperbolic and parabolic fluxes. A space–time fully-discrete scheme is derived using a finite volume framework. A hydrostatic reconstruction strategy, an upwinding of the sources at the interfaces, and a careful choice of the central discretisation of the parabolic fluxes are used to achieve the well-balancing property for hydrostatic steady states. Results of several numerical case studies are presented to substantiate the theoretical claims and to verify the robustness of the scheme.

A cell-centered implicit-explicit Lagrangian scheme for a unified model of nonlinear continuum mechanics on unstructured meshes

A cell-centered implicit-explicit updated Lagrangian finite volume scheme on unstructured grids is proposed for a unified first-order hyperbolic formulation of continuum fluid and solid mechanics, namely the Godunov-Peshkov-Romenski (GPR) model. The scheme provably respects the stiff relaxation limits of the continuous model at the fully discrete level, thus it is asymptotic preserving. Furthermore, the GCL is satisfied by a compatible discretization that makes use of a nodal solver to compute vertex-based fluxes that are used both for the motion of the computational mesh as well as for the time evolution of the governing PDEs. Second order of accuracy in space is achieved using a TVD piecewise linear reconstruction, while an implicit-explicit (IMEX) Runge-Kutta time discretization allows the scheme to obtain higher accuracy also in time. Particular care is devoted to the design of a stiff ODE solver, based on approximate analytical solutions of the governing equations, that plays a crucial role when the visco-plastic limit of the model is approached. We demonstrate the accuracy and robustness of the scheme on a wide spectrum of material responses covered by the unified continuum model that includes inviscid hydrodynamics, viscous heat conducting fluids, elastic and elasto-plastic solids in multidimensional settings.

High order modal Discontinuous Galerkin Implicit–Explicit Runge Kutta and Linear Multistep schemes for the Boltzmann model on general polygonal meshes

Deterministic solutions of the Boltzmann equation represent a real challenge due to the enormous computational effort which is required to produce such simulations and often stochastic methods such as Direct Simulation Monte Carlo (DSMC) are used instead due to their lower computational cost. In this work, we show that combining different technologies for the discretization of the velocity space and of the physical space coupled with suitable time integration techniques, it is possible to compute very precise deterministic approximate solutions of the Boltzmann model in different regimes, from extremely rarefied to dense fluids, with CFL conditions only driven by the hyperbolic transport term. To that aim, we develop modal Discontinuous Galerkin (DG) Implicit–Explicit Runge Kutta schemes (DG-IMEX-RK) and Implicit–Explicit Linear Multistep Methods based on Backward-Finite-Differences (DG-IMEX-BDF) for solving the Boltzmann model on multidimensional unstructured meshes. The solution of the Boltzmann collision operator is obtained through fast spectral methods, while the transport term in the governing equations is discretized relying on an explicit shock-capturing DG method on polygonal tessellations in the physical space. A novel class of WENO-type limiters, based on a shifting of the moments of inertia for each zone of the mesh, is used to control spurious oscillations of the DG solution across discontinuities. The use of Linear Multistep Methods (LMM) allows the Boltzmann solutions to be consistent not only with the compressible Euler limit but also with the Navier–Stokes asymptotic regime. In addition, as numerically proven, they also permit to strongly reduce the computational effort compared to Runge–Kutta approaches while maintaining the same or even larger accuracy. The performances of these different time discretization techniques are measured comparing both precision and efficiency. At the same time, comparisons against simpler relaxation type kinetic models such as the BGK model are proposed. The order of convergence is numerically measured for different regimes and found to agree with the theoretical findings. The new methods are validated considering two-dimensional benchmark test cases typically used in the fluid dynamics community. A prototype engineering problem consisting of a supersonic flow around a NACA 0012 airfoil with space–time-dependent boundary conditions is also presented for which the pressure coefficients are measured.

A compact Fourth-Order Implicit-Explicit Runge-Kutta Type Method for Solving Diffusive Lotka–Volterra System

This paper aims to developed a high-order and accurate method for the solution of one-dimensional Lotka-Volterra-diffusion with Numman boundary conditions. A fourth-order compact finite difference scheme for spatial part combined with implicit-explicit Runge Kutta scheme in temporal are proposed. Furthermore, boundary points are discretized by using a compact finite difference scheme in terms of fourth order accuracy. A key idea for proposed scheme is to take full advantage of method of line (MOL), this is consequently enabling us to use implicit-explicit Runge Kutta method, that are of fourth order in time. We constructed fourth order accuracy in both space and time and is unconditionally stable. This is consequently leading to a reduction in the computational cost of the scheme. Numerical experiments show that the combination of the compact finite difference with IMEX- RK methods give an accurate and reliable for solving the Lotka-Volterra-diffusion.

New optimized implicit-explicit Runge-Kutta methods with applications to the hyperbolic conservation laws

This paper discusses a new class of optimized implicit-explicit Runge-Kutta methods for the numerical solution of the dispersive and non-dispersive hyperbolic systems. Optimized implicit-explicit methods are formulated for the better stability and dispersion properties. Moreover, for present methods inversion of the coefficient matrices is not necessary, which makes these methods very attractive in terms of computational cost and complexity. To validate the efficiency of the developed methods, we have solved the one- and two-dimensional dispersive rotating shallow water equations and benchmark problems from acoustics. Computed solutions are also compared with the exact and experimental results available in the literature. The present methods compete well with the existing multi-stage time-integration methods in terms of accurately resolving the physical characteristics for the chosen problems. Furthermore, the computational costs of the proposed methods are significantly lower as compared to the four-stage, fourth-order explicit Runge-Kutta (RK4) method.

A Low Mach Number IMEX Flux Splitting for the Level Set Ghost Fluid Method

Considering droplet phenomena at low Mach numbers, large differences in the magnitude of the occurring characteristic waves are presented. As acoustic phenomena often play a minor role in such applications, classical explicit schemes which resolve these waves suffer from a very restrictive timestep restriction. In this work, a novel scheme based on a specific level set ghost fluid method and an implicit-explicit (IMEX) flux splitting is proposed to overcome this timestep restriction. A fully implicit narrow band around the sharp phase interface is combined with a splitting of the convective and acoustic phenomena away from the interface. In this part of the domain, the IMEX Runge-Kutta time discretization and the high order discontinuous Galerkin spectral element method are applied to achieve high accuracies in the bulk phases. It is shown that for low Mach numbers a significant gain in computational time can be achieved compared to a fully explicit method. Applications to typical droplet dynamic phenomena validate the proposed method and illustrate its capabilities.

Error estimates of local discontinuous Galerkin method with implicit-explicit Runge Kutta for two-phase miscible flow in porous media

Two-phase miscible flow in porous media is described by a coupled system of nonlinear partial differential equations. In this paper, we develop and analyze the fully discrete local discontinuous Galerkin method coupled with the implicit-explicit Runge Kutta (IMEX-RK) method for the problems. Second-order time discretization scheme in two adjacent time layer can be obtained with the IMEX-RK method. For the temporal-spatial discretization schemes, we give the optimal convergence analysis for both the pressure and concentration in L2-norm. We present several numerical examples for two-phase miscible flow problem to demonstrate the effectiveness and robustness of our method.

A Discontinuous Galerkin Method for Three-Dimensional Poroelastic Wave Propagation: Forward and Adjoint Problems

We develop a numerical solver for three-dimensional wave propagation in coupled poroelastic-elastic media, based on a high-order discontinuous Galerkin (DG) method, with the Biot poroelastic wave equation formulated as a first order conservative velocity/strain hyperbolic system. To derive an upwind numerical flux, we find an exact solution to the Riemann problem, including the poroelastic-elastic interface; we also consider attenuation mechanisms both in Biot's low- and high-frequency regimes. Using either a low-storage explicit or implicit-explicit (IMEX) Runge-Kutta scheme, according to the stiffness of the problem, we study the convergence properties of the proposed DG scheme and verify its numerical accuracy. In the Biot low frequency case, the wave can be highly dissipative for small permeabilities; here, numerical errors associated with the dissipation terms appear to dominate those arising from discretisation of the main hyperbolic system. We then implement the adjoint method for this formulation of Biot's equation. In contrast with the usual second order formulation of the Biot equation, we are not dealing with a self-adjoint system but, with an appropriate inner product, the adjoint may be identified with a non-conservative velocity/stress formulation of the Biot equation. We derive dual fluxes for the adjoint and present a simple but illuminating example of the application of the adjoint method.

High order semi-implicit weighted compact nonlinear scheme for viscous Burgers’ equations

This paper develops a high order semi-implicit weighted compact nonlinear scheme (WCNS) for one- and two-dimensional viscous Burgers’ equations. This semi-implicit WCNS combines a fifth-order WCNS with a third-order implicit-explicit Runge–Kutta (IMEX-RK) time-stepping method. The advection terms of viscous Burgers’ equations are treated explicitly, while the viscosity terms are treated implicitly. Stability analysis shows that the CFL condition of the semi-implicit WCNS is controlled only by the advection terms. Compared to the explicit time-stepping method, the semi-implicit method has the advantage in terms of computational efficiency. Numerical results validate the accuracy and efficiency of the semi-implicit WCNS. This method can also solve viscous Burgers’ equations with large Reynolds numbers and has high-resolution shock-capturing ability.

On the preserving of the maximum principle and energy stability of high-order implicit-explicit Runge-Kutta schemes for the space-fractional Allen-Cahn equation

We put forward and analyze the high-order (up to fourth) strong stability-preserving implicit-explicit Runge-Kutta schemes for the time integration of the space-fractional Allen-Cahn equation, which inherits the maximum principle preserving and energy stability. The space-fractional Allen-Cahn equation with homogeneous Dirichlet boundary condition is first discretized in the spatial direction by using a second-order fractional centered difference scheme that preserves the semi-discrete maximum principle. It is subsequently integrated in the temporal direction by a class of strong stability-preserving implicit-explicit Runge-Kutta schemes that are specifically designed to preserve the maximum principle to the optimal time step size. The convergence order in the discrete $L^{\infty }$ norm and energy boundedness are provided by using the established maximum principle. Finally, a series of numerical experiments are carried out to demonstrate the high-order convergence, maximum principle preserving, and energy stability of the proposed schemes.

Modified shallow water model for viscous fluids and positivity preserving numerical approximation

Shallow water equations are widely used in the simulation of those geophysical flows for which the flow horizontal length scale is much greater than the vertical one. Inspired by the example of lava flows, we consider here a modified model with an additional transport equation for a scalar quantity (e.g., temperature), and the derivation of the shallow water equations from depth-averaging the Navier-Stokes equations is presented. The assumption of constant vertical profiles for some of the model variables is relaxed allowing the presence of vertical profiles, and it follows that the non-linearity of the flux terms results in the introduction of appropriate shape coefficients.The space discretization of the resulting system of hyperbolic partial differential equations is obtained with a modified version of the finite volume central-upwind scheme introduced by Kurganov and Petrova in 2007. The time discretization is based on an implicit-explicit Runge-Kutta method which couples properly the hyperbolic part and the stiff source terms, avoiding the use of a very small time step; the use of complex arithmetic increases accuracy in the implicit treatment of stiff terms. The whole scheme is proved to preserve the positivity of flow thickness and the stationary steady-states.Some numerical experiments are performed to validate the proposed method and to show the incidence on the numerical solutions of shape coefficients introduced in the model.

A linear, high-order, and unconditionally energy stable scheme for the epitaxial thin film growth model without slope selection

The epitaxial thin film growth model without slope selection is the L2-gradient flow of energy with a logarithmic potential in terms of the gradient of a height function. A challenge to numerically solving the model is how to treat the nonlinear term to preserve energy stability without compromising accuracy and efficiency. To resolve this problem, we present a high-order energy stable scheme by placing the linear and nonlinear terms in the convex and concave parts, respectively, and employing the specially designed implicit–explicit Runge–Kutta method. As a result, our scheme is linear, high-order accurate in time, and unconditionally energy stable. We show analytically that the scheme is unconditionally uniquely solvable and energy stable. Numerical experiments are presented to demonstrate the accuracy, efficiency, and energy stability of the proposed scheme.