arbitrary and the cored isothermal sphere (CIS, Sec. III.2), with a finite-
density core of size $x_{c}$. For each lens we presented the physical scales,
GO structure (images, caustics, magnifications), bGO corrections, and WO
predictions, including the convergence to the SIS limit. We mostly focused on
the strong-lensing regime with mild magnification, although we showed examples
of WO signatures in the single-image regime (weak lensing) and for sources
near the caustics (Figs. 8, 9).

A distinctive feature of these lenses is the existence of a central image
associated with the maximum of the Fermat potential. Central images produce a
secondary modulation in the amplification factor, encoding additional
information about the lens. For sources near a caustic, the additional
modulation of the amplification becomes very broad in frequency space and
pushes convergence to GO towards very high dimensionless frequency $w\gg
1/\text{min}\,|\phi_{I}-\phi_{J}|$. The central image vanishes for SIS and
narrow gSIS profiles ($k\leq 1$): its detection is thus a direct probe of the
inner structure of halos. We discussed qualitative aspects and observational
prospects (Sec. III.3), including parameter estimation in the $w\gg 1$ limit,
the prospect of detecting central images and mismatches between GWs lensed by
gSIS and CIS, relative to SIS.

We then addressed the potential for GW detectors to probe lens features (Sec.
IV). We used a Fisher-matrix analysis to forecast the ability of LISA and LIGO
to reconstruct the parameters of the lens (mass, impact parameter, slope/core
size) and of the source (signal amplitude and initial phase), as well as their
correlations. We focused on equal-mass, non-spinning sources and kept their
parameters fixed (except the amplitude and overall phase). While our analysis
is simplified, it produces reasonable agreement with more detailed treatments
of the source and the detector [36]. The large degeneracy between lens
parameters makes the Fisher matrix very sensitive to numerical errors in the
computation. To circumvent this issue, we used approximations to the Fisher
matrix at low/high frequencies and focused on typical strong-lensing
situations (away from caustic) where convergence to GO occurs at relatively
low dimensionless frequencies. We also addressed the validity of the linear
signal approximation (LSA) in different cases.

LISA and LIGO can identify the CIS lens parameters with precision
$\sim\mathcal{O}(1/\text{SNR})$ for large $w_{\rm chirp}\propto M_{Lz}/M_{\rm
BBH}$ (GO limit). The precision is somewhat higher for the core size and lower
for the lens mass. The marginalised constraints saturate in the GO limit
because of the large degeneracies between lens parameters, although the 2D
posteriors become thinner as $M_{Lz}$ increases. This behaviour is generic in
strongly-lensed systems with a central image (three images total), allowing to
recover up to four parameters in GO (two time delays with high accuracy $\sim
1/M_{Lz}$ and two magnification ratios $\sim 1/\text{SNR}$). The constraints
degrade at lower lens masses (such that the merger occurs below the onset of
magnification), due to both parameter degeneracies and diffraction suppressing
the amplification factor.

Wave optics can substantially improve the constraints in regimes where GO does
not provide enough information. We saw this explicitly for the narrow gSIS
lens, which only produces two GO images. In this case, the most stringent
constraints are found when the onset of magnification occurs close to the
merger, $w_{\rm chirp}\propto M_{Lz}/M_{\rm BBH}\sim 1-10$. The constraints
are driven by WO and degrade for increasing $M_{Lz}$, as the bGO corrections
become negligible. The study of the gSIS revealed that the LSA overestimates
the sensitivity at high $M_{Lz}$ due to the contribution from the lens’ cusp.
For this reason, in addition to full WO we quoted bGO-only results, which
neglect the cusp contribution and are thus conservative for $w_{\rm chirp}\gg
1$. Our Fisher matrix results are supported by the mismatch analysis (Fig.
11), which does not rely on the LSA. Given the limitations we uncovered for
the Fisher analysis on lensed signals, it would be interesting to perform more
sophisticated forecasts by sampling the full likelihood and better asses past
results obtained in the literature (e.g. [4]).

Identifying lens features at percent and sub-percent precision (for LIGO and
LISA respectively) opens up novel applications of GW lensing. While GW lensing
measures projected quantities, limits on the lens mass and physical scales can
be obtained under reasonable assumptions regarding the cosmology and the halo
mass function (Sec. V.1 and Fig. 21). In particular $M_{\rm vir}$, $(\xi_{0}$)
has a hard lower (upper) limit, and their 95% credible regions can be obtained
within a factor of a few, except for very high redshifts. These limitations
may be lifted if the lens redshift can be identified, e.g. via EM follow-up
[112, 117].

The possibility of simultaneously constraining the lens mass and core size
enables tests of Dark Matter scenarios. We addressed the prospects of probing
self-interacting and ultra-light DM models (Sec. V.2 and V.3, respectively),
which predict cored halos whose size depends on the properties of DM and the
halo mass. Deriving simple expressions for the projected core size in these
models, we translated the sensitivity from our forecasts into the parameter
space accessible to strongly lensed GWs. GWs can probe SIDM cross sections as
low as $10^{-4}-10^{-5}{\rm cm}^{2}/{\rm g}$, with the maximum sensitivity
achieved around the onset of magnification (Fig. 22). For ULDM, the core sizes
of DM halos can be tested if $m_{\phi}\lesssim 10^{-17}\,{\rm eV}$ (cf. Fig.
23). These bounds assume the largest possible projected core size for a given
source and are thus conservative within our model. Nonetheless, the assumption
of a $\rho\propto r^{-2}$ profile might give overly optimistic constraints for
SIDM, compared to a shallower central profile such as NFW.

The capacity to detect very small cores (recall $x_{c}$ is normalized by the
Einstein radius) leads to high sensitivity, particularly for lighter halos,
and in some cases well beyond current limits. The capacity of GWs to constrain
$M_{Lz}$ and $x_{c}$ independently may allow observations to discern DM-
dominated halos ($M_{\rm vir}\lesssim 10^{9}\,M_{\odot}$), and thus whether
the core is likely to represent a signature of DM properties, rather than
baryonic effects. We note that turning these prospects into actual constraints
requires further modelling of lens features, astrophysical effects and lensing
cross sections, as well as robust detection of lensed GW signals. To
illustrate the need for more realistic modelling, we discuss how our analyses,
based on axisymmetric lenses, can be generalized to elliptical lenses (Sec.
V.4). Despite the challenges ahead, our results warrant further investigations
into the potential of DM tests with strongly-lensed GWs.

Our analysis highlights several ways in which GW lensing is highly
complementary to electromagnetic observations:

  1. 1.

The low frequency of GW sources enables observation of WO effects in lenses
with mass $M_{Lz}\gtrsim 10^{2}$, $10^{5}\,M_{\odot}$ (point-lens) and $M_{\rm
vir}\gtrsim 10^{4}$, $10^{7}\,M_{\odot}$ (extended) for LVK and LISA,
respectively. WO provides information about a region of the lens plane of size
$\sim w^{-1/2}$.

  2. 2.

WO effects provide additional information through the onset of magnification
and mild frequency-dependence at low/high dimensionless frequencies (beyond
GO), as well as the contributions from cusps. These effects persist in the
single-image (weak-lensing) regime at moderate $w$ (Fig. 9).

  3. 3.

GW lensing is a unique probe of the inner structure of halos. Central,
demagnified images are brighter for GWs than for EM sources
($\sqrt{|\mu_{H}|}\gg|\mu_{H}|$), cf. Fig. 10. Unlike EM signals, central
images of GW are not blocked or outshined by matter in the inner regions of
halos.

  4. 4.

Central images cause additional frequency-dependent modulation of the
amplification pattern and carry the strongest frequency dependence via beyond
GO corrections (Figs. 4, 7). In the vicinity of a caustic, both the onset of
this modulation and convergence to GO is pushed to very high frequencies (Fig.
8).

  5. 5.

GW lensing effects enable measurements of individual lens parameters with
precision $\sim\mathcal{O}(1/{\rm SNR})$, despite strongly correlated
posteriors, Figs. 17. Constraints on strongly-lensed signals ($y<y_{rc}$)
improve for large impact parameters, due to larger time delays, Figs. 15, 19.
This feature favours observational prospects, since differential probabilities
scale as $\propto y^{2}$.

  6. 6.

While lensed waveforms are sensitive to projected parameters ($M_{Lz}$, $y$,
$x_{c}$), limits on the physical scales and virial mass can be placed under
reasonable assumptions, Fig. 21. In particular, $M_{\rm vir}$ can be obtained
within a factor of a few.

  7. 7.

Central images may enable novel tests of astrophysics and fundamental physics,
such as probing DM through lens cores (Figs. 22, 23). Variations of these
tests may rely on measuring the slope of the inner halo cusp or the effect of
microlenses and supermassive black holes close to the central image.

Based on these prospects, we envision further extensions of our work regarding
lens models, source properties, statistical analysis and towards other lensing
regimes.

While we focused on simple symmetric lenses, we expect many of our conclusions
to hold in more realistic situations. For instance, the central GO image
probes the inner region of the lens, which may differ from the outer, averaged
lens distribution. While this dependence is exact only for GO and symmetric
lenses (the lens equation depends on the integrated projected mass), it should
approximate more general settings. Further characterization may enable
analyses based on central images complementary to those proposed for EM
sources, e.g. [16, 17, 20, 21, 22]. Still, more complex lenses may mimic or
obscure some of the signatures of cores and cusps. Studying other lenses
[105], adding large-scale lens features (variable slope, ellipticity, external
convergence & shear) and small-scale substructure (DM subhalos, stellar
population) is thus relevant. Future work needs to address these modelling
challenges and opportunities.

Another interesting direction is a more detailed exploration of DM scenarios.
By modelling DM halos as CIS, our work has shown the promise of GW lensing
under rather simplistic assumptions. As discussed in Sec. V, it will be
necessary to model lenses more accurately, along the lines mentioned above.
Our work can, in principle, be extended to other DM models. Warm dark matter
halos have been shown to have small cores, $r_{c}/r_{200}\lesssim 10^{-3}$
[136]. Another promising signature is the inner halo slope $\rho\propto
r^{-3/2}$ [137, 138]. This feature is also present in other models such as
compact DM structures [63], with a sharp slope $\rho\propto r^{-9/4}$
(although see Ref. [139]) and constant density core, whose size is set by the
DM velocity dispersion or annihilation cross-section [140]. Ultimately,
deriving constraints will require a detailed assessment of lensing
probabilities, including false positives due to baryonic or complex
structures.

A key advantage of GWs is the existence of accurate source models. While we
considered the source parameters as fixed (except for distance and phase) and
simple sources (equal mass, non-spinning binaries), further understanding the
interplay between source models and lensing effects is necessary. Besides
varying the source parameters, other directions include the study of different
sources (e.g. Extreme-mass-ratio inspirals), and more realistic waveform
models (higher harmonics, spin-precession, eccentricity, see Ref. [36]).
Ultimately, full Bayesian parameter estimation of lensed signal injections
will be necessary to derive accurate posteriors and address the shortcomings
of the Fisher matrix. These studies will allow an understanding of both
parameter degeneracies and potential systematics, e.g. if a lensed event is
analyzed with an unlensed waveform. Eventually, WO amplification factors need
to be included in GW searches and parameter estimations, perhaps taking
advantage of machine learning algorithms, e.g. [141].

Our results also warrant investigating other regimes of GW lensing. Here we
have focused on the strong-lensing regime, at mild magnifications. An
interesting extension is highly magnified configurations in which the source
is close to a caustic [39, 142], and the maximum magnification is limited by
wavelength rather than source size [143, 144]. This requires pushing our
numerical methods to new limits, especially in computing forecasts. Another
future direction is the study of line-of-sight halos [145, 146, 147], and
explore the interplay between multi-plane lensing and wave optics [148] to
enable ray-tracing in cosmological simulations. Finally, wave optics in the
single-image regime is also promising, as the frequency dependence may be
detectable up to large impact parameters [4, 33, 36, 6, 149]. Since the
lensing probabilities scale as the square of the impact parameter, weak-
lensing effects might be the first detectable signature, or even contribute as
a systematic in the waveform reconstruction.

As the detection rate increases, GW lensing is bound to become a reality. It
is a matter of chance whether LVK searches return confident identification of
lensed signals in the coming years. However, the capacity of next-generation
detectors to observe binary mergers in the entire Universe [11] will make GW
lensing into a reality: even the least likely lens configurations will stand a
chance to be observed. Other proposed missions will offer new opportunities,
such as WO effects for galactic scale lenses from $\mu{\rm Hz}$ space-borne
detectors [103]. As GW astronomy matures, wave optics effects may serve both
to identify lensed sources and to study the matter distribution in the
Universe. The complementarity between lensing of GW and of EM sources may
enable novel probes of the structure of halos, perhaps helping us identify the
nature of DM and opening new insights into astrophysics, fundamental physics,
and the history and evolution of the Universe.

###### Acknowledgements.

It is a pleasure to thank G. Brando, JM. Ezquiaga, A. Kumar Mehta, J. Poon, S.
Savastano and H. Villarrubia-Rojo for useful discussions and comments. L.D.
acknowledges the research grant support from the Alfred P. Sloan Foundation
(Award Number FG-2021-16495). M.H.-Y.C. is a Croucher Scholar supported by the
Croucher Foundation. M.H.-Y.C. is also supported by NSF Grants No.
AST-2006538, PHY-2207502, PHY-090003 and PHY-20043, and NASA Grants No.
19-ATP19-0051, 20-LPS20- 0011 and 21-ATP21-0010.

##  Appendix A Lensing Recap

Here we introduce briefly wave-optics computations in the low-frequency limit
(A.1) and via numerical methods (A.2), as well as the singular-isothermal
sphere lens (A.3), whose generalizations we used throughout the paper.

###  A.1 Wave-Optics: Low-frequency expansion

Let us summarise the results obtained in [5] regarding the low-$w$ expansion
of $F(w)$. At low frequencies, the GO expansion fails and one needs a
different approach to obtain approximations for the amplification factor. In
the limit $w\ll 1$, the wavelength of the wave becomes larger than the
characteristic scale of the lens, and the signal remains unperturbed:
$F(w)\sim 1$. Corrections to this result, at small impact parameters and for
axially-symmetric lensing potentials, are computed as a series expansion in
powers of $\psi(x)$ [5]:

$\displaystyle F(w)$ $\displaystyle\simeq
e^{-iw\phi_{m}}\bigg{[}1-\int_{0}^{\infty}\,{\rm d}z\,z\,e^{-z^{2}/2}\times$
$\displaystyle\left(iw\psi(q)+\frac{w^{2}}{2}\psi(q)^{2}+\mathcal{O}(w^{3}\psi(q)^{3})\right)\bigg{]}\;,$
(77)

where $q\equiv e^{i\pi/4}z/\sqrt{w}$. This expansion is well defined provided
that $\psi(q)$ grows less rapidly than the quadratic part of the Fermat
potential as $z\to\infty$.

For the lenses we are going to discuss (see Tab. 1), the integrals in Eq.
(A.1) can be evaluated explicitly. By suppressing the overall factor
$e^{-iw\phi_{m}}$, for the point lens $\psi(x)=\log x$ we have

$F^{\rm pl}(w)\simeq
1+\frac{w}{4}\left(\pi+2i\gamma_{\mathrm{E}}+2i\log\frac{w}{2}\right)\;,$ (78)

where $\gamma_{\mathrm{E}}$ is the Euler’s constant. Similarly, for SIS

$F^{\rm SIS}(w)\simeq 1-(-1)^{3/4}\sqrt{\frac{\pi w}{2}}-iw\;.$ (79)

While for gSIS we have

$F^{\rm gSIS}(w)\simeq 1+(-w/2)^{k/2}\Gamma\left((2-k)/2\right)\;.$ (80)

Finally, for the CIS:

$\displaystyle F^{\mathrm{CIS}}(w)$ $\displaystyle\simeq
1-\left[(-1)^{3/4}\sqrt{\frac{\pi
w}{2}}+\frac{wx_{c}}{4}\left(2i\gamma_{\mathrm{E}}\right.\right.$
$\displaystyle\left.\left.+\pi+2i\log\frac{w}{2}+4i\log\left(2x_{c}\right)+\frac{4i}{x_{c}}\right)\right]\;.$
(81)

We will use these expressions to compute $F(w)$ and its derivatives at low
frequencies.

###  A.2 Wave Optics: Numerical Methods

Numerical methods are necessary to compute the amplification factor at
intermediate frequencies, bridging the gap between the low- and high-frequency
expansions described above. In Ref. [5] we developed and validated two
numerical methods, regularized contour flow [3] and complex deformation [150].
Both yield the same results at high accuracy, but the regularized contour flow
is substantially faster and will be our method of choice. We will give an
overview now; further details can be found in Ref. [5].

The starting point is the Fourier transform of the integral in (14),
$\tilde{I}(\tau)\equiv\int{\rm d}^{2}\boldsymbol{x}\int\frac{{\rm
d}w}{2\pi}\exp\left(iw\left(\phi(\boldsymbol{x},\boldsymbol{y})-\tau\right)\right)$,
or

$\tilde{I}(\tau)=\int{\rm
d}^{2}\boldsymbol{x}\,\delta_{D}\left(\phi(\boldsymbol{x},\boldsymbol{y})-\tau\right)\,,$
(82)

where $\delta_{D}(x)$ is the Dirac-delta function. The above integral is
computed on contours of equal Fermat potential
$\phi(\boldsymbol{x},\boldsymbol{y})=\tau$. As contours end on critical
points, the first step is to solve the lens equation to find the initial/final
conditions for the contours. The integral is then sampled on the nodes that
define each contour. Then the nodes are evolved to the next value of $\tau$,
adding or removing more nodes depending on the curvature of the contour. The
next value of $\tau$ is chosen adaptively and depends on the derivatives of
$\tilde{I}(\tau)$.

The amplification factor (14) in the frequency domain is then obtained via a
fast-Fourier transform (FFT). This is computed by combining all contours,
interpolated on a homogeneous grid, as defined by the frequency range of the
FFT. In order to improve the computation, we split $\tilde{I}(\tau)$ into a
regular and a singular part. The singular part is the time-domain counterpart
to the GO amplification factor (21). It can be computed analytically and
consists of discontinuities (for maxima/minima of $\phi$) and logarithmic
divergences (for saddle points of $\phi$), which if not dealt with separately
contribute spurious high-frequency noise to the FFT.

This algorithm is efficient and accurate at intermediate frequencies. At high
frequency, we use the bGO approximation described in Eq. (24). The transition
between the two regimes, $w_{\rm high}$, is chosen such that the relative
deviation $|F_{\rm bGO}/F_{\rm WO}-1|$ is sub-percent. At the same time,
$w_{\rm high}$ is chosen relatively low to avoid numerical aliasing from the
FFT, usually $w_{\rm high}\sim\mathcal{O}(100)$. The choice for $w_{\rm high}$
is more delicate when computing derivatives of $F(w)$ with respect to the lens
parameters since numerical noise is higher in these cases: Appendix B
discusses this issue in the context of Fisher-matrix forecasts.

###  A.3 Singular Isothermal Sphere

The SIS is characterised by the following density distribution

$\rho(r)=\frac{\sigma_{v}^{2}}{2\pi Gr^{2}}\;,$ (83)

where $\sigma_{v}$ is the line of sight velocity dispersion. This lens is
often used to model dark matter halos. Because of the falloff as $\sim
1/r^{2}$ of the density profile, the halo mass is formally infinite, and the
model thus requires a physical cutoff for $r$, after which the density is
imagined to drop to zero quickly. In the context of lensing, this is not an
issue as long as the cutoff is much larger than the scale $\xi_{0}$ specified
below.

The projected mass density of this lens is (Eq. (5))

$\Sigma(\boldsymbol{\xi})=\frac{\sigma_{v}^{2}}{2G|\boldsymbol{\xi}|}\,.$ (84)

This simple lens model falls within the class of axisymmetric lenses since the
projected mass only depends on $|\boldsymbol{\xi}|$. The lensing potential is
then obtained from Eq. (11) as $\psi(x)=\sigma_{v}^{2}x/(\xi_{0}G\,\Sigma_{\rm
cr})$. In our work we use a standard choice for the normalisation of the lens,
$\xi_{0}=\sigma_{v}^{2}/(G\,\Sigma_{\rm cr})$. With this choice the lensing
potential becomes particularly simple

$\psi(x)=x\;.$ (85)

In the GO limit we can have one or two images depending on whether the impact
parameter is outside or inside the caustic $y_{\rm cr}=1$ (see Fig. 3, Panel
C). For $y<y_{\rm cr}$, the two images have magnifications $\mu_{\pm}=1/y\pm
1$, time delays $\phi_{\pm}=\mp y-1/2$ and Morse phases $n_{+}=0$ (minimum),
$n_{-}=1/2$ (saddle). For $y>y_{\rm cr}$ only the minimum image survives.
Consistently with the prescription of Eq. (4), the minimum time delay needs to
be set to zero by adding the appropriate $\phi_{m}(y)=-y-1/2$.

In order to properly interpret our results, we now discuss the relationship
between the quantities above and the physical properties of the SIS. To define
the total mass of the SIS it is necessary to truncate the mass density at a
finite radius. We will follow the standard definition of considering the
virial radius $r_{\rm vir}$ as the radius at which the density reaches
$200\rho_{c}$ , with $\rho_{c}=3H_{0}^{2}/(8\pi G)$ being the critical density
of the Universe today ($H_{0}=100\,h\,\rm{km}\,\rm{s}^{-1}\,\rm{Mpc}^{-1}$ is
the Hubble constant), and choose this as the truncation radius. Then, the
virial mass $M_{\rm vir}$ is defined as the mass within the virial radius and
coincides with the mass of the halo. Using Eq. (83) and the expression for
$r_{\rm vir}$, $M_{\rm vir}$ is obtained as

$\displaystyle M_{\rm vir}$ $\displaystyle=4\pi\int_{0}^{r_{\rm vir}}{\rm
d}r\,r^{2}\rho(r)=\frac{2\sigma_{v}^{3}}{5\sqrt{6}GH_{0}}$ (86)
$\displaystyle=1.8\cdot 10^{11}\,M_{\odot}\left(\frac{\sigma_{v}}{70\,{\rm
km/s}}\right)^{3}\;.$

Alternatively, the expression above can be seen as a relation between the
velocity dispersion and the halo mass.

With our normalisations, the effective lens mass is instead given by
$M_{Lz}=4\pi^{2}(1+z_{L})^{2}\sigma_{v}^{4}d_{\rm eff}/G$ (see Eq. (16) and
our choice of $\xi_{0}$). Using Eq. (86) we can then relate $M_{Lz}$ with the
virial mass:

$\displaystyle M_{Lz}$ $\displaystyle=\frac{4\pi^{2}}{G}(1+z_{L})^{2}d_{\rm
eff}\left(\frac{5\sqrt{6}}{2}GH_{0}M_{\rm vir}\right)^{4/3}$ (87)
$\displaystyle=2.3\cdot 10^{6}\,M_{\odot}(1+z_{L})^{2}\left(\frac{d_{\rm
eff}}{1\,{\rm Gpc}}\right)\left(\frac{M_{\rm
vir}}{10^{9}\,M_{\odot}}\right)^{4/3}\;.$

This relation, together with Eq. (19), suggests that for LISA frequencies and
virial masses around $M_{\rm vir}\sim 10^{9}\,M_{\odot}$, lensing is best
described by wave optics.

##  Appendix B Forecast implementation and tests

Figure 24: Convergence of the amplification factor and the Fisher derivatives
to bGO and rGO respectively for a CIS lens (left) and gSIS (right). The
bGO/rGO values are used for high frequencies (shaded region) to avoid aliasing
noise.

In this appendix, we are going to discuss the methods and validation
procedures we employed in evaluating Fisher-matrix elements. The numerical
derivatives of $F(w)$ with respect to the lens parameters $\theta_{l}$ are
particularly delicate. Indeed, as we are going to mention soon, the large
degeneracies in the Fisher matrix can amplify numerical errors when computing
its inverse. It is then also important to develop approximation methods at low
and high frequencies in order to validate our results in these regions.

At high frequencies we compute derivatives of the amplification factor using
the beyond GO expression, Eq. (24), which can be arranged as a sum over GO
images $I$,

$\partial_{l}F(w)=\sum_{I}\left(iw\psi_{I,l}+\frac{\mu_{I,l}}{2\mu_{I}}+\frac{i}{w}\Delta_{I,l}\right)F_{I}(w)\,.$
(88)

Here uppercase Latin indices $I,J,...$ refer to GO images, lowercase indices
$l,m,...$ to lens parameters and commas denote partial derivatives with
respect to the lens parameters ($f_{,l}\equiv\partial f/\partial\theta_{l}$).
Additionally, $F_{I}(w)$ is the GO amplification of the $I$-th image, Eq.
(21). Note that the Morse phases $n_{I}$ are assumed to be constant in the
derivatives.

Parametrizing the effective lens mass by its logarithm provides a convenient
rescaling of the Fisher matrix. Since $M_{Lz}$ only enters through the
dimensionless frequency $w$, Eq. (15), the associated derivative reads

$\frac{\partial F(w)}{\partial\log M_{Lz}}=\frac{\partial F(w)}{\partial\log
w}=\sum_{I}\left(iw\psi_{I}-\frac{i}{w}\Delta_{I}\right)F_{I}(w)\,.$ (89)

Here the first equality is fully general and will be used in the intermediate
and low-frequency regimes. The second expression is used in the bGO regime.
Similar expressions are used for the derivatives of the cusp terms appearing
in the gSIS lens, Eq. (26).

At low frequencies we compute $F(w)$ and its derivatives using the low-$w$
expansion of Sec. A.1 (one needs to re-introduce the minimum time delay
$\phi_{m}$ in these expressions, as it depends on the lensing parameters and
so it affects the derivatives). Moreover, we checked these results do not
change when using the full WO result at low frequency (computed using Eq. (16)
of Ref. [5]).

Derivatives of $F(w)$ in the intermediate region are computed in the time
domain via finite differences and then Fourier transformed. This process
significantly worsens the numerical noise at high frequencies, in a manner
analogous to the introduction of terms $\propto w$ in the bGO derivatives, Eq.
(88). Given the sensitivity of the inverse Fisher matrix to numerical errors,
it is necessary to achieve sufficient accuracy. To this end, we use a
supplementary regularization procedure, described in App. A of Ref. [5]. In
addition, we choose $w_{\rm high}$, the transition between WO and bGO, such
that the errors between both computations are smaller than $1\%$ for $F(w)$
and $3\%$ for the derivatives with respect to the lens parameters. Figure 24
shows how these tests are implemented for the CIS and gSIS. We have found that
achieving our target precision is hardest for the core size and the density
slope.

The large degeneracy between lens parameters makes the inversion of the Fisher
matrix (i.e. the posterior) very sensitive to numerical errors. This can be
seen in the condition numbers, which are typically $\mathcal{O}(10^{4})$ at
best and much larger at very low frequencies $w\ll 1$, as the lack of
magnification prevents the recovery of lens parameters (see Fig. 17 in the
main text and Fig. 27 for the dependence of the eigenvalues with the lens mass
in a specific example). The conditioning number also grows as $M_{Lz}^{2}$
when $w\gg 1$ due to the degeneracies in lens parameters in the geometric
optics limit. These degeneracies are responsible for the saturation of the
precision at high lens masses (for fixed SNR), e.g. Figs. 13, 15 and 15.

As a further test of the stability of the computation, we have checked the
convergence between our full calculation, GODA and bGODA (introduced in App.
C) at sufficiently high lens masses (see Fig. 25 in App. C). In the gSIS case
with $k>1$ (only two images) the Fisher matrix for GODA is not invertible:
with only two images one can extract at most two lens parameters. The
inclusion of bGO terms allows instead to extract more information, although
the precision falls with $M_{Lz}$, since GODA is recovered in this limit.
Finally, the rGODA terms (containing also the cusp contribution) allow to
formally obtain constant precision at high $M_{Lz}$, cf. Fig. 26 in App. C. We
explain this (paradoxical) result in Sec. IV.3, motivating that it cannot be
trusted, as the LSA becomes invalid in this limit.

##  Appendix C Geometric Optics Diagonal Approximation

Figure 25: 68% c.l. marginalised constraints on the lens parameters as a
function of the lens mass. $y=0.3$, $x_{c}=10^{-2}$ and the source
$10^{6}\,M_{\odot}$ at $z=3$. Results are normalized to a lensed ${\rm
SNR}=10^{3}$. The solid black line shows the marginalised errors using the
full calculation, including wave-optics effects. The red-dotted and magenta
dash-dotted correspond to Geometric-optics diagonal approximation and its
beyond GO extension. The upper horizontal axes show the dimensionless
frequency corresponding to the innermost stable circular orbit (ISCO) given
the source and lens mass.

Figure 26: 68% c.l. marginalised constraints on the lens parameters for gSIS
as a function of the lens mass. The fiducial parameters are $y=0.3$,
$k=1+10^{-4}$ and the source mass $10^{6}M_{\odot}$ at $z=3$. Results are
normalized to a lensed ${\rm SNR}=10^{3}$. The solid black line shows the
marginalised errors using the full calculation, including wave-optics effects.
The red-dotted and magenta dash-dotted correspond to beyond Geometric-optics
diagonal and the rGO diagonal approximations. The upper horizontal axes show
the dimensionless frequency corresponding to the innermost stable circular
orbit (ISCO) given the source and lens mass. See discussions in Sec. B and
IV.3.

We will introduce an approximation valid at high frequencies, based on the bGO
Fisher derivatives (88). Here we will give expressions for GO, but bGO (and
rGO) generalizations can be obtained straightforwardly. In the GO limit the
Fisher-matrix elements involving lens parameters $l,m$ read

$\mathcal{F}_{lm}^{\rm GO}=\sum_{IJ}\sqrt{|\mu_{I}\mu_{J}|}\int{\rm
d}f\frac{|\tilde{h}_{0}(f)|^{2}}{S_{n}(f)}\mathcal{I}_{Il,Jm}\;.$ (90)

The integrand associated with each term is

$\displaystyle\mathcal{I}_{Il,Jm}$
$\displaystyle=\left(\frac{\mu_{I,l}\mu_{J,m}}{\mu_{I}\mu_{J}}+4w^{2}\phi_{I,l}\phi_{J,m}\right)\cos(\Delta_{IJ})$
$\displaystyle\quad+2w\left(\phi_{I,l}\frac{\mu_{J,m}}{\mu_{J}}-\frac{\mu_{I,l}}{\mu_{I}}\phi_{J,m}\right)\sin(\Delta_{IJ})\;.$
(91)

The oscillatory pieces are given by the phase difference between images

$\Delta_{IJ}=w(\phi_{I}-\phi_{J})+n_{I}-n_{J}\;.$ (92)

At high frequencies, the GO Fisher matrix will in general be dominated by the
diagonal terms $I=J$, for which the oscillatory contribution cancels. We,
therefore, introduce the Geometric Optics Diagonal Approximation (GODA) for
the Fisher matrix, defined by neglecting the cross-image terms as follows:

$\displaystyle\mathcal{F}_{lm}^{\rm GD}$
$\displaystyle\equiv\sum_{I}\bigg{[}\frac{\mu_{I,l}\mu_{I,m}}{4|\mu_{I}|}(h_{0}|h_{0})$
$\displaystyle\quad\quad\quad\quad+|\mu_{I}|\phi_{I,l}\phi_{I,m}(wh_{0}|wh_{0})\bigg{]}\;.$
(93)

The first term is proportional to the unlensed signal’s $\text{SNR}^{2}$. The
second term involves the signal weighted by the dimensionless frequency
$(wh_{0}|wh_{0})\propto M_{Lz}^{2}(fh_{0}|fh_{0})$ and typically dominates in
the high-mass limit. Note that GODA is equivalent to GO in the weak-lensing
regime, as there is a single image.

Figure 25 shows the convergence of the full calculation to GO and bGO for an
example lensed LISA source. Good convergence is achieved above lens masses
such that the dimensionless frequency corresponding to the innermost stable
circular orbit (ISCO) is large, $w_{\rm ISCO}>10$. For lower lens masses the
GODA overestimates the sensitivity: it traces the shape of the error in the
mass, but completely overestimates other lens parameters. This is because
information about the mass is only contained in the phase/time delays, and
hence lost in the limit $w(\phi_{I}-\phi_{J})\lesssim 1$. Other lens
parameters affect the GO magnification factors, and hence appear to be
measurable in the GODA limit even when $w\ll 1$, as GODA fails to account for
the lack of magnification in the low-$w$ regime. WO has an additional source
of information on the lens mass, namely the onset of magnification around
$w\sim 1$.

The neglected cross-image terms ($I\neq J$) lead to oscillatory integrals and
are in general subdominant for large $w$. This can be seen as follows,
assuming that the coefficient $|h_{0}|^{2}/S_{n}$ depends weakly on frequency.
The integral of a term proportional to $w^{\alpha}$ in Eq. (C) would be
proportional to $w^{\alpha+1}$ for $I=J$, while to
$\Delta_{IJ}^{-1}w^{\alpha}\cos(\Delta_{IJ}w)\leq w^{\alpha+1}$ for $I\neq J$.
Cross-image terms may be important in some cases, e.g. if the source is
observed over a narrow range of frequencies. In addition, if an image is very
faint compared to other images (e.g. $\mu_{I}^{2}\ll|\mu_{I}\mu_{J}|$), its
GODA contribution can be smaller than the $I\neq J$ cross terms. Finally, non-
diagonal terms might be important if the unlensed waveform $\tilde{h}_{0}$ has
a modulation with period $w\sim 1/\Delta_{IJ}$.

##  Appendix D LSA validity for CIS

At this point we would like to check the validity of the Fisher-matrix
formalism for our analysis, starting from the CIS lens. Let us write the full
likelihood $\mathcal{L}$ as

$\mathcal{L}=e^{-\Lambda}=e^{-\frac{1}{2}(s-h(\boldsymbol{\theta})|s-h(\boldsymbol{\theta}))}\;,$
(94)

We then require that the $1\sigma$ posteriors of the full likelihood (94) fall
within the 2$\sigma$ contours obtained through the Fisher analysis.141414This
criterion is less conservative than the one proposed in [108] to check the
validity of the LSA: ours can be satisfied even in cases where the full
likelihood and the one from the Fisher analysis differ substantially.
Nonetheless, we believe it to be informative since it bounds the size of the
correlations of the full likelihood. If this condition is satisfied, it means
we are not missing important degeneracies that could hinder our results. In
practice, we check a simplified version of this condition, namely that the
full log-likelihood $\Lambda$ evaluated at the $2\sigma$ contours (obtained
from the Fisher) is larger than the Fisher log-likelihood evaluated at
$1\sigma$ (from our definition of $\Lambda$ in Eq. (94) this means the
probability density is smaller).

In order to evaluate this likelihood for a large number of points, we work in
the regime of high $w$ (high $M_{Lz}$), where we can obtain lensed signals
quickly using the GO/bGO amplification factor. This is also the region where
the forecast gives the best sensitivity and, as we have motivated in the main
text (see Sec. IV.1), the breakdown of the LSA condition is expected to
appear. A point on the $2\sigma$ contours is parameterized as
$\Delta\theta^{i}=\sum_{k=1}^{5}u_{(k)}^{i}c_{k}(\lambda_{(k)}/2)^{-1/2}$,
where $\lambda_{(k)}$ and $u_{(k)}^{i}$ are the eigenvalues and eigenvectors
defined in Eq. (67) and $c_{k}$ is a unit-norm vector on the 5D parameter-
space sphere. To generate more points, we draw the values of $c_{k}$ from a
normal distribution [108]. We then evaluate the full log-likelihood
$\Lambda=(h(\boldsymbol{\bar{\theta}})-h(\boldsymbol{\bar{\theta}}+\Delta\boldsymbol{\theta})|h(\boldsymbol{\bar{\theta}})-h(\boldsymbol{\bar{\theta}}+\Delta\boldsymbol{\theta}))/2$.
Again, $h(\boldsymbol{\theta})$ is computed for simplicity using bGO. Our
criterion is satisfied if $\Lambda$ is larger than the Fisher log-likelihood
for most of the points $\Delta\boldsymbol{\theta}$.

This criterion is found to be satisfied $99.7\,\%$ of the times in our typical
fiducial values $y=0.3$, $x_{c}=10^{-2}$, $M_{Lz}=10^{10}\,M_{\odot}$, $M_{\rm
BBH}=10^{6}\,M_{\odot}$ and $z=3$ (for $10^{3}$ points).

Even when the condition is satisfied, we find that the log-likelihood
$\Lambda$ is not matching the Fisher one. This does not necessarily call for
concern, in the sense we just explained, but in any case points towards a
breakdown of the LSA in the high-$M_{Lz}$ limit. These conclusions motivate
the explanation given in the main text, Sec. IV.1 and IV.2.

Notice that this issue is common to all lens models in the strong-lensing
regime, not only to the CIS. To our knowledge, this issue was not investigated
previously in the literature. It would be interesting to revisit the Fisher
results in the literature and characterize the sensitivities with a more
sophisticated analysis.

As we are soon going to mention, this issue becomes more problematic in the
gSIS case, where the Fisher analysis yields very optimistic results that
contradict the number of parameters that can be reconstructed in the GO limit,
Eq. (57). We give a concrete analytical example of this issue in App. E. As we
comment at the end, the analysis of this example also supports the conclusion
that a violation of the LSA condition for the CIS should not call for concern.

Figure 27: Top: Normalized Fisher matrix for a LISA source with $M_{\rm
BBH}=10^{6}\,M_{\odot}$ and CIS with $M_{Lz}=10^{4}\,M_{\odot}$,
$10^{7}\,M_{\odot}$ and $10^{10}\,M_{\odot}$, $y=0.3$, $x_{c}=10^{-2}$.
Middle: Normalized inverse Fisher matrix. Bottom: Principal components of the
Fisher matrix (left to right, normalized to the largest component). The
precision is given on top of each eigenvector and it is rescaled to a lensed
${\rm SNR}=10^{3}$.

##  Appendix E LSA approximation for cusps

In the main text we argued that the linearisation of the likelihood used in
the Fisher matrix formalism might not be always appropriate, even for
moderately high $\rm SNR$. In turn, the linearisation procedure can
overestimate the sensitivity on parameters, yielding unrealistic, or even
paradoxical, results. Specifically, this seems to be the case for the
parameter $k$ in the gSIS forecast. For a fiducial value $k=1$ the Fisher
matrix predicts $\Delta k\propto 1/{\rm SNR}$ in the GO limit (large lens
masses $M_{Lz}$). On the other hand, in this limit the characteristic feature
due to $k$ (the signal from the cusp) quickly decays. Thus, no meaningful
constraint on $k$ should be present. Here we provide a toy example
highlighting this feature: by studying the full probability distribution we
can see that indeed the error on $k$ grows for large $M_{Lz}$.

Let us consider a toy lensing signal with a single generic lensing parameter
$\theta$ (that could represent $k$ in the case of the gSIS). The lensed GW
waveform $h_{L}(f)$ is taken to be

$h_{L}(f)=h_{0}(f)\left(1+\frac{e^{iw\theta}}{w}\right)\;,$ (95)

where $h_{0}(f)$ is the unlensed waveform, $w=8\pi GM_{Lz}$ and $M_{Lz}$ is
taken as fixed. This signal is indeed analogous to the cusp obtained for the
gSIS in the limit of $k\to 1$, see Eq. (26). One can easily generalise what
follows to cases with different $w$ dependencies.

Let us write $\theta=\bar{\theta}+\Delta\theta$, where $\bar{\theta}$ is the
fiducial value. We want to compute the probability distribution for
$\Delta\theta$, assuming some analytic behaviour for $h_{0}(f)$ and the noise.

Let us call $\Delta h\equiv h(\bar{\theta})-h(\theta)$. Then, the exponent
$\Lambda$ of the likelihood is defined as

$\Lambda=\frac{1}{2}\left(\Delta h|\Delta h\right)\;,$ (96)

where we are neglecting the noise term. In this way the likelihood is written
as $\mathcal{L}\propto e^{-\Lambda}$. We can write the full $\Lambda$ without
expanding in $\Delta\theta$. It is easy to get

$\Lambda=8\int\frac{{\rm
d}f}{S_{n}(f)}\frac{|h_{0}(f)|^{2}}{w^{2}}\left[1-\cos(w\Delta\theta)\right]\;.$
(97)

If we were to expand the square bracket for small $\Delta\theta$, we would
recover the Fisher formalism. Notice that if we expand and trust the result
for $w\Delta\theta>1$, then we are _overestimating_ $\Lambda$. In the full
result, we see indeed that the square bracket remains bounded for large
$\Delta\theta$.

At this point, we can understand the behaviour of the integral by taking a
simplified form for $h_{0}(f)$ and $S_{n}(f)$. The lower band of the LISA
noise is characterised by $S_{n}(f)\sim s_{n}f^{-4}$, where $s_{n}$ is a
dimensionful constant that depends on the noise properties of LISA. This is a
good approximation for events with large masses $M_{\rm BBH}\gtrsim
10^{7}\,M_{\odot}$ since the merger happens before the noise curve reaches its
minimum. For the waveform, we take simply $h_{0}(f)=Af^{-7/6}e^{i\Psi(f)}$,
where $A$ is the amplitude and $\Psi(f)$ the phase evolution [151]. We perform
the integration over $f$ from some low $f_{0}$, where the signal enters the
LISA band, and we cut the signal at the ISCO frequency, which we indicate as
$f_{1}$.

The SNR for this signal is then

${\rm SNR}^{2}=(h_{0}|h_{0})\simeq\frac{3A^{2}}{2s_{n}}f_{1}^{8/3}\;,$ (98)

where we assumed $f_{1}\gg f_{0}$.

We can use the expression above to simplify $\Lambda$, and get

$\Lambda\simeq\frac{8\,{\rm SNR}^{2}}{3f_{1}^{8/3}}\int_{f_{0}}^{f_{1}}{\rm
d}f\,\frac{f^{5/3}}{w^{2}}\left[1-\cos(w\Delta\theta)\right]\;.$ (99)

This integral can be solved explicitly in terms of exponential integrals.
However, for our purposes, we can focus on the two limits of $w\Delta\theta\ll
1$ and $w\Delta\theta\gg 1$, where $w$ represents the dimensionless frequency
contributing mainly to the integral (i.e. $w$ at the ISCO).

For small $w\Delta\theta$ we can expand the square bracket, recovering a
quadratic function of $\Delta\theta$: $\Lambda\simeq{\rm
SNR}^{2}\Delta\theta^{2}$. This corresponds with the Fisher-matrix result.
Naively, it would predict $\Delta\theta\sim 1/{\rm SNR}$, even though the
correction due to $\theta$ in Eq. (95) becomes negligible at high $M_{Lz}$.
From this discussion it is now clear that this result is only applicable for
$w_{\rm ISCO}\Delta\theta<1$ and therefore breaks down at high $M_{Lz}$. This
is analogous to the LSA condition in Eq. (66) for this simplified case.

Instead, in the opposite limit, we can estimate Eq. (99) by noting that the
cosine function becomes rapidly oscillating. This means its contribution to
the integral is negligible: the log-likelihood then becomes independent of
$\Delta\theta$:

$\Lambda\simeq\frac{4\,{\rm SNR}^{2}}{w_{\rm ISCO}^{2}}\;.$ (100)

This result shows $\theta$ cannot be efficiently constrained when $w_{\rm
ISCO}\sim{\rm SNR}$, since its probability distribution stops being localized
around $\Delta\theta=0$ and becomes flat, with $\mathcal{L}\propto
e^{-\Lambda}\sim 1$. In the more realistic situation where, on top of the cusp
contribution, there are more GO images the situation is of course more
involved. Indeed since the cusp is very faint, its contribution to the
likelihood might be dominated by the product with other images. The general
logic is expected to apply nonetheless.

In conclusion, we expect the Fisher approximation to stop being accurate
around $\Delta\theta\sim w_{\rm ISCO}^{-1}$, with the constraint on $\theta$
becoming weaker and weaker as we increase $M_{Lz}$. Eventually $\theta$
becomes unconstrained when we reach $w_{\rm ISCO}\sim{\rm SNR}$.

We can also comment on the case in which the lensing effects are not
suppressed by $1/w$, as in Eq. (95), but are due to a GO image. (Notice also
that bGO terms do not have the form of Eq (95), since the phase factor would
be factorised outside of the bracket.) It is straightforward to obtain the
analogous of Eq. (99) in this case (the factor $w^{-2}$ disappears).
Linearisation becomes unreliable when $\Delta\theta\sim w_{\rm ISCO}^{-1}$, as
in the previous case. However, the likelihood is still highly suppressed for
large $w_{\rm ISCO}$, since $\Lambda\sim{\rm SNR}^{2}$ (the tails of the
distribution differ from a Gaussian, but the associated probability remains
small). Applied to the CIS case, this suggests the Fisher analysis still
reliably captures the size of the sensitivities even after the LSA condition
is violated.

## References

  * Bartelmann [2010] M. Bartelmann, Class. Quant. Grav. 27, 233001 (2010), arXiv:1010.3829 [astro-ph.CO] .
  * Bertone and Tait [2018] G. Bertone and T. Tait, M. P., Nature 562, 51 (2018), arXiv:1810.01668 [astro-ph.CO] .
  * Ulmer and Goodman [1995] A. Ulmer and J. Goodman, Astrophys. J. 442, 67 (1995), arXiv:astro-ph/9406042 .
  * Takahashi and Nakamura [2003] R. Takahashi and T. Nakamura, Astrophys. J. 595, 1039 (2003), arXiv:astro-ph/0305055 .
  * Tambalo  _et al._ [2022] G. Tambalo, M. Zumalacárregui, L. Dai, and M. H.-Y. Cheung, (2022), arXiv:2210.05658 [gr-qc] .
  * Savastano  _et al._ [2023] S. Savastano, G. Tambalo, H. Villarrubia-Rojo, and M. Zumalacarregui, (2023), arXiv:2306.05282 [gr-qc] .
  * Oguri [2019] M. Oguri, Rept. Prog. Phys. 82, 126901 (2019), arXiv:1907.06830 [astro-ph.CO] .
  * Liao  _et al._ [2022] K. Liao, M. Biesiada, and Z.-H. Zhu, (2022), arXiv:2207.13489 [astro-ph.HE] .
  * Katz  _et al._ [2018] A. Katz, J. Kopp, S. Sibiryakov, and W. Xue, JCAP 12, 005 (2018), arXiv:1807.11495 [astro-ph.CO] .
  * Katz  _et al._ [2020] A. Katz, J. Kopp, S. Sibiryakov, and W. Xue, Mon. Not. Roy. Astron. Soc. 496, 564 (2020), arXiv:1912.07620 [astro-ph.CO] .
  * Kalogera  _et al._ [2021] V. Kalogera _et al._ , (2021), arXiv:2111.06990 [gr-qc] .
  * Schneider  _et al._ [1992] P. Schneider, J. Ehlers, and E. E. Falco,  _Gravitational Lenses_ (1992).
  * Shajib  _et al._ [2022] A. J. Shajib, G. Vernardos, T. E. Collett, V. Motta, D. Sluse, L. L. R. Williams, P. Saha, S. Birrer, C. Spiniello, and T. Treu, (2022), arXiv:2210.10790 [astro-ph.GA] .
  * Mao  _et al._ [2001] S. Mao, H. J. Witt, and L. V. E. Koopmans, Mon. Not. Roy. Astron. Soc. 323, 301 (2001), arXiv:astro-ph/0007011 .
  * Rusin and Ma [2001] D. Rusin and C.-P. Ma, Astrophys. J. Lett. 549, L33 (2001), arXiv:astro-ph/0009079 .
  * Keeton and Madau [2001] C. R. Keeton and P. Madau, Astrophys. J. Lett. 549, L25 (2001), arXiv:astro-ph/0101058 .
  * Evans and Hunter [2002] N. W. Evans and C. Hunter, Astrophys. J. 575, 68 (2002), arXiv:astro-ph/0204206 .
  * Keeton [2003] C. R. Keeton, Astrophys. J. 582, 17 (2003), arXiv:astro-ph/0206243 .
  * Rusin  _et al._ [2005] D. J. Rusin, C. R. Keeton, and J. N. Winn, Astrophys. J. Lett. 627, L93 (2005), arXiv:astro-ph/0412034 .
  * Li  _et al._ [2012a] N. Li, S. Mao, L. Gao, A. Loeb, and R. Di Stefano, Mon. Not. Roy. Astron. Soc. 419, 2424 (2012a), arXiv:1110.0887 [astro-ph.CO] .
  * Hezaveh  _et al._ [2015] Y. D. Hezaveh, P. J. Marshall, and R. D. Blandford, Astrophys. J. Lett. 799, L22 (2015), arXiv:1501.01757 [astro-ph.GA] .
  * Quinn  _et al._ [2016] J. Quinn _et al._ , Mon. Not. Roy. Astron. Soc. 459, 2394 (2016), arXiv:1604.00173 [astro-ph.CO] .
  * Inada  _et al._ [2008] N. Inada, M. Oguri, E. E. Falco, T. J. Broadhurst, E. O. Ofek, C. S. Kochanek, K. Sharon, and G. P. Smith, Publ. Astron. Soc. Jap. 60, L27 (2008), arXiv:0808.1769 [astro-ph] .
  * Blanchet [2014] L. Blanchet, Living Rev. Rel. 17, 2 (2014), arXiv:1310.1528 [gr-qc] .
  * Buonanno and Sathyaprakash [2014] A. Buonanno and B. S. Sathyaprakash, “Sources of Gravitational Waves: Theory and Observations,” (2014) arXiv:1410.7832 [gr-qc] .
  * Schmidt [2020] P. Schmidt, Front. Astron. Space Sci. 7, 28 (2020).
  * Venumadhav  _et al._ [2017] T. Venumadhav, L. Dai, and J. Miralda-Escudé, Astrophys. J. 850, 49 (2017), arXiv:1707.00003 [astro-ph.CO] .
  * Diego  _et al._ [2018] J. M. Diego _et al._ , Astrophys. J. 857, 25 (2018), arXiv:1706.10281 [astro-ph.CO] .
  * Zumalacarregui and Seljak [2018] M. Zumalacarregui and U. Seljak, Phys. Rev. Lett. 121, 141101 (2018), arXiv:1712.02240 [astro-ph.CO] .
  * Birrer  _et al._ [2022] S. Birrer, M. Millon, D. Sluse, A. J. Shajib, F. Courbin, L. V. E. Koopmans, S. H. Suyu, and T. Treu, (2022), arXiv:2210.10833 [astro-ph.CO] .
  * Suyu  _et al._ [2023] S. H. Suyu, A. Goobar, T. Collett, A. More, and G. Vernardos, (2023), arXiv:2301.07729 [astro-ph.CO] .
  * Sereno  _et al._ [2010] M. Sereno, A. Sesana, A. Bleuler, P. Jetzer, M. Volonteri, and M. C. Begelman, Phys. Rev. Lett. 105, 251101 (2010), arXiv:1011.5238 [astro-ph.CO] .
  * Gao  _et al._ [2022] Z. Gao, X. Chen, Y.-M. Hu, J.-D. Zhang, and S.-J. Huang, Mon. Not. Roy. Astron. Soc. 512, 1 (2022), arXiv:2102.10295 [astro-ph.CO] .
  * Gais  _et al._ [2022] J. Gais, K. Ng, E. Seo, K. W. K. Wong, and T. G. F. Li, (2022), arXiv:2201.01817 [gr-qc] .
  * Savastano  _et al._ [2022] S. Savastano, F. Vernizzi, and M. Zumalacárregui, (2022), arXiv:2212.14697 [gr-qc] .
  * Çalışkan  _et al._ [2022a] M. Çalışkan, L. Ji, R. Cotesta, E. Berti, M. Kamionkowski, and S. Marsat, (2022a), arXiv:2206.02803 [astro-ph.CO] .
  * Lai  _et al._ [2018] K.-H. Lai, O. A. Hannuksela, A. Herrera-Martín, J. M. Diego, T. Broadhurst, and T. G. F. Li, Phys. Rev. D 98, 083005 (2018), arXiv:1801.07840 [gr-qc] .
  * Christian  _et al._ [2018] P. Christian, S. Vitale, and A. Loeb, Phys. Rev. D 98, 103022 (2018), arXiv:1802.02586 [astro-ph.HE] .
  * Diego  _et al._ [2019] J. M. Diego, O. A. Hannuksela, P. L. Kelly, T. Broadhurst, K. Kim, T. G. F. Li, G. F. Smoot, and G. Pagano, Astron. Astrophys. 627, A130 (2019), arXiv:1903.04513 [astro-ph.CO] .
  * Diego [2020] J. M. Diego, Phys. Rev. D 101, 123512 (2020), arXiv:1911.05736 [astro-ph.CO] .
  * Meena and Bagla [2020] A. K. Meena and J. S. Bagla, Mon. Not. Roy. Astron. Soc. 492, 1127 (2020), arXiv:1903.11809 [astro-ph.CO] .
  * Cheung  _et al._ [2021] M. H. Y. Cheung, J. Gais, O. A. Hannuksela, and T. G. F. Li, Mon. Not. Roy. Astron. Soc. 503, 3326 (2021), arXiv:2012.07800 [astro-ph.HE] .
  * Mishra  _et al._ [2021] A. Mishra, A. K. Meena, A. More, S. Bose, and J. S. Bagla, Mon. Not. Roy. Astron. Soc. 508, 4869 (2021), arXiv:2102.03946 [astro-ph.CO] .
  * Basak  _et al._ [2022] S. Basak, A. Ganguly, K. Haris, S. Kapadia, A. K. Mehta, and P. Ajith, Astrophys. J. Lett. 926, L28 (2022), arXiv:2109.06456 [gr-qc] .
  * Biesiada and Harikumar [2021] M. Biesiada and S. Harikumar, Universe 7, 502 (2021), arXiv:2111.05963 [gr-qc] .
  * Seo  _et al._ [2021] E. Seo, O. A. Hannuksela, and T. G. F. Li, in  _17th International Conference on Topics in Astroparticle and Underground Physics_ (2021) arXiv:2110.03308 [astro-ph.HE] .
  * Yeung  _et al._ [2021] S. M. C. Yeung, M. H. Y. Cheung, J. A. J. Gais, O. A. Hannuksela, and T. G. F. Li, (2021), arXiv:2112.07635 [gr-qc] .
  * Suvorov [2021] A. G. Suvorov, (2021), arXiv:2112.01670 [astro-ph.HE] .
  * Meena  _et al._ [2022] A. K. Meena, A. Mishra, A. More, S. Bose, and J. S. Bagla, Mon. Not. Roy. Astron. Soc. 517, 872 (2022), arXiv:2205.05409 [astro-ph.GA] .
  * Dai  _et al._ [2018] L. Dai, S.-S. Li, B. Zackay, S. Mao, and Y. Lu, Physical Review D 98, 104029 (2018).
  * Gil Choi  _et al._ [2021] H. Gil Choi, C. Park, and S. Jung, (2021), arXiv:2103.08618 [astro-ph.CO] .
  * Cremonese  _et al._ [2021] P. Cremonese, D. F. Mota, and V. Salzano, (2021), arXiv:2111.01163 [astro-ph.CO] .
  * Guo and Lu [2022] X. Guo and Y. Lu, (2022), arXiv:2207.00325 [astro-ph.CO] .
  * Buckley and Peter [2018] M. R. Buckley and A. H. G. Peter, Phys. Rept. 761, 1 (2018), arXiv:1712.06615 [astro-ph.CO] .
  * Jung and Shin [2019] S. Jung and C. S. Shin, Phys. Rev. Lett. 122, 041103 (2019), arXiv:1712.01396 [astro-ph.CO] .
  * Oguri and Takahashi [2020] M. Oguri and R. Takahashi, Astrophys. J. 901, 58 (2020), arXiv:2007.01936 [astro-ph.CO] .
  * Urrutia and Vaskonen [2021] J. Urrutia and V. Vaskonen, Mon. Not. Roy. Astron. Soc. 509, 1358 (2021), arXiv:2109.03213 [astro-ph.CO] .
  * Fairbairn  _et al._ [2022] M. Fairbairn, J. Urrutia, and V. Vaskonen, (2022), arXiv:2210.13436 [astro-ph.CO] .
  * Oguri and Takahashi [2022] M. Oguri and R. Takahashi, (2022), arXiv:2204.00814 [astro-ph.CO] .
  * Zhou  _et al._ [2022] H. Zhou, Z. Li, K. Liao, and Z. Huang, Mon. Not. Roy. Astron. Soc. 518, 149 (2022), arXiv:2206.13128 [astro-ph.CO] .
  * Buckley and DiFranzo [2018] M. R. Buckley and A. DiFranzo, Phys. Rev. Lett. 120, 051102 (2018), arXiv:1707.03829 [hep-ph] .
  * Savastano  _et al._ [2019] S. Savastano, L. Amendola, J. Rubio, and C. Wetterich, Phys. Rev. D 100, 083518 (2019), arXiv:1906.05300 [astro-ph.CO] .
  * Arvanitaki  _et al._ [2020] A. Arvanitaki, S. Dimopoulos, M. Galanis, L. Lehner, J. O. Thompson, and K. Van Tilburg, Phys. Rev. D 101, 083014 (2020), arXiv:1909.11665 [astro-ph.CO] .
  * Chakraborty  _et al._ [2022] A. Chakraborty, P. K. Chanda, K. L. Pandey, and S. Das, Astrophys. J. 932, 119 (2022), arXiv:2204.09628 [astro-ph.CO] .
  * Hu  _et al._ [2000] W. Hu, R. Barkana, and A. Gruzinov, Phys. Rev. Lett. 85, 1158 (2000), arXiv:astro-ph/0003365 .
  * Hui  _et al._ [2017] L. Hui, J. P. Ostriker, S. Tremaine, and E. Witten, Phys. Rev. D 95, 043541 (2017), arXiv:1610.08297 [astro-ph.CO] .
  * Ferreira [2021] E. G. M. Ferreira, Astron. Astrophys. Rev. 29, 7 (2021), arXiv:2005.03254 [astro-ph.CO] .
  * Hui [2021] L. Hui, Ann. Rev. Astron. Astrophys. 59, 247 (2021), arXiv:2101.11735 [astro-ph.CO] .
  * Spergel and Steinhardt [2000] D. N. Spergel and P. J. Steinhardt, Phys. Rev. Lett. 84, 3760 (2000), arXiv:astro-ph/9909386 .
  * Tulin and Yu [2018] S. Tulin and H.-B. Yu, Phys. Rept. 730, 1 (2018), arXiv:1705.02358 [hep-ph] .
  * Bullock and Boylan-Kolchin [2017] J. S. Bullock and M. Boylan-Kolchin, Ann. Rev. Astron. Astrophys. 55, 343 (2017), arXiv:1707.04256 [astro-ph.CO] .
  * Haris  _et al._ [2018] K. Haris, A. K. Mehta, S. Kumar, T. Venumadhav, and P. Ajith, (2018), arXiv:1807.07062 [gr-qc] .
  * Li  _et al._ [2019] A. K. Y. Li, R. K. L. Lo, S. Sachdev, C. L. Chan, E. T. Lin, T. G. F. Li, and A. J. Weinstein, (2019), arXiv:1904.06020 [gr-qc] .
  * McIsaac  _et al._ [2020] C. McIsaac, D. Keitel, T. Collett, I. Harry, S. Mozzon, O. Edy, and D. Bacon, Phys. Rev. D 102, 084031 (2020), arXiv:1912.05389 [gr-qc] .
  * Liu  _et al._ [2021] X. Liu, I. M. Hernandez, and J. Creighton, Astrophys. J. 908, 97 (2021), arXiv:2009.06539 [astro-ph.HE] .
  * Lo and Magaña Hernandez [2021] R. K. L. Lo and I. Magaña Hernandez, (2021), arXiv:2104.09339 [gr-qc] .
  * Janquart  _et al._ [2022a] J. Janquart, O. A. Hannuksela, K. Haris, and C. Van den Broeck, in  _56th Rencontres de Moriond on Gravitation_ (2022) arXiv:2203.06444 [gr-qc] .
  * Ali  _et al._ [2022] S. Ali, E. Stoikos, E. Meade, M. Kesden, and L. King, (2022), arXiv:2210.01873 [gr-qc] .
  * Janquart  _et al._ [2022b] J. Janquart, A. More, and C. Van Den Broeck, (2022b), arXiv:2205.11499 [gr-qc] .
  * Hannuksela  _et al._ [2019] O. A. Hannuksela, K. Haris, K. K. Y. Ng, S. Kumar, A. K. Mehta, D. Keitel, T. G. F. Li, and P. Ajith, Astrophys. J. Lett. 874, L2 (2019), arXiv:1901.02674 [gr-qc] .
  * Kim  _et al._ [2022] K. Kim, J. Lee, O. A. Hannuksela, and T. G. F. Li,  (2022), Arxiv:2206.08234v1 .
  * Dai and Venumadhav [2017] L. Dai and T. Venumadhav, (2017), arXiv:1702.04724 [gr-qc] .
  * Ezquiaga  _et al._ [2021] J. M. Ezquiaga, D. E. Holz, W. Hu, M. Lagos, and R. M. Wald, Phys. Rev. D 103, 064047 (2021), arXiv:2008.12814 [gr-qc] .
  * Wang  _et al._ [2021] Y. Wang, R. K. L. Lo, A. K. Y. Li, and Y. Chen, Phys. Rev. D 103, 104055 (2021), arXiv:2101.08264 [gr-qc] .
  * Vijaykumar  _et al._ [2022] A. Vijaykumar, A. K. Mehta, and A. Ganguly, (2022), arXiv:2202.06334 [gr-qc] .
  * Abbott  _et al._ [2021] R. Abbott _et al._ (LIGO Scientific, VIRGO), Astrophys. J. 923, 14 (2021), arXiv:2105.06384 [gr-qc] .
  * Dai  _et al._ [2020] L. Dai, B. Zackay, T. Venumadhav, J. Roulet, and M. Zaldarriaga, (2020), arXiv:2007.12709 [astro-ph.HE] .
  * Diego  _et al._ [2021] J. M. Diego, T. Broadhurst, and G. Smoot, Phys. Rev. D 104, 103529 (2021), arXiv:2106.06545 [gr-qc] .
  * Ng  _et al._ [2018] K. K. Y. Ng, K. W. K. Wong, T. Broadhurst, and T. G. F. Li, Phys. Rev. D 97, 023012 (2018), arXiv:1703.06319 [astro-ph.CO] .
  * Oguri [2018] M. Oguri, Mon. Not. Roy. Astron. Soc. 480, 3842 (2018), arXiv:1807.02584 [astro-ph.CO] .
  * Mukherjee  _et al._ [2021] S. Mukherjee, T. Broadhurst, J. M. Diego, J. Silk, and G. F. Smoot, Mon. Not. Roy. Astron. Soc. 506, 3751 (2021), arXiv:2106.00392 [gr-qc] .
  * Wierda  _et al._ [2021] A. R. A. C. Wierda, E. Wempe, O. A. Hannuksela, L. e. V. E. Koopmans, and C. Van Den Broeck, Astrophys. J. 921, 154 (2021), arXiv:2106.06303 [astro-ph.HE] .
  * Cusin  _et al._ [2021] G. Cusin, R. Durrer, and I. Dvorkin, Universe 8, 19 (2021).
  * Xu  _et al._ [2022] F. Xu, J. M. Ezquiaga, and D. E. Holz, Astrophys. J. 929, 9 (2022), arXiv:2105.14390 [astro-ph.CO] .
  * Çalışkan  _et al._ [2022b] M. Çalışkan, J. M. Ezquiaga, O. A. Hannuksela, and D. E. Holz, (2022b), arXiv:2201.04619 [astro-ph.CO] .
  * Dai  _et al._ [2017] L. Dai, T. Venumadhav, and K. Sigurdson, Phys. Rev. D 95, 044011 (2017), arXiv:1605.09398 [astro-ph.CO] .
  * Cusin and Tamanini [2021] G. Cusin and N. Tamanini, Mon. Not. Roy. Astron. Soc. 504, 3610 (2021), arXiv:2011.15109 [astro-ph.CO] .
  * Auclair  _et al._ [2022] P. Auclair _et al._ (LISA Cosmology Working Group), (2022), arXiv:2204.05434 [astro-ph.CO] .
  * Harte [2019] A. I. Harte, Gen. Rel. Grav. 51, 14 (2019), arXiv:1808.06203 [gr-qc] .
  * Andersson  _et al._ [2021] L. Andersson, J. Joudioux, M. A. Oancea, and A. Raj, Phys. Rev. D 103, 044053 (2021), arXiv:2012.08363 [gr-qc] .
  * Oancea  _et al._ [2022] M. A. Oancea, R. Stiskalek, and M. Zumalacárregui, (2022), arXiv:2209.06459 [gr-qc] .
  * Suyama  _et al._ [2005] T. Suyama, R. Takahashi, and S. Michikoshi, Phys. Rev. D 72, 043001 (2005), arXiv:astro-ph/0505023 .
  * Sesana  _et al._ [2021] A. Sesana _et al._ , Exper. Astron. 51, 1333 (2021), arXiv:1908.11391 [astro-ph.IM] .
  * Takahashi [2004] R. Takahashi, Astron. Astrophys. 423, 787 (2004), arXiv:astro-ph/0402165 .
  * Keeton [2001] C. R. Keeton, (2001), arXiv:astro-ph/0102341 .
  * Hinshaw and Krauss [1987] G. Hinshaw and L. M. Krauss, Astrophys. J. 320, 468 (1987).
  * Flores and Primack [1996] R. A. Flores and J. R. Primack, Astrophys. J. Lett. 457, L5 (1996), arXiv:astro-ph/9512063 .
  * Vallisneri [2008] M. Vallisneri, Phys. Rev. D 77, 042001 (2008), arXiv:gr-qc/0703086 .
  * Marsat  _et al._ [2021] S. Marsat, J. G. Baker, and T. Dal Canton, Phys. Rev. D 103, 083011 (2021), arXiv:2003.00357 [gr-qc] .
  * Husa  _et al._ [2016] S. Husa, S. Khan, M. Hannam, M. Pürrer, F. Ohme, X. Jiménez Forteza, and A. Bohé, Phys. Rev. D 93, 044006 (2016), arXiv:1508.07250 [gr-qc] .
  * Khan  _et al._ [2016] S. Khan, S. Husa, M. Hannam, F. Ohme, M. Pürrer, X. Jiménez Forteza, and A. Bohé, Phys. Rev. D 93, 044007 (2016), arXiv:1508.07253 [gr-qc] .
  * Hannuksela  _et al._ [2020] O. A. Hannuksela, T. E. Collett, M. Çalışkan, and T. G. F. Li, Mon. Not. Roy. Astron. Soc. 498, 3395 (2020), arXiv:2004.13811 [astro-ph.HE] .
  * [113] “Updated advanced ligo sensitivity design curve,” .
  * Tinker  _et al._ [2008] J. L. Tinker, A. V. Kravtsov, A. Klypin, K. Abazajian, M. S. Warren, G. Yepes, S. Gottlober, and D. E. Holz, Astrophys. J. 688, 709 (2008), arXiv:0803.2706 [astro-ph] .
  * Diemer [2018] B. Diemer, Astrophys. J. Suppl. 239, 35 (2018), arXiv:1712.04512 [astro-ph.CO] .
  * Aghanim  _et al._ [2020] N. Aghanim _et al._ (Planck), Astron. Astrophys. 641, A6 (2020), [Erratum: Astron.Astrophys. 652, C4 (2021)], arXiv:1807.06209 [astro-ph.CO] .
  * Wempe  _et al._ [2022] E. Wempe, L. V. E. Koopmans, A. R. A. C. Wierda, O. Akseli Hannuksela, and C. van den Broeck, arXiv e-prints , arXiv:2204.08732 (2022), arXiv:2204.08732 [astro-ph.HE] .
  * Rocha  _et al._ [2013] M. Rocha, A. H. G. Peter, J. S. Bullock, M. Kaplinghat, S. Garrison-Kimmel, J. Onorbe, and L. A. Moustakas, Mon. Not. Roy. Astron. Soc. 430, 81 (2013), arXiv:1208.3025 [astro-ph.CO] .
  * Vogelsberger  _et al._ [2012] M. Vogelsberger, J. Zavala, and A. Loeb, Mon. Not. Roy. Astron. Soc. 423, 3740 (2012), arXiv:1201.5892 [astro-ph.CO] .
  * Elbert  _et al._ [2015] O. D. Elbert, J. S. Bullock, S. Garrison-Kimmel, M. Rocha, J. Oñorbe, and A. H. G. Peter, Mon. Not. Roy. Astron. Soc. 453, 29 (2015), arXiv:1412.1477 [astro-ph.GA] .
  * Brinckmann  _et al._ [2018] T. Brinckmann, J. Zavala, D. Rapetti, S. H. Hansen, and M. Vogelsberger, Mon. Not. Roy. Astron. Soc. 474, 746 (2018), arXiv:1705.00623 [astro-ph.CO] .
  * Nadler  _et al._ [2020] E. O. Nadler, A. Banerjee, S. Adhikari, Y.-Y. Mao, and R. H. Wechsler, Astrophys. J. 896, 112 (2020), arXiv:2001.08754 [astro-ph.CO] .
  * Kaplinghat  _et al._ [2016] M. Kaplinghat, S. Tulin, and H.-B. Yu, Phys. Rev. Lett. 116, 041302 (2016), arXiv:1508.03339 [astro-ph.CO] .
  * Madau  _et al._ [2014] P. Madau, S. Shen, and F. Governato, Astrophys. J. Lett. 789, L17 (2014), arXiv:1405.2577 [astro-ph.GA] .
  * Oñorbe  _et al._ [2015] J. Oñorbe, M. Boylan-Kolchin, J. S. Bullock, P. F. Hopkins, D. Kerěs, C.-A. Faucher-Giguère, E. Quataert, and N. Murray, Mon. Not. Roy. Astron. Soc. 454, 2092 (2015), arXiv:1502.02036 [astro-ph.GA] .
  * Chan  _et al._ [2015] T. K. Chan, D. Kereš, J. Oñorbe, P. F. Hopkins, A. L. Muratov, C. A. Faucher-Giguère, and E. Quataert, Mon. Not. Roy. Astron. Soc. 454, 2981 (2015), arXiv:1507.02282 [astro-ph.GA] .
  * Tollet  _et al._ [2016] E. Tollet _et al._ , Mon. Not. Roy. Astron. Soc. 456, 3542 (2016), [Erratum: Mon.Not.Roy.Astron.Soc. 487, 1764 (2019)], arXiv:1507.03590 [astro-ph.GA] .
  * Graham  _et al._ [2015] P. W. Graham, I. G. Irastorza, S. K. Lamoreaux, A. Lindner, and K. A. van Bibber, Ann. Rev. Nucl. Part. Sci. 65, 485 (2015), arXiv:1602.00039 [hep-ex] .
  * Di Luzio  _et al._ [2020] L. Di Luzio, M. Giannotti, E. Nardi, and L. Visinelli, Phys. Rept. 870, 1 (2020), arXiv:2003.01100 [hep-ph] .
  * Schive  _et al._ [2014] H.-Y. Schive, M.-H. Liao, T.-P. Woo, S.-K. Wong, T. Chiueh, T. Broadhurst, and W. Y. P. Hwang, Phys. Rev. Lett. 113, 261302 (2014), arXiv:1407.7762 [astro-ph.GA] .
  * Dalal and Kravtsov [2022] N. Dalal and A. Kravtsov, Phys. Rev. D 106, 063517 (2022), arXiv:2203.05750 [astro-ph.CO] .
  * May and Springel [2021] S. May and V. Springel, Mon. Not. Roy. Astron. Soc. 506, 2603 (2021), 2101.01828 [astro-ph.CO] .
  * Mocz  _et al._ [2020] P. Mocz _et al._ , Mon. Not. Roy. Astron. Soc. 494, 2027 (2020), arXiv:1911.05746 [astro-ph.CO] .
  * Kormann  _et al._ [1994] R. Kormann, P. Schneider, and M. Bartelmann, Astronomy and Astrophysics (ISSN 0004-6361), vol. 284, no. 1, p. 285-299 284, 285 (1994).
  * Keeton  _et al._ [1997] C. R. Keeton, C. S. Kochanek, and U. Seljak, Astrophys. J. 482, 604 (1997), arXiv:astro-ph/9610163 .
  * Villaescusa-Navarro and Dalal [2011] F. Villaescusa-Navarro and N. Dalal, JCAP 03, 024 (2011), arXiv:1010.3008 [astro-ph.CO] .
  * Ishiyama [2014] T. Ishiyama, Astrophys. J. 788, 27 (2014), arXiv:1404.1650 [astro-ph.CO] .
  * Delos and White [2022] M. S. Delos and S. D. M. White,  (2022), 10.1093/mnras/stac3373, arXiv:2207.05082 [astro-ph.CO] .
  * Berezinsky  _et al._ [2010] V. Berezinsky, V. Dokuchaev, Y. Eroshenko, M. Kachelriess, and M. A. Solberg, Phys. Rev. D 81, 103529 (2010), arXiv:1002.3444 [astro-ph.CO] .
  * Li  _et al._ [2012b] F. Li, A. L. Erickcek, and N. M. Law, Phys. Rev. D 86, 043519 (2012b), arXiv:1202.1284 [astro-ph.CO] .
  * Dax  _et al._ [2021] M. Dax, S. R. Green, J. Gair, J. H. Macke, A. Buonanno, and B. Schölkopf, Phys. Rev. Lett. 127, 241103 (2021), arXiv:2106.12594 [gr-qc] .
  * Diego [2019] J. M. Diego, Astron. Astrophys. 625, A84 (2019), arXiv:1806.04668 [astro-ph.GA] .
  * Ohanian [1974] H. C. Ohanian, Int. J. Theor. Phys. 9, 425 (1974).
  * Matsunaga and Yamamoto [2006] N. Matsunaga and K. Yamamoto, JCAP 01, 023 (2006), arXiv:astro-ph/0601701 .
  * Çaǧan Şengül  _et al._ [2020] A. Çaǧan Şengül, A. Tsang, A. Diaz Rivero, C. Dvorkin, H.-M. Zhu, and U. Seljak, Phys. Rev. D 102, 063502 (2020), arXiv:2006.07383 [astro-ph.CO] .
  * Şengül  _et al._ [2021] A. c. Şengül, C. Dvorkin, B. Ostdiek, and A. Tsang, (2021), arXiv:2112.00749 [astro-ph.CO] .
  * Fleury  _et al._ [2021] P. Fleury, J. Larena, and J.-P. Uzan, JCAP 08, 024 (2021), arXiv:2104.08883 [astro-ph.CO] .
  * Feldbrugge [2020] J. Feldbrugge, (2020), arXiv:2010.03089 [astro-ph.CO] .
  * Çalışkan  _et al._ [2023] M. Çalışkan, N. Anil Kumar, L. Ji, J. M. Ezquiaga, R. Cotesta, E. Berti, and M. Kamionkowski, (2023), arXiv:2307.06990 [astro-ph.CO] .
  * Feldbrugge  _et al._ [2019] J. Feldbrugge, U.-L. Pen, and N. Turok, (2019), arXiv:1909.04632 [astro-ph.HE] .
  * Maggiore [2007] M. Maggiore,  _Gravitational waves: Volume 1: Theory and experiments_ (OUP Oxford, 2007).