Logo
ScienceStack
Table of Contents
Fetching ...
Sign In

Shock with Confidence: Formal Proofs of Correctness for Hyperbolic Partial Differential Equation Solvers

Jonathan Gorard, Ammar Hakim

TL;DR

The paper develops a formal verification pipeline in Racket to produce provably correct first-order hyperbolic PDE solvers and automatically generate verifiable C implementations, addressing the numerical challenges posed by shocks and discontinuities. By encoding PDEs, flux functions, and CFL constraints in a domain-specific language, and applying floating-point-aware symbolic theorem-proving and automatic differentiation, the work proves key properties such as hyperbolicity preservation, $L^1$/$L^2$/$L^\infty$ stability, and flux continuity for selected schemes like Lax-Friedrichs and Roe, including second-order extensions with flux limiters. The results demonstrate full proofs for scalar equations and several vector systems (notably Maxwell and isothermal Euler) with partial proofs in some isothermal Euler components, highlighting both the potential and current limitations of the approach. The work provides a solid proof-of-concept toward widespread formal verification of numerical solvers, with plans to extend to ODE integrators, DG methods, and broader solver families, thereby increasing reliability in scientific computing pipelines such as Gkeyll.

Abstract

First-order systems of hyperbolic partial differential equations (PDEs) occur ubiquitously throughout computational physics, commonly used in simulations of fluid turbulence, shock waves, electromagnetic interactions, and even general relativistic phenomena. Such equations are often challenging to solve numerically in the non-linear case, due to their tendency to form discontinuities even for smooth initial data, which can cause numerical algorithms to become unstable, violate conservation laws, or converge to physically incorrect solutions. In this paper, we introduce a new formal verification pipeline for such algorithms in Racket, which allows a user to construct a bespoke hyperbolic PDE solver for a specified equation system, generate low-level C code which verifiably implements that solver, and then produce formal proofs of various mathematical and physical correctness properties of the resulting implementation, including L^2 stability, flux conservation, and physical validity. We outline how these correctness proofs are generated, using a custom-built theorem-proving and automatic differentiation framework that fully respects the algebraic structure of floating-point arithmetic, and show how the resulting C code may either be used to run standalone simulations, or integrated into a larger computational multiphysics framework such as Gkeyll.

Shock with Confidence: Formal Proofs of Correctness for Hyperbolic Partial Differential Equation Solvers

TL;DR

The paper develops a formal verification pipeline in Racket to produce provably correct first-order hyperbolic PDE solvers and automatically generate verifiable C implementations, addressing the numerical challenges posed by shocks and discontinuities. By encoding PDEs, flux functions, and CFL constraints in a domain-specific language, and applying floating-point-aware symbolic theorem-proving and automatic differentiation, the work proves key properties such as hyperbolicity preservation, // stability, and flux continuity for selected schemes like Lax-Friedrichs and Roe, including second-order extensions with flux limiters. The results demonstrate full proofs for scalar equations and several vector systems (notably Maxwell and isothermal Euler) with partial proofs in some isothermal Euler components, highlighting both the potential and current limitations of the approach. The work provides a solid proof-of-concept toward widespread formal verification of numerical solvers, with plans to extend to ODE integrators, DG methods, and broader solver families, thereby increasing reliability in scientific computing pipelines such as Gkeyll.

Abstract

First-order systems of hyperbolic partial differential equations (PDEs) occur ubiquitously throughout computational physics, commonly used in simulations of fluid turbulence, shock waves, electromagnetic interactions, and even general relativistic phenomena. Such equations are often challenging to solve numerically in the non-linear case, due to their tendency to form discontinuities even for smooth initial data, which can cause numerical algorithms to become unstable, violate conservation laws, or converge to physically incorrect solutions. In this paper, we introduce a new formal verification pipeline for such algorithms in Racket, which allows a user to construct a bespoke hyperbolic PDE solver for a specified equation system, generate low-level C code which verifiably implements that solver, and then produce formal proofs of various mathematical and physical correctness properties of the resulting implementation, including L^2 stability, flux conservation, and physical validity. We outline how these correctness proofs are generated, using a custom-built theorem-proving and automatic differentiation framework that fully respects the algebraic structure of floating-point arithmetic, and show how the resulting C code may either be used to run standalone simulations, or integrated into a larger computational multiphysics framework such as Gkeyll.

Paper Structure

This paper contains 14 sections, 8 theorems, 47 equations, 5 tables.

Table of Contents

  1. Introduction
  2. Preliminaries
  3. Hyperbolic PDE Solvers
  4. The Lax-Friedrichs Flux
  5. The Roe Flux
  6. Flux Extrapolation and Limiters
  7. Methodology and Results
  8. Automatic Code Generation
  9. Symbolic Theorem-Proving
  10. Automatic Differentiation
  11. Stability, Hyperbolicity and Convexity
  12. Symbolic Limits and Symmetry
  13. Results
  14. Conclusions and Future Work

Key Result

theorem 1

A Lax-Friedrichs solver preserves hyperbolicityhartman_lemma_1960 (i.e. preserves well-posedness of the Cauchy problem along all non-characteristic hypersurfaces) for an arbitrary hyperbolic PDE system if the flux Jacobian with respect to the conserved variable vector ${\mathbf{U}}$: is diagonalizable, with purely real eigenvalues.

Theorems & Definitions (8)

  • theorem 1
  • theorem 2
  • theorem 3
  • theorem 4
  • theorem 5
  • theorem 6
  • theorem 7
  • theorem 8