Priorities in gravitational waveforms for future space-borne detectors: vacuum accuracy or environment?

ABSTRACT

In preparation for future space-borne gravitational-wave (GW) detectors, should the modelling effort focus on high-precision vacuum templates or on the astrophysical environment of the sources? We perform a systematic comparison of the phase contributions caused by (1) known environmental effects in both gaseous and stellar matter backgrounds, or (2) high-order post-Newtonian (PN) terms in the evolution of mHz GW sources during the inspiral stage of massive binaries. We use the accuracy of currently available analytical waveform models as a benchmark value, finding the following trends: the largest unmodelled phase contributions are likely environmental rather than PN for binaries lighter than ∼10⁷/(1 + z)² M_⊙, where z is the redshift. Binaries heavier than ∼10⁸/(1 + z) M_⊙ do not require more accurate inspiral waveforms due to low signal-to-noise ratios (SNRs). For high-SNR sources, environmental phase contributions are relevant at low redshift, while high-order vacuum templates are required at z ≳ 4. Led by these findings, we argue that including environmental effects in waveform models should be prioritized in order to maximize the science yield of future mHz detectors.

1 CHALLENGES AND OPPORTUNITIES IN THE MILLIHERTZ BAND

From the early post-Newtonian (PN) results (see e.g. Damour & Deruelle 1985; Damour, Soffel & Xu 1991; Jaranowski & Schäfer 1998, 1999), formalisms such as the effective-one-body (see e.g. Buonanno & Damour 1999; Damour 2001; Barack, Damour & Sago 2010; Akcay et al. 2012; Pan et al. 2014; Ossokine et al. 2020) and self-force perturbation theory (see e.g. Teukolsky 1973; Quinn & Wald 1997; Burko 2004; Gralla & Wald 2011; Lackeos & Burko 2012), combined with advances in numerical relativity (see e.g. Pretorius 2005; Ajith et al. 2011; Pürrer 2014; Szilágyi et al. 2015), have pushed the validity of analytical waveform templates to higher and higher orders. In combination with Bayesian inference techniques, they have made it possible to extract an impressive amount of information from gravitational-wave (GW) events detected by the Laser Interferometer GW Observatory (LIGO; see e.g. Abbott et al. 2019; Chatziioannou et al. 2019; Romero-Shaw et al. 2020; Islam et al. 2021). Now, the promise of space-borne mHz detectors in the early 2030s is breathing new life into the waveform modelling effort. Sources of GWs in the mHz band present ulterior challenges with respect to the ones routinely detected by ground-based observatories. First, crucial parameters such as mass and mass ratio can vary by several orders of magnitude. Secondly, signal-to-noise ratios (SNRs) of several thousands are expected for appropriate sources (see e.g. Amaro-Seoane et al. 2017). These effects combined greatly increase the accuracy and breadth demands of the required waveform templates, motivating several research programmes with the aim of developing more sophisticated analytical approximations as well as expanding the parameter space in which numerical relativity is viable (see e.g. Lousto & Healy 2022; Nagar et al. 2022, for some recent work). As an interesting recent example, template formalisms based on the PN and post-Minkowskian series have benefited from a recent influx of particle physics techniques (see e.g. Bern et al. 2019; Mogull, Plefka & Steinhoff 2021; Buonanno et al. 2022), which are rapidly succeeding in computing high-order terms.

A further, crucial difference between Hz and mHz binary sources is that the latter are more likely to be affected by their astrophysical environment. The presence of gas and massive third bodies can influence the source’s evolution within the mHz band, confounding expectations based on vacuum templates, limiting the effectiveness of parameter estimation (and thus tests of general relativity), and introducing spurious biases (see e.g. Barausse & Rezzolla 2008; Gair et al. 2013; Barausse, Cardoso & Pani 2014; Caputo et al. 2020; Chen, Xuan & Peng 2020). Conversely, detecting these deviations represents a unique opportunity to measure properties of the source’s environment (from constraining accretion disc physics to detecting dark matter and exoplanets; see e.g. Chakrabarti 1993; Inayoshi et al. 2017; Derdzinski et al. 2019, 2021; Tamanini & Danielski 2019; Cardoso & Maselli ; Coogan et al. 2022; Speri et al. 2022; Zwick et al. 2022). A question is often posed in the context of the scientific groundwork required for missions such as the Laser Interferometer Space Antenna (LISA; Barack et al. 2019; Thorpe et al. 2019; Amaro-Seoane et al. 2022) and TianQin (Luo et al. 2016; Mei et al. 2021): what is the relative importance of the environment with respect to the vacuum evolution of a GW source? In this work, we aim to to provide a simple but general framework to estimate the relative importance of environmental effects without needing to specify any particular waveform model. We apply this framework to analyse massive black hole (BH) binaries of varying mass and mass ratio, since they are often considered to be largely unaffected by their environments while radiating GWs in the mHz band.

2 COMPARING THE PHASE OF GRAVITATIONAL WAVES

2.1 Vacuum waveforms

The ability of a waveform template to accurately match the phase of a real signal is a crucial benchmark with regards to its accuracy. Over the course of an observation, the total phase ϕ of a GW is determined by the source’s frequency evolution, |$\dot{f}$|⁠. It reads (see e.g. Cutler & Flanagan 1994)

where f is the observed GW frequency, f_r = f(1 + z) is the rest frame frequency at a redshift z, f_min is the frequency at which the GW source enters a detector band, and f_max is the maximum frequency reached within the observation window. Sophisticated waveform models expand upon the historical result derived in Peters & Mathews (1963) by including contributions from higher-order GW modes as well as PN corrections in the frequency evolution (see e.g. Blanchet 2014):

where M is the binary’s total mass, q its mass ratio, v its characteristic orbital velocity, c the speed of light in vacuum, and G Newton’s constant, and we assume quasi-circular orbits. Corrections to the leading order expression proportional to powers of (v/c)²ⁿ correspond to the nth relative PN order. For the purposes of our phenomenological analysis of environmental effects, we make the two following assumptions regarding vacuum waveforms:

• The phase accuracy of a waveform template is a proxy for its capacity to recover source parameters without bias. We thus neglect the merger and ring-down contributions to the SNR of a GW signal.¹

• The phase accuracy of any waveform template can be translated into an equivalent accuracy in terms of PN orders, regardless of the original formalism used to construct the template.

Both of these assumptions are motivated throughout the inspiral phase of the binary source, where environmental effects are most likely to be significant (see e.g. Blanchet 2014). Therefore, we deem a waveform template to be accurate to the nth PN order if it can match the phase of a true signal up to an error smaller than the next PN contribution ∼(v/c)^{2n + 1}, as defined, e.g. by equations (2) and (3). In this language, current state-of-the-art waveform models range in their phasing accuracy depending on several simplifying assumptions such as spin alignment, mass ratio, or lack of eccentricity, going as high as ∼20-PN (see e.g. Fujita 2015; Munna 2020). In the case of a generic LISA source with arbitrary spins, moderate mass ratios (1 ≲ q ≲ 0.01), and small eccentricities, purely analytical methods have achieved 4-PN accuracy, while 5-PN accuracy has been achieved by calibrating against numerical relativity simulations. A general overview of recent literature leads us to set a reasonable benchmark of 5-PN phasing precision to be the current standard for waveform templates (based on several works; see e.g. Messina et al. 2019; Accettulli Huber et al. 2021; Cho, Pardo & Porto 2021; Khalil 2021; Blümlein et al. 2022; Chattaraj et al. 2022; Cho, Porto & Yang 2022; Nagar et al. 2022). While such a benchmark is arbitrary (and is bound to change in the following years), it will serve as a useful comparison tool to assess the current state of the field.

In order to mimic vacuum waveform templates of arbitrary precision, we model the phase evolution of a source using a phenomenological form for equation (2), based on the PN series as well as dimensional arguments (see also Zwick et al. 2021). It reads

where the exponent j denotes the (n/2)th effective PN order and we arbitrarily set all dimensionless coefficients A_j = 1, thus only preserving the physical (v/c) scaling information (see e.g. Cardoso & Maselli 2020b). Note that, while equation (3) is only a crude simplification, there are several fundamental uncertainties in the modelling of environmental effects which will overshadow any loss of precision due to setting A_j = 1. Note also that we are assuming that the phase evolution can be well described by orbit-averaged equations, losing information regarding the initial true anomaly of the source, an important parameter in many waveform models.

Even in the case of vacuum sources, the accuracy required to extract a maximal amount of information is bounded by the inherent SNR limitations of GW detectors. A simple approximation of the SNR of a GW event can be found by using the following formula (see e.g. Klein et al. 2016):

where S_t is the detector’s power spectral density and h_c the source’s characteristic strain. For the former, we take the LISA specifications as reported in Robson, Cornish & Liu (2019). The latter reads

where D_l is the source’s luminosity distance and N_c is the number of cycles it completes at a given rest frame frequency. For an observation window of 4 yr, the latter reads

In order to define whether a small phase shift δϕ is distinguishable from noise, we consider the following criterion, commonly used in more qualitative treatments of environmental effects on GWs (see e.g. Kocsis, Yunes & Loeb 2011):

where δSNR is a detectability threshold, customarily chosen to be equal to 8. With the aid of equation (4) and assuming a constant phase shift, equation (7) can be rearranged into a useful form,

which essentially states that the phase of a GW signal can be reconstructed with an accuracy of ∼2π/SNR (also see e.g. Katz et al. 2021). Phenomenologically, waveform templates are therefore only required to achieve an accuracy comparable to the limit imposed by equation (8), since any smaller contribution to the phase of the real signal would be washed away by noise. Adopting a more accurate criterion would likely decrease the sensitivity to account for the phase’s accumulation rate along with degeneracies and other subtleties of more accurate GW data analysis. Note also, that equations (4)–(6) are technically only valid at Newtonian order (Mangiagli et al. 2019). However, they suffice for the purposes of this work as they only serve to produce a reference SNR value, which we apply as a detectability criterion equally for both PN and environmental phase shifts: a crude but fair comparison in line with the phenomenological form of equation (3) and the intrinsic uncertainties of astrophysical environmental effects.

2.2 A sample of environmental effects

Environmental influences can introduce additional terms that modify the frequency evolution of a source of GWs. In this work, we consider a minimal model for three types of environmental effects that are considered typical in the astrophysical setting of LISA sources, i.e. we only consider a simple phase contribution. Richer environmental signatures can also be produced (see e.g. Zwick et al. 2022 and Cardoso et al. 2022 for gas-embedded binaries, or Torres-Orjuela et al. 2019; Torres-Orjuela, Chen & Amaro Seoane 2021 for sources with a peculiar velocity), but are beyond the scope of this work. The influence of additional environments such as dark matter or other baryonic fields has also been explored in the literature (see e.g. Eda et al. 2013; Macedo et al. 2013; Cardoso et al. 2016; Baumann et al. 2022; Cole et al. 2022).

Gas torques can act on a GW source by transferring energy and angular momentum between the binary and the surrounding gas. The presence of gas is likely in the case of supermassive BH (SMBH) LISA sources, since the galaxy mergers responsible for the binary itself can trigger large inflows of gas towards the central regions, because of both tidal torques (due to gravitational forces; see e.g. Barnes & Hernquist 1996; Hopkins & Quataert 2010; Capelo et al. 2015) and hydrodynamical torques (e.g. large-scale ram-pressure shocks; Barnes 2002; Capelo & Dotti 2017; Blumenthal & Barnes 2018). Furthermore, gas is thought to be one of the key agents that can aid binary hardening below pc scales (see e.g. Souza Lima et al. 2020). For sources embedded in a circumbinary disc, gas torques are well described as resulting from viscous forces (see e.g. D’Orazio & Duffell 2021):

where |$\dot{M}$| is the gas accretion rate on to the binary and A is a dimensionless pre-factor that depends on disc properties. Following D’Orazio & Duffell (2021), we adopt a value of A = 3, which is appropriate for binaries with 10⁻² < q < 1. We scale the mass accretion rate with the Eddington (1916) limit:

where f_Edd is the Eddington ratio, m_P the proton mass, σ_T the Thomson cross section, and η = 0.1 the radiative efficiency (see e.g. Marconi et al. 2004). Typical values for f_Edd in active galactic nuclei at z ∼ 2 range from 10⁻³ to 1, peaked around the commonly assumed value of ∼10⁻¹ (see e.g. Suh et al. 2015). Larger Eddington ratios are thought to be common at higher redshifts (Willott et al. 2010).

The gravitational influence of a third body can induce tidal fields and/or linear accelerations that affect the orbit of the inner binary and can produce subtle modifications to the source’s GW emission (see e.g. Bonetti et al. 2017; Randall & Xianyu 2019; Torres-Orjuela et al. 2019, 2021; Liu et al. 2022; Xuan, Naoz & Chen 2022). As in the rest of this work, we will limit our analysis to the lowest-order effect on the inner binary’s frequency evolution and also assume that the latter is the sole source of GWs. In a hierarchical triplet, the average effect of the third body is to induce a tidal term E_t which modifies the inner binary’s binding energy E_b (see e.g. Will 2021, for a recent PN treatment of the three-body problem). Dimensionally, the ratio between the energies reads

where a is the inner binary’s separation, r ≫ a is the distance to the third body, m₃ its mass, and B is a dimensionless pre-factor of order unity that depends on the configuration of the system, which we set equal to one. As shown in Zwick et al. (2021), the tidal term affects the inspiral rate of the inner binary by a factor proportional to E_t/E_b, modifying its frequency evolution:

A third body of mass m₃ at a distance r also induces an acceleration of the inner binary’s centre of mass, which can produce a time-dependent change in the peculiar velocity of the GW source. While a constant Doppler shift is degenerate with redshift, a peculiar acceleration along the line of sight can cause a time variation in the Doppler shift which, if integrated over an observation time, causes a non-degenerate shift in the source’s GW phase. We model this effect by Doppler shifting the frequency by an amount proportional to the time-dependent line-of-sight velocity V_p(t) which is solely caused by the third body:

where V_p(t) = (Gm₃/r²)t. We can then integrate the evolution equations as usual to obtain an accumulated, non-degenerate dephasing caused by the centre of mass acceleration. We parametrize this effect by considering the line-of-sight velocity, σ_v, reached after a time T of acceleration:

In our calculations, σ_v therefore represents the maximum value of V_p(t) reached within an observation window. Typical values for σ_v are

Both unmodelled velocities of the order of km s⁻¹ and tidal deformations of the order E_t/E_b ∼ 10⁻⁶ can be produced in an astrophysical context: the presence of a heavy star or stellar-mass BH in the innermost ∼10⁻² to ∼10⁻⁴ pc of a nuclear cluster is a likely consequence of relaxation and mass segregation (see e.g. Bahcall & Wolf 1977; Alexander & Hopman 2009; Linial & Sari 2022) or the presence of a large-scale accretion disc (see e.g. Goodman & Tan 2004; Levin 2007; Bellovary et al. 2016). Similar values can be produced by pc-scale SMBH triplets, another channel that can produce hard SMBH binaries despite possible angular momentum barriers or other delays (Bonetti et al. 2019).

2.3 Methodology

Having set up our models for both vacuum waveform templates and environmental effects, we devise a simple strategy to compare the importance of typical environmental perturbations to the accuracy of an arbitrarily precise vacuum template. We compare contributions to the total phase of the GW signal produced by additional PN orders and environmental effects by means of equation (1) and the various forms of equation (2) discussed in the previous section. Every GW source is integrated from the time it enters the LISA band for a 4-yr period or, if it occurs first, until it reaches a separation of 12GM/c². An example of our computations can be seen in Fig. 1, in which we plot PN contributions to the GW phase for a range of possible sources, as well as the detectability criterion defined by equation (8).

We term an environmental perturbation to be of ‘comparative nth PN order’ if it produces a phase shift comparable to the nth PN correction in equation (3). In Fig. 2, we show the comparative PN order of selected environmental effects for two particular sources of GWs at z = 2 with a binary’s total mass of 10⁵ M_⊙ and a mass ratio q = 1/3 and 1/10, as a function of the dimensionless variables f_Edd, σ_v/c, and E_t/E_b. The comparative PN order of environmental effects strongly depends on the value of these parameters. For realistic values of f_Edd ∼ 0.1 and σ_v/c ∼ E_t/E_b ∼ 10⁻⁶, it tends to be between the 4th and the 6th PN order, comparable to the benchmark 5-PN precision of available waveform templates. Note how very recent work on PN waveforms (Owen et al. 2023) has determined that the truncation of 5PN terms will lead to a systematic bias in the parameter estimation of generic vacuum sources. By analogy, mis-modelling environmental effects of comparative 5PN order is expected to lead to similar bias.

3 VACUUM OR ENVIRONMENT?

3.1 Comparable-mass mergers

We apply our methodology on BH binaries with total masses between 10³ and 10⁸ M_⊙, mass ratios between 10⁻² and 1, and redshift between 0 and 5. Sources with an SNR < 8 according to equation (4) are automatically discarded. The results of our analysis are visualized in Fig. 3, in which we assume representative values of f_Edd = 0.1, σ_v = 1 km s⁻¹, and E_t/E_b = 10⁻⁶. With the aid of equation (7), we show the detectability regions of environmental effects, represented by the brown (detectable) and blue (undetectable) areas. The red contours denote the maximum detectable PN order consistent with SNR limitations, also according to equation (7). The regions highlighted in gray are defined by enforcing the following criteria: (1) environmental effects are not detectable; and (2) phase contributions due to PN orders higher than the 5 PN benchmark are detectable. The significant trends in the figure can be roughly summarized as follows:

• In unequal-mass binaries with total masses of the order of ≲10⁷/(1 + z)² M_⊙, environmental effects are likely to produce a larger phase contribution than the benchmark 5 PN precision.

• The accumulated phase of heavy binaries (≳10⁸/(1 + z) M_⊙) can be adequately modelled with available waveforms due to their lower SNRs and total accumulated phase.

• Sources with the highest SNRs land around a total mass of ∼10⁶ M_⊙ and mass ratio of ≲0.2. In this range, environmental phase contributions are relevant at low redshift (z ≲ 2), while high-order vacuum phase contributions dominate at z ≳ 4.

Varying the choice of f_Edd, σ_v, or E_t/E_b strongly affects the results, as suggested by Fig. 2. For example, a larger Eddington fraction implies detectable phase contributions even at higher redshift (see Garg et al. 2022, for a detailed analysis).

The results of Fig. 3 must be interpreted in light of the expected merger rates of BH binaries. Several estimates suggest that lighter sources (≲few 10⁵ M_⊙) at relatively low redshift (≲3) are expected to dominate event rates, although significant variation can be caused by changing seeding prescriptions (see e.g. Rhook & Wyithe 2005; Sesana et al. 2005; Sesana, Volonteri & Haardt 2007; Volonteri et al. 2020). If these estimates are accurate, they imply that a majority of LISA sources will likely fall in the brown regions of Fig. 3, strongly suggesting that the study of environmental effects should also become priority when it comes to LISA SMBH binaries.

3.2 A comment on extreme mass-ratio inspirals

A second class of sources requiring complex waveforms are extreme mass-ratio inspirals (EMRIs). In the standard dynamical (dry) formation channel, stellar-mass BHs are scattered on to low peri-apsis orbits around SMBHs by relaxation processes (see e.g. Alexander & Hopman 2003; Amaro-Seoane, Freitag & Spurzem 2004; Merritt 2013; Vázquez-Aceves et al. 2022) and are able to complete ∼10⁵ orbits only a few Schwarzschild radii above the event horizon of the primary BH. Because of their high initial eccentricity and rapid circularization, EMRIs can ‘skip’ over low-frequency GW emission: they enter the LISA band directly in a regime where environmental perturbations are completely negligible.

However, an alternative formation pathway can potentially produce similar expected EMRI rates. In the so-called active galactic nucleus (wet) channel, compact objects align or form within the central SMBH’s accretion disc (see e.g. Rauch 1995; Pan, Lyu & Yang 2021; Derdzinski & Mayer 2022). They are dragged into the LISA band by means of gas torques, producing low-eccentricity EMRIs. They will likely enter the LISA band at lower frequencies, damping out high-order PN effects and reducing the overall SNR. While this might be a blow to precision general-relativity measurements, gas-embedded EMRIs could be powerful probes of accretion disc physics, enabling measurements inaccessible to electromagnetic instruments (see e.g. Chakrabarti (see e.g. Chakrabarti 1993; Ryan 1995; Levin 2007; Barausse & Rezzolla 2008; Kocsis et al. 2011; Barausse et al. 2014; Derdzinski et al. 2019, 2021; Destounis et al. 2022; Polcar, Lukes-Gerakopoulos & Witzany 2022; Speri et al. 2022; Zwick et al. 2022).

Because of the high uncertainties in the event rates, a choice between prioritising vacuum templates is premature. In any case, the detection of EMRIs is likely to provide a wealth of scientific observations, be it in fundamental physics (see e.g. Babak et al. 2017; Arun et al. 2022) or astrophysics.

4 SETTING THE PRIORITIES

In short, our analysis suggests that the majority of massive binary mHz GW sources are likely to populate the regions of parameter space in which (1) currently available waveforms can already precisely match the phase of realistic signals if SNR limitations are considered, or (2) environmental phase contributions are likely dominant over vacuum corrections of order higher than the current benchmark value of 5 PN.

The aim of developing high PN order templates is to perform the precision measurements required to test fundamental physics, in particular modifications to general relativity. Unless the environment is properly modelled, these measurements will most likely require special sources (denoted by the grey areas in Fig. 3) in which the SNR is sufficiently high and the environment is negligible. Additionally, most of the astrophysical information regarding the BH binary itself (mass, mass ratio, component spin, and orientation) can, in principle, be recovered with lower PN order phenomenological templates, since even complex dynamics such as spin-orbit and spin-spin couplings produce phase shifts at much less than 5th PN order. On the other hand, further development on the environmental side would (1) reduce the possibility of unknown biases and (2) allow to constrain the astrophysical surrounding of GW sources. The latter is especially promising if richer environmental signatures are also taken into account.

Led by these considerations, we argue that systematically including environmental effects in waveform templates should take priority with respect to further increasing the accuracy of inspiral vacuum templates. If the goal is to maximize the science yield of future missions, the community could be better served by shifting the focus from the source of GWs to its surroundings.

Acknowledgement

The authors acknowledge support from the Swiss National Science Foundation under the Grant 200020_192092.

DATA AVAILABILITY

The data underlying this article will be shared on reasonable request to the authors.

Footnotes

References