PaperlayerView exact revision on arXiv

Source-linked AI summary

Constraint Preserving AFD-WENO Schemes for Relativistic Hydrodynamics with General Equations of State

Pramodit Mishra, Shubham Upadhyay, Rakesh Kumar, Biswarup Biswas

arXiv:2608.29654v1math.NA

TL;DR

Existing PCP AFD-WENO methods do not yet address relativistic hydrodynamics with general equations of state, despite the need for robust methods under complex thermodynamic closures. This paper develops a state-and-flux-limited PCP AFD-WENO framework and rigorously proves constraint preservation, with benchmarks showing robustness, accuracy, and effective resolution across one- and two-dimensional tests and varied EOS.

  • Problem

    The applicability of PCP AFD-WENO methods to relativistic hydrodynamics with general equations of state remains unexplored under complex thermodynamic closures.

  • Method

    The paper develops a PCP AFD-WENO framework using state and flux limiting alongside efficient state-variable WENO interpolation for relativistic hydrodynamics with a general EOS.

  • Results

    The proposed scheme rigorously preserves physical constraints and demonstrates robustness, accuracy, and effectiveness across extensive one- and two-dimensional tests with varied EOS.

  • Takeaways & Limitations

    The framework provides a robust high-order PCP approach for relativistic hydrodynamics across different equations of state and challenging wave interactions.

  • Takeaways & Limitations

    The formulation requires an equation-of-state closure, while the ideal-gas EOS can poorly approximate semi-relativistic or two-component fluids.

Abstract

from arXiv · show

We develop a high-order physical-constraint-preserving (PCP) alternative finite difference weighted essentially non-oscillatory (AFD-WENO) scheme for the special relativistic hydrodynamics equations with general equations of state. The proposed scheme comprises two key limiters: a state limiter, which acts after the WENO state interpolation step, and a flux limiter, which acts on the final high-order fluxes. The state limiter ensures that the interpolated states are physically admissible, while the flux limiter ensures that the numerical fluxes are physically admissible. The resulting scheme is rigorously proved to satisfy the physical constraints. Incorporating multiple WENO interpolation techniques, including an improved adaptive-order formulation (WENO-AOI), the method is validated through extensive one- and two-dimensional numerical benchmarks with various equations of state. The numerical results demonstrate high-order accuracy, sharp resolution of discontinuities, and robust stability in extreme relativistic regimes.

1 Introduction

Relativistic hydrodynamics requires robust numerical methods because nonlinear conserved–primitive coupling complicates analysis, while physical solutions must preserve density, pressure, and subluminal-velocity constraints. This work addresses the unexplored combination of PCP AFD-WENO methods and general equations of state using state and flux limiting plus efficient state-variable interpolation.

  • Motivation: Highly nonlinear RHD equations involve implicit Lorentz-factor dependence and intricate conserved–primitive coupling, making robust numerical methods indispensable.Relevant applications include relativistic jets, gamma-ray bursts, pulsar wind nebulae, supernovae, and neutron star mergers [35] [37] [18] [22] [2].
  • Physical constraints: Physical RHD solutions require positive rest-mass density, positive pressure, and subluminal velocity; violating these conditions can cause instability and simulation breakdown.The risk is especially pronounced for strong shocks and low-density regions, where early high-order schemes often failed to maintain admissibility.
  • Equation of state: General equations of state are essential for realistic relativistic temperatures, variable composition, and dense nuclear matter, beyond the ideal-gas EOS’s algebraic simplicity.The EOS determines the qualitative and quantitative character of the resulting solution.
  • Research gap: Although PCP frameworks exist for classical WENO and discontinuous Galerkin methods with general EOS and for AFD-WENO with an ideal-gas EOS [48] [3], their combination had remained unexplored.This gap motivates extending physical-constraint preservation to AFD-WENO under general thermodynamic closures.
  • Contribution: The work develops a physical-constraint-preserving AFD-WENO framework for relativistic hydrodynamics with a general EOS, combining flux and state limiting with efficient state-variable WENO interpolation.This targets the previously unexplored intersection of PCP AFD-WENO methods and general thermodynamic closures.

2 Preliminaries

This section formulates the two-dimensional relativistic hydrodynamics system, specifies the general and benchmark equations of state, and describes conservative-to-primitive conversion. It also defines the admissible physical-state constraints whose discrete preservation is essential for robust and stable numerical schemes.

  • Governing equations: The two-dimensional RHD equations are written in conservative form with conserved density, momentum, and total energy variables and fluxes in the x and y directions.The formulation uses c = 1; D is conserved mass density, m_x and m_y are momentum densities, and E is total energy density.
  • Equations of state: The system requires an EOS satisfying relativistic hyperbolicity conditions, including sound speed 0 < c_s < 1 and the stated specific-enthalpy inequality.The conserved and primitive variables are related through Lorentz-factor-dependent transformations, and the general EOS formulation simplifies their conversion.
  • Equations of state: Because the ideal-gas EOS can poorly approximate semi-relativistic or two-component flows, the study also considers TM-EOS [33], IP-EOS, and RC-EOS [40].The ideal-gas model uses Γ ∈(1, 2], while the additional EOS choices provide better relativistic approximations.
  • Variable conversion: Conservative-to-primitive inversion is nonlinear, so the method uses the robust and efficient procedure from, with the RC-EOS-specific method from [4].Both procedures are described as provably robust and efficient across the relevant equations of state.
  • Physical constraints: Physically admissible states require positive density, positive pressure, and subluminal velocity, making discrete constraint preservation crucial for avoiding instability and conversion failure.The admissible conservative-state set is denoted G_p = {u : ρ(u) > 0, p(u) > 0, |v(u)| < 1}.

3 AFD-WENO schemes

This section formulates AFD-WENO schemes by reconstructing conserved-variable point values and combining a low-order Riemann flux with a high-order correction. It also reviews WENO-JS, WENO-Z, and WENO-AO interpolations and introduces WENO-AOI(5,3) as an improved high-order, non-oscillatory variant.

  • 3.1 AFD-WENO framework: AFD-WENO reconstructs conserved-variable point values and evaluates each interface flux as a low-order Riemann flux plus a high-order correction.The correction is built from physical flux values and provides the target order of accuracy.
  • 3.1 AFD-WENO framework: The fifth-order AFD-WENO scheme uses r = 3 and a specified correction stencil to achieve fifth-order accuracy.The correction coefficients are given for different orders, with the work employing the fifth-order case.
  • 3.1 AFD-WENO framework: Component-wise conserved-variable WENO interpolation is computationally efficient but may oscillate near strong discontinuities, motivating local characteristic-space reconstruction.In the LCD approach, stencil values are projected onto interface-local characteristic variables, interpolated component-wise, and transformed back to physical space.
  • 3.2 WENO interpolation procedures: The interpolation component reviews WENO-JS5, WENO-Z5, and WENO-AO(5,3) using Legendre-basis polynomial reconstructions and smoothness indicators.These procedures combine candidate stencil polynomials through nonlinear weights designed from stencil smoothness information.
  • 3.2 WENO interpolation procedures: WENO-AOI(5,3) preserves the WENO-AO reconstruction structure and high-order, non-oscillatory properties while using a different global smoothness indicator to enhance performance.Its nonlinear weights are otherwise computed as in WENO-AO.

4 Extension to Two Dimensions

The AFD-WENO scheme extends to two dimensions on a uniform Cartesian mesh by applying one-dimensional reconstruction and flux evaluation independently in each coordinate direction. The resulting semi-discrete formulation provides the basis for the subsequent physical-constraint-preserving analysis.

  • Computational grid: The two-dimensional domain is partitioned into rectangular cells on a uniform Cartesian mesh with spacings Δx and Δy.The grid approximation u_i,j represents u(t, x_i, y_j) at each grid point.
  • Dimension-by-dimension extension: The two-dimensional AFD-WENO extension applies one-dimensional reconstruction and flux evaluation independently in the x- and y-directions.The resulting semi-discrete finite-difference scheme uses numerical fluxes obtained from the one-dimensional AFD-WENO reconstruction.
  • Constraint-preserving formulation: The semi-discrete formulation, including reconstructed left and right interface states and directional numerical fluxes, underpins the subsequent physical-constraint-preserving analysis.The numerical fluxes are defined separately in the x- and y-directions.

5 Physical Constraint Preservation by AFD-WENO Schemes · 5.1 Admissible Set in Conservative Variables · 5.2 Wave-Speed Estimates and Grid Ratios

The PCP framework reformulates physical admissibility in conservative variables, establishes convexity, and specifies wave-speed, CFL, and grid-ratio conditions for constraint-preserving AFD-WENO updates. These ingredients support admissible state limiting and LLF flux construction, with local coefficients rigorously bounded to preserve physical admissibility.

  • 5.1 Admissible Set in Conservative Variables: The conservative admissible set G is defined by D > 0 and q(u) > 0, and is equivalent to the primitive-variable physical set Gp.Here, D > 0 ensures positive conserved mass density, while q(u) > 0 encodes pressure positivity and subluminal velocity.
  • 5.1 Admissible Set in Conservative Variables: The q-function is a single scalar certificate for pressure positivity and the subluminal velocity constraint, making it the quantity monitored by the numerical scheme.The conservative reformulation therefore directly targets both physical restrictions through q(u).
  • 5.1 Admissible Set in Conservative Variables: The admissible set G is convex because the q-function is concave and D varies linearly under convex combinations.This convexity underpins the state limiter and the multi-dimensional constraint-preservation proof.
  • 5.2 Wave-Speed Estimates and Grid Ratios: The framework uses spectral radii of the coordinate-direction flux Jacobians to define global wave-speed estimates and a uniform CFL-controlled time step.The global estimates also define the scaling parameters used in the update.
  • 5.2 Wave-Speed Estimates and Grid Ratios: Directional splitting weights βx and βy are formed from wave speeds and grid spacings, satisfy βx + βy = 1, and accompany mesh-ratio scaling parameters Λx and Λy.The construction explicitly incorporates the directional grid ratios into the constraint-preserving update.
  • 5.2 Wave-Speed Estimates and Grid Ratios: Global wave speeds determine the uniform time step, while local wave speeds are used in the LLF fluxes.This separates global CFL control from local flux dissipation estimates.

5.3 The PCP Algorithm … Step 3: Low-Order Local Lax–Friedrichs Flux

The PCP algorithm combines component-wise or characteristic WENO interface reconstruction with admissibility-preserving state limiting and low-order LLF flux evaluation at each time step. State limiting uses admissible anchors and separate D and q constraints before subsequent flux calculations.

  • 5.3 The PCP Algorithm: At each time step, the numerical procedure operates on grid points (i,j) and time level t_n, with y-direction operations symmetric to x-direction operations.
  • Step 1: High-Order Interface Reconstruction: WENO interpolation reconstructs left and right interface states component-wise or through local characteristic decomposition.
  • Step 2: State-Limiting Procedure: Because high-order reconstruction may violate admissibility, each reconstructed state is limited toward a provably admissible anchor before flux evaluation.
  • Step 2: State-Limiting Procedure: The anchor state is the arithmetic average of neighboring grid-point values and belongs to G whenever those neighboring values are admissible.
  • Step 2: State-Limiting Procedure: State limiting first enforces D(u_D) ≥ ε_D through a convex blend, then enforces q(u_PCP) ≥ ε_q through a second blend toward the anchor.
  • Step 2: State-Limiting Procedure: The D-stage parameter is computed analytically, whereas the nonlinear q-constraint parameter is obtained by bisection on [0,1].The tolerances ε_D and ε_q are typically 10^-13; limited states replace reconstructed states in all subsequent steps.
  • Step 3: Low-Order Local Lax–Friedrichs Flux: Low-order local Lax–Friedrichs numerical fluxes are evaluated from grid-point values using local wave-speed estimates at each interface.

Step 4: High-Order AFD-WENO Flux

The high-order AFD-WENO fluxes combine LLF fluxes evaluated on PCP-limited interface states with a high-order correction computed from grid-point physical fluxes. Because the correction has no independent admissibility guarantee, the resulting forward-Euler candidate states may be inadmissible, motivating subsequent flux limiting.

  • Flux construction: LLF fluxes are evaluated on PCP-limited reconstructed interface states before forming the high-order AFD-WENO fluxes.The interface states use PCP-limited values rather than grid-point values in the LLF flux evaluation.
  • Flux construction: A high-order correction term from Section 3 is added to the LLF fluxes to obtain the final AFD-WENO fluxes.
  • Flux limiting motivation: The correction terms use only grid-point physical fluxes and lack an independent admissibility guarantee, so flux limiting is required afterward.Without flux limiting, the high-order fluxes may produce inadmissible forward-Euler candidate states.

Step 5: Two-Stage Flux Limiter

The two-stage flux limiter restores admissibility by minimally blending the high-order AFD-WENO flux toward an LLF flux whose forward-Euler candidate states are admissible under the CFL condition. It first enforces the D constraint and then the q constraint, with a symmetric construction in the y-direction.

  • Step 5: Two-Stage Flux Limiter: The limiter minimally blends each interface’s AFD-WENO flux toward the LLF flux until the one-sided forward-Euler candidate states are admissible.The LLF candidate states lie in G under CFL condition (9), providing the admissible reference flux.
  • Step 5: Two-Stage Flux Limiter: The same limiting procedure is applied symmetrically in the y-direction to obtain the PCP flux ĝ_PCP.The x-direction construction is the representative case.
  • Step 5: Two-Stage Flux Limiter: Stage II blends the intermediate flux with the LLF flux to enforce the q constraint and produce the final PCP flux.Small positive thresholds ε̃_D and ε̃_q are chosen from the strictly positive LLF values; 10^-13 is used in simulations.
  • Step 5: Two-Stage Flux Limiter: Stage I limits only the D component, preserving the high-order AFD-WENO momentum and energy components.The resulting intermediate flux is denoted by f̂^D.

Step 6: Conservative Update

The semi-discrete scheme is advanced with an admissibility-preserving time-stepping method, presented using forward Euler; Figure 1 summarizes the complete algorithm.

  • Step 6: Conservative Update: The update advances grid-point values to time level t_n+1 using forward Euler, a time-stepping method that preserves admissibility.This discretizes the semi-discrete scheme (5) in time while maintaining admissibility.
  • Step 6: Conservative Update: Figure 1 provides an overview of the complete algorithm.

5.4 Theoretical Analysis

The analysis proves that the two-stage PCP flux limiter preserves the admissible set G. LLF positivity, staged D/q enforcement, and a final convex-combination argument establish constraint preservation at every time step.

  • Lemma 5.2: Under CFL ≤ 1/2 and bounded local wave speed, LLF forward-Euler candidate states belong to G.The proof uses nonnegative coefficients summing to one, so each candidate is a convex combination of admissible states.
  • Lemma 5.3: Stage I enforces D(u±,D) ≥ ε̃_D for both intermediate candidate states.Because the first correction modifies only the D-component, the limiter either leaves the state unchanged or sets the component exactly to ε̃_D.
  • Lemma 5.4: Stage II enforces q(u±,PCP) ≥ ε̃_q while retaining D(u±,PCP) ≥ ε̃_D.The q constraint follows from concavity and the selected limiter parameter, while D preservation follows from convexly combining states already satisfying the D bound.
  • Theorem 5.5: Theorem 5.5 proves that the updated grid-point values produced by the two-stage PCP flux limiter remain in G.The proof applies Lemma 5.4 to four one-sided forward-Euler states and rewrites the conservative update as their convex combination.
  • Theorem 5.5: The final update is admissible because coefficients β_x/2, β_x/2, β_y/2, and β_y/2 are nonnegative and sum to one.With all four one-sided states in G, the convex-combination representation and Lemma 5.1 imply the updated state lies in G.

6 Numerical Results

Across one- and two-dimensional benchmarks, the schemes accurately capture exact solutions, shocks, contacts, and complex relativistic flow structures across four equations of state. WENO-AOI generally provides the sharpest resolution while maintaining physical-constraint preservation and stability in extreme regimes.

  • Convergence tests: Case-I and the low-density, low-pressure Case-II converge to the exact solution, with Case-II specifically testing robustness under challenging physical conditions.The convergence results are reported in Tables 1 and 2 for Case-I and Case-II, respectively.
  • One-dimensional comparisons: WENO-AOI provides slightly sharper resolution than WENO-JS, WENO-Z, and WENO-AO in one-dimensional tests, including rarefaction peaks and shock-contact structures across all four EOS.WENO-AO and WENO-AOI are sharper than WENO-JS and WENO-Z, with WENO-AOI performing best among the compared schemes.
  • Overall benchmark performance: Across the benchmark suite, all schemes accurately resolve essential solution features and strong discontinuities across the four equations of state, while maintaining PCP in challenging interactions.The results include close agreement with reference solutions, accurate shocks and contacts, and robustness for strong discontinuities and complex wave interactions.
  • Two-dimensional wave interactions: WENO-AOI consistently captures complex two-dimensional wave interactions and sharp structures across EOS, while ID-EOS can produce visibly different low-density cores or detailed vortex structures.Other EOS often exhibit similar overall flow patterns, but EOS choice affects detailed density, vortex, and wave structures.
  • Shock-interaction benchmarks: For shock–bubble and double-Mach-reflection tests, WENO-AOI resolves transmitted and reflected waves, deformed interfaces, vortices, and characteristic shock patterns across all four EOS.EOS choices substantially change post-shock density distributions and bubble dynamics, while the principal structures remain consistently resolved.
  • Relativistic flow tests: In shear-layer and relativistic-jet tests, WENO-AOI produces comparable primary-vortex structures across EOS and remains stable and well resolved for increasingly collimated, near-light-speed jets.The jet tests show stronger axial compression and sharper internal structures in cold configurations as beam velocity approaches the speed of light.

7 Conclusion

The paper develops a high-order accurate, physical-constraint-preserving AFD-WENO scheme for relativistic hydrodynamics with a general equation of state. Rigorous analysis and one- and two-dimensional experiments establish its physical admissibility, broad EOS applicability, robustness, accuracy, and effectiveness.

  • 7 Conclusion: The proposed scheme is a high-order accurate, physical-constraint-preserving AFD-WENO method for RHD equations with a general EOS.
  • 7 Conclusion: The scheme’s PCP property is rigorously proven, ensuring numerical solutions remain within the physically admissible set.
  • 7 Conclusion: Extensive one- and two-dimensional numerical experiments demonstrate robustness, accuracy, and effectiveness across a wide range of EOS.
Loading 2608.29654v1…