Modeling Relative Peak Times of Gravitational Wave Harmonics

TL;DR

This work tackles the problem of predicting the relative peak times of gravitational-wave harmonics in binary black-hole mergers. It introduces two semi-analytical approaches—the Backwards-One-Body (BOB) frequency-evolution method and an equatorial geodesic model in the remnant Kerr spacetime—that use only remnant properties $M_f$, $\chi_f$, and mode peak frequencies from Numerical Relativity (NR) to predict mode-peak timings. Across quasi-circular, non-precessing binaries with $l=|m|\le 8$, both methods achieve high accuracy (mean/median errors $\lesssim 1 \; M_i$) and outperform leading IMR models in cases with large timing differences, highlighting their potential to reduce NR reliance in strong-field predictions. The results suggest that relative peak timings are governed by linear dynamics on the Kerr background and may be connected to effective-potential minima, with implications for improving subdominant-mode modeling for next-generation detectors. Limitations include the equatorial, non-precessing focus, and future work will extend to precession and additional modes.

Abstract

Accurate modeling of gravitational waves from binary black hole mergers is essential for extracting their rich physics. A key detail for understanding the physics of mergers is predicting the precise time when the amplitude of the gravitational wave strain peaks, which can differ significantly among the different harmonic modes. We propose two semi-analytical methods to predict these differences using the same three inputs from Numerical Relativity (NR): the remnant mass and spin and the instantaneous frequency of each mode at its peak amplitude. The first method uses the frequency evolution predicted by the Backwards-One-Body model, while the second models the motion of an equatorial timelike geodesic in the remnant black hole spacetime. We compare our models to the SXS waveform catalog for quasi-circular, non-precessing systems and find excellent agreement for $l = |m|$ modes up to $l=8$, with mean and median differences from NR below 1$M$ in nearly all cases across the parameter space. We compare our results to the differences predicted by leading Effective-One-Body and NR surrogate waveform models and find that in cases corresponding to the largest timing differences, our models can provide significant increases in accuracy.

Figures (5)

• Figure 1: Characteristic radii $r_{22}$ calculated using Eq. (\ref{['eq:radius']}) as a function of the radius at which the minimum of the scalar effective potential $W$ occurs for the $(2,2)$ mode. For reference, an exact agreement between the two quantities is shown by a green line. Equal mass configurations are shown with filled red triangles, all other configurations are shown with filled black circles. The kink seen in the upper right may be due to poor NR frequency data, resulting in deviations from the linear trend.
• Figure 3: Comparison of predicted time differences between the peak amplitude of the $(2, 2)$ and corresponding $(l, m)$ modes for equal-mass binaries as a function of the remnant mass. Results from our Geodesic model (brown diamonds), BOB model (green triangles), NRHybSur3dq8 (blue stars), and SEOBNRv5HM (red squares) are compared against SXS NR data (black circles). For $l\,\hbox{$<$}\hbox{$-$}\,5$, we selected systems with initial dimensionless spin $|\chi_i| \,\hbox{$<$}\hbox{$-$} \,0.8$ for consistency with NRHybSur3dq8's training limits. SEOBNRv5HM and NRHybSur3dq8 do not model $l>5$ modes so comparisons to NR are not possible for $l>5$ modes. The error on the SXS, Geodesic, and BOB data points, calculated through comparisons using lower resolution data, are not shown but are all under $1M_i$.
• Figure 4: Difference between NR and the BOB model for the time between the peak amplitude of the $(2, 2)$ and corresponding $(l, m)$ modes across all available configurations as a function of the effective spin $\chi_{\mathrm{eff}}$. The colorbar represents the initial mass ratio $q$ of each configuration.The data is in units of $M_i$ and [$l,l'$] indicates that the difference is being calculated between the $(l,m=l)$ and $(l',m'=l')$ modes.
• Figure 5: Difference between NR and the Geodesic model for the time between the peak amplitude of the $(2, 2)$ and corresponding $(l, m)$ modes across all available configurations as a function of the effective spin $\chi_{\mathrm{eff}}$. The colorbar represents the initial mass ratio $q$ of each configuration. The data is in units of $M_i$ and [$l,l'$] indicates that the difference is being calculated between the $(l,m=l)$ and $(l',m'=l')$ modes.
• Figure 6: Differences between NR and the SEOBNRv5HM (upper) and NRHybSur3dq8 (lower) IMR models for the time between the peak amplitude of the $(2, 2)$ and corresponding $(l, m)$ modes across many configurations as a function of the effective spin $\chi_{\mathrm{eff}}$. Systems are restricted to $|\chi_i|\,\hbox{$<$}\hbox{$-$}\,0.8$ and $q\,\hbox{$<$}\hbox{$-$}\,8$ for consistency with NRHybSur3dq8's training set. The colorbar represents the initial mass ratio $q$ of each configuration.