Well-balanced high order finite difference WENO schemes for a first-order Z4 formulation of the Einstein field equations

TL;DR

This work targets the need for very high-order, robust numerical methods in numerical relativity by developing well-balanced finite-difference WENO schemes for the first-order non-conservative Z4 (FO-Z4) formulation of the Einstein equations. It introduces two families of schemes, FD-WENO (reconstruction-based) and AFD-WENO (interpolation-based), tailored to handle non-conservative hyperbolic systems with constraint damping and auxiliary variables, while enabling well-balancing to preserve stationary equilibria such as Kerr–Schild spacetimes. The methods are validated on a suite of standard NR tests, including gauge waves, robust stability, Gowdy waves, stationary black holes, and head-on BH collisions, demonstrating accurate high-order convergence, long-time stability, and the ability to maintain exact equilibria with well-balancing activated. The results suggest that memory-efficient FD-WENO/AFD-WENO approaches are viable competitors for production NR codes, offering high accuracy, strong shock-capturing capabilities when needed, and a flexible ecosystem of enhancements such as divergence-curl preservation and PCP properties.

Abstract

In this work we aim at developing a new class of high order accurate well-balanced finite difference (FD) Weighted Essentially Non-Oscillatory (WENO) methods for numerical general relativity, which can be applied to any first-order reduction of the Einstein field equations, even if non-conservative terms are present. We choose the first-order non-conservative Z4 formulation of the Einstein equations, which has a built-in cleaning procedure that accounts for the Einstein constraints and that has already shown its ability in keeping stationary solutions stable over long timescales. Upon the introduction of auxiliary variables, the vacuum Einstein equations in first order form constitute a ...

Figures (16)

• Figure 1: The schematic structure of a fifth order FD-WENO scheme in conservation form. Please focus on the purple zone boundary at $i+1/2$. The numerical flux at that zone boundary is made up of two split flux contributions $\mathbf{F}^+$ and $\mathbf{F}^-$. The stencils contributing to the reconstruction of the right-going flux $\mathbf{F}^+$ at the purple zone boundary are shown. The stencils contributing to the reconstruction of the left-going flux $\mathbf{F}^-$ at the purple zone boundary are also shown. .
• Figure 2: Panel shows part of the mesh around zone “i”. The mesh functions are collocated at the zone centers, as shown by the thick dots. The zone boundaries are shown by the vertical lines. The figure also shows the stencils associated with the zone “$i$” for the fifth order WENO--AO reconstruction/interpolation. We have three smaller third order stencils and a large fifth order stencil. The reconstructed/interpolated variables at the zone boundaries are shown with a caret. The variables with a superscript star are resolved states obtained by the pointwise application of a simple HLL or LLF Riemann solver at the zone boundaries.
• Figure 3: Panel shows part of the mesh around zone boundary “$i+1/2$”. The fluxes are evaluated pointwise at the zone centers, as shown by the thick dots. The zone boundaries are shown by the vertical lines. The figure also shows the stencils associated with the zone boundary “$i+1/2$” for the third and fifth order AFD-WENO schemes. We have two smaller third order stencils and a large sixth order stencil. For a third order AFD-WENO scheme, the two smaller stencils can be non-linearly hybridized. In that case, the second derivatives of the flux can be obtained at the zone boundary when the smoothness in the solution warrants it. For fifth order AFD-WENO, the two smaller stencils can be non-linearly hybridized along with the larger stencil. In that case, the second and fourth derivatives of the flux can be obtained at the zone boundary when the smoothness in the solution warrants it. The process described here can be done for Adaptive Order and Multiresolution WENO interpolation.
• Figure 4: \ref{['sec:gaugewave']}: The left panel (Figs. \ref{['fig:FDgauge']}a, \ref{['fig:FDgauge']}c, \ref{['fig:FDgauge']}e, \ref{['fig:FDgauge']}g) shows the difference among the numerical solution and the exact one for the variable $\alpha$ (lapse) at time $t=1000$ for various order accurate FD--WENO schemes. The right panel (Figs. \ref{['fig:FDgauge']}b, \ref{['fig:FDgauge']}d, \ref{['fig:FDgauge']}f, \ref{['fig:FDgauge']}h) shows the $L_2$-error evolution of the ADM constraints for the corresponding order scheme.
• Figure 5: \ref{['sec:gaugewave']}: The left panel (Figs. \ref{['fig:AFDgauge']}a, \ref{['fig:AFDgauge']}c, \ref{['fig:AFDgauge']}e, \ref{['fig:AFDgauge']}g) shows the difference among the numerical solution and the exact one for the variable $\alpha$ (lapse) at time $t=1000$ for various order accurate FD--WENO schemes. The right panel (Figs. \ref{['fig:AFDgauge']}b, \ref{['fig:AFDgauge']}d, \ref{['fig:AFDgauge']}f, \ref{['fig:AFDgauge']}h) shows the $L_2$-error evolution of the ADM constraints for the corresponding order scheme.