Source-linked AI summary
An Efficient and Accurate Two-Stage Fourth-order Gas-kinetic Scheme for the Navier-Stokes Equations
Liang Pan, Kun Xu, Qibing Li, Jiequan Li
TL;DR
Higher-order CFD schemes can offer greater accuracy but are often less robust and more complicated. This paper constructs a two-stage fourth-order gas-kinetic scheme from a second-order flux function, reporting robustness comparable to the second-order method and an extension path to fifth-order time accuracy.
Problem
Higher-order numerical methods can provide more accurate solutions than second-order schemes, but are less robust and more complicated.
Source: §1 · p. 2Method
Results
The fourth-order scheme is reported as robust as the second-order method, with favorable shock capturing and a CFL number on the order of 0.5.
Source: §1 · p. 2, §5 · p. 25Takeaways & Limitations
The two-stage framework can be extended to fifth-order temporal accuracy by including first- and second-order time derivatives of the flux function.
Source: §1 · p. 3In the double Mach reflection problem, slip-line oscillations are not as large as with many other higher-order schemes because physical dissipation stabilizes the shear instability.
Source: §4.5 · p. 18
Abstract
from arXiv · showhide
For computational fluid dynamics (CFD), the generalized Riemann problem (GRP) solver and the gas-kinetic kinetic scheme (GKS) provide a time-accurate flux function starting from a discontinuous piecewise linear flow distributions around each cell interface. With the use of time derivative of the flux function, a two-stage Lax-Wendroff-type (L-W for short) time stepping method has been recently proposed in the design of a fourth-order time accurate method [18]. In this paper, based on the same time-stepping method and the second-order GKS flux function [34], a fourth-order gas-kinetic scheme is constructed for the Euler and Navier-Stokes equations. In comparison with the formal one-stage time-stepping third-order gas-kinetic solver [21], the current fourth-order method not only reduces the complexity of the flux function, but also improves the accuracy of the scheme, even though the third- and fourth-order schemes have similar computation cost. Most importantly, the robustness of the fourth-order GKS is as good as the second-order one. Perfect numerical solutions can be obtained from the high Reynolds number boundary layer solutions to the hypersonic viscous heat conducting flow computations. Many numerical tests, including many difficult ones for the Navier-Stokes solvers, have been used to validate the current fourth-order method. Following the two-stage time-stepping framework, the one-stage third-order GKS can be easily extended to a fifth-order method with the usage of both first-order and second-order time derivatives of the flux function. The use of time-accurate flux function may have great impact on the development of higher-order CFD methods.
1. Introduction
The paper develops a two-stage, fourth-order gas-kinetic scheme for the Euler and Navier–Stokes equations using a time-accurate flux and its time derivative. It targets higher-order accuracy while retaining the robustness of the second-order scheme.
1. Introduction
- The two-stage Lax–Wendroff-type discretization uses the gas-kinetic flux function and its time derivative to construct a fourth-order scheme for Euler and Navier–Stokes equations.The gas-kinetic scheme supplies a time-accurate flux through its kinetic-to-hydrodynamic relaxation process.Source: §1 · p. 2
- The fourth-order scheme is reported to be as robust as the second-order method, with favorable shock capturing and a CFL number on the order of 0.5.The introduction describes its performance from subsonic to hypersonic viscous heat-conducting flows.Source: §1 · p. 2
- The two-stage framework can be extended to fifth-order temporal accuracy by including first- and second-order time derivatives of the flux function.The paper says the extension process can continue to even higher orders.Source: §1 · p. 3
2. Fourth-order temporal discretization
The temporal discretization applies a two-stage Lax–Wendroff-type update to a finite-volume formulation. Its fourth-order accuracy follows from using time derivatives of the spatial flux operator at the initial and intermediate states.
2. Fourth-order temporal discretization
- The initial condition is specified at t_n, and the time derivatives of the spatial flux operator are obtained with the Cauchy–Kovalevskaya method.The operator L denotes the spatial derivative of flux.Source: §2 · p. 3
- The intermediate-stage state is used to obtain corresponding time derivatives before the conservative state is updated.The supplied passages describe the intermediate state and subsequent update as separate steps.Source: §2 · p. 3
- The two-stage update using Eq.(3) and Eq.(4) is stated to produce a fourth-order time-accurate solution at t_n + Δt.The accuracy follows from a time-accurate solution for ∂L/∂t.Source: §2 · pp. 3–4
- The method applies the two-stage Lax–Wendroff framework to a finite-volume equation for conservative variables, with interface fluxes updating cell averages.The formulation is presented for one dimension, with a similar finite-volume form in two and three dimensions.Source: §2 · p. 4
3. A fourth-order gas-kinetic scheme
The scheme constructs a fourth-order gas-kinetic method for Euler and Navier–Stokes solutions by combining a second-order flux function with two-stage time discretization. Its flux solver, spatial reconstruction, and flux-derivative evaluation supply the components of that construction.
3. A fourth-order gas-kinetic scheme
- The fourth-order scheme for Euler and Navier–Stokes solutions is built from a second-order gas-kinetic flux function within a two-stage temporal framework.The paper relates the GRP solver and GKS through a generalized Lax–Wendroff-type procedure for discontinuous cases.Source: §3 · p. 4
3.1. Second-order gas-kinetic flux solver
- The flux solver evolves the gas distribution at a cell interface from initial data, with Euler and Navier–Stokes behavior obtained through zeroth- and first-order truncations in collision time.The Navier–Stokes formulation uses the first-order truncation; the Euler formulation uses the zeroth-order truncation.Source: §3.1 · p. 5
3.2. Spatial reconstruction
- Fifth-order WENO reconstruction supplies high-order macroscopic flow data, while fifth-order Gauss quadrature provides tangential reconstruction values in two-dimensional computations.The one-dimensional reconstruction also determines interface states and their spatial derivatives.Source: §3.2 · pp. 6–7
3.3. Two-stage gas-kinetic scheme
- Because gas evolution makes the flux a complicated function of time, the method approximates it linearly over a time interval to obtain flux time derivatives for the two-stage algorithm.The initial flux and its time derivative are computed from flux values at the half-step and full-step times.Source: §3.3 · p. 7
4. Numerical tests
The numerical tests assess the fourth-order GKS on one- and two-dimensional accuracy problems, using fifth-order WENO reconstruction. They report fifth-order accuracy in smooth advection tests, with temporal error limiting accuracy when the CFL number is reduced.
4. Numerical tests
- For viscous-flow tests, the collision time reduces to τ = µ/p in smooth regions; an added term enlarges kinetic-scale physics near discontinuities.The added dissipation supports construction of a numerical shock structure through particle free transport and inadequate particle collision.Source: §4 · pp. 8–9
- For strong discontinuities, WENO reconstruction can use characteristic variables, obtained by projecting conservative variables with the right eigenvectors of the flux Jacobian.The reconstruction is applied to the characteristic variables, then the conservative variables are recovered by inverse projection.Source: §4 · p. 9
4.1. Accuracy tests
- The density-perturbation advection test uses a one-dimensional sinusoidal initial density profile and periodic boundaries.The exact solution propagates with velocity (1, 1) in the two-dimensional accuracy test.Source: §4.1 · p. 9
- Fifth-order accuracy is retained through mesh size N = 640 for density-perturbation advection, while the second-order GKS achieves only second-order accuracy with the same spatial reconstruction and CFL = 0.4.For the second-order GKS, a smaller CFL = 0.1 initially restores fifth-order accuracy, but refinement makes temporal error dominant and reduces accuracy to second order.Source: §4.1 · p. 9
- The two-dimensional isentropic-vortex test achieves the expected accuracy, with L1 and L∞ errors and orders reported at t = 10.The results are presented for uniform N × N meshes in Table 3.Source: §4.1 · p. 10, Table 3 · p. 10
4.2. One dimensional Riemann problems
The one-dimensional Riemann tests cover Sod, the Woodward-Colella blast wave, and Shu-Osher shock-acoustic interaction. The reported profiles agree well with exact solutions for Sod and blast-wave cases, and the scheme is basically stable under conventional CFL conditions.
4.2. One dimensional Riemann problems
- The Sod test uses a one-dimensional domain [0, 1] with 100 uniform mesh points and non-reflecting boundaries.Its initial states are (1, 0, 1) left of x = 0.5 and (0.125, 0, 0.1) to the right.Source: §4.2 · p. 10
- The Woodward-Colella blast-wave test uses 200 and 400 uniform mesh points on [0, 100], with reflected boundaries.The initial pressure is 1000 for x < 10, 0.01 for 10 ≤ x < 90, and 100 for x ≥ 90.Source: §4.2 · p. 12
- The Sod and blast-wave profiles agree well with exact solutions and resolve wave profiles, particularly local extreme values.Sod is shown at t = 0.2 and blast wave at t = 3.8; the figure includes 100 Sod cells and 200 and 400 blast-wave cells.Source: Figure 1 · p. 11, §4.2 · p. 12
- The Shu-Osher shock-acoustic interaction is shown as a computed density profile at t = 1.8 with 400 mesh points.The initial flow combines a left state with a sinusoidal density field to its right.Source: §4.2 · p. 12
- The scheme is basically stable for CFL numbers from 0.2 to 0.7, with wave profiles well resolved around CFL = 0.5.The following numerical tests use CFL = 0.4 unless otherwise stated.Source: §4.2 · p. 12
4.3. Two-dimensional Riemann problems
The two Riemann tests show the scheme capturing intricate shock and vortex-sheet structures. The authors report a cost trade-off favoring fourth order when matching accuracy, alongside robustness comparable to the second-order scheme.
4.3. Two-dimensional Riemann problems
- The rarefaction-wave and vortex-sheet interaction produces a roll-up that the scheme captures.The density distribution is shown for 600 × 600 and 1000 × 1000 mesh points.Source: §4.3 · p. 13, Figure 5 · p. 15
- Four initial shocks develop a complex pattern: S−23 bifurcates into a reflected shock, Mach stem, and slip line, while small-scale structures are captured.The density distribution is shown for two mesh resolutions.Source: §4.3 · p. 13, Figure 4 · p. 14
- About 3× the second-order scheme’s cost is reported for fourth order, while matching its accuracy with second order requires refinement that raises two-dimensional cost 8×.The authors conclude that the fourth-order scheme is more efficient than the second-order method.Source: §4.3 · pp. 14–15, Table 4 · p. 14
- The fourth-order scheme is reported as accurate and as robust as the second-order scheme; its robustness is attributed to a first-order relaxation model with a fully hyperbolic system and local source term.Source: §4.3 · p. 16
4.4. Shock vortex interaction
The section tests a small vortex interacting with a stationary Mach 1.1 shock. The reported computations capture shock reflection and show accurate centerline density profiles at three mesh resolutions.
4.4. Shock vortex interaction
- A small vortex perturbs the flow around a stationary Mach 1.1 shock positioned at x = 0.5, normal to the x-axis.The computational domain is [0, 2] × [0, 1].Source: §4.4 · p. 16
- By t = 0.8, a shock-bifurcation branch reaches and reflects from the top boundary, and the reflection is captured; centerline density profiles at three mesh spacings demonstrate accuracy.Source: §4.4 · p. 16
- The vortex perturbation uses κ = 0.3, μ = 0.204, and r_c = 0.05, with κ setting strength, μ decay rate, and r_c the radius of maximum strength.Source: §4.4 · p. 16
4.5. Double Mach reflection problem
The tests show the fourth-order GKS resolving the triple Mach stem, capturing hypersonic cylinder flows, and matching boundary-layer reference solutions. The reported comparisons also show agreement with experiment and accurate profiles on a coarse mesh.
4.5. Double Mach reflection problem
- The density contours resolve the flow structure beneath the triple Mach stem, while physical dissipation stabilizes shear instability and limits slip-line oscillations.The contours are shown on 720 × 240 and 1440 × 480 meshes.Source: §4.5 · pp. 17–18, Figure 8 · p. 18, Figure 9 · p. 18
4.6. Hypersonic flow past a cylinder
- For hypersonic inviscid flow past a cylinder, Mach-number distributions at Ma = 5, 10, and 20 show that the scheme captures strong shocks very well without carbuncle phenomenon.Source: §4.6 · p. 19, Figure 10 · p. 19
- For viscous, heat-conducting cylinder flow at Ma∞ = 8.03 and Re = 1.835 × 10^5, computed surface pressure and heat flux agree well with experimental data.A non-uniform 60 × 160 mesh uses near-wall cell width 1/2000 to resolve the boundary layer.Source: §4.6 · p. 19, Figure 12 · p. 20
4.7. Laminar boundary layer
- Flat-plate boundary-layer velocity profiles match the exact Blasius solution, with the upstream profile captured using four grid points within the layer.The test uses a 260 × 90 mesh.Source: §4.7 · p. 19, Figure 11 · p. 20
4.8. Lid-driven cavity flow
The scheme matches benchmark data for lid-driven cavity flow at two Reynolds numbers and gives mesh-convergent solutions for reflecting shock-boundary layer interactions. At Re = 1000, mesh convergence is basically obtained with mesh size decreasing to 1/1500.
4.8. Lid-driven cavity flow
- At Re = 1000 and 3200, cavity simulations match benchmark data; the Re = 3200 case demonstrates higher-order accuracy with 65 mesh points per direction.The results include streamlines and centerline velocity profiles.Source: §4.8 · pp. 20–21, Figure 14 · p. 21, Figure 15 · p. 22
4.9. Viscous shock tube problems
- For the reflecting shock-boundary layer interaction at Re = 200, density solutions are mesh-convergent and the primary-vortex height agrees with reference data.The paper shows density distributions and lower-wall profiles at different mesh resolutions.Source: Table 5 · p. 22, Figure 16 · p. 23, Figure 17 · p. 23
- At Re = 1000, the flow structure is more complicated, and mesh convergence is basically obtained as mesh size decreases to 1/1500; results agree well with reference data.Density distributions and lower-wall profiles are shown for mesh sizes 1/1500 and 1/2000.Source: §4.9 · p. 24, Figure 18 · p. 24, Figure 19 · p. 24
5. Conclusion
The conclusion presents a robust fourth-order gas-kinetic scheme for inviscid and viscous flows, with fifth-order spatial and fourth-order temporal accuracy in the reported formulation. It identifies higher-order time derivatives of the flux as a route to still higher temporal accuracy.
5. Conclusion
- The fourth-order GKS uses two stages for fourth-order time accuracy, compared with four stages in conventional methods based on first-order Riemann solvers.The scheme can use a CFL number on the order of 0.5.Source: §5 · p. 25
- With fifth-order WENO reconstruction, the scheme has fifth-order spatial and fourth-order temporal accuracy, favorable shock capturing, and robustness comparable to the second-order scheme.The conclusion reports expected accuracy for smooth flows and favorable shock capturing for discontinuous solutions.Source: §5 · p. 25
- Including the second-order time derivative of the flux function is presented as a way to develop even higher-order schemes, including a fifth-order scheme in the Appendix.The conclusion also identifies first-order Riemann solvers as the bottleneck hindering higher-order accurate and robust compressible-flow schemes.Source: §5 · p. 25
Appendix: Extension to higher order
The appendix extends the two-stage time-stepping framework to fifth-order temporal accuracy. It uses a quadratic approximation of the time-dependent flux and states that the update equations provide the claimed accuracy.
Appendix: Extension to higher order
- Third-order GRP and GKS flux solvers, which include first- and second-order time derivatives, can be used with two-stage discretization to develop fifth-order temporal accuracy.The appendix frames the time-dependent flux as the key ingredient for the two-stage fourth-order method and its extension.Source: Appendix: Extension to higher order · p. 25
- The appendix states that Eqs. (17) and (18) give a fifth-order temporally accurate approximation at t = tn + ∆t with the specified coefficients.The derivation uses a Taylor expansion of the operator and solves for coefficients by comparing terms.Source: Appendix: Extension to higher order · p. 26
- For fifth-order temporal accuracy, the time-dependent flux is approximated by a quadratic function.The appendix introduces an intermediate state and an update involving the second time derivative of the operator at that state.Source: Appendix: Extension to higher order · pp. 26–27