Convergence, gauge, and memory-mode systematics in AthenaK simulations of the binary black hole merger GW150914

The first AthenaK (GPU-accelerated) binary black hole simulations targeted to GW150914 produced remnant properties and gravitational-wave multipoles extracted at future null infinity, including nonlinear memory effects, but the main results were obtained at a single baseline grid resolution ($N=128$), with only a two-point comparison against one higher-resolution run. The open problem, stated explicitly by the authors, is to complete the validation of AthenaK binary black hole waveforms: perform a systematic multi-resolution convergence study quantifying discretization errors in the waveform multipoles (especially near merger and in memory-dominated modes such as $(\ell,m)=(4,0)$), measured with the superrest-frame convention held fixed across resolutions, and identify and mitigate gauge effects — including the determination of the binary's superrest frame — in comparisons against the SXS and RIT simulations of the same event. No published work to date has delivered this convergence and gauge study for AthenaK binary black hole simulations.

Background

GW150914, the first gravitational-wave detection, was the merger of two black holes observed by LIGO in 2015; accurate waveform modeling of such mergers requires numerical relativity (NR), the direct numerical solution of the full Einstein equations. Different NR codes use different formulations and techniques: the SXS collaboration's SpEC code evolves the generalized-harmonic form of the equations with spectral methods and singularity excision, the RIT code LazEV uses the BSSN formulation (the Baumgarte–Shapiro–Shibata–Nakamura rewriting of the Einstein equations into a standard first-order evolution system) with the moving-puncture gauge (the coordinate treatment that lets finite-difference grids evolve the puncture singularities of merging black holes without excising them) on finite-difference grids, and AthenaK evolves the Z4c formulation (a modification of BSSN with improved constraint propagation) with the moving-puncture gauge on a GPU-accelerated, block-octree adaptive mesh refinement (AMR) framework. Gravitational-wave signals are decomposed into spin-weighted spherical-harmonic multipoles $h_{\ell m}$, with the quadrupolar $(2,2)$ mode dominant. Two extraction approaches exist: finite-radius extraction of the Weyl scalar $\Psi_4$ on coordinate spheres, and Cauchy-characteristic extraction (CCE), which propagates the metric data on a world tube through Einstein's equations to future null infinity $\mathcal{I}^{+}$, yielding strain waveforms free of finite-radius extraction artifacts and containing nonlinear memory effects — permanent, non-decaying offsets in the waveform, visible chiefly in the $m=0$ modes, that would produce a lasting displacement in a detector. Asymptotic waveforms carry residual gauge freedom associated with the BMS (Bondi–Metzner–Sachs) group of asymptotic symmetries; fixing the translations, boosts, and supertranslations determines the superrest frame of the binary, and the choice of frame affects the distribution of signal power among multipoles, particularly the memory-dominated modes. Radice, Gamba, Zhu, and Rashti (2025) simulated GW150914 with AthenaK, computed the remnant mass, spin, and recoil with the isolated-horizon formalism (a quasi-local method for measuring black-hole mass and spin directly from simulation data), extracted waveforms at $\mathcal{I}^{+}$ via CCE, and compared against the highest-resolution SXS (SXS:BBH:305) and RIT (RIT:BBH:0062) simulations of the same event, finding excellent agreement in remnant properties and a maximum $(2,2)$ dephasing of $\Delta\phi\simeq 0.35$ rad with $\Delta A/A\simeq 0.4\%$), and used the waveform to re-analyze the GW150914 detector data with broadly consistent parameter posteriors. However, the baseline simulation used 128 points per direction on the base grid (finest resolution $\Delta=0.0078125\,M$ across 12 AMR levels), a single grid resolution for the main results; the paper's appendix reports one higher-resolution run ($N=192$) whose comparison with the baseline ($\Delta\phi_{22}=0.38$ rad at merger between the two resolutions) showed good agreement for most modes but revealed that the resolution change alters the determined superrest frame, so that memory-dominated modes such as $(\ell,m)=(4,0)$ appear rather different between resolutions. A two-point comparison cannot establish the observed convergence order or provide a per-multipole discretization-error budget, and the authors identify the single resolution and the unquantified gauge effects in the cross-code comparison as the work's main limitations, with the necessary studies stated as ongoing future work.

Problem Statement

Complete the accuracy validation of AthenaK binary black hole simulations targeted to GW150914, as explicitly called for by Radice, Gamba, Zhu, and Rashti (2025), by: (1) performing a systematic convergence study in which simulations at three or more grid resolutions — differing from one another only in grid spacing, and optionally including the existing $N=128$ baseline and $N=192$ appendix run — are evolved through inspiral, merger, and ringdown, so that the observed convergence order and the discretization-error budget are quantified for the waveform multipoles extracted at $\mathcal{I}^{+}$, especially near merger; and (2) investigating and mitigating possible gauge effects in the comparison of AthenaK waveforms with the SXS (SXS:BBH:305) and RIT (RIT:BBH:0062) simulations of the same event, including how the procedure used to determine each binary's superrest frame — which the authors showed changes with resolution and drives the differences in memory-dominated modes such as $(\ell,m)=(4,0)$ — contributes to the residual cross-code dephasing. The common objective is a quantified error budget for AthenaK binary black hole waveforms, including their memory content, sufficient to establish their reliability for gravitational-wave data analysis and for validation against other numerical-relativity codes. This objective comprises the measured convergence order and a per-multipole discretization-error budget over the required resolution set, together with the gauge and superrest-frame analysis of the cross-code comparison.

References

1. D. Radice, R. Gamba, H. Zhu, A. Rashti, 'AthenaK simulations of the binary black hole merger GW150914', Classical and Quantum Gravity 42, 185003 (2025). DOI: 10.1088/1361-6382/adfffa, arXiv:2506.06838, https://arxiv.org/abs/2506.06838
2. C. Musolino, C. Ecker, K. Topolski, M. Cassing, K. Miler, H. Ho-Yin Ng, K. Pierre, E. R. Most, L. Rezzolla, 'GRACE: An Open-Source Framework for GPU-Accelerated Numerical Relativity' (2026). arXiv:2607.09854, https://arxiv.org/abs/2607.09854
3. M.-Z. Han, K. Kiuchi, M. Shibata, 'SACRA-K: A Performance-Portable Numerical Relativity Code with Kokkos' (2026). arXiv:2607.08743, https://arxiv.org/abs/2607.08743
4. H. Zhu, F. Pretorius, J. M. Stone, 'Trapping, Irregular Waveforms, and Efficient Radiation in Ultra-relativistic Black Hole Encounters' (2026). arXiv:2604.26253, https://arxiv.org/abs/2604.26253
5. J. Fields, D. Radice, 'AthenaK Simulations of Magnetized Binary Neutron Star Mergers' (2025). arXiv:2511.06058, https://arxiv.org/abs/2511.06058
6. A. Rashti, R. Gamba, K. Chandra, D. Radice, B. Daszuta, W. Cook, S. Bernuzzi, 'Binary Black Hole Waveforms from High-Resolution GR-Athena++ Simulations' (2024, v3 2025). arXiv:2411.11989, https://arxiv.org/abs/2411.11989
7. E. Shukla, A. Rashti, R. Gamba, D. Radice, K. Chandra, 'GR-Athena++ Simulations of Spinning Binary Black Hole Mergers' (2025, v2 2026). arXiv:2510.13963, https://arxiv.org/abs/2510.13963
8. J. Fontbuté, S. Bernuzzi, S. Albanesi, D. Radice, A. Rashti, W. Cook, B. Daszuta, A. Nagar, 'Covariant and Gauge-invariant Metric-based Gravitational-waves Extraction in Numerical Relativity' (2025). arXiv:2508.03799, https://arxiv.org/abs/2508.03799