BSSN Evolution Systems

BSSN rewrites the spatial metric and extrinsic curvature into conformal variables whose numerical evolution is substantially more robust than a direct ADM update in many black-hole calculations. That empirical success depends on gauge, constraint additions, derivative stencils, dissipation, and boundary treatment; “using BSSN” does not by itself define an algorithm. This page specifies the evolved fields, numerical controls, standard test sequence, and invariant checks needed to support a finite evolution claim. Continuum formulation and hyperbolicity remain with Formulations of General Relativity.

Helpful background. Formulations of General Relativity develops the continuum evolution systems, Ordinary Differential Equations supplies time integration and stability, Computational Relativity supplies the volume-wide evidence pattern, and Black-Hole Dynamics and Collapse supplies the physical evolution problems.

The evolved conformal variables

Use $K_{ij}=-\tfrac12\mathcal L_n\gamma_{ij}$. In Cartesian spatial coordinates, a common BSSN choice is

$$\gamma_{ij}=e^{4\phi}\tilde\gamma_{ij},\qquad K_{ij}=e^{4\phi}\left(\tilde A_{ij}+\frac13\tilde\gamma_{ij}K\right),$$

with $\det\tilde\gamma=1$, $\tilde\gamma^{ij}\tilde A_{ij}=0$, and conformal connection functions $\tilde\Gamma^i=\tilde\gamma^{jk}\tilde\Gamma^i{}_{jk}=-\partial_j\tilde\gamma^{ij}$. The final equality uses unit conformal determinant; curvilinear implementations need a consistent reference-metric or tensor-density formulation. Some codes evolve $\chi=e^{-4\phi}$ or $W=e^{-2\phi}$ near punctures. The variable choice, enforcement of algebraic constraints, placement of momentum-constraint substitutions, and advective derivative convention must be recorded.

The conformal curvature here differs from the tensor used in Conformal Decompositions of Initial Data: BSSN uses $\tilde A_{ij}=\psi^{-4}A_{ij}$ with $\psi=e^\phi$, whereas the initial-data page’s CTT tensor obeys $A^{ij}=\psi^{-10}\tilde A_{\rm CTT}^{ij}$. Raising the BSSN tensor with the conformal metric gives $\tilde A_{\rm BSSN}^{ij}=\psi^{-6}\tilde A_{\rm CTT}^{ij}$. Passing a CTT array directly into a BSSN variable without this conversion changes the physical extrinsic curvature.

Derive one update and test its meaning

The metric equation $\partial_t\gamma_{ij}=-2\alpha K_{ij}+\mathcal L_\beta\gamma_{ij}$ implies

$$\partial_t\ln\gamma=-2\alpha K+\beta^k\partial_k\ln\gamma +2\partial_k\beta^k.$$

Since $\phi=\ln\gamma/12$, the conformal-factor update is

$$(\partial_t-\beta^k\partial_k)\phi =-\frac16\alpha K+\frac16\partial_k\beta^k.$$

The shift-divergence term is essential even in flat spacetime: it accounts for a changing coordinate volume. In zero shift, positive $K$ contracts that volume and decreases $\phi$. Deriving this scalar equation is a useful first implementation check before the longer curvature and connection equations; see Baumgarte and Shapiro 2010, §11.5, Eqs. (11.32)–(11.40), pp. 387–388. This calculation explains one component of BSSN, not a complete evolution system.

Spatial derivatives may use centered finite differences plus upwinded shift advection, finite elements, or high-order dissipation. A method-of-lines integrator supplies time evolution under a CFL condition. Baumgarte and Shapiro 1999, §§II–IV and Figs. 1–3 compares ADM and conformal-traceless evolution and identifies the numerical role of the reformulation. Its results do not replace a test of a modern implementation.

Constraints are evolved diagnostics

Monitor the Hamiltonian $\mathcal H$, momentum $\mathcal M_i$, connection constraint $\mathcal G^i=\tilde\Gamma^i+\partial_j\tilde\gamma^{ij}$, determinant, and trace-free condition. Recompute them with an independent derivative operator. Algebraically resetting $\det\tilde\gamma$ and the trace can control drift, but the reset magnitude must be output; otherwise it conceals how strongly the unconstrained update violates the chosen variable manifold.

Use volume-weighted $L_2$ and $L_\infty$ norms, refinement-scaled plots, and localized slices through difficult regions. Constraint reduction alone is not sufficient: lapse collapse can freeze coordinate fields while horizon mass or curvature invariants drift. Track invariant or quasi-local observables alongside coordinate components.

A staged benchmark suite

Begin with robust stability: small random perturbations of flat spacetime test high-frequency growth. Continue with a linear wave, gauge wave, shifted gauge wave, and expanding Gowdy spacetime. These tests exercise propagation, gauge, advection, and strong-field dynamics in increasing order. Daverio, Dirian, and Mitsou 2018, §§4–8 and Figs. 21–23 compares conformal formulations with the Apples-with-Apples tests and common measures.

For a smooth exact test, use at least three resolutions and demonstrate the designed convergence order before boundary contamination. For robust-stability noise tests, retain perturbations at the grid scale and follow the chosen benchmark’s resolution-dependent amplitude and norm definitions; otherwise the newly resolved high-frequency modes are never tested. Record seeds and compare ensembles so a different random draw does not masquerade as convergence. A fixed physical-spectrum perturbation is a useful separate control, not a replacement for the grid-scale stress test. Vary dissipation and CFL independently; stability that exists only under excessive damping is not accuracy.

Then evolve a stationary Schwarzschild or Kerr black hole with declared slicing and gauge. Monitor horizon area, spin, constraints, and any extracted radiation expected to vanish. Move the outer boundary and distinguish startup gauge settling from persistent drift. This strong-field benchmark connects the testbed to the intended application without pretending to validate a binary merger.

The grid hierarchy must resolve both the near-hole gradients and the shortest scientifically relevant wavelength in the test. Repeat the stationary case with shifted refinement boundaries and, where possible, a second coordinate representation. Agreement of horizon mass under that change distinguishes a geometric result from fortuitous cancellation on one mesh. Archive time series rather than only final norms, because transient exponential growth can disappear after a coarse output cadence or a late gauge freeze.

Environment and rerun record

Pin the exact BSSN equations, variable choice, gauge and parameters, constraint substitutions, algebraic projection, stencil and advection bias, dissipation, Runge–Kutta scheme, CFL factor, precision, grid hierarchy, boundary closure, test-data revision, and code commit. A suitable interface is evolve-bssn --config gauge-wave.toml --output run/, emitting fields, all constraints, projection corrections, invariants, convergence differences, timings, and environment metadata. This is the record an implementation should provide, not a claim that code was run for this prose.

Reject an evolution with resolution-accelerating constraint growth, loss of convergence before the declared interval, unexplained dependence on dissipation, boundary reflections, or stable coordinate fields accompanied by drifting horizon invariants. Further failure modes include inconsistent advective signs, under-resolved puncture gradients, applying a periodic benchmark with nonperiodic boundary closures, and comparing grids after arbitrary coordinate alignment.

Diagnostic

A gauge-wave test remains bounded for $1000$ crossings when strong dissipation is used. The wave amplitude decays by half, and this decay does not decrease under refinement. Has BSSN passed the test?

Solution

No. Boundedness does not establish convergence to the intended solution. Check the dissipation operator, its powers of grid spacing, and its dimensionless coefficient; correctly scaled high-order artificial dissipation should vanish on a fixed smooth solution as the grid is refined. Repeat at several dissipation strengths and report waveform and constraint errors. A fixed dimensionless dissipation coefficient is not itself an error; the full operator scaling determines the continuum effect.

References

• Baumgarte, T. W., and Shapiro, S. L. (1999). “On the numerical integration of Einstein’s field equations.” Physical Review D 59, 024007. arXiv.
• Baumgarte, T. W., and Shapiro, S. L. (2010). Numerical Relativity: Solving Einstein’s Equations on the Computer. Cambridge University Press. DOI.
• Daverio, D., Dirian, Y., and Mitsou, E. (2018). “Apples with Apples comparison of 3+1 conformal numerical relativity schemes.” arXiv:1810.12346. arXiv.