Source-linked AI summary
The BR1 Scheme is Stable for the Compressible Navier-Stokes Equations
Gregor J. Gassner, Andrew R. Winters, Florian J. Hindenlang, David A. Kopriva
TL;DR
BR1’s stability for compressible Navier–Stokes equations had not been established, despite its use and the need for ad hoc stabilization. The paper analyzes a BR1-type DGSEM on three-dimensional curvilinear meshes under specified discretization conditions. It shows that BR1 interface fluxes add no artificial dissipation at interior faces and that the resulting scheme is energy-stable for linear equations and entropy-stable for nonlinear equations.
Problem
BR1’s stability for compressible Navier–Stokes equations had not been shown, although ad hoc stabilization is often needed for DG approximations.
Source: §1 · p. 2Method
The paper analyzes a BR1-type DGSEM on three-dimensional curvilinear meshes with discrete metric identities, averaged metric terms, and entropy-variable gradients for viscous terms.
Source: §1 · p. 3, §5 · p. 35Results
BR1 fluxes add no artificial dissipation at interior faces, and the scheme is energy-stable for linear or entropy-stable for nonlinear compressible Navier–Stokes equations.
Source: §5 · p. 35Takeaways & Limitations
The work presents a full proof of entropy stability for DGSEM discretizations of compressible Euler terms on three-dimensional curvilinear meshes.
Source: §5 · p. 35Entropy stability requires positive density and pressure and does not guarantee positive density or temperature in computed solutions.
Source: §4.2 · p. 29
Abstract
from arXiv · showhide
We show how to modify the original Bassi and Rebay scheme (BR1) [F. Bassi and S. Rebay, A High Order Accurate Discontinuous Finite Element Method for the Numerical Solution of the Compressible Navier-Stokes Equations, Journal of Computational Physics, 131:267--279, 1997] to get a provably stable discontinuous Galerkin collocation spectral element method (DGSEM) with Gauss-Lobatto (GL) nodes for the compressible Navier-Stokes equations (NSE) on three dimensional curvilinear meshes. Specifically, we show that the BR1 scheme can be provably stable if the metric identities are discretely satisfied, a two-point average for the metric terms is used for the contravariant fluxes in the volume, an entropy conserving split form is used for the advective volume integrals, the auxiliary gradients for the viscous terms are computed from gradients of entropy variables, and the BR1 scheme is used for the interface fluxes. Our analysis shows that even with three dimensional curvilinear grids, the BR1 fluxes do not add artificial dissipation at the interior element faces. Thus, the BR1 interface fluxes preserve the stability of the discretization of the advection terms and we get either energy stability or entropy-stability for the linear or nonlinear compressible NSE, respectively.
1 Introduction
The paper addresses the unproven stability of BR1 for compressible Navier–Stokes equations and develops a stable version for three-dimensional curvilinear grids. Its analysis treats linear energy stability and nonlinear entropy stability under specific discretization choices.
1 Introduction
- BR1 rewrites second-order equations as a first-order system with gradients as new unknowns and uses arithmetic means for viscous fluxes and solution values at interfaces.This makes it a simple variant for discretizing second-order terms.Source: §1 · p. 2
- Despite BR1’s popularity, its stability for compressible Navier–Stokes equations had not been shown, while ad hoc stabilization is often needed in practice.Earlier concerns about BR1 for elliptic problems include instability and a widened stencil.Source: §1 · p. 2
- The paper provides energy analysis for the linearized equations and entropy analysis for the nonlinear equations, extending the argument to three-dimensional unstructured curvilinear grids.The authors introduce the basic steps through one-dimensional scalar examples before the general three-dimensional analysis.Source: §1 · p. 3
- BR1-type discretization is provably stable when discrete metric identities, averaged metric terms, split-form advection, entropy-variable gradients, and BR1 interface fluxes are used.The analysis covers compressible Navier–Stokes equations on three-dimensional curvilinear grids.Source: §1 · p. 3
2 The DGSEM for Compressible Viscous Flows
This section formulates a DGSEM for compressible viscous flows on mapped hexahedral elements and specifies its polynomial, quadrature, and flux approximations. It also explains why the usual approximation may lack linear stability and identifies discrete metric identities as a key geometric condition.
2 The DGSEM for Compressible Viscous Flows
- The spectral element method maps conforming hexahedral elements, including curved elements, to a reference cube and transforms gradients, divergence, and fluxes using geometric metrics.The transformed Navier–Stokes system is written on the reference element.Source: §2 · pp. 5–6
2.1 The Spectral Element Approximation
- The DGSEM approximates states and fluxes with degree-N polynomials at nodal points, uses Gauss–Lobatto quadrature, and couples elements through numerical fluxes.Gauss–Lobatto nodes include boundary points, avoiding interpolation of volume data to element surfaces when constructing surface terms.Source: §2.1 · pp. 7–9
- The discrete metric identities are enforced through the metric construction, ensuring discrete free-stream preservation and supporting stability.The metric terms are approximated with degree-N polynomials.Source: §2.1 · p. 7
- Gauss and Gauss–Lobatto quadrature yield a summation-by-parts formula, which extends to all space dimensions and underlies the discrete weak forms.The weak forms couple elements through numerical fluxes.Source: §2.1 · pp. 8–9
- The usual DGSEM approximation is not necessarily linearly stable; the analysis constructs stable advection and diffusion approximations while retaining BR1 surface viscous fluxes.The section then examines contributors to instability in underresolved flows.Source: §2.1 · p. 9
3 Linear Stability Analysis
This section analyzes stability for split-form DGSEM discretizations with BR1 diffusion, from a 1D advection–diffusion setting to the 3D linearized compressible Navier–Stokes equations. It examines the role of metric identities, interface dissipation, and stable treatment of physical boundaries in the resulting energy bound.
3.1 BR1 is Stable: Linear Scalar Advection-Diffusion in 1D
- The scalar analysis uses a variable-coefficient advection–diffusion problem as a blueprint, retaining split form even though constant advection speed makes it equivalent to conservative form.The setup takes constant positive advection speed and positive variable diffusion coefficient.Source: §3.1 · pp. 9–10
3.1.1 Continuous Energy Analysis in 1D
- The continuous 1D energy analysis derives an energy bound that serves as the target for stability of the discrete approximation.The boundary value problem is described as well-posed, and the bound follows from integrating the energy relation over time and applying the initial condition.Source: §3.1.1 · p. 10
3.1.2 An Energy Stable DGSEM in 1D
- The 1D DGSEM is stable for advection–diffusion with split-form advection and BR1 diffusion.The upwind advective interface contribution is non-positive, and the boundary conditions are dissipative.Source: §3.1.2 · p. 13
3.2 BR1 is Stable: 3D Linearized Compressible Navier-Stokes Equations
- The linearized 3D Navier–Stokes analysis uses symmetrization and an energy method, assuming stable approximations of advection and physical boundary conditions.The treatment focuses on those two steps, leaving the approximation of physical boundary conditions to future work.Source: §3.2 · p. 14
3.2.1 Continuous Energy Analysis in 3D
- For the continuous 3D linearized system, the flux-volume terms cancel, so L2 energy growth is determined by the boundary integral.The energy analysis tests with a symmetrized state and considers homogeneous boundary data.Source: §3.2.1 · pp. 15–17
3.2.2 An Energy Stable DGSEM in 3D
- In 3D, the discrete metric identities make the transformed coefficient-matrix divergence exactly zero; without them, aliasing could disrupt stability.The DGSEM uses Gauss–Lobatto quadrature and polynomial interpolants; its advective volume terms cancel, and the physical boundary terms are assumed dissipative.Source: §3.2.2 · pp. 18–20
- With upwind or central advective fluxes, interface contributions are dissipative, while the BR1 viscous surface contribution vanishes.The advective interface terms contribute nonpositively to the energy; BR1 itself adds no energy dissipation.Source: §3.2.2 · pp. 20–21
- With stably approximated physical boundary terms, the 3D DGSEM satisfies an energy stability bound, with physical dissipation and, for σ > 0, jump-dependent surface dissipation.The BR1 scheme itself has no effect on the energy; stability requires properly posed and stably implemented physical boundary conditions.Source: §3.2.2 · p. 21
3.3 Why is Instability Seen in Practice?
- For smooth flows, stability depends on the advective approximation, not BR1, which neither causes instability nor adds stabilizing dissipation.In the standard-form DGSEM, aliasing creates a product-rule error of indeterminate sign; under-resolved solutions or high Reynolds numbers can leave insufficient physical or interface dissipation to counter it.Source: §3.3 · pp. 21–22
- Stable approximation of the advective terms is ultimately key; Riemann-solver dissipation can counter aliasing, with Lax-Friedrichs more dissipative than exact upwind.For well-resolved solutions, aliasing is spectrally small and physical or interface dissipation may suffice; changing BR1 alone is not the central stability requirement.Source: §3.3 · p. 22
4 Nonlinear Stability Analysis
This section develops a nonlinear entropy-stability argument, using 1D viscous Burgers to motivate an entropy estimate for the 3D compressible Navier–Stokes DGSEM. The estimate uses entropy-conserving advection with averaged metric terms and BR1 viscous fluxes, under positivity and boundary-condition assumptions, with viscous contributions at interior faces canceling exactly.
4.1.1 Continuous Entropy Analysis in 1D
- In the continuous Burgers analysis, positive viscosity makes the viscous entropy contribution non-positive, bounding total entropy by its initial value.The analysis tests the equation with the entropy variable and uses the viscous flux relation.Source: §4.1.1 · pp. 23–24
4.1.2 An Entropy-Stable DGSEM in 1D
- The standard DGSEM is unstable regardless of surface-flux choice because aliasing in the advective volume terms creates an unbounded discrete entropy source or sink.Under coarse resolution, entropy can accumulate and the solution can blow up.Source: §4.1.2 · p. 24
- For 1D viscous Burgers, the entropy-conserving advection discretization paired with BR1 viscous fluxes yields a discrete entropy estimate, with additional dissipation only at physical boundaries.The viscous contribution is non-positive; the estimate mimics the continuous one.Source: §4.1.2 · pp. 26–28
4.2 BR1 is Entropy-Stable: 3D Nonlinear Compressible Navier-Stokes Equations
- For compressible Navier–Stokes, entropy-stability is stronger than linear stability but is not full nonlinear stability: it requires positive density and pressure and can still permit invalid states.The positivity requirement makes the mapping between conservative and entropy variables one-to-one; additional strategies are needed to guarantee positivity.Source: §4.2 · p. 29
4.2.1 Continuous Entropy Analysis in 3D
- The continuous 3D Navier–Stokes entropy estimate follows by testing with entropy variables and viscous flux, rewriting viscous terms using entropy-variable gradients, and assuming non-positive boundary contributions.Under that boundary assumption, entropy does not increase in time.Source: §4.2 · p. 29, §4.2.1 · p. 30
4.2.2 An Entropy-Stable DGSEM in 3D
- The 3D DGSEM analysis begins by transforming the Navier–Stokes equations to reference space and introducing an alternative expression for the viscous flux.Source: §4.2.2 · p. 31
- The advective volume discretization uses an entropy-conservative two-point flux and averaged metric terms to obtain entropy stability and address metric aliasing.The standard Euler DGSEM is not entropy-stable because nonlinear fluxes generate aliasing errors; the entropy-stable form separates the averaged metric terms from the specially averaged flux.Source: §4.2.2 · pp. 31–32
- The entropy estimate assumes a non-negative interpolant of the element mapping Jacobian determinant at Gauss-Lobatto nodes and boundary fluxes chosen to make boundary effects dissipative.The advective surface flux is selected so the surface-quadrature argument vanishes at each interior-face quadrature point.Source: §4.2.2 · pp. 33–35
- BR1 viscous numerical flux contributions at interior element faces vanish exactly.The viscous numerical surface fluxes are chosen as averages.Source: §4.2.2 · p. 34
- The final discrete entropy statement bounds the total entropy rate by physical-boundary contributions and added dissipation from numerical surface fluxes.Stability requires choosing physical-boundary terms so their effect is dissipative.Source: §4.2.2 · pp. 34–35
5 Conclusions
The analysis establishes stability for a BR1-type approximation on three-dimensional curvilinear meshes under stated discretization conditions, and extends the stability estimate to general advection-diffusion problems when their advection DGSEM is stable. It also reports stable, virtually dissipation-free behavior in an under-resolved Taylor-Green vortex test.
5 Conclusions
- Stability requires discrete metric identities, two-point metric averages in the volume, entropy-variable gradients for viscous terms, and BR1 interface fluxes.The estimate extends to general advection-diffusion problems when stable DGSEM approximations of the advection terms are known.Source: §5 · p. 35
- BR1 fluxes add no artificial dissipation at interior faces, preserving advection stability and yielding energy stability for linear or entropy-stability for nonlinear compressible NSE.The result holds on three-dimensional curvilinear grids.Source: §5 · p. 35
- Numerical investigations report robustness against standard DGSEM and ad hoc stabilization for under-resolved compressible turbulence simulations.The paper also presents a full proof of entropy-stability for DGSEM compressible Euler terms on three-dimensional curvilinear meshes.Source: §5 · p. 35
- At 643 degrees of freedom, the modified scheme does not crash without added interface dissipation, while standard DGSEM immediately crashes, even with Riemann-solver dissipation.The Taylor-Green test uses N = 7 and 83 elements; at later times its numerical Reynolds number approaches the physical Re = 1600.Source: §5 · pp. 35–36, Figure 1 · p. 36
A An Entropy Conserving Euler Flux
The appendix presents an entropy-conserving Euler flux and an entropy-stable variant with matrix dissipation. It also specifies metric averaging for contravariant fluxes and notes that the chosen dissipation resolves stationary contact discontinuities.
A An Entropy Conserving Euler Flux
- The paper explicitly presents a kinetic-energy-preserving and entropy-conservative flux for the compressible Euler equations, first derived by Chandrashekar [8].Contravariant entropy-conservative fluxes incorporate an average of the metric terms.Source: A An Entropy Conserving Euler Flux · pp. 37–38
- A numerically stable computation of the logarithmic mean is available when its two arguments are approximately equal.The logarithmic mean is identified in the flux definition.Source: A An Entropy Conserving Euler Flux · p. 37
- Adding matrix dissipation gives an entropy-stable approximation, and the selected discrete dissipation operator exactly resolves stationary contact discontinuities.The stated proof follows the structure presented by Chandrashekar [8].Source: A An Entropy Conserving Euler Flux · p. 37
B.1 Proof of Entropy Conservation for the One-Dimensional Volume Integral
The proof rewrites the one-dimensional volume integral using the summation-by-parts matrix and the entropy-conservation condition. Symmetry and reindexing then cancel terms; the stated result depends on summation by parts rather than on DG specifically.
B.1 Proof of Entropy Conservation for the One-Dimensional Volume Integral
- The proof expresses the volume integral as a sum and uses the summation-by-parts matrix Q = MD, with entries Q_im = ω_i D_im.The boundary evaluation matrix is B = diag([-1, 0, ..., 0, 1]), and Q satisfies a property because the derivative of a constant is exact.Source: §B.1 · pp. 38–39
- Reindexing the sums and using the symmetry of the entropy-conserving flux allows terms to be recombined and canceled to obtain the stated result.The derivation also uses the entropy potential and the flux consistency condition at equal states.Source: §B.1 · pp. 39–40
- The stated result depends only on the summation-by-parts property, not specifically on the discontinuous Galerkin approximation.Source: §B.1 · p. 40
B.2 Proof of Entropy Conservation at Interelement Interfaces for Curvilinear Elements
The interface proof assumes a watertight mesh and uses discrete metric identities to ensure continuity of the surface metric. It then moves the metric into the jump and establishes pointwise vanishing of the integrand at interior-face quadrature points.
B.2 Proof of Entropy Conservation at Interelement Interfaces for Curvilinear Elements
- For a conforming hexahedral mesh, discrete metric identities and the stated construction of normals and surface elements ensure discrete surface-metric continuity under the watertight-mesh assumption.The assumption is that the normal vector and surface element are continuous across element interfaces.Source: §B.2 · p. 40
- Surface-metric continuity permits using the master-side metric and moving it inside the jump in the interface proof.The proof writes a normal contravariant surface flux using Cartesian fluxes and the metric.Source: §B.2 · pp. 40–41
- The resulting integrand vanishes at each quadrature point of an interior face individually.The proof treats the approximate surface integral at its Gauss-Lobatto quadrature points.Source: §B.2 · pp. 40–41
B.3 Proof of Entropy Conservation in 3D Curvilinear Coordinates
The 3D curvilinear-coordinate proof establishes the entropy condition for the volume approximation when the discrete metric identities hold. It extends the one-dimensional argument across all three coordinate directions using entropy potentials and summation-by-parts manipulations.
B.3 Proof of Entropy Conservation in 3D Curvilinear Coordinates
- The proof introduces one entropy potential for each Cartesian direction to rewrite the two-point flux entropy-conservation condition.Source: §B.3 · p. 41
- For the ξ-direction volume term, the proof follows the one-dimensional argument using summation-by-parts, re-indexing, and symmetry properties.It uses Q_im = ω_iD_im and the stated properties of the metric term and two-point flux.Source: §B.3 · pp. 42–44
- The same result holds for the other two coordinate-direction terms, completing the three-dimensional volume approximation argument.Source: §B.3 · p. 44
- The regrouped boundary and volume terms satisfy the desired entropy condition provided the metric identities hold discretely.The condition applies at all Gauss-Lobatto points within an element.Source: §B.3 · p. 45