Source-linked AI summary
A Posteriori Subcell Limiting of the Discontinuous Galerkin Finite Element Method for Hyperbolic Conservation Laws
Michael Dumbser, Olindo Zanotti, Raphael Loubere, Steven Diot
TL;DR
DG schemes are high-order and stable but require limiting near shocks without losing their subcell resolution. This paper extends MOOD with a posteriori detection and ADER-WENO recomputation on subcells, demonstrating high-order multidimensional simulations and large-scale parallel execution.
Problem
DG methods can generate oscillations and fail near strong shocks or steep gradients, despite nonlinear L2 stability, while limiting should preserve high-order accuracy and subcell resolution.
Method
The method computes an unlimited one-step ADER-DG candidate, tests positivity and a relaxed discrete maximum principle, and recomputes troubled cells on a 2N + 1-subcell grid with ADER-WENO.
Results
The scheme is demonstrated in two and three space dimensions with DG orders up to tenth order and a 3D simulation using 10 billion space-time degrees of freedom per time step.
Takeaways & Limitations
The a posteriori subcell strategy provides a high-order DG limiter that retains the method’s subcell resolution while locally recomputing invalid cells.
Takeaways & Limitations
The authors identify future work for unstructured meshes, broader PNPM schemes, and additional conservation-law systems, including magnetohydrodynamics.
Abstract
from arXiv · showhide
The purpose of this work is to propose a novel a posteriori finite volume subcell limiter technique for the Discontinuous Galerkin finite element method for nonlinear systems of hyperbolic conservation laws in multiple space dimensions that works well for arbitrary high order of accuracy in space and time and that does not destroy the natural subcell resolution properties of the DG method. High order time discretization is achieved via a one-step ADER approach that uses a local space-time discontinuous Galerkin predictor method. Our new limiting strategy is based on the so-called MOOD paradigm, which aposteriori verifies the validity of a discrete candidate solution against physical and numerical detection criteria. Within the DG scheme on the main grid, the discrete solution is represented by piecewise polynomials of degree N. For those troubled cells that need limiting, our new limiter approach recomputes the discrete solution by scattering the DG polynomials at the previous time step onto a set of N_s=2N+1 finite volume subcells per space dimension. A robust but accurate ADER-WENO finite volume scheme then updates the subcell averages of the conservative variables within the detected troubled cells. The choice of N_s=2N+1 subcells is optimal since it allows to match the maximum admissible time step of the finite volume scheme on the subgrid with the maximum admissible time step of the DG scheme on the main grid. We illustrate the performance of the new scheme via the simulation of numerous test cases in two and three space dimensions, using DG schemes of up to tenth order of accuracy in space and time (N=9). The method is also able to run on massively parallel large scale supercomputing infrastructure, which is shown via one 3D test problem that uses 10 billion space-time degrees of freedom per time step.
1. Introduction
DG schemes provide high-order, locally entropy-stable approximations for nonlinear hyperbolic conservation laws, but require nonlinear limiting near shocks. The paper extends the a posteriori MOOD paradigm to DG using a fine subgrid and ADER-WENO recomputation to preserve subcell resolution.
- DG schemes satisfy a local cell entropy inequality for any polynomial degree N, yielding nonlinear L2 stability at arbitrary high order.
- Explicit DG time discretizations impose a restrictive CFL condition whose maximum Courant number scales approximately as 1/(2N +1).
- MOOD generates an unlimited high-order finite-volume candidate, tests admissibility, and locally recomputes troubled cells with lower-order reconstructions until validity is restored.
- Classical DG limiting detects and corrects troubled cells, but differs from MOOD because it typically postprocesses the current candidate at one time level.
- The proposed limiter recomputes troubled DG cells on a finer subgrid with an accurate ADER-WENO finite-volume scheme and transfers data through subcell projection and reconstruction.
- The subgrid size is selected so its finite-volume time-step limit matches the DG limit, producing Ns = 2N + 1 subcells per space dimension.
2. The ADER Discontinuous Galerkin Method
The ADER-DG method advances piecewise polynomial solutions with a one-step predictor-corrector formulation. A local space-time DG predictor supplies high-order flux information, while the corrector couples neighboring cells through numerical fluxes under a restrictive order-dependent CFL condition.
- The method solves nonlinear hyperbolic conservation laws with pure PNPN discontinuous Galerkin schemes, providing high-order accuracy in space and time.
- Each main-grid cell stores the state as a piecewise polynomial of maximum degree N in a basis of Lagrange polynomials at Gauss-Legendre points.
- A local space-time Galerkin predictor evolves each cell independently while neglecting neighbor influence before the corrector restores intercell coupling.
- The PDE is mapped to a space-time reference element [0; 1]d+1, with time represented by t = tn + ∆t τ.
- The corrector integrates the governing equation over each space-time control volume and replaces boundary fluxes with a numerical Riemann-solver flux.
- The explicit DG CFL number decreases with polynomial degree, roughly following 1/(2N + 1), while high-degree polynomials provide multiple subcell features within one cell.
- For sufficiently smooth solutions, the unlimited ADER-DG scheme achieves effective order N + 1 but cannot prevent oscillations near steep gradients or shocks.
3. Nonlinear stability via a priori limiting
DG solutions can become strongly oscillatory near discontinuities despite nonlinear L2 stability, motivating nonlinear limiting. Existing a priori approaches add dissipation or modify polynomial representations, but must preserve accuracy and subcell features.
- Smooth initial data can evolve into discontinuous solutions, and strong shocks may cause DG oscillations that ultimately fail without nonlinear limiting.
- A troubled-cell indicator should localize limiting, while the correction should suppress oscillations near discontinuities without damaging high-order accuracy or subcell resolution.
- The main a priori DG limiter families are artificial-viscosity, slope or moment, and WENO or HWENO methods.
- Artificial-viscosity approaches stabilize shocks by adding dissipative terms, with models differing in how shock regions and viscosity magnitudes are determined.
- Slope /moment reduction: Slope and moment limiters constrain or remove high-order polynomial components, but must avoid degrading accuracy at smooth extrema.
- Some DG limiters dissipate subscale features, making simultaneous preservation of subcell accuracy and shock stability difficult.
- Existing subcell limiters use finite-volume subcells or switches to first-order subgrid methods, whereas the proposed approach uses a posteriori detection and higher-order subgrid finite volumes.
4. Nonlinear stability via a posteriori sub-cell limiting (SCL)
The method extends the MOOD a posteriori limiting paradigm to DG by recomputing troubled cells on a fine finite-volume subgrid while preserving DG subcell information and high-order accuracy.
- A posteriori detection: The MOOD procedure computes an unlimited candidate, tests physical admissibility and a relaxed discrete maximum principle, then flags invalid cells for recomputation.Physical criteria include positivity of relevant quantities such as density and pressure for Euler equations.
- Subcell representation: DG subcell averages are obtained by projecting degree-N cell polynomials onto (N_s)^d piecewise-constant subcells.The projection computes integral averages on equidistant Cartesian subcells.
- Subcell representation: For N_s ≥ N + 1, reconstruction recovers the original DG polynomial without loss of accuracy, including when the reconstruction system is overdetermined.The overdetermined system is solved with constrained least squares.
- Subcell representation: The subgrid representation is nominally equivalent in accuracy to a degree-N DG representation on the main cell.This equivalence preserves the information needed for subcell limiting.
- Subcell recomputation: Troubled cells are recomputed with a robust ADER-WENO finite-volume scheme on the subgrid and gathered back into degree-N DG polynomials.Unlike postprocessing limiters, the approach replaces the troubled-cell update with a more robust scheme.
- Subgrid size: N_s = 2N + 1 matches the finite-volume subgrid time step to the DG time step and is optimal for accuracy and computational efficiency.Coarser or finer subgrids are respectively more dissipative or more restrictive in time stepping.
5. Numerical results
The numerical study applies the MPI-parallel ADER-DG method with the subcell limiter to compressible Euler tests spanning smooth flows, shocks, multidimensional interactions, and a large-scale 3D problem.
- Implementation: The implementation uses an MPI-parallel 3D Cartesian-grid code for hyperbolic systems and supplements ADER-DG-P5 or ADER-DG-P9 with the subcell limiter.The subcell solver is a third-order ADER-WENO finite-volume method.
- Implementation: The visualization uses N_s = N + 1 for unlimited cells and N_s = 2N + 1 for troubled cells to inspect subscale solution structure.This supports checking physical validity within large cells and discontinuity quality at shocks.
- Shock resolution: The limiter is intended to deposit entropy and dissipation on subcell scales rather than spreading shocks across excessively large coarse-cell regions.The numerical section uses the tests to provide evidence for this property.
- Test cases: The test suite covers smooth-vortex convergence, oscillatory shock tubes, shock reflections, supersonic flow over a step, two-dimensional Riemann problems, and shock-vortex interaction.These tests target accuracy, small-scale feature capture, complex shock patterns, and multidimensional wave interactions.
- Test cases: The 3D explosion problem evaluates separate spherical waves and demonstrates operation in a massively parallel MPI framework.The reported configuration uses 8000 CPU cores and 10 billion space-time degrees of freedom per time step and conserved variable.
5.1. Sod and Lax shock tube
The Sod and Lax shock tubes test one-dimensional wave resolution on a very coarse two-dimensional mesh. ADER-DG-P9 with an a posteriori WENO3 subcell limiter agrees closely with exact solutions while resolving key discontinuities sharply.
- Problem setup: The Sod and Lax problems use a 20 × 5 structured mesh, γ = 1.4, and a discontinuity at x = 0.5.The computational domain is Ω = [0, 1] × [−0.5, 0.5], with Dirichlet boundaries in x and periodic boundaries in y.
- Results: The shock is resolved in one cell, while the Sod contact discontinuity is resolved within one cell in an almost S-type shape.The sharp contact resolution is attributed to the degree-N = 9 DG polynomials.
- Numerical method: ADER-DG-P9 with the a posteriori WENO3 subcell limiter is compared against exact density, velocity, and pressure profiles.Figure 3 samples each numerical solution along x and compares it with the exact solution.
- Results: Troubled cells are limited with the subcell method, whereas blue cells retain the unlimited ADER-DG-P9 update on the main grid.Figure 4 displays the density variable and marks limited cells red and unlimited cells blue.
5.3. Shu-Osher oscillatory shock tube
The Shu-Osher oscillatory shock tube tests whether the method captures small scales generated when a Mach 3 shock interacts with a sinusoidal density fluctuation. ADER-DG-P9 with the a posteriori WENO3 subcell limiter produces an excellent result on a very coarse mesh.
- Problem setup: The test uses a Mach 3 shock interacting with a sinusoidal density fluctuation of wavelength λ = 5 and amplitude ε = 0.2.The initial shock is at x = −4, and the final time is tfinal = 0.18.
- Numerical method: The computation uses a very coarse mesh of 40 cells in x and 5 cells in y with ADER-DG-P9 and an a posteriori WENO3 subcell limiter.An Osher-type numerical flux is used.
- Results: The result is reported as excellent for the 40-cell mesh against a reference solution from a third-order ADER-WENO scheme on a very fine mesh.The DG polynomials are represented by 10 sample points per cell in the one-dimensional view.
5.4. Double Mach reflection problem
The double Mach reflection tests assess the limiter on a Mach 10 shock interacting with a ramp, using ADER-DG schemes of degrees 1, 2, 5, and 9. The method detects shocks while preserving smooth vortex structures and high-order accuracy.
- Problem setup: The test uses a Mach 10 shock striking a 30° ramp on a 350 × 100 grid with h = 1/100.ADER-DG-PN schemes with N = 1, 2, 5, and 9 are supplemented by an a posteriori ADER-WENO3 subcell limiter.
- Detection behavior: All polynomial degrees detect the shock waves, while ADER-DG-P5 and ADER-DG-P9 require little limiting in the vortex-rich interaction zone.The higher-order schemes treat the vortex structures as smooth subscale features, whereas ADER-DG-P2 flags more troubled cells.
- Resolution: Increasing the polynomial degree captures more vortex structures and reduces shock smearing on the same coarse grid.ADER-DG-P9 is less dissipative than the lower-order schemes, whose dissipation destroys small-scale vortices.
- Accuracy: The subcell limiter provides valid troubled-cell detection and preserves the overall accuracy of the high-order DG scheme.The ADER-WENO3 update on the subgrid is reported to maintain the main-grid DG accuracy.
5.5. Forward facing step
The forward facing step test examines ADER-DG-P5 with a posteriori ADER-WENO3 limiting in a Mach 3 wind tunnel. The simulation resolves shock waves, shear-layer instabilities, and small-amplitude acoustic waves while limiting mainly near shocks.
- Problem setup: The forward facing step problem models a Mach 3 wind tunnel with reflective walls and a stepped computational domain.The initial gas has density ρ = γ, pressure p = 1, velocity (u,v) = (3,0), and γ = 1.4.
- Flow structures: The calculation resolves Kelvin-Helmholtz instabilities along the shear wave and acoustic waves generated by unsteady vortices.These features remain visible as they cross shock waves, although more dissipative schemes often smear the acoustic waves.
5.6. 2D Riemann problem
The two-dimensional Riemann tests evaluate ADER-DG-P5 with an a posteriori WENO3 subcell limiter across five configurations. Limiting remains concentrated near strong discontinuities, while smooth regions retain the unlimited scheme's precision and the main structures are captured.
- Problem setup: Five two-dimensional Riemann configurations are simulated on a 100 × 100 main grid using ADER-DG-P5 with an a posteriori WENO3 subcell limiter.The configurations correspond to entries listed in Table 5 and previously studied problem families.
- Visualization: The results show density isolines alongside meshes marking limited cells red and unlimited cells blue.Figures 9 and 10 present RP1–RP3 and RP4–RP5, respectively.
- Limiter behavior: The limiter is active mainly along strong discontinuities and inactive in smooth regions, retaining the optimal precision of unlimited ADER-DG-P5.Some parasitical cells and start-up error are reported for RP3.
- Comparison: RP3 is compared with a sixth-order ADER-WENO finite-volume scheme using adaptive mesh refinement at equivalent resolution.The comparison is illustrated in Figure 11, with the DG result on the left and the AMR finite-volume result on the right.
- Conclusion: The coupled method captures discontinuous waves without spurious oscillations and smooth flow without excessive numerical dissipation.This conclusion is stated for the classical test cases considered in the section.
5.7. Shock-vortex interaction
The shock-vortex interaction test probes ADER-DG-P5 on a flow containing both a stationary shock and a vortex. The reported visualizations compare density structures and limited-cell locations for multidimensional flow features.
- Problem setup: The test places a stationary normal shock at x = 0.5 and a vortex at (x_c,y_c) = (0.25,0.5) in the domain [0,2] × [0,1].The vortex–shock benchmark combines smooth features with discontinuous waves.
- Related visualizations: The section presents density isolines and limited-cell maps for two-dimensional Riemann problems RP1–RP5, with RP4 and RP5 shown in Figure 10.Figure 10 uses the same caption structure as Figure 9.
- Comparison: A comparison for RP3 places the ADER-DG-P5 subcell-limited result beside a sixth-order ADER-WENO6 calculation with adaptive mesh refinement.The finite-volume comparison is shown in Figure 11 at equivalent resolution.
- Parameters: The vortex parameters include γ = 1.4, a = 0.075, b = 0.175, M_S = 1.5, M_V = 0.7, p_0 = 1, and ρ_0 = 1.The downstream state is computed using Rankine–Hugoniot conditions.
- Numerical experiment: The simulation uses ADER-DG-P5 on a grid with characteristic length h = 1/100 and reports density at tfinal = 0.7.The accompanying visualization marks troubled cells red and unlimited cells separately.
5.8. 3D explosion problem
The three-dimensional explosion test evaluates ADER-DG-P9 with the a posteriori ADER-WENO subcell limiter on coarse and fine Cartesian grids. The results compare density profiles with a one-dimensional reference solution and show the method’s subcell resolution and large-scale computational setting.
- Problem setup: The explosion problem extends the Sod problem to three dimensions on Ω = [−1; 1]3, with a spherical initial discontinuity of radius R = 0.5.The ideal-gas equation of state uses γ = 1.4, and the spherically symmetric reference solution comes from a radial one-dimensional problem with geometric source terms.
- Numerical method: The simulation uses ADER-DG-P9 with the new ADER-WENO subcell limiter on 25^3 and 100^3 element grids.For N = 9, each element’s space-time predictor contains (N + 1)^4 = 10^4 space-time degrees of freedom.
- Results: Even on the 25^3 grid, the high-degree DG polynomial captures the outgoing shock wave in sampled density profiles.The profiles are sampled at 125 equidistant points and compared with the one-dimensional reference solution at t_final = 0.2.
- Large-scale computation: The 100^3-element three-dimensional explosion simulation corresponds to 10 billion space-time degrees of freedom per time step.Figure 13 identifies limited cells in red, unlimited DG cells in blue, and the density distribution on the plane z = 0.
- Results: At t_final = 0.2, Figure 14 compares the one-dimensional reference solution with numerical solutions from 25^3 and 100^3 element grids.The comparison is made using density profiles along the x-axis.
6. Conclusion and Perspectives
The paper presents an a posteriori subcell limiter that detects invalid ADER-DG candidates and recomputes troubled cells with an ADER-WENO finite volume scheme. It also identifies possible extensions for greater robustness, alternative time discretizations, meshes, approximation families, and physical systems.
- Conclusion: The limiter first computes an unlimited high-order ADER-DG candidate and accepts each cell only if positivity and a relaxed discrete maximum principle are satisfied.For troubled cells, the method projects the previous-time DG polynomial onto piecewise-constant subcell averages, evolves them with ADER-WENO, and reconstructs DG polynomials.
- Conclusion: The a posteriori subcell strategy preserves DG’s subcell resolution by recomputing troubled cells on a finer subgrid instead of simply decrementing polynomial order.The subgrid update uses a robust and high-order ADER-WENO finite volume scheme.
- Perspectives: The limiter could also be combined with explicit RK-DG schemes or pointwise WENO schemes using TVD-Runge-Kutta time discretization.These alternatives are proposed for updating cell averages on the subgrid inside troubled zones.
- Perspectives: The WENO3 stage may still encounter positivity problems for extremely complex hyperbolic systems, motivating a more dissipative TVD or Godunov subgrid layer.The proposed cascade would apply the MOOD criteria again after the subcell WENO solver.
- Perspectives: Future work includes extending the limiter to unstructured meshes, PNPM schemes, and systems such as compressible multimaterial flows and magnetohydrodynamics.The paper states that these extensions remain to be investigated.
a priori a posteriori
The figure adds a fallback layer beneath the subcell WENO procedure when the latter fails the MOOD admissibility checks. This layer uses a lower-order TVD finite volume scheme.
- a priori a posteriori: If the WENO procedure fails the MOOD detection criteria, a standard first- or second-order TVD finite volume scheme provides an additional fallback layer.The extra layer follows the WENO procedure in the limiting cascade.