On detectability of non-linear modes in ringdown through Parameter Estimation

In general relativity, a perturbed black hole (BH) produced in a binary merger relaxes to the Kerr solution by emitting GWs in a discrete set of characteristic QNMs whose complex frequencies depend solely on the mass $\mathcal{M}$ and angular momentum $\mathcal{J}$ of the remnant BH. The standard picture, where the ringdown is modeled as a linear superposition of exponentially decaying sinusoids, has been very successful; however, given the intrinsically nonlinear nature of Einstein’s equations, it is interesting that the transition from inspiral to merger and ringdown appears so smooth.

Recent studies have provided evidence that nonlinear (quadratic) modes may be excited during merger, even in quasicircular binaries and high-energy head-on collisions. These quadratic modes emerge due to mode-mode coupling in second-order perturbation theory. In essence, the nonlinearity of Einstein’s gravity guarantees that products of first-order (linear) perturbations act as effective sources for second-order modes. Such effects have been predicted by analytical studies and corroborated by numerical simulations. In this project, we develop the tools to inject and recover quadratic modes in ringdown signals, aiming to assess their detectability and explore how their inclusion might enrich black hole spectroscopy.

The gravitational wave signal emitted by a perturbed Kerr BH can be decomposed into its two polarizations using the Teukolsky formalism. The plus and cross polarizations are expressed in terms of the radial Teukolsky solution $R_{\ell m\omega}$ as $h_{+}-ih_{\times}=-\frac{2}{r^4}\int_{-\infty}^{\infty}\frac{d\omega}{\omega^2}\sum\limits_{\ell,m} \ ^{-2}S_{\ell m}(\iota,\beta) R_{\ell m\omega}(r)$, where the asymptotic behavior, $R_{\ell m\omega}\rightarrow r^3\mathcal{Z}_{\ell mn}^{\text{out}}\exp(-i\omega r)$ defines the complex amplitude $\mathcal{Z}_{\ell mn}^{\text{out}}$. Here, $_{-2}S_{\ell m}(\iota,\beta)$ are the spin-2 weighted spheroidal harmonics, which naturally arise as eigenfunctions of the angular part of the Teukolsky equation in the Kerr background. Their dependence on the angular coordinates $(\iota,\beta)$ encapsulates the intrinsic spin effects of the BH, and in the nonspinning limit, these functions reduce to the familiar spin-weighted spherical harmonics, $_{-2}Y_{\ell m}(\theta,\phi)$. In practice, numerical relativity waveform extraction often employs the spherical basis for convenience; however, for spinning BHs, the spheroidal basis more accurately represents the true eigenfunctions of the perturbative operator.

Assuming that the ringdown signal is a linear superposition of exponentially damped sinusoids, the waveform can be written as $h_{+}-ih_{\times}=\frac{1}{r}\sum_{\ell m n} \mathcal{Z}_{\ell mn}^{\text{out}}\exp\bigl(i\omega_{\ell mn}t\bigr)\exp(-t/\tau_{\ell mn}) _{-2}S_{\ell m}(\iota, \beta)$, where $n$ labels the overtone and $t$ represents the retarded time ($t-r$). By decomposing the complex amplitude as $\mathcal{Z}_{\ell mn}^{\text{out}}=M\mathcal{A}_{\ell mn}\exp(i\Phi_{\ell mn})$, we obtain $h_{+}-ih_{\times}=\frac{M}{r}\sum_{\ell m n} \mathcal{A}_{\ell mn}\exp\bigl[i\omega_{\ell mn}t+\Phi_{\ell mn}\bigr]\exp(-t/\tau_{\ell mn}) _{-2}S_{\ell m}(\iota, \beta)$.

The complex frequency $\omega_{\ell mn}$ encodes both the oscillation frequency (real part) and the damping rate (inverse of the imaginary part). The retarded time $t-r$ ensures causality, as the QNMs are only defined in the exterior region of the BH. The use of spin-weighted spheroidal harmonics accurately captures the angular structure of perturbations, particularly for rotating BHs.

Beyond the linear approximation, the nonlinear nature of Einstein’s equations where first-order perturbations will act as effective sources for second-order (quadratic) modes. At this level, the total waveform can be expressed as a sum of linear and quadratic contributions. Schematically, one writes $h_{+}-ih_{\times}=\frac{M}{r}\sum_{k\equiv\ell m n)}\mathcal{A}_{k}\exp\bigl[i\omega_{k}t+\Phi_{k}\bigr] Y_{k}$, where the index $k$ runs over both linear modes and their nonlinear (quadratic) counterparts.

For quadratic modes generated by the coupling of two linear modes denoted $k_1\equiv(\ell_1 m_1 n_1)$ and $k_2\equiv(\ell_2 m_2 n_2)$, the following relations are expected $\omega_{\ell m n} = \omega_{\ell_1 m_1 n_1}+\omega_{\ell_2 m_2 n_2}$, $\frac{1}{\tau_{\ell m n}} = \frac{1}{\tau_{\ell_1 m_1 n_1}}+\frac{1}{\tau_{\ell_2 m_2 n_2}}$, $\mathcal{A}_{\ell m n} \propto \mathcal{A}_{\ell_1 m_1 n_1}\mathcal{A}_{\ell_2 m_2 n_2}$, $\Phi_{\ell m n} = \Phi_{\ell_1 m_1 n_1}+\Phi_{\ell_2 m_2 n_2}+\text{const.}$.

The additive nature of the frequencies reflects that the quadratic mode is driven by a bilinear source term such that its oscillatory behaviour is the combined response of the two parent modes. The damping rate of a quadratic mode is governed by the sum of the inverse damping times of the parent modes, indicating the nonlinearity of the complex frequency. The amplitude scales as the product of the parent amplitudes, which is a direct consequence of the second-order source being quadratic in the first-order perturbations. The phase is the sum of the parent phases (up to an additive constant), ensuring coherence between the coupled modes.

Furthermore, angular momentum conservation imposes strict selection rules on the allowed mode couplings. Only those combinations that satisfy $|\ell_1 - \ell_2| < \ell < \ell_1 + \ell_2 \quad \text{and} \quad |m|=|m_1|+|m_2|$ will contribute. In practice, the precise nonlinear coupling coefficients are determined by Clebsch-Gordan coefficients, which enforce these selection rules and ensure that the overall angular momentum is conserved during the coupling process. Furthermore, the nonlinear source term involves the product $Y_{\ell_1 m_1 n_1}(\theta,\phi) Y_{\ell_2 m_2 n_2}(\theta,\phi)$, which can be decomposed via the Clebsch–Gordan series as $Y_{\ell_1 m_1 n_1} Y_{\ell_2 m_2 n_2} = \sum_{\ell,m,n} \langle \ell_1, m_1, n_1; \ell_2, m_2, n_2 \mid \ell, m, n \rangle Y_{\ell m n}$. The corresponding Clebsch–Gordan coefficient is obtained from the overlap integral $\int d\Omega \ Y_{\ell_1m_1n_1}(\theta,\phi)Y_{\ell_2m_2n_2}(\theta,\phi)Y^{*}_{\ell_1+\ell_2, m_1+m_2, n_1+n_2}(\theta,\phi)$ which quantifies the coupling strength between the parent modes.

Additionally, when solving the inhomogeneous second-order perturbation equation, the Green’s function enters the integration over the radial coordinate, providing the dynamical response of the system. Both the angular overlap factor encoded by the Clebsch–Gordan coefficient and the radial integration via Green’s function are absorbed into the overall coefficient of proportionality. This coefficient thereby relates the quadratic mode amplitude to the product of the linear mode amplitudes, ensuring that the nonlinear coupling is fully captured by our model.

Drafting is currently in progress. Sorry for any inconvenience. Please check back soon!