Source-linked AI summary

A pressure-based semi-implicit space-time discontinuous Galerkin method on staggered unstructured meshes for the solution of the compressible Navier-Stokes equations at all Mach numbers

Maurizio Tavelli, Michael Dumbser

arXiv:1612.08866v1math.NAcs.NA

TL;DR

Explicit time stepping is severely restricted for low-Mach flows, and shocks can make unlimited high-order DG solutions oscillatory. The paper develops a pressure-based semi-implicit staggered DG method and pairs artificial viscosity with a posteriori MOOD limiting; reported tests reproduce shock positions without significant oscillations and match a Taylor–Green reference solution closely.

  • Problem

    Explicit schemes have severe low-Mach time-step restrictions because their CFL condition includes sound speed, while unlimited high-order schemes can oscillate near shocks and produce unphysical pressure or density.

    Source: §1 · p. 2
  • Method

    The method places pressure on the primal grid and density and velocity on the staggered dual grid, solves a pressure-only linear system, and combines artificial viscosity with a posteriori MOOD detection.

    Source: §1 · pp. 2–3, §3.3 · p. 13, §3.5 · p. 14
  • Results

    The Sod test reproduces the shock position without significant oscillations, while the Taylor–Green results fit the reference solution well, with slightly higher dissipation at Re = 800.

    Source: §4.1.10 · p. 26, §4.2.3 · p. 31
  • Takeaways & Limitations

    For piecewise constant time polynomials, the pressure system experimentally retains the good properties observed in the incompressible case, and the formulation permits implicit viscous-stress and heat-flux treatment.

    Source: §5 · p. 32

    For large viscosities, especially when artificial viscosity is used, an explicit viscous-term treatment can impose too restrictive a time-step condition, making implicit treatment necessary.

    Source: §3.4 · p. 14

Abstract

from arXiv · show

We propose a new arbitrary high order accurate semi-implicit space-time discontinuous Galerkin (DG) method for the solution of the two and three dimensional compressible Euler and Navier-Stokes equations on staggered unstructured curved meshes. The method is pressure-based and semi-implicit and is able to deal with all Mach number flows. In our scheme, the discrete pressure is defined on the primal grid, while the discrete velocity field and the density are defined on a face-based staggered dual grid. All convective terms are discretized explicitly, while the pressure terms appearing in the momentum and energy equation are discretized implicitly. Substitution of the momentum equation into the energy equation yields a linear system for the scalar pressure as the only unknown. The enthalpy and the kinetic energy are taken explicitly and are then updated using a simple Picard procedure. Thanks to the use of a staggered grid, the final pressure system is a very sparse block five-point system for three dimensional problems. The viscous terms and the heat flux are also discretized making use of the staggered grid by defining the viscous stress tensor and the heat flux vector on the dual grid, which corresponds to the use of a lifting operator on the dual grid. The time step of our new numerical method is limited by a CFL condition based only on the fluid velocity and not on the sound speed. This makes the method particularly interesting for low Mach number flows. Finally, a very simple combination of artificial viscosity and the a posteriori MOOD technique allows to deal with shock waves and thus permits also to simulate high Mach number flows. We show computational results for a large set of two and three-dimensional benchmark problems, including both low and high Mach number flows and using polynomial approximation degrees up to p=4.

1. Introduction

Explicit schemes face severe sound-speed time-step restrictions for low-Mach flows, while high-order methods need limiting near shocks. This paper extends semi-implicit pressure-based methods to a high-order staggered DG scheme and combines artificial viscosity with a posteriori MOOD limiting.

1. Introduction

  • Explicit schemes face severe low-Mach time-step restrictions because their CFL condition includes sound speed, motivating semi-implicit pressure-based approaches.The compressible Navier–Stokes equations cover both low-Mach and high-Mach flow regimes.Source: §1 · p. 2
  • The proposed scheme extends prior semi-implicit finite-volume ideas to high-order staggered DG for compressible Euler and Navier–Stokes equations, with a velocity-based CFL restriction and local/global conservation.Pressure is defined on the primal grid and density and velocity on the staggered dual grid; substituting momentum into energy yields a pressure-only linear system.Source: §1 · p. 2
  • Near shocks, unlimited high-order schemes can oscillate and produce unphysical negative pressure or density, so limiting is needed.The paper combines artificial viscosity with a posteriori MOOD detection for staggered-grid semi-implicit methods.Source: §1 · pp. 2–3

2. Staggered space-time DG scheme for the compressible Navier-Stokes equations

The scheme formulates compressible flow on primal and staggered dual meshes, with governing equations closed by an equation of state. Spatial and temporal polynomial bases are constructed on the corresponding space-time elements.

2.1. Governing equations

  • The compressible equations express mass, momentum, and energy conservation; the energy flux uses enthalpy and kinetic energy, and the EOS is constrained to be linear in pressure.The framework allows a general e(p, ρ) relation with this constraint, including ideal-gas and van der Waals equations of state.Source: §2.1 · pp. 3–4

2.2. Staggered unstructured grid

  • The 2D primal mesh uses triangles with edge-based quadrilateral dual cells, while the 3D primal tetrahedra have face-based five-point hexahedral dual elements.In 3D, each dual element is associated with the shared face of two tetrahedra.Source: §2.2 · p. 4

2.3. Space-time basis functions

  • Spatial bases use nodal polynomials on 2D triangles and quadrilaterals and on 3D tetrahedra, with a modal physical-space basis on dual hexahedra; time bases use Lagrange polynomials.The space-time bases are formed by tensor products of spatial and temporal bases.Source: §2.3 · pp. 5–6

3. Semi-implicit space-time DG scheme

This section presents a semi-implicit space-time DG scheme on staggered primal/dual grids, combining implicit pressure treatment and explicit convection with Picard iteration for a linear pressure solve. It also describes its velocity-based CFL restriction and troubled-cell limiter, which uses artificial viscosity and a detection-cost measure to address low- and high-Mach flow requirements.

3. Semi-implicit space-time DG scheme

  • Pressure and total energy are represented on the primal grid, while velocity, density, enthalpy, and momentum are represented on the staggered dual grid.The discretization represents the solution with piecewise space-time polynomials in primal and dual control volumes.Source: §3 · p. 6

3.1. Auxiliary variables

  • High-order averaging operators transfer quantities between the staggered dual and primal grids, while derived dual-grid quantities use pointwise evaluation in 2D and L2 projection in 3D.The averaging operators map in both directions between the grids.Source: §3.1 · pp. 7–8

3.2. Convective terms

  • The convective terms are discretized explicitly using a space-time weak formulation, with an upwind flux in time and a Rusanov-type spatial flux.Numerical Riemann fluxes resolve jumps at element interfaces.Source: §3.2 · p. 8

3.3. Derivation of the semi-implicit space-time DG scheme

  • Each Picard sweep updates density, evaluates convective and viscous terms, solves a linear system with pressure as its only unknown, then updates momentum and energy; the flux-form scheme conserves mass, momentum, and total energy.Pressure is treated implicitly, while convective terms, enthalpy, and kinetic energy are taken from the previous Picard iteration; the iteration linearizes the otherwise nonlinear pressure system.Source: §3.3 · pp. 13–14

3.4. Time step restriction

  • The CFL restriction depends on local flow velocity rather than sound speed, but explicit viscous stresses and heat flux impose a severe quadratic time-step restriction at large viscosity.The explicit viscous treatment remains computationally efficient for small viscosities without the artificial-viscosity limiter; large viscosities or use of that limiter require implicit viscous and heat-flux treatment.Source: §3.4 · p. 14

3.5. A posteriori DG limiter with artificial viscosity

  • For shock-containing high-Mach flows, the limiter first detects troubled cells a posteriori and then applies artificial viscosity only on those cells.The method combines artificial viscosity with MOOD detection because sub-cell limiting is not readily generalized to staggered unstructured meshes.Source: §3.5 · p. 14

3.5.1. Detection criteria

  • A candidate solution is rejected if density or pressure is nonpositive or if any conserved variable violates a relaxed discrete maximum principle over the element and its Voronoi neighbors.The admissibility checks are applied a posteriori at each time step.Source: §3.5.1 · p. 15
  • To reduce the cost of evaluating elementwise solution extrema, the detector checks only the known degrees of freedom stored in each element.The passage also specifies detector parameters δ0 = 2 · 10−3 and ϵ = 10−3.Source: §3.5.1 · p. 15

3.5.2. Limiter

  • When troubled cells are detected, the candidate solution is discarded and added viscosity is used there; otherwise, the candidate is accepted unchanged.The added viscosity is chosen to make the troubled-cell mesh Reynolds number satisfy the limiter criterion.Source: §3.5.2 · p. 15
  • The artificial heat-conduction coefficient is set so that the effective Prandtl number is 1, and the added viscosity is extended to all dual elements.The passage specifies Pri = µiγcv/λi = 1.Source: §3.5.2 · p. 15
  • Because the added viscosity can be high when the maximum wave speed includes sound speed, implicit viscous and heat-flux treatment becomes advantageous.After computing artificial viscosity across the mesh, the time step is recomputed with it active in troubled cells.Source: §3.5.2 · p. 15

4.1. Two dimensional numerical tests

The 2D tests assess smooth-flow accuracy on an isentropic vortex and acoustic wave. The reported results show polynomial-degree-dependent convergence for the vortex and good agreement with a reference solution for the acoustic wave.

4.1. Two dimensional numerical tests

  • The two-dimensional compressible Navier–Stokes tests compare numerical solutions with exact or numerical references in the literature; the isentropic-vortex test uses successively refined grids.The acoustic-wave test uses a nonorthogonal unstructured grid at very low effective Mach number, with the time step set by a flow-velocity-based CFL condition rather than a sound-speed-based one.Source: §4.1 · p. 15, §4.1.1 · p. 16, §4.1.2 · p. 16

4.1.1. Isentropic vortex

  • For the isentropic vortex, odd polynomial degrees yield optimal convergence, while even degrees yield suboptimal convergence.For this smooth unsteady problem, high-order polynomials in time matter because time error significantly affects total error.Source: §4.1.1 · p. 16, Table 1 · p. 16

4.1.2. Smooth acoustic wave propagation in 2D

  • The acoustic wave’s position is well reproduced, with good agreement against the reference solution along y = 0.The numerical solution also remains angularly symmetric despite the fully unstructured, nonsymmetric grid.Source: §4.1.2 · p. 17, Figure 4 · p. 17

4.1.3. Heat conduction

The heat-conduction test agrees well with its reference solution, and the viscous-shock test agrees very well with the exact solution, including derivative-based quantities.

4.1.3. Heat conduction

  • The heat-conduction test shows good agreement with the reference for temperature, heat flux, and density.The comparison uses a cut along y = 0 at t = 1.Source: §4.1.3 · p. 18, Figure 5 · p. 18

4.1.4. Viscous shock

  • The viscous-shock results agree very well with the exact solution for density, velocity, pressure, viscous stress, and heat flux.The benchmark includes viscosity, heat conduction, nonlinear convection, and pressure forces.Source: §4.1.4 · p. 19, Figure 7 · p. 19

4.1.5. Compressible mixing layer

The proposed scheme reproduces the compressible mixing-layer vorticity pattern in good qualitative agreement with reference solutions.

4.1.5. Compressible mixing layer

  • The final-time vorticity pattern shows good qualitative agreement with three reference solutions.The comparison includes a compact finite-difference scheme, a sixth-order PNPM scheme, and an ADER-WENO finite-volume scheme for the GPR model.Source: §4.1.5 · p. 20

4.1.6. Two dimensional lid-driven cavity flow

For the lid-driven cavity, the computed velocity profiles match incompressible reference results at Reynolds numbers 100 and 1000, while streamlines show main and secondary corner vortices.

4.1.6. Two dimensional lid-driven cavity flow

  • The velocity profiles fit Ghia’s incompressible reference solution very well at Re = 100 and Re = 1000.Figures 10 and 11 show the velocity cuts and streamline plots.Source: §4.1.6 · p. 22, Figure 10 · p. 22, Figure 11 · p. 22
  • Streamlines show the main recirculation zone and secondary corner vortices.The cavity uses a moving upper boundary and no-slip conditions on the remaining boundaries.Source: §4.1.6 · pp. 20, 22

4.1.7. Low Reynolds number flow past a circular cylinder

The cylinder calculation captures a von Karman street and a Strouhal number matching a reported value; the two viscous Mach-reflection tests show agreement with experimental and numerical results, including wedge boundary layers.

4.1.7. Low Reynolds number flow past a circular cylinder

  • The computed Strouhal number is St = 0.1834, matching Müller’s value and agreeing with the 0.182 reported in.The vortex street forms behind a circular cylinder at M = 0.2 and Re = 150.Source: §4.1.7 · pp. 22–23

4.1.8. Viscous single Mach reflection problem

  • The viscous single Mach reflection at Ms = 1.7 agrees with experimental results and other numerical experiments.A zoom shows the boundary layer generated by the viscous no-slip wall.Source: §4.1.8 · p. 23, Figure 14 · p. 24

4.1.9. Viscous double Mach reflection problem

  • For the viscous double Mach reflection at Ms = 10.0, the density contours resemble those in, and the upper boundary shows a boundary layer.Vorticity, Mach number, and pressure contours are also reported.Source: §4.1.9 · p. 24, Figure 15 · p. 25, Figure 16 · p. 25

4.1.10. Riemann problems

The Riemann tests examine whether the limiter yields physically admissible, non-oscillatory solutions across shock and rarefaction problems. The results show both its necessity in difficult cases and localized errors or oscillations in some tests.

Riemann problems

  • RP1 reproduces the shock position without significant oscillations, with limiting confined to the shock wave.The initial discontinuity is smoothed before projection onto the polynomial space, and the plotted profile is compared with the exact solution.Source: §4.1.10 · p. 26, Figure 17 · p. 26
  • The RP2 solution matches the exact profile but retains small oscillations; without the limiter, pressure becomes negative after only a few time steps.The limiter acts on only a small number of elements, and the authors use MOOD detection with artificial viscosity set to a mesh Reynolds number of unity without tuning its coefficients.Source: §4.1.10 · p. 27, Figure 18 · p. 27
  • In RP3, limiting leaves an unphysical density at the origin, whereas disabling it gives better agreement there because no artificial viscosity is added.No cells are limited at the final time, although steep initial gradients activate the limiter early in the simulation.Source: §4.1.10 · p. 27, Figure 20 · p. 28
  • RP4 shows small oscillations and an unexpected density behavior at the origin, consistent with a local wall-cooling error also observed in the second-order scheme.The solution otherwise has a relatively good shape compared with the exact reference.Source: §4.1.10 · pp. 27–28, Figure 19 · p. 28
  • RP5 agrees well with the exact solution, has small oscillations between the shock and rarefaction, and shows no sonic glitch inside the rarefaction fan.The initial states and solution times for the Riemann tests are listed in Table 2.Source: Table 2 · p. 26, §4.1.10 · p. 27

4.1.11. Circular explosion problem

The two-dimensional explosion test compares the staggered DG solution with a radial reference and examines where the limiter acts. The computed flow agrees well with the reference, with limiting concentrated at the discontinuity.

Circular explosion problem

  • The test covers [-1, 1]^2 with 5,616 triangles, uses (p, pγ) = (3, 0) and transmissive boundaries everywhere, and runs the simulation to t_end = 0.25.The radial reference reduces the 2D Euler equations in polar coordinates to a 1D radial system with geometric source terms, solved using a second-order TVD scheme with Osher-type flux.Source: §4.1.11 · p. 28
  • The explosion solution shows good agreement with the radial reference, and the limiter acts only at the discontinuity.The comparison uses a one-dimensional cut and scatter plots of the main flow variables against a reference obtained from the radially reduced equations.Source: §4.1.11 · p. 29, Figure 23 · p. 30

4.2. Three dimensional numerical tests

The three-dimensional tests cover smooth acoustic propagation, a lid-driven cavity, a Taylor–Green vortex, and flow around a sphere. Across these cases, the reported solutions agree well with reference data or experiments, with specific dissipation and wake-frequency comparisons.

4.2. Three dimensional numerical tests

  • The 3D test suite spans low to moderate Mach flows, including smooth acoustic propagation, cavity flow, vortex evolution, and a sphere wake.The sphere case uses curved elements to describe the geometry.Source: §4.2 · p. 30, §4.2.4 · p. 31

4.2.1. Smooth acoustic wave propagation in 2D

  • For the smooth acoustic wave, the solution is symmetric on an unstructured mesh and has the correct propagation speed, with some dissipation near velocity and pressure peaks.The solution also agrees reasonably well with a reduced one-dimensional reference.Source: §4.2.1 · pp. 30–31

4.2.2. 3D Lid-driven cavity

  • The 3D cavity results at Re=400 and Re=1000 agree very well with literature reference data and reproduce secondary Taylor–Görtler-like vortex structures.The calculations use 1,380 tetrahedra at Re=400 and 11,040 at Re=1000.Source: §4.2.2 · p. 31

4.2.3. 3D Taylor-Green vortex

  • At Re=400, the Taylor–Green vortex energy decay agrees well with DNS, including at t>7; at Re=800, the computed dissipation rate is slightly higher.The test develops small-scale structures from a simple initial state and compares kinetic-energy dissipation with Brachet et al.'s DNS.Source: §4.2.3 · p. 31

4.2.4. Flow around a sphere

  • For flow around a sphere at M∞=0.5, the Strouhal number is 0.139, close to the reported values 0.137 in and by Johnson and Patel.The velocity and density patterns are also reported as similar to experimental observations.Source: §4.2.4 · pp. 31–32

5. Conclusions

The paper concludes that its semi-implicit staggered DG method applies to compressible flows across a broad Mach-number range. The authors identify complex-geometry, geophysical, meteorological, natural-convection, and multiphase or magnetohydrodynamic applications as future directions.

Conclusions

  • The method is tested on two- and three-dimensional problems from very low to very high Mach numbers using the same discretization, supporting its characterization as an all-Mach solver.The formulation is in principle arbitrary high order in space and time, and an outer Picard iteration makes the nonlinear terms explicit.Source: §5 · p. 32
  • Implicit treatment of viscous stress and heat flux is particularly useful with artificial viscosity, whose explicit treatment would impose very small time steps.The conclusion also reports the pressure system's good properties for piecewise constant polynomials in time.Source: §5 · p. 32
  • Future applications include geophysical and meteorological flows, natural convection, complex-geometry industrial and biological flows, and extensions to MHD and Baer–Nunziato multiphase equations.The authors state that high-order semi-implicit discretization would allow larger time steps than purely explicit DG for natural convection.Source: §5 · p. 32
Loading 1612.08866v1…