Source-linked AI summary
A Two-Stage Fourth Order Time-Accurate Discretization for Lax-Wendroff Type Flow Solvers. I. Hyperbolic Conservation Laws
Jiequan Li, Zhifang Du
TL;DR
High-order time discretization for hyperbolic conservation laws traditionally relies on multistage Runge-Kutta methods. This paper develops a two-stage fourth-order Lax-Wendroff approach using generalized Riemann problem solvers, saving computational cost while retaining high-order accuracy.
Problem
The paper seeks to extend Lax-Wendroff-type schemes to higher-order temporal accuracy as an alternative to classical Runge-Kutta methods.
Method
The method uses a second-order Lax-Wendroff solver, with a generalized Riemann problem solver as its building block, in a two-stage procedure for fourth-order temporal accuracy.
Results
At least about 20% computational cost is saved for 1-D problems and 30% for 2-D problems relative to the standard Runge-Kutta approach, while numerical examples validate accuracy and performance.
Takeaways & Limitations
The Lax-Wendroff approach extends from second to higher-order accuracy without successive replacement of temporal derivatives by spatial derivatives in a one-stage method.
Takeaways & Limitations
The study applies the approach only to hyperbolic conservation laws in a finite-volume framework on rectangular meshes, with numerical stability remaining an open theoretical problem.
Abstract
from arXiv · showhide
In this paper we develop a novel two-stage fourth order time-accurate discretization for time-dependent flow problems, particularly for hyperbolic conservation laws. Different from the classical Runge-Kutta (R-K) temporal discretization for first order Riemann solvers as building blocks, the current approach is solely associated with Lax-Wendroff (L-W) type schemes as the building blocks. As a result, a two-stage procedure can be constructed to achieve a fourth order temporal accuracy, rather than using well-developed four stages for R-K methods. The generalized Riemann problem (GRP) solver is taken as a representative of L-W type schemes for the construction of a two-stage fourth order scheme.
1. Introduction
The paper develops a two-stage fourth-order temporal discretization that extends Lax-Wendroff-type solvers for hyperbolic conservation laws. Using a GRP solver as a building block, it reduces the stages and reconstruction work required by conventional Runge-Kutta approaches while retaining fourth-order temporal accuracy.
- Motivation: High-order CFD methods commonly use Runge-Kutta iterations with first-order numerical fluxes, requiring four stages for fourth-order temporal accuracy.
- Method: The proposed method extends second-order Lax-Wendroff methodology to fourth-order temporal accuracy through two stages using solution and time-derivative information.The stages consist of a Lax-Wendroff step at t = tn and a solution-advancing step at the intermediate time.
- Advantages: Compared with traditional Runge-Kutta schemes, the two-stage method removes two stages of data reconstruction and saves at least about 20% of 1-D and 30% of 2-D computational cost.
- Advantages: The governing equation and solution time derivatives are explicitly incorporated into the Lax-Wendroff approximation, which the authors associate with sharper discontinuity capture.
- Scope: The framework is applicable beyond hyperbolic conservation laws and can be used with DG or finite-difference formulations when corresponding Lax-Wendroff or Cauchy-Kovalevskaya tools are available.
- Method: The approach uses a generalized Riemann problem solver as its representative Lax-Wendroff building block for time-accurate flow computation.The paper also identifies ADER and gas kinetic solvers as alternative building blocks.
2. A high order temporal discretization for time-dependent problems
The section constructs a two-stage fourth-order temporal discretization for Lax-Wendroff-type solvers by introducing an intermediate-time approximation and matching Taylor expansions. The resulting parameters are uniquely fixed by the fourth-order accuracy conditions.
- Construction: The method introduces an intermediate value at t = tn + A∆t and uses it to formulate the middle stage.The intermediate value is constructed with third-order accuracy before determining the final update.
- Accuracy conditions: Taylor expansions determine the coefficients A, B0, B1, C0, and C1 required by the accuracy conditions.The coefficient equations are obtained by comparing the scheme expansion with the direct Taylor expansion of the time integral.
- Accuracy conditions: The iterations (2.5)–(2.7) provide a fourth order accurate approximation at t = tn + ∆t.The parameters are uniquely determined for the fourth-order accuracy requirement.
- Algorithm: The resulting algorithm advances the solution through intermediate values and a final update using the explicitly defined stages.The section presents the procedure as Step 1, defining intermediate values, followed by Step 2, advancing the solution.
- Interpretation: The midpoint stage at tn + A∆t can be interpreted as an Hermite-type approximation to the time integral.The midpoint interpretation distinguishes the construction from classical Runge-Kutta iteration.
3. Fourth order accurate temporal discretization for hyperbolic conservation laws
The paper extends the two-stage Lax-Wendroff construction to hyperbolic conservation laws using instantaneous solutions and temporal derivatives from the governing equation. GRP-based flux evaluation supplies the needed treatment of discontinuities in one and two dimensions.
- One-dimensional method: The hyperbolic-conservation-law extension uses Lax-Wendroff-type solvers together with instantaneous solutions and their temporal derivatives.The governing equation supplies the relation needed to compute temporal variation from spatial information.
- GRP solver: The GRP solver provides analytical or approximate instantaneous solutions, with ADER identified as an approximate alternative.The temporal derivative can be obtained from the GRP solution or related Lax-Wendroff formulations.
- One-dimensional method: The one-dimensional algorithm reconstructs initial data, computes instantaneous values, forms an intermediate state, and advances to the next time level.The procedure uses the Lax-Wendroff solver twice, at tn and at the midpoint time.
- Flux evaluation: The numerical flux at cell interfaces is a central component of the finite-volume update and must be accurately approximated.The discussion identifies interface flux evaluation as crucial to the method.
- Discontinuity treatment: Strong discontinuities can produce large flux discrepancies because jump-based flux errors are not proportional to mesh size or time-step length.The method addresses this issue by analytically solving the generalized Riemann problem and deriving the temporal derivative.
- Accuracy verification: Numerical verification indicates that the two-stage method provides fourth-order temporal accuracy.For linear equations, the quadrature formula is exact; numerical examples further verify fourth-order temporal accuracy.
- Two-dimensional extension: The same two-stage framework extends to two-dimensional problems using finite volumes, rectangular meshes, Gauss quadrature, and quasi-one-dimensional GRP problems at Gauss points.The two-dimensional algorithm repeats reconstruction and instantaneous-value computation before the final update.
4. Numerical Examples
The numerical examples test GRP4-HWENO5 against RK4-WENO5 on scalar conservation laws, one-dimensional Euler equations, and two-dimensional Euler problems. The reported results verify expected fourth-order accuracy and demonstrate performance on smooth flows, shocks, rarefactions, vortices, and blast-wave interactions.
- Experiment design: The experiments cover linear advection, Burgers flow, one-dimensional Euler equations, and two-dimensional Euler equations, using CFL number 0.5 except for the large-density-ratio problem.GRP4-HWENO5 denotes the GRP solver with fifth-order HWENO reconstruction, while RK4-WENO5 uses fifth-order WENO reconstruction.
- Scalar conservation laws: The linear advection results achieve the expected accuracy, with L1, L∞ errors and convergence orders compared for GRP4-HWENO5 and RK4-WENO5.The comparison uses m cells and reports results at t = 10.
- One-dimensional Euler equations: The smooth Euler test verifies the expected accuracy order, while the shock-turbulence density profile is comparable with results from other schemes.The smooth test uses periodic boundary conditions; the shock-turbulence comparison uses 400 cells and a reference solution.
- 2-D Examples: The isentropic vortex test achieves the expected fourth-order accuracy, and two-dimensional Riemann tests resolve small-scale vortices with performance comparable to an adaptive moving-mesh GRP method.The two-dimensional Riemann examples include shock, vortex-sheet, and rarefaction interactions.
- 2-D Examples: The double Mach reflection results are reported with excellent performance, using a high-resolution two-dimensional calculation.The figure uses GRP4-HWENO5 with 960 × 240 cells and shows the solution at t = 0.2.
5. Discussions and Prospectives
The paper develops a two-stage fourth-order temporal discretization based on Lax-Wendroff flow solvers, applies it with HWENO5 spatial reconstruction, and validates it on hyperbolic conservation laws. The approach differs from classical Runge-Kutta discretization and remains limited to the paper’s finite-volume rectangular-mesh setting, with theoretical stability questions unresolved.
- Contributions: The proposed method combines two-stage fourth-order time accuracy with fifth-order HWENO spatial accuracy for hyperbolic conservation laws.Numerical examples assess its accuracy and computational performance on complex flow problems.
- Methodological distinction: The discretization extends Lax-Wendroff evolution beyond second order without successive one-stage replacement of temporal derivatives by spatial derivatives.The authors describe the construction as Hermite-type rather than the classical Runge-Kutta Simpson-type approach.
- Scope: The approach is applied here to hyperbolic conservation laws in a finite-volume framework over rectangular meshes.The authors state that broader applications require Lax-Wendroff-type solvers on other formulations or meshes.
- Performance: Numerical experiments show that the fourth-order scheme can use a CFL number as large as the second-order GRP scheme without losing accuracy.The CFL number can exceed 1/2 when the computed waves are not very strong.
- Open problems: Numerical stability and other theoretical problems remain for further study.The authors characterize the work as a starting point for designing high-order accurate methods.
Appendix A. The GRP solver
The appendix documents the GRP solver used in the implementation and specifies the source-term treatment underlying the one- and two-dimensional formulations.
- Solver documentation: The appendix includes the GRP solver used in coding for completeness and reader convenience.Further details are referred to separate treatments of Euler equations and general hyperbolic systems.
- Source-term setting: The paper focuses on the homogeneous one-dimensional case with g(u, x) ≡ 0.For two-dimensional problems, tangential effects at cell interfaces can be treated as source terms.
A.1. 1-D GRP.
The one-dimensional GRP solver starts from two polynomial pieces of initial data, uses fifth-order HWENO reconstruction in this study, and provides acoustic and genuinely nonlinear versions.
- Initial data: The one-dimensional GRP solver assumes that the initial data consist of two pieces of polynomials.
- Initial data: The study uses HWENO to construct fifth-order polynomial initial data for the GRP solver.
- Solver variants: The GRP solver has acoustic and genuinely nonlinear versions.These are identified as the two solver variants used for different flow regimes.
- Acoustic GRP: The acoustic GRP is designed for weak discontinuities or smooth flows.
A.1.1. Acoustic GRP solver.
The acoustic GRP computes instantaneous time derivatives by linearizing around a reference state and using characteristic information, while strong discontinuities require the genuinely nonlinear GRP solver.
- Strong discontinuities: The jump across the initial data is not necessarily small in the genuinely nonlinear treatment.
- Acoustic GRP: The acoustic GRP linearizes the governing relation around u0 before computing the instantaneous time derivative.
- Characteristic calculation: Its characteristic formulation uses eigenvalues, an eigenmatrix, and positive and negative eigenvalue parts of the system matrix.These quantities define the characteristic decomposition used by the derivative calculation.
- Acoustic GRP: The acoustic GRP is the G1 scheme and is consistent with Toro’s ADER solver.
- Strong discontinuities: Large jumps make the acoustic GRP insufficient for strong discontinuities, so the associated generalized Riemann problem is solved analytically by the genuinely nonlinear GRP.The paper names this solver G∞GRP.
A.1.2. Nonlinear GRP solver.
The Euler-equation GRP resolution obtains the instantaneous state from a Riemann solver and the instantaneous time derivative by solving algebraic equations. The density derivative is then determined using the equation of state and quantified entropy variation.
- For the Euler equations, the instantaneous value u0 comes from the Riemann solver, while (∂u/∂t)0 comes from a pair of algebraic equations.
- The algebraic-equation coefficients ai, bi, and di are given explicitly in terms of the initial data.
- The instantaneous density derivative is obtained from the equation of state p = p(ρ, s) after quantifying the entropy variation.
A.2. Quasi-1-D GRP solver.
The quasi-1-D GRP solver treats multidimensional flow near an interface by fixing a transverse coordinate and incorporating tangential effects as a source term. It supports local wave propagation normal to the interface and includes multidimensional effects.
- The two-dimensional problem is handled with a quasi-1-D GRP solver.
- Polynomial states on neighboring cells define the left and right data, while tangential effects are treated as a source when constructing the normal interface flux.
- Fixing y reduces the interface calculation to a 1-D GRP at (0, ˜y), with g(u)y accounting for local wave propagation.
- The quasi-1-D solver has two versions, differing from the 1-D solver by inclusion of multidimensional effects.
- The quasi-1-D acoustic case applies when neighboring states are approximately equal but their gradients differ in norm.
A.2.1. Quasi-1-D acoustic case.
The solver distinguishes acoustic and genuinely nonlinear configurations using neighboring-state differences, then constructs the local GRP solution through eigen-decomposition or a two-step analytical procedure.
- Quasi-1-D acoustic case: In the acoustic case, the flux Jacobian A(u0) is decomposed as A(u0) = RΛR−1 with Λ containing the characteristic eigenvalues.
- Genuinely nonlinear case: When the neighboring-state difference is large, the solver treats the configuration as genuinely nonlinear and solves the quasi-1-D GRP analytically.
- Genuinely nonlinear case: The genuinely nonlinear construction begins by solving the local 1-D planar Riemann problem to obtain the local Riemann solution u0.
- Genuinely nonlinear case: The construction then solves the quasi-1-D GRP to determine the interface time derivative, with further details referenced externally.
- The paper acknowledges contributors for improving the English presentation and lists institutional and grant support.