A dual-pairing summation-by-parts finite difference framework for nonlinear conservation laws

Robust and convergent high-order numerical methods for solving partial differential equations are highly attractive due to their efficiency on modern and next-generation hardware architectures. However, designing such methods for nonlinear hyperbolic conservation laws remains a significant challenge. In this work, we introduce a framework based on dual-pairing (DP) and upwind summation-by-parts (SBP) finite difference (FD) and discontinuous Galerkin (DG) finite element methods, aimed at achieving accurate and robust numerical approximations of nonlinear conservation laws. The framework ensures entropy consistency and features an intrinsic high-order accurate filter designed to detect and resolve regions where the solution is poorly captured or discontinuities are present. The DP SBP FD/DG operators form a dual pair of discrete derivative operators that collectively preserve the SBP property. Furthermore, these operators are constructed to be upwind, allowing them to incorporate dissipation within the elements themselves.This contrasts with traditional SBP and collocated DG spectral element methods, which typically induce dissipation solely through numerical fluxes at element interfaces. Our framework facilitates the systematic combination of DP SBP FD/DG operators with skew-symmetric and upwind flux splitting techniques. This integration enables the development of robust, high-order accurate schemes for nonlinear hyperbolic conservation laws.

💡 Research Summary

This paper introduces a novel high‑order numerical framework for solving nonlinear hyperbolic conservation laws, built on the combination of dual‑pairing (DP) operators and upwind summation‑by‑parts (SBP) techniques for both finite‑difference (FD) and discontinuous‑Galerkin (DG) discretizations. The authors identify two fundamental shortcomings of existing high‑order SBP/DG‑SEM schemes when applied to nonlinear problems with shocks or turbulence: (1) traditional central‑difference SBP operators provide dissipation only at element interfaces, leaving interior regions vulnerable to spurious high‑frequency modes, and (2) the lack of a discrete chain‑rule or product‑rule forces the continuous equations to be rewritten in skew‑symmetric form, which complicates the construction of entropy‑conserving schemes.

To overcome these issues, the paper proposes a dual‑pairing SBP construction. Two discrete derivative matrices, (D) and (\tilde D), together with their respective mass matrices (P) and (\tilde P), satisfy the relation
(D^{!T}P + \tilde P,\tilde D = B),
where (B) contains only boundary contributions. This dual pair preserves the SBP property while allowing each operator to embed an upwind (dissipative) component directly inside the element. Consequently, numerical viscosity is no longer confined to interfaces; it is distributed throughout the element, improving robustness for discontinuous solutions.

Entropy consistency is achieved by employing a skew‑symmetric flux splitting. The continuous conservation law (\partial_t u + \nabla!\cdot f(u)=0) is rewritten as a convex combination of the flux form and the quasi‑linear form, parameterized by a scalar (\alpha). This reformulation yields a differential operator (F(u,f,\partial_x)) that satisfies
((g, F) = (1, \nabla!\cdot q))
with entropy variables (g = \partial_u e(u)) and entropy flux (q). Because only integration‑by‑parts is required, the discrete entropy analysis mirrors the continuous one, avoiding the need for a discrete chain rule.

The semi‑discrete scheme obtained from the DP‑SBP operators and the skew‑symmetric split is shown to be (i) conservative, (ii) entropy‑stable (i.e., the discrete total entropy does not increase), and (iii) provably high‑order accurate. A high‑order “filter” is also introduced: a polynomial‑based smoother that activates automatically in regions where the solution exhibits large gradients or discontinuities, damping spurious high‑frequency components without degrading the formal order of accuracy.

Numerical experiments cover three benchmark problems: (a) inviscid Burgers’ equation, (b) the nonlinear shallow‑water equations, and (c) the two‑dimensional compressible Euler equations. In all cases, the method attains the expected convergence rates, preserves the designed entropy dissipation at shocks, and suppresses non‑physical oscillations that plague traditional SBP/DG‑SEM schemes. Multi‑block configurations demonstrate that interior upwind dissipation works independently within each block while inter‑block interfaces are treated with standard numerical fluxes, ensuring scalability on modern parallel architectures.

Overall, the paper delivers a comprehensive framework that simultaneously addresses high‑order accuracy, robustness, and physical fidelity for nonlinear conservation laws. By embedding upwind dissipation directly into the dual‑pairing SBP operators and coupling them with skew‑symmetric flux splitting, the authors achieve provable entropy stability without sacrificing spectral accuracy. The work opens avenues for extending the approach to more complex systems (e.g., magnetohydrodynamics, multiphase flows), unstructured meshes, and hardware‑accelerated implementations on GPUs and next‑generation exascale platforms.