# Constraining runaway dilaton models using joint gravitational-wave and
electromagnetic observations

Arnab Dhani<EMAIL_ADDRESS>Institute for Gravitation and the Cosmos,
Department of Physics, Pennsylvania State University, University Park, PA
16802, USA  Anuradha Gupta<EMAIL_ADDRESS>Department of Physics and
Astronomy, The University of Mississippi, University, Mississippi 38677, USA
B. S. Sathyaprakash<EMAIL_ADDRESS>Institute for Gravitation and the Cosmos,
Department of Physics, Pennsylvania State University, University Park, PA
16802, USA  Department of Astronomy & Astrophysics, Pennsylvania State
University, University Park, PA 16802, USA  School of Physics and Astronomy,
Cardiff University, Cardiff, UK, CF24 3AA

###### Abstract

With the advent of gravitational-wave astronomy it has now been possible to
constrain modified theories of gravity that were invoked to explain the dark
energy. In a class of dilaton models, distances to cosmic sources inferred
from electromagnetic and gravitational wave observations would differ due to
the presence of a friction term. In such theories, the ratio of the Newton’s
constant to the fine structure constant varies with time. In this paper we
explore the degree to which it will be possible to test such models. If
collocated sources (e.g. supernovae and binary neutron star mergers), but not
necessarily multimessengers, can be identified by electromagnetic telescopes
and gravitational-wave detectors one can probe if light and gravitational
radiation are subject to the same laws of propagation over cosmological
distances. This helps in constraining the variation of Newton’s constant
relative to fine-structure constant. The next generation of gravitational wave
detectors, such as the Cosmic Explorer and Einstein Telescope, in tandem with
the Vera Rubin Observatory and gamma ray observatories such as the Fermi Space
Observatory will be able to detect or constrain such variations at the level
of a few parts in 100. We apply this method to GW170817 with distances
inferred by the LIGO and Virgo detectors and the observed Kilonova.

Gravitational waves, electromagnetic waves, fine structure constant,
gravitational constant

##  I Introduction

Gravitational waves (GWs) and electromagnetic (EM) waves follow the same
propagation equations in General Relativity (GR) Maggiore (2007).
Consequently, the various distance measures in cosmology (e.g., luminosity
distance, angular diameter distance, comoving distance, etc.) are identical
for both GW and EM. Several alternative theories of gravity with additional
scalar degrees of freedom Fujii and Maeda (2007); Bertolami  _et al._ (2007,
2008); Bertolami and Paramos (2008); Sotiriou and Faraoni (2008); De Felice
and Tsujikawa (2010); Harko  _et al._ (2013); Das and Banerjee (2008); Bisabr
(2012); Moffat and Toth (2012); Shiralilou  _et al._ (2022) modify the
propagation of either or both by altering the friction term in the wave
equations due to the evolution of the scalar field. We will, however, restrict
to a class of scalar-tensor theories in which the dispersion relation remains
unchanged. Hence, in these scalar-tensor theories distance to an astronomical
source inferred from gravitational-wave observation will be different from
that inferred using electromagnetic radiation.

The presence of a scalar field is also motivated by the low energy effective
field theories of Loop Quantum Gravity Rovelli and Smolin (1994); Domagala
_et al._ (2010) and String Theory Green  _et al._ (1988); Uzan (2011); Damour
and Polyakov (1994a, b); Gasperini  _et al._ (2002); Minazzoli and Hees
(2013). Furthermore, dark energy Ratra and Peebles (1988); Caldwell  _et al._
(1998); Peebles and Ratra (2003), inflation Guth (1981); Linde (1982);
Albrecht and Steinhardt (1982); Linde (2008), and variations of the
fundamental constants are often modeled using a scalar field Bekenstein
(1982); Sandvik  _et al._ (2002); Dvali and Zaldarriaga (2002); Olive and
Pospelov (2008); Damour (2012). In fact, it has been claimed that the
requirement of gauge and diffeomorphism invariances would invariably lead to
scalar-tensor theories with minimal/non-minimal coupling to the matter sector
Armendariz-Picon (2002). The coupling of the scalar field to the gravitational
sector in such theories has been tightly constrained in the weak-field limit
using solar system tests Adelberger  _et al._ (2003, 2007, 2009); Kapner  _et
al._ (2007); Will (2006). If the scalar field couples non-minimally to the
matter sector, the Einstein equivalence principle is broken. The equivalence
principle, likewise, has been tested to a very high accuracy within the solar
system Rosenband  _et al._ (2008); Will (2006); Adelberger  _et al._ (2009);
Williams  _et al._ (2012). A variety of decoupling Tseytlin and Vafa (1992);
Damour and Vilenkin (1996); Damour  _et al._ (2002); Jarv  _et al._ (2008);
Damour and Nordtvedt (1993); Damour and Polyakov (1994b); Minazzoli and Hees
(2013) or screening Khoury (2010); Khoury and Weltman (2004a, b); Hees and
Fuzfa (2012); Hinterbichler and Khoury (2010); Hinterbichler  _et al._ (2011)
mechanisms have, therefore, been proposed to keep these theories viable for
cosmological evolution.

The propagation of waves on a modified background allows one to test for the
presence of a scalar field on cosmological scales. High redshift quasar
absorption spectra Webb  _et al._ (2001); King  _et al._ (2012); Webb  _et
al._ (2011), galaxy clustering data Holanda  _et al._ (2016), and 21cm neutral
hydrogen intensity mapping Khatri and Wandelt (2007) have been used to place
limits on the spatio-temporal evolution of the fine structure constant which
can be modeled using a scalar field. Type Ia supernova (SNeIa) data is used to
fit the EM luminosity distance-redshift relation and constrain models of
dynamical dark energy Riess  _et al._ (2016). Other studies use the EM
luminosity distance estimates from SNeIa in parallel with the EM angular
diameter distance measurements from X-ray and Sunyaev-Zel’dovich observations
of galaxy clusters to directly constrain the violation of the distance-duality
relation in the EM sector Cao and Liang (2011); Hees  _et al._ (2014).

Gravitational wave astronomy has opened a new means of revealing the presence
of a scalar degree of freedom. Coincident measurements of the luminosity
distance from GW observations and the redshift from follow-up EM observations
of “bright” sirens, such as the first observation of gravitational waves from
a binary neutron star merger, GW170817 Abbott  _et al._ (2017a, b), have been
used to put limits on the modified friction term in $f(R)$ and scalar-tensor
theories of gravity Fanizza  _et al._ (2020); Finke  _et al._ (2021). The
luminosity distance-redshift relation has also been constrained for “dark”
sirens (GW observations without an EM counterpart) by cross-correlations with
galaxy catalogs Finke  _et al._ (2021); Mukherjee  _et al._ (2020). In these
methods, the modified friction term is constrained together with the standard
cosmological parameters, however, the two sets of parameters are strongly
correlated with each other. Mukherjee  _et al._ (2020) propose the use of
Baryon acoustic oscillation (BAO) data together with luminosity-distance
measurements from GW observations and redshifts from galaxy catalog cross-
correlations to directly constrain the ratio between the GW luminosity-
distance and EM luminosity distance in terms of the modified friction
parameter.

In this study, we propose the direct use of the EM luminosity distance from
SNeIa/kilonova concurrently with the GW luminosity distance from “bright”
sirens and the redshift obtained from photometric/spectroscopic studies of the
identified galaxy or galaxy cluster to directly constrain the ratio of the two
luminosity distances for a class of scalar-tensor theories with a non-minimal
multiplicative coupling between the scalar field and the matter sector. The
crucial distinction with the method described in Mukherjee  _et al._ (2020) is
their use of the BAO data to convert the angular diameter distance to EM
luminosity distance via the distance-duality relation, which is broken for us
due to the non-minimal coupling of the scalar field to the matter sector. In
other words, their procedure is valid for alternative theories of gravity in
which gravity is minimally coupled to the matter sector whereas our method
applies to more general theories. Furthermore, they infer the redshifts to GW
sources using galaxy correlation and as a result also measure some
cosmological parameters. We restrict ourselves to “bright” sirens and,
therefore, have a direct measurement of the redshift. In this way, our
parameter constraints do not suffer from degeneracies with the other
cosmological parameters.

The class of scalar-tensor theories considered in this study arise as low
energy action of string theories and satisfy the solar system tests for both
the modifications to the gravitational sector and the breakage of the
equivalence principle. This class of theories, known as the runaway dilaton
models Gasperini  _et al._ (2002); Damour  _et al._ (2002); Minazzoli and Hees
(2014); Hees  _et al._ (2014), has a Brans-Dicke type gravitational
interaction and a universal multiplicative coupling between the scalar field
and the matter sector which breaks the equivalence principle. The unequal
coupling of the scalar field to the metric and the matter sector leads to
(distinct) modified propagation equations for gravitational waves and
electromagnetic waves.

We parameterize the ratio of the electromagnetic and gravitational-wave
luminosity distance using a parameter $\eta_{0}$. We find that the planned
upgrades to the second-generation of advanced gravitational-wave detector
networks (e.g. the A+ upgrade Reitze  _et al._ (2019); Abbott  _et al._ (2018)
of Advanced LIGO and similar upgrades to Advanced Virgo Acernese  _et al._
(2015a), KAGRA Aso  _et al._ (2013); Somiya (2012) and LIGO-India Unnikrishnan
(2013); Saleem  _et al._ (2022)) constrains $\eta_{0}$ to $|\eta_{0}|<0.2,$
while the proposed improvement of the network to Voyager sensitivity Adhikari
_et al._ (2020a) refines the constraint to $|\eta_{0}|<0.05$. The proposed
third-generation of ground-based gravitational-wave detector network (Cosmic
Explorer Evans  _et al._ (2021); Reitze  _et al._ (2019) and Einstein
Telescope Punturo  _et al._ (2010a, b); Hild  _et al._ (2011)) will place the
best limits on $\eta_{0}$ at $|\eta_{0}|<0.01$.

In Sec. II, we briefly describe runaway dilaton models and their EM and GW
propagation equations. We also describe how the ratio of the luminosity
distance of each sector can be related to the variation of the fundamental
constants. In Sec. III, we discuss the gravitational-wave detectors considered
in this study, the simulations we performed, and the electromagnetic data that
we used. We outline our main results and forecasts in Sec. IV and the
constraints that can be placed using GW170817 in Sec. V. Sec. VI concludes the
paper.

##  II Background

In this section, we briefly review the equations of motion for runaway dilaton
models, derive the propagation equations for electromagnetic and gravitational
waves on a homogeneous and isotropic background, parameterize the ratio of the
luminosity distances as a function of redshift, and relate it to the redshift
variation of the fundamental constants.

###  II.1 Runaway dilaton models

The action for runaway dilaton models Gasperini  _et al._ (2002); Damour  _et
al._ (2002); Minazzoli and Hees (2014); Hees  _et al._ (2014), a class of
scalar-tensor theories with a generic multiplicative coupling $h(\phi)$
between a scalar field $\phi$ and the matter Lagrangian
$\mathcal{L}_{m}[g_{\mu\nu,\Psi}]$, is given by

$S=\int d^{4}x\sqrt{-g}\Bigg{[}\frac{1}{2\kappa}\left(\phi
R-\frac{\omega(\phi)}{\phi}\nabla_{\mu}\phi\nabla^{\mu}\phi-V(\phi)\right)\\\
+h(\phi)\mathcal{L}_{m}[g_{\mu\nu,\Psi}]\Bigg{]}\,,$ (1)

where $\kappa=8\pi G$ with $G$ being the gravitational coupling constant, $R$
is the Ricci scalar, and $\Psi$ consists of all the Standard Model fields.

The gravitational equations of motion, given by the variation of the action
with respect to the metric, takes the form,

$R_{\mu\nu}-\frac{1}{2}g_{\mu\nu}R=\kappa\frac{h(\phi)}{\phi}T_{\mu\nu}+\frac{1}{\phi}(\nabla_{\mu}\nabla_{\nu}-g_{\mu\nu}\Box)\phi\\\
+\frac{\omega(\phi)}{\phi^{2}}(\nabla_{\mu}\phi\nabla_{\nu}\phi-\frac{1}{2}g_{\mu\nu}\nabla_{\alpha}\phi\nabla^{\alpha}\phi)-g_{\mu\nu}\frac{V(\phi)}{2\phi}.$
(2)

Similarly, one can obtain the equation of motion for the scalar field by
varying the action with respect to it. Upon replacing the Ricci scalar in the
resulting equation with the trace of the gravitational equations of motion Eq.
(2), one finds that the Klein-Gordon equation for the scalar field is given
by,

$\frac{2\omega(\phi)+3}{\phi}\Box\phi=\kappa\left(\frac{h(\phi)}{\phi}T^{\alpha}_{\alpha}-2h^{\prime}(\phi)\mathcal{L}_{m}\right)\\\
-\frac{\omega^{\prime}(\phi)}{\phi}\nabla_{\alpha}\phi\nabla^{\alpha}\phi+V^{\prime}(\phi)-2\frac{V(\phi)}{\phi},$
(3)

where a prime denotes the derivative with respect to $\phi$. The stress-energy
tensor $T_{\mu\nu}$ in the above equations can be defined by the variation of
the matter Lagrangian with the metric $g_{\mu\nu}$,

$T_{\mu\nu}\equiv\frac{-2}{\sqrt{-g}}\frac{\delta}{\delta
g^{\mu\nu}}(\sqrt{-g}\mathcal{L}_{m}).$ (4)

We will consider the background to be a homogeneous, isotropic, and spatially
flat Universe described by the Friedmann-Lemaître-Robertson-Walker (FLRW)
metric,

$ds^{2}=-dt^{2}+a(t)^{2}\delta_{ij}dx^{i}dx^{j},$ (5)

where the size of the homogeneous, isotropic, and spatially flat 3-surface is
given by the scale factor $a(t)$. The background spacetime is considered to be
sourced by a perfect fluid with its stress-energy tensor given by

$T^{\mu\nu}=(\rho+p)u^{\mu}u^{\nu}+pg^{\mu\nu},$ (6)

where $\rho$ is the total energy density and $p$ is the pressure of the fluid
in its rest-frame, and $u^{\mu}$ is the 4-velocity of the fluid with respect
to an observer.

The Friedmann equations that describe the evolution of the background
spacetime are obtained by substituting Eq. (5) and Eq. (6) in to the
gravitational field equations (2):

$H^{2}=\kappa\frac{h(\phi)}{3\phi}\rho+\frac{V(\phi)}{6\phi}+\frac{\omega(\phi)}{6}\left(\frac{\dot{\phi}}{\phi}\right)^{2}-H\frac{\dot{\phi}}{\phi},$
(7)
$2\dot{H}+3H^{2}=-\kappa\frac{h(\phi)}{\phi}p+\frac{V(\phi)}{2\phi}-2H\frac{\dot{\phi}}{\phi}-\frac{\omega(\phi)}{2}\left(\frac{\dot{\phi}}{\phi}\right)^{2}-\frac{\ddot{\phi}}{\phi},$
(8)

where $H(t)$ is the Hubble parameter defined as $H(t)\equiv\dot{a}(t)/a(t)$
and dots denote derivatives with respect to the time coordinate $t$. In Eqs.
(7) and (8), if the scalar field $\phi$ is a constant, only the first two
terms on the right-hand side are non-zero and we recover the standard
Friedmann equations for $\Lambda$CDM cosmology up to the normalization of the
field $\phi$.

###  II.2 Propagation of gravitational waves

In this and the following subsection, we will derive the equations describing
the propagation of gravitational and electromagnetic waves, respectively, on
the background spacetime.

Gravitational waves are propagating tensor perturbations of the background
spacetime. To get the equations of motion for tensor perturbations, we perturb
our FLRW metric as

$ds^{2}=-dt^{2}+a(t)^{2}(\delta_{ij}+h_{ij})dx^{i}dx^{j},$ (9)

where $h_{ij}$ is a small perturbation of the background geometry and is
transverse ($\partial_{i}h^{ij}=0$) and trace-less ($h_{i}^{i}=0$) in the
chosen coordinate system. Note that this is not a generic perturbation of the
background. A generic perturbation can be decomposed in to scalar, vector, and
tensor components that do not mix under diffeomorphisms. Furthermore, at the
leading order, the equations of motion for these components are decoupled.
Here, since we are only interested in GWs, it is sufficient to perturb the
background with the tensor component which is transverse and trace-less.

The equations of motion for gravitational-wave propagation are then given by

$\ddot{h}_{ij}+\left(3H+\frac{\dot{\phi}}{\phi}\right)\dot{h}_{ij}-\frac{\nabla^{2}h_{ij}}{a(t)^{2}}=0,$
(10)

where $\dot{\phi}/\phi$ is the modified friction term that would change the
observed GW amplitude and as a result the luminosity distance with respect to
GR. We note that the luminosity distance is additionally modified since the
Friedmann equations get altered due to the presence of the scalar field. In
other words, the evolution of the scale factor $a(t)$ is different from that
in GR. Note, however, that the dispersion relation is unchanged with respect
to GR and, hence, GWs travel at the speed of light.

Throughout this study, we are interested in solutions of the wave equations
under geometric optics approximation. This is because the length scales of the
signals of interest to us (stellar-mass compact binary mergers and SNeIa) are
much smaller than the Hubble scale. In this limit, the metric perturbations
can be written as

$h_{ij}=\mathcal{R}\\{(b_{ij}+\epsilon
c_{ij}+\mathcal{O}(\epsilon^{2}))\,e^{i\theta/\epsilon}\\},$ (11)

where $\theta$ is the phase of the plane wave and $\epsilon$ is an order-
counting parameter, which can be set to 1 at the end of the calculation.
Substituting Eq. (11) into Eq. (10) and collecting terms of the same order, we
get

$\displaystyle k_{\mu}k^{\mu}$ $\displaystyle=0,$ (12) $\displaystyle
k^{\nu}\nabla_{\nu}k^{\mu}$ $\displaystyle=0,$

at the $\mathcal{O}(\epsilon^{2})$ and

$\nabla_{\mu}(b^{2}k^{\mu})=-b^{2}k^{\mu}\nabla_{\mu}\ln\phi,$ (13)

at the $\mathcal{O}(\epsilon)$, where $k_{\mu}=\nabla_{\mu}\theta$, the wave
vector, is null and follows null geodesics and $b=||b_{ij}||$ is the Euclidean
norm of the leading-order amplitude. The latter equation is the one which is
modified with respect to GR and represents the non-conservation of the
graviton number as it propagates on the background spacetime.

The luminosity distance can then be calculated following Minazzoli and Hees
(2014)111The derivation of the luminosity distance is carried out for EM waves
but the procedure is the same for GWs. and is given by

$d_{L}^{\rm GW}=(1+z)\sqrt{\frac{\phi_{0}}{\phi}}\int_{0}^{z}\frac{dz}{H(z)},$
(14)

where where $\phi_{0}$ is the value of the field in the present epoch and
$H(z)$ is the modified Hubble relation [Eq. (7)].

###  II.3 Propagation of electromagnetic waves

The field equations that govern the propagation of electromagnetic waves can
be obtained by the variation of the action Eq. (1) with respect to the EM
4-potential $A^{\mu}$ and are given by

$\nabla_{\nu}(h(\phi)F^{\mu\nu})=0,$ (15)

where $F_{\mu\nu}=\nabla_{\mu}A_{\nu}-\nabla_{\nu}A_{\mu}$ is the
electromagnetic field tensor.

Solving the above equation in the geometric optics limit, as for metric
perturbations, yields similar equations for photons, namely, they travel on
null geodesics. The equation for photon number ‘non-conservation’ is given by:

$\nabla_{\mu}(\overline{b}^{2}k^{\mu})=-\overline{b}^{2}k^{\mu}\nabla_{\mu}\ln
h(\phi),$ (16)

where the electromagnetic potential in the geometric optics limit is given by

$A^{\mu}=\mathcal{R}\\{(b^{\mu}+\epsilon
c^{\mu}+\mathcal{O}(\epsilon^{2}))\,e^{i\theta/\epsilon}\\},$ (17)

with $\overline{b}=||b^{\mu}||$ being the norm taken with respect to the
background FLRW metric.

The luminosity distance is then given by Minazzoli and Hees (2014)

$d_{L}^{\rm
EM}=(1+z)\sqrt{\frac{h(\phi_{0})}{h(\phi)}}\int_{0}^{z}\frac{dz}{H(z)}\,,$
(18)

where $\phi_{0}$ is the value of the field at the present epoch and $H(z)$ is
given by the modified Friedman equations. Of note here is that the two
luminosity distances differ only if $h(\phi)\neq\phi$ which is the premise we
are working under.

###  II.4 Parameterizing the modified luminosity distances

Now that we have obtained the luminosity distance-redshift relation for both
electromagnetic and gravitational waves, we can parameterize the scalar field
dependence of their ratio. We choose to do this parameterization in terms of
the violation of the distance-duality relation for both the sectors. This
choice helps connect our results to those of the numerous experiments in the
electromagnetic sector that parameterize this deviation Hees  _et al._ (2014).
The distance-duality relation connects the luminosity distance to the angular
diameter distance. The latter is defined by:

$d_{A}(z)=\frac{1}{1+z}\int_{0}^{z}\frac{dz}{H(z)}.$ (19)

This is a geometric quantity that can be derived by integrating the geodesic
equation. For the class of theories considered in this study, both
gravitational and electromagnetic waves travel on null geodesics and,
therefore, their angular diameter distances are unchanged from that of GR,
apart from a modification to the Friedmann equations. We can, then, write the
parameterization as

$\eta(z)=\frac{d_{L}(z)}{d_{A}(z)(1+z)^{2}}.$ (20)

In GR, the distance-duality relation implies $\eta(z)=1$.

We consider $\eta(z)$ to have the functional form

$\displaystyle\eta^{\rm
GW}(z)=\sqrt{\frac{\phi_{0}}{\phi}}=1+\eta_{1}\frac{z}{1+z},$ (21)
$\displaystyle\eta^{\rm{}EM}(z)=\sqrt{\frac{h(\phi_{0})}{h(\phi)}}=1+\eta_{2}\frac{z}{1+z},$

parameterizing the deviation from the gravitational (electromagnetic)
distance-duality relation by $\eta_{1}$ ($\eta_{2}$). Additionally, we
parameterize the ratio of the luminosity distances as

$\frac{d_{L}^{\rm GW}}{d_{L}^{\rm EM}}=\frac{\eta^{\rm GW}(z)}{\eta^{\rm
EM}(z)}=\sqrt{\frac{\phi_{0}/h(\phi_{0})}{\phi/h(\phi)}}=1+\eta_{0}\frac{z}{1+z}$
(22)

from which one can deduce that $\eta_{0}\approx\eta_{1}-\eta_{2}$. The above
form of the parameterization was introduced in Holanda  _et al._ (2012) with
the advantage being that it avoids divergence at large redshifts which the
linear expansions suffer from. Given a simultaneous measurement of the GW and
EM luminosity distances, either from the same source or the same galaxy or the
same galaxy cluster (see Sec. III.2.3 and VI for a discussion), one can place
constraints on the parameter $\eta_{0}$.

At this point, we note that studies in Mukherjee  _et al._ (2020); Fanizza
_et al._ (2020); Finke  _et al._ (2021) have constrained the ratio of the two
luminosity distances, albeit in the context of modifying the frictional term
in the gravitational sector alone, through the parameterization,

$\frac{d_{L}^{\rm GW}(z)}{d_{L}^{\rm
EM}(z)}=\Xi_{0}+\frac{1-\Xi_{0}}{(1+z)^{n}}.$ (23)

where $\Xi_{0}=1$ in GR and $n$ gives the rate at which the ratio saturates to
its asymptotic value $\Xi_{0}$. Our parameter $\eta_{0}$ is related to
$(\Xi_{0},n)$ as

$\eta_{0}=1-\Xi_{0}\quad{\rm for}\quad n=1,$ (24)

i.e., our parameterization is a subclass of the $(\Xi_{0},n)$ parameterization
for a fixed saturation rate.

The errors on $\eta_{0}$ can be calculated from the errors on the GW and EM
luminosity distances, assuming the redshift to the source is known, using the
standard error propagation formula for independent variables ($d_{L}^{\rm GW}$
and $d_{L}^{\rm EM}$ in this case) as,

$\sigma_{\eta_{0}}^{2}=\left(\frac{\partial\eta_{0}}{\partial d_{L}^{\rm
GW}}\right)^{2}\sigma_{d_{L}^{\rm
GW}}^{2}+\left(\frac{\partial\eta_{0}}{\partial d_{L}^{\rm
EM}}\right)^{2}\sigma_{d_{L}^{\rm EM}}^{2},$ (25)

where $\sigma_{X}$ denotes $1$-$\sigma$ error in the quantity $X$. Simplifying
the above equation by evaluating the derivative expressions leads to

$\sigma_{\eta_{0}}=\frac{1+z}{z}\frac{d_{L}^{\rm GW}}{d_{L}^{\rm
EM}}\sqrt{\left(\frac{\sigma_{d_{L}^{\rm GW}}}{d_{L}^{\rm
GW}}\right)^{2}+\left(\frac{\sigma_{d_{L}^{\rm EM}}}{d_{L}^{\rm
EM}}\right)^{2}}.$ (26)

###  II.5 Redshift variation of fundamental constants

The non-minimal coupling of the scalar field to the matter and gravitational
sectors leads to the dependence of the fundamental constants on the scalar
field and they, therefore, evolve with the evolution of the scalar field Yunes
_et al._ (2010, 2016); Vijaykumar  _et al._ (2021). From the action given by
Eq. (1), it can be read out that the fine structure constant $\alpha$ and the
gravitational constant $G$ depend on the scalar field via

$\displaystyle\alpha\sim h^{-1}(\phi),$ (27) $\displaystyle G\sim\phi^{-1},$

and, hence, their redshift variation can be written as

$\displaystyle\frac{\Delta\alpha(z)}{\alpha_{0}}\equiv\frac{\alpha(z)-\alpha_{0}}{\alpha_{0}}=\frac{h(\phi_{0})}{h(\phi)}-1=\eta^{\rm
EM}(z)^{2}-1$ (28) $\displaystyle\frac{\Delta
G(z)}{G_{0}}=\frac{G(z)-G_{0}}{G_{0}}=\frac{\phi_{0}}{\phi}-1=\eta^{\rm
GW}(z)^{2}-1\,,$

where $\alpha_{0}$ and $G_{0}$ are the values of $\alpha$ and $G$ at the
current epoch, respectively.

Given the experimental constraints on the ratio of the two luminosity
distances, one can constrain the temporal variation of $\frac{G}{\alpha}(z)$
in the current epoch as

$\beta\equiv\left.\frac{\tfrac{d}{dt}(G/\alpha)}{G/\alpha}\right|_{0}=-2H_{0}\left.\frac{d\eta}{dz}\right|_{0}=-2H_{0}\eta_{0}\,,$
(29)

where $\eta(z)=\eta^{\rm GW}(z)/\eta^{\rm EM}(z)$ and $H_{0}$ is the present
value of the Hubble parameter. If one uses constraints from other
electromagnetic probes Hees  _et al._ (2014), the temporal variations of both
$\alpha$ and $G$ can be separately constrained.

We point out here that $G_{0}$ is not the effective gravitational constant
$G_{\rm eff}$ that enters the Poisson equation at the Newtonian order and
should not be interpreted as the strength of the gravitational force between
two test masses separated by a unit distance. The two are related by

$G_{\rm
eff}=G_{0}\left(1+\frac{1-2\phi_{0}\frac{h^{\prime}(\phi_{0})}{h(\phi_{0})}}{2\omega(\phi_{0})+3}\right)\frac{h(\phi_{0})}{\phi_{0}}.$
(30)

In the absence of the scalar field, $G_{0}$ and $G_{\rm eff}$ coincide, as
expected.

##  III Method

In this section, we describe the different gravitational wave detector
networks considered in this study, outline our procedure for simulating
gravitational-wave sources, calculate the rate of spatially coincident EM and
GW signals, and estimate the distribution of luminosity distance errors for
the coincidentally observed population of sources.

###  III.1 Gravitational wave detector networks

We consider three GW detector networks across three technology generations.
The 2G+ network consists of the five second generation GW detectors with three
LIGO detectors Aasi  _et al._ (2015) (LIGO-Hanford, LIGO-Livingston, LIGO-
India) operating at A+ sensitivity, the Virgo Acernese  _et al._ (2015b) and
the KAGRA Akutsu  _et al._ (2019) detectors at AdV+ and KAGRA+ sensitivities,
respectively. The Voy+ network includes the same five second generation
detectors but with the LIGO detectors upgraded to a proposed ‘Voyager’
Adhikari  _et al._ (2020b) technology. The final network, ECC, consists of
three proposed third generation detectors, specifically, two Cosmic Explorer
Reitze  _et al._ (2019) detectors and an Einstein Telescope Punturo  _et al._
(2010b). We show the noise power spectral densities (PSDs) for the individual
detectors in Fig. 1. The locations of these detectors and the technologies
used in a network are given in Table 1.

Figure 1: The noise power spectral density (PSD) estimates used for the individual detectors considered in this study. We use a low frequency cutoff of 5Hz for all but the advanced Virgo detector for which the PSD starts at 10Hz. Network label |  Detector location (technology)
---|---
2G+ |  Hanford WA (A+), Livingston LA (A+), Cascina Italy (AdV+), Kamioka Japan (KAGRA+), Hingoli India (A+)
Voy+ |  Hanford WA (Voyager), Livingston LA (Voyager), Cascina Italy (AdV+), Kamioka Japan (KAGRA+), Hingoli India (Voyager)
ECC |  Cascina Italy (ET-D), fiducial US site (CE1_cb), fiducial Australian site (CE1_cb)
Table 1: An overview of the three networks used in the study. The location
determines the detector antenna patterns, while the technology indicates the
used power spectral density. The Voyager and Cosmic Explorer power spectral
densities are chosen to be low-frequency optimized and in the case of the
latter for a detector arm length of $40\,\text{km}$.

###  III.2 Rates

####  III.2.1 Binary neutron star merger rates

We simulate a population of binary neutron star (BNS) merger events up to a
redshift of $z=1.5$ assuming a uniform mass distribution between 1$M_{\odot}$
and 2.5$M_{\odot}$ for the individual NSs Abbott  _et al._ (2020). The other
parameters, cosine of the inclination angle $\cos\iota$, location of the
source on the plane of the sky $\Omega$ (cosine of the declination angle
$\cos\delta$ and right ascension $\alpha$), polarization angle $\psi$, and the
phase of the signal at coalescence $\phi_{0}$, of the fiducial BNS population
are drawn from a uniform distribution across their domains. We assume 10 years
of observing time for each network with an 80% duty cycle for each detector
Belgacem  _et al._ (2019). The redshift distribution for our BNS population is
given by the following probability distribution,

$p(z)=\frac{R_{z}(z)}{\int_{0}^{10}R_{z}(z)dz},$ (31)

where an upper limit of $z=10$ is justified since the contribution to the
integral from redshifts larger than 10 is negligible. $R_{z}(z)$, the merger
rate density in the observer frame, can be expressed as

$R_{z}(z)=\frac{R_{m}(z)}{1+z}\frac{dV(z)}{dz}.$ (32)

Here $R_{m}(z)$ is the merger rate per comoving volume in the source frame and
$dV/dz$ is the comoving volume element. The former is given by

$R_{m}(z)=\int_{t_{\rm min}}^{t_{\rm max}}R_{f}[t(z)+t_{d}]P(t_{d})dt_{d},$
(33)

where $R_{f}(t)$ is the binary star formation rate (SFR) which we assume
follows the Vangioni cosmic SFR Vangioni  _et al._ (2015). The delay time (the
time it takes for a binary to coalesce after formation) distribution is taken
to be $P(t_{d})\propto 1/t_{d}$ with $t_{\rm min}=20\,\rm Myr$ and $t_{\rm
max}$ set to the Hubble time $1/H_{0}$. The value of $R_{m}$ at $z=0$ is
estimated from the population properties of the third LIGO–Virgo
Gravitational-Wave Transient Catalog, GWTC-3 Abbott  _et al._ (2021) to be
between

$R_{m}(z=0)=\mbox{13--1900}\;\rm Gpc^{-3}yr^{-1}.$ (34)

We present results for both the optimsitic and pessimistic local merger rates.

####  III.2.2 Electromagnetic counterpart

The Voy+ and the ECC network of GW detectors have a reach beyond the horizon
distance of the current and future EM telescopes for kilonova which can be
observed up to a redshift of about $z=0.5$ (see, e.g. Table 2.2 in Ref.
Kalogera  _et al._ (2021)). Hence, BNS events beyond a redshift of $z=0.5$ is
electromagnetically observable only through short gamma ray burst events.
Therefore, in this study, we assume that 10% of the BNS events up to a
redshift of 0.5 will have a dedicated EM follow-up search to detect their
kilonova emissions (we assume this to be in addition to possible GRB
detection, which do not need a dedicated search owing to the near all-sky
sensitivity of GRB detectors) and for BNS observations beyond a redshift of
0.5, we assume a coincident electromagnetic detection to consist of only GRBs.

Figure 2: Number of GW (left panel) and GW+GRB (right panel) detections as a
function of redshift for different detector networks considered in this study.
The optimistic case (solid lines) represents the upper limit on the local BNS
merger rate and the pessimistic case (dashed lines) the lower limit. The range
of values for the local BNS merger rate is given in Eq. (34). The lifetime of
a network is assumed to be 10 years with an 80% duty cycle for each detector.

We calculate the rate of a coincident GRB detection following the procedure
outlined in Belgacem  _et al._ (2019) and sketched it out here for
completeness. We assume a Gaussian structured jet profile Howell  _et al._
(2018) for a GRB burst and the luminosity $L(\theta_{V})$ is given by

$L(\theta_{V})=L_{p}\exp\left(-\frac{\theta_{V}^{2}}{2\theta_{c}^{2}}\right),$
(35)

where $\theta_{V}$ is the viewing angle and $\theta_{c}=4.7\degree$ represents
the variation in the GRB jet opening angle. $L_{p}$ is the peak luminosity of
each burst assuming isotropic emission in the rest frame in the $1-10^{4}$ keV
energy range and can be sampled from the probability distribution

$\Phi(L_{p})\propto\begin{cases}(L_{p}/L_{*})^{\alpha},\qquad L_{p}<L_{*},\\\
(L_{p}/L_{*})^{\beta},\qquad L_{p}\geq L_{*},\end{cases}$ (36)

where the parameters of the broken power-law distribution are $L_{*}=2\times
10^{52}\,\rm erg/s$, $\alpha=-1.95$, and $\beta=-3$ Wanderman and Piran
(2015). A GRB is assumed to be detected if the observed peak flux
$F_{P}(\theta_{V})=L(\theta_{V})/4\pi d_{L}^{2}$, given the GW luminosity
distance and inclination angle, is greater than the flux limit of
$1.1\rm\,ph\,s^{-1}\,cm^{-2}$ Belgacem  _et al._ (2019) in the $50$–300 keV
band for Fermi-GBM. The total time-averaged observable sky fraction for the
Fermi-GBM is taken to be 0.6 Burns  _et al._ (2016).

####  III.2.3 Rates for spatially coincident SNeIa

Following the arguments presented in Sec. 3 of Gupta  _et al._ (2019), we now
estimate the rates for a spatial coincidence of SNeIa and BNSs in a galaxy
cluster. Gupta  _et al._ (2019) concluded that the rate of spatial coincidence
of SNeIa and BNS mergers in a galaxy given their rates Li  _et al._ (2011);
Abbott  _et al._ (2019) is extremely small. Moreover, in their Sec. 5 Gupta
_et al._ (2019) showed that there is ${\cal O}(1\%)$ error in the distance
estimation of SNeIa if calibrated through a BNS in the same galaxy cluster
instead of the same galaxy. Therefore, coincident of SNeIa with a BNS in the
same galaxy cluster is sufficient to obtain the redshift information of BNSs.

The current volumetric merger rate of BNSs is $13$–1900 $\rm Gpc^{-3}yr^{-1}$
Abbott  _et al._ (2021) and that of local SNeIa is $3.0^{+0.6}_{-0.6}\times
10^{4}\,\rm Gpc^{-3}yr^{-1}$ Li  _et al._ (2011). Considering the median of
local SNeIa rates, it implies that there will be roughly 15 to 2300 SNeIa for
a BNS merger in a galaxy. As in Gupta  _et al._ (2019), we assume that the
ratio of the SNeIa and BNS rates will be similar in rich galaxy clusters as
well since both types of populations involve compact object mergers. For
sources up to the redshift of $z=1.5$, we use SNeIa rate to be
$0.65^{+0.61}_{-0.49}\times 10^{-12}\,\rm L_{B,\odot}^{-1}yr^{-1}$ in rich
galaxy clusters Friedmann and Maoz (2018), where $L_{\odot}$ is the bolometric
luminosity in solar units. Consequently, these numbers suggest that every year
there will be $\approx 3$ SNeIa in a Coma-like cluster with bolometric
luminosity of $L_{B}\approx 5.0\times 10^{12}L_{B,\odot}$ Girardi  _et al._
(2002) which is sufficient to confirm the association with BNSs and derive
their redshifts.

###  III.3 Luminosity distance errors

####  III.3.1 Errors from gravitational-wave observations

We simulate a population of neutron star binaries using the procedure outlined
in Sec. III.2. The redshift of a source is converted to its luminosity
distance, the GW observable, using Planck18 Aghanim  _et al._ (2020)
cosmology. For a BNS merger to be detectable, we require a network signal-to-
noise ratio (SNR) threshold of 12 for each binary in the population but do not
demand a minimum SNR for individual detectors. Note that the probability of
having just one detector online in a 5 detector network with a duty cycle of
80% for each detector is less than a percent Belgacem  _et al._ (2019). We
calculate the errors in the estimation of the binary parameters using the
publicly available code, gwbench Borhanian (2020), which implements a Fisher-
matrix formalism Cutler and Flanagan (1994); Poisson and Will (1995) for error
calculation. We use the IMRPhenomPv2_NRTidal waveform model in our Fisher
analysis, with a fixed effective tidal parameter $\tilde{\Lambda}=100$. We do
not compute the error on the $\tilde{\Lambda}$ measurement because this
parameter is not expected to appreciably affect the luminosity distance
estimate. We assume that the electromagnetic counterpart accurately provides
the sky location, so we do not compute an error on it. We further take the
chirp mass to be given because it is well estimated and mostly not degenerate
with the luminosity distance.

Network | GW events | GW events ($z<0.5$) | GW + GRB events | GW + GRB events ($z>0.5$) | GW + EM counterpart
---|---|---|---|---|---
2G+ | 10,259 (83) | 10,259 (83) | 304 (2) | 0 (0) | 1330 (11)
Voy+ | 83,697 (589) | 81,415 (571) | 825 (5) | 17 (0) | 8,967 (62)
ECC | 5,286,423 (36,001) | 505,073 (3454) | 2,810 (19) | 1,657 (15) | 53,317 (364)
Table 2: Number of GW events detected by the three networks in 10 years,
together with the coincident GRB detection rate and the same for sources with
$z>0.5$, assuming the detector characteristics of Fermi-GBM.

This slightly underestimates the distance errors but it would not affect our
results significantly. From a technical perspective, this renders some of the
otherwise ill-conditioned Fisher matrices of the 3G network to behave well. We
are then left with a seven dimensional Fisher matrix consisting of the
symmetric mass ratio $\eta$, the luminosity distance $d_{L}^{\rm GW}$, the
inclination angle $\iota$, the polarization angle $\psi$, the time of
coalescence $t_{c}$, and the phase of coalescence $\phi_{c}$. We subsequently
extract the errors in the measurement of the luminosity distance which is the
parameter of interest here.

Figure 2 shows the redshift distribution of the detected GW events in our
population (left panel), together with those that have an observable GRB
counterpart in the Fermi-GBM detector (right panel) for the three different
networks considered. The distribution is shown for both the optimistic case
(solid lines) and the pessimistic case (dashed lines) corresponding to the
range of the local BNS merger rates given in Eq. (34). We note that only the
3G network can observe BNS coalescences from the furthest redshifts
considered. In Tab. 2, we quote the figures for the expected number of GW
events, the corresponding number whose redshift is less than 0.5, the total
number of events with a GRB counterpart, the number of events with a GRB
counterpart above redshift $z>0.5$, and the cumulative number of events
expected to contribute to the measurement of $\eta_{0}$ according to our
assumptions in Sec. III.2.2. The numbers in parenthesis correspond to the
pessimistic case.

We see that the 2G+ network has a horizon reach of less than $z<0.3$ and the
total number of coincident electromagnetic detections for the optimistic case
are $\sim 1330$ (304 + 10% of 10,259). The corresponding numbers for the Voy+
and ECC networks are $\sim 8,970$ and $\sim 53,320$, respectively, as given in
the last column of Tab. 2.

Figure 3: Fractional error in the measurement of gravitational luminosity
distance as a function of redshift for simulated the BNS population. We see
that only the third generation detector network detects sources from the
highest redshift ($z=1.5$) considered in this study.

In Fig. 3, we show the fractional error in the measurement of GW luminosity
distance as a function of redshift for the three detector networks for our
detected population. To get the average behavior, we distribute the sources
into redshift bins and calculate the median of the fractional errors of the
sources in each redshift bin. We model the fractional luminosity distance
errors as a function of redshift as a series of Heaviside step functions,
which entails taking the errors in each redshift bin to be a constant.

We note that a fit for the fractional error in luminosity distance as a
function of redshift can be found in Belgacem  _et al._ (2019). The reason we
do not directly use their fits in our study is because the only parameter in
their Fisher matrix is the luminosity distance and, therefore, their errors
are unrealistic. Crucially, they ignore the correlations between the
luminosity distance and inclination angle, which is known to increase the
errors significantly Marković (1993); Cutler and Flanagan (1994).

As can be seen from Fig. 3, the horizon distance for the 2G+ network is $z\sim
0.3$ and hence we consider the full detectable population to have a possible
kilonova counterpart detection. Another point of note from the figure is that
the largest redshift considered in this study ($z=1.5$) is within the horizon
distance for the ECC network. We do not consider higher redshift sources
because we are limited by the farthest observed SNeIa in the Union2 data-set
(see Sec. III.3.2).

####  III.3.2 Errors from electromagnetic observations

We model the EM luminosity distance errors using the Union2 Amanullah  _et
al._ (2010) SNeIa compilation. Supernovae distances are measured in units of
distance modulus $\mu,$ which is related to luminosity distance by

$\mu=5\log_{10}\left(\frac{10^{5}d_{L}}{Mpc}\right).$ (37)

It can be easily seen from the above expression that the fractional error in
the luminosity distance is given by

$\frac{\Delta d_{L}}{d_{L}}=\frac{5}{\log_{e}(10)}\Delta\mu.$ (38)

Figure 4 shows the fractional errors in the EM luminosity distance in the
Union2 data-set as a function of redshift. We do not see any functional
behaviour in the luminosity distance errors across redshift bins (with the
errors approximately constant across bins) and, therefore, do not attempt at a
fit, instead treating the redshift behaviour of the errors as a piece-wise
step function.

We also note that for the low redshift events of 2G+ and Voy+ networks, we are
limited by the SNeIa luminosity distance errors and the median errors for the
ECC network is always less than their SNeIa counterpart except around the
redshift limit of $z=1.5$.

The median errors in the distance modulus for the Union2 data-set is 0.19. The
Rubin Observatory Legacy Survey of Space and Time is expected to observe about
half a million supernovae in its survey life cycle of 10 years with a large
fraction of them expected to have distance modulus errors of order 0.12 which
is a $\sim 40\%$ improvement over the Union2 data-set LSS , which would
further improve our estimates.

Figure 4: Fractional errors on SNe luminosity distances as a function of
redshift. The red steps denote the median error for the corresponding redshift
bin.

##  IV Results

We calculate the errors on $\eta_{0}$ for our simulation as follows. From our
sub-population of observed sources of GWs and their EM counterparts, we
randomly select $N$ binaries. Given that we know the redshift to each of our
sources, we get the median fractional error in the GW and EM luminosity
distances for each of the $N$ detections from our modeling of the same as
described in Sec. III.3. Assuming that the central value for both the
luminosity distances are the same, we use Eq. (26) to calculate the errors on
$\eta_{0}$ for each source. The combined error for $N$ independent
observations is given by

$\frac{1}{\sigma^{2}}=\sum_{i=1}^{N}\frac{1}{\sigma_{i}^{2}},$ (39)

where $\sigma_{i}$ is the error for each event. We show the resultant errors
on $\eta_{0}$ in Fig. 5 as a function of the number of observed events, in
increments of 5 up to the expected number of observations for the optimistic
case, for the three networks under consideration. The dashed vertical lines
show the expected number of observations for the pessimistic case rounded to
the nearest multiple of 5 for easy reading of the associated error. Also
depicted in the figure on the right axis is the same error converted to the
temporal variation of the ratio of the gravitational and fine structure
constant at the current epoch [see Eq. (29)] where $H_{0}=73.04\pm 1.04\,\rm
km/s/Mpc$ Riess  _et al._ (2021). We fit the $1/\sqrt{N}$ asymptotic behavior
of the errors and quote the typical error for an observation in each network
in Table 3.

Figure 5: The expected error on the measurement of $\eta_{0}$ ($\sigma_{\eta_{0}}$) as a function of the number of observations for the three detector networks examined in this study. The axis on the right enumerates the same errors in terms of the temporal variation of the ratio of the gravitational and fine structure constant at the present epoch, $\beta$. The maximum number of events for each network denotes the expected number of total observations for the optimistic case. The corresponding number for the pessimistic case (rounded to the nearest multiple of 5) is shown by the vertical dashed lines. Network |  | $\sigma_{\eta_{0}}$ |  |  $\sigma_{\beta}$ $[\times 10^{-9}\rm yr^{-1}]$
---|---|---|---|---
2G+ |  | 7.3 |  | 1.1
Voy+ |  | 4.1 |  | 0.61
ECC |  | 1.6 |  | 0.24
Table 3: Typical value of the error in the measurement of $\eta_{0}$ and the
same in terms of the temporal variation of $G/\alpha(z=0)$ for a single GW
detection for the three networks considered in this study.

We note that the constraints from the violation of the distance-duality
relation directly translates to constraints on the variation of fundamental
constants. We now compare our results to complimentary EM experiments that
look for the variation of the fine structure constant from cosmological data.
We stress that this comparison can only be done for a restricted class of
models that do not modify the gravitational sector as the relevant EM
experiments are oblivious to these modifications. Hees  _et al._ (2014)
briefly reviewed such EM probes and the constraints from various probes are
quoted in Table i@ of their paper. Holanda  _et al._ (2012) and Cao and Liang
(2011) use the same parameterization of $\eta(z)$ as ours and quote average
$1\sigma$ errors of 0.12 and 0.22 on $\eta_{0}$, respectively. The latter is
on par with the capability of the 2G+ network at the end of its observing
cycle in the optimistic case.

Other studies use different parameterizations but one can deduce that although
the constraints using 2G+ network’s forecast to be of the same order or
slightly more than the electromagnetic ones, Voy+ network would be able to
place limits that are a few times better than most of the EM
experiments—except constraints from high redshift quasar absorption spectra
Webb  _et al._ (2001); King  _et al._ (2012); Webb  _et al._ (2011)—for the
optimistic case. The ECC network improves the Voy+ network constraints by an
additional factor of 5. The pessimistic case yields constraints that are an
order of magnitude poorer than the optimistic case for all the three networks
studied here, which is in line with the $1/\sqrt{N}$ behaviour of errors – the
optimistic case has two orders of magnitude more events than the pessimistic
case. Furthermore, observations of quasars from high redshifts suggest a
possible spatial variation of the fine structure constant King  _et al._
(2012); Webb  _et al._ (2011). The large number of gravitational wave
observations, albeit from smaller redshifts, would allow for local constraints
on the spatial variation of the fine structure constant too.

##  V Constraints based on GW170817

In the previous sections, we focused on the detection of a spatially
coincident SNeIa to provide the EM luminosity distance. This is because a
SNeIa is a standard candle and, hence, has constant absolute luminosity in the
source frame. In addition, the systematic uncertainties of modeling SNeIa as a
standard candles is well understood and, therefore, provides an unbiased
estimate for the luminosity distance to the source.

Recently, there have been efforts to model the kilonova emissions following a
BNS merger as a standard candle Kashyap  _et al._ (2019); Coughlin  _et al._
(2020). This would provide another independent measure of the EM luminosity
distance for low redshift sources with the added benefit of not having to
search for a spatially coincident supernova. We, however, do not forecast the
constraints that can be placed on $\eta_{0}$ for a population of joint GW-
kilonova sources using standardised kilonova emissions since these models are
at a very nascent stage of development with large systematic uncertainties.
Nevertheless, we use the EM luminosity distance estimates of Coughlin  _et
al._ (2020) for GW170817 Abbott  _et al._ (2017a) to place the first
constraints on our deviation parameter $\eta_{0}$. The gravitational-wave
luminosity distance for GW170817 was estimated to be $d_{L}^{\rm
GW}=43.8^{+2.9}_{-6.9}\,\text{Mpc}$ Abbott  _et al._ (2017c). Coughlin  _et
al._ (2020) give three measurements of the EM luminosity distance. The first
value, $d_{L}^{\rm EM}=31^{+17}_{-11}\,\text{Mpc}$, is a direct measurement
from the lightcurve based on the analysis of Kasen  _et al._ (2017). The other
two are inferred from ejecta parameters based on the analyses of Kasen  _et
al._ (2017) and Bulla (2019) and are given by $d_{L}^{\rm
EM}=37^{+8}_{-7}\,\text{Mpc}$ and $d_{L}^{\rm EM}=40^{+9}_{-8}\,\text{Mpc}$,
respectively. These three luminosity distance measurements correspond to
$\eta_{0}=42^{+79}_{-55}$, $18^{+27}_{-29}$, and $10^{+26}_{-28}$,
respectively. Unsurprisingly, the estimate of $\eta_{0}$ is consistent with 0.

The above measurements of $\eta_{0}$ give the temporal variation of
$G/\alpha(z)$ [see Eq. (29)] – in units of $[\times 10^{-9}\rm yr^{-1}]$ – at
the current epoch to be

$\beta=-6^{+8}_{-12},\;-3^{+4}_{-4},\;-1^{+4}_{-4},$ (40)

respectively. A Hubble constant value of $H_{0}=73.04\pm 1.04\,\rm km/s/Mpc$
as reported by Riess  _et al._ (2021) is used for this calculation.

##  VI Conclusion

In this paper, we focused on constraining the ratio of gravitational-wave and
electromagnetic-wave luminosity distances and, consequently, the variation in
the ratio of the gravitational and fine structure constant, using coincident
gravitational and electromagnetic wave observations for a class of scalar-
tensor theories known as runaway dilaton models. These theories have
multiplicative couplings of generic scalar fields to gravitational and
electromagnetic sectors. To constrain the modified propagation in such
theories without having to fit for other cosmological parameters, as is done
while using a luminosity distance-redshift relation, a second distance measure
is necessary. We used a spatially coincident supernova as the EM probe to
provide the complimentary EM luminosity distance estimate.

We find that the planned upgrade to the current second-generation ground-based
detector network (2G+) can constrain the parameter modeling the ratio of the
EM and GW luminosity distance to below $|\eta_{0}|\lesssim 0.2,$ while the
proposed improvement of the 2G+ network to Voy+ sensitivity would be able to
place an upper limit of $|\eta_{0}|\lesssim 0.05$, at the end of an 8-year
effective observing cycle if no deviation from the GR value of $\eta_{0}=0$ is
measured. The proposed next-generation ground-based detector network (ECC) can
further improve the constraints to $|\eta_{0}|\lesssim 0.01$. We see that the
constraints using this method for the sub-class of theories that modify the EM
sector alone would be competitive with most of the current EM probes in the
literature Hees  _et al._ (2014) for the 2G+ network. The Voy+ network would
improve these estimates by a factor of 4 and the ECC network improves the Voy+
network constraints by a factor of 5. We also showed how these numbers
translate in to the temporal variation of the fundamental constants.

We, further, make use of recent progress in kilonova light-curve modeling and,
consequently, the EM luminosity distance estimates from them to place the
first constraints on our $\eta_{0}$ parameter for GW170817. As expected, we
find consistency with GR.

We expect a number of BNS merger observations with counterpart in the fourth
observing run of aLIGO/aVirgo/KAGRA and plan to use them to constrain this
class of theories.

###### Acknowledgements.

AD was supported by NSF grant PHY-2012083 and BSS was supported in part by NSF
grants PHY-1836779, PHY-2012083 and AST-2006384.

## References

  * Maggiore (2007) M. Maggiore,  _Gravitational Waves. Vol. 1: Theory and Experiments_ , Oxford Master Series in Physics (Oxford University Press, 2007).
  * Fujii and Maeda (2007) Y. Fujii and K. Maeda,  _The scalar-tensor theory of gravitation_, Cambridge Monographs on Mathematical Physics (Cambridge University Press, 2007).
  * Bertolami  _et al._ (2007) O. Bertolami, C. G. Boehmer, T. Harko, and F. S. N. Lobo, Phys. Rev. D 75, 104016 (2007), arXiv:0704.1733 [gr-qc] .
  * Bertolami  _et al._ (2008) O. Bertolami, F. S. N. Lobo, and J. Paramos, Phys. Rev. D 78, 064036 (2008), arXiv:0806.4434 [gr-qc] .
  * Bertolami and Paramos (2008) O. Bertolami and J. Paramos, Class. Quant. Grav. 25, 245017 (2008), arXiv:0805.1241 [gr-qc] .
  * Sotiriou and Faraoni (2008) T. P. Sotiriou and V. Faraoni, Class. Quant. Grav. 25, 205002 (2008), arXiv:0805.1249 [gr-qc] .
  * De Felice and Tsujikawa (2010) A. De Felice and S. Tsujikawa, Living Rev. Rel. 13, 3 (2010), arXiv:1002.4928 [gr-qc] .
  * Harko  _et al._ (2013) T. Harko, F. S. N. Lobo, and O. Minazzoli, Phys. Rev. D 87, 047501 (2013), arXiv:1210.4218 [gr-qc] .
  * Das and Banerjee (2008) S. Das and N. Banerjee, Phys. Rev. D 78, 043512 (2008), arXiv:0803.3936 [gr-qc] .
  * Bisabr (2012) Y. Bisabr, Phys. Rev. D 86, 127503 (2012), arXiv:1212.2709 [gr-qc] .
  * Moffat and Toth (2012) J. W. Moffat and V. T. Toth, Int. J. Mod. Phys. D 21, 1250084 (2012), arXiv:1001.1564 [gr-qc] .
  * Shiralilou  _et al._ (2022) B. Shiralilou, T. Hinderer, S. M. Nissanke, N. Ortiz, and H. Witek, Class. Quant. Grav. 39, 035002 (2022), arXiv:2105.13972 [gr-qc] .
  * Rovelli and Smolin (1994) C. Rovelli and L. Smolin, Phys. Rev. Lett. 72, 446 (1994), arXiv:gr-qc/9308002 .
  * Domagala  _et al._ (2010) M. Domagala, K. Giesel, W. Kaminski, and J. Lewandowski, Phys. Rev. D 82, 104038 (2010), arXiv:1009.2445 [gr-qc] .
  * Green  _et al._ (1988) M. B. Green, J. H. Schwarz, and E. Witten,  _SUPERSTRING THEORY. VOL. 2: LOOP AMPLITUDES, ANOMALIES AND PHENOMENOLOGY_ (Cambridge University Press, 1988).
  * Uzan (2011) J.-P. Uzan, Living Rev. Rel. 14, 2 (2011), arXiv:1009.5514 [astro-ph.CO] .
  * Damour and Polyakov (1994a) T. Damour and A. M. Polyakov, Gen. Rel. Grav. 26, 1171 (1994a), arXiv:gr-qc/9411069 .
  * Damour and Polyakov (1994b) T. Damour and A. M. Polyakov, Nucl. Phys. B 423, 532 (1994b), arXiv:hep-th/9401069 .
  * Gasperini  _et al._ (2002) M. Gasperini, F. Piazza, and G. Veneziano, Phys. Rev. D 65, 023508 (2002), arXiv:gr-qc/0108016 .
  * Minazzoli and Hees (2013) O. Minazzoli and A. Hees, Phys. Rev. D 88, 041504 (2013), arXiv:1308.2770 [gr-qc] .
  * Ratra and Peebles (1988) B. Ratra and P. J. E. Peebles, Phys. Rev. D 37, 3406 (1988).
  * Caldwell  _et al._ (1998) R. R. Caldwell, R. Dave, and P. J. Steinhardt, Phys. Rev. Lett. 80, 1582 (1998), arXiv:astro-ph/9708069 .
  * Peebles and Ratra (2003) P. J. E. Peebles and B. Ratra, Rev. Mod. Phys. 75, 559 (2003), arXiv:astro-ph/0207347 .
  * Guth (1981) A. H. Guth, Phys. Rev. D 23, 347 (1981).
  * Linde (1982) A. D. Linde, Phys. Lett. B 108, 389 (1982).
  * Albrecht and Steinhardt (1982) A. Albrecht and P. J. Steinhardt, Phys. Rev. Lett. 48, 1220 (1982).
  * Linde (2008) A. D. Linde, Lect. Notes Phys. 738, 1 (2008), arXiv:0705.0164 [hep-th] .
  * Bekenstein (1982) J. D. Bekenstein, Phys. Rev. D 25, 1527 (1982).
  * Sandvik  _et al._ (2002) H. B. Sandvik, J. D. Barrow, and J. Magueijo, Phys. Rev. Lett. 88, 031302 (2002), arXiv:astro-ph/0107512 .
  * Dvali and Zaldarriaga (2002) G. R. Dvali and M. Zaldarriaga, Phys. Rev. Lett. 88, 091303 (2002), arXiv:hep-ph/0108217 .
  * Olive and Pospelov (2008) K. A. Olive and M. Pospelov, Phys. Rev. D 77, 043524 (2008), arXiv:0709.3825 [hep-ph] .
  * Damour (2012) T. Damour, Class. Quant. Grav. 29, 184001 (2012), arXiv:1202.6311 [gr-qc] .
  * Armendariz-Picon (2002) C. Armendariz-Picon, Phys. Rev. D 66, 064008 (2002), arXiv:astro-ph/0205187 .
  * Adelberger  _et al._ (2003) E. G. Adelberger, B. R. Heckel, and A. E. Nelson, Ann. Rev. Nucl. Part. Sci. 53, 77 (2003), arXiv:hep-ph/0307284 .
  * Adelberger  _et al._ (2007) E. G. Adelberger, B. R. Heckel, S. A. Hoedl, C. D. Hoyle, D. J. Kapner, and A. Upadhye, Phys. Rev. Lett. 98, 131104 (2007), arXiv:hep-ph/0611223 .
  * Adelberger  _et al._ (2009) E. G. Adelberger, J. H. Gundlach, B. R. Heckel, S. Hoedl, and S. Schlamminger, Prog. Part. Nucl. Phys. 62, 102 (2009).
  * Kapner  _et al._ (2007) D. J. Kapner, T. S. Cook, E. G. Adelberger, J. H. Gundlach, B. R. Heckel, C. D. Hoyle, and H. E. Swanson, Phys. Rev. Lett. 98, 021101 (2007), arXiv:hep-ph/0611184 .
  * Will (2006) C. M. Will, Living Rev. Rel. 9, 3 (2006), arXiv:gr-qc/0510072 .
  * Rosenband  _et al._ (2008) T. Rosenband, D. B. Hume, P. O. Schmidt, C. W. Chou, A. Brusch, L. Lorini, W. H. Oskay, R. E. Drullinger, T. M. Fortier, J. E. Stalnaker, S. A. Diddams, W. C. Swann, N. R. Newbury, W. M. Itano, D. J. Wineland, and J. C. Bergquist, Science 319, 1808 (2008), https://science.sciencemag.org/content/319/5871/1808.full.pdf .
  * Williams  _et al._ (2012) J. G. Williams, S. G. Turyshev, and D. Boggs, Class. Quant. Grav. 29, 184004 (2012), arXiv:1203.2150 [gr-qc] .
  * Tseytlin and Vafa (1992) A. A. Tseytlin and C. Vafa, Nucl. Phys. B 372, 443 (1992), arXiv:hep-th/9109048 .
  * Damour and Vilenkin (1996) T. Damour and A. Vilenkin, Phys. Rev. D 53, 2981 (1996), arXiv:hep-th/9503149 .
  * Damour  _et al._ (2002) T. Damour, F. Piazza, and G. Veneziano, Phys. Rev. Lett. 89, 081601 (2002), arXiv:gr-qc/0204094 .
  * Jarv  _et al._ (2008) L. Jarv, P. Kuusk, and M. Saal, Phys. Rev. D 78, 083530 (2008), arXiv:0807.2159 [gr-qc] .
  * Damour and Nordtvedt (1993) T. Damour and K. Nordtvedt, Phys. Rev. Lett. 70, 2217 (1993).
  * Khoury (2010) J. Khoury, (2010), arXiv:1011.5909 [astro-ph.CO] .
  * Khoury and Weltman (2004a) J. Khoury and A. Weltman, Phys. Rev. Lett. 93, 171104 (2004a), arXiv:astro-ph/0309300 .
  * Khoury and Weltman (2004b) J. Khoury and A. Weltman, Phys. Rev. D 69, 044026 (2004b), arXiv:astro-ph/0309411 .
  * Hees and Fuzfa (2012) A. Hees and A. Fuzfa, Phys. Rev. D 85, 103005 (2012), arXiv:1111.4784 [gr-qc] .
  * Hinterbichler and Khoury (2010) K. Hinterbichler and J. Khoury, Phys. Rev. Lett. 104, 231301 (2010), arXiv:1001.4525 [hep-th] .
  * Hinterbichler  _et al._ (2011) K. Hinterbichler, J. Khoury, A. Levy, and A. Matas, Phys. Rev. D 84, 103521 (2011), arXiv:1107.2112 [astro-ph.CO] .
  * Webb  _et al._ (2001) J. K. Webb, M. T. Murphy, V. V. Flambaum, V. A. Dzuba, J. D. Barrow, C. W. Churchill, J. X. Prochaska, and A. M. Wolfe, Phys. Rev. Lett. 87, 091301 (2001), arXiv:astro-ph/0012539 .
  * King  _et al._ (2012) J. A. King, J. K. Webb, M. T. Murphy, V. V. Flambaum, R. F. Carswell, M. B. Bainbridge, M. R. Wilczynska, and F. E. Koch, Mon. Not. Roy. Astron. Soc. 422, 3370 (2012), arXiv:1202.4758 [astro-ph.CO] .
  * Webb  _et al._ (2011) J. K. Webb, J. A. King, M. T. Murphy, V. V. Flambaum, R. F. Carswell, and M. B. Bainbridge, Phys. Rev. Lett. 107, 191101 (2011), arXiv:1008.3907 [astro-ph.CO] .
  * Holanda  _et al._ (2016) R. F. L. Holanda, S. J. Landau, J. S. Alcaniz, I. E. Sanchez G., and V. C. Busti, JCAP 05, 047 (2016), arXiv:1510.07240 [astro-ph.CO] .
  * Khatri and Wandelt (2007) R. Khatri and B. D. Wandelt, Phys. Rev. Lett. 98, 111301 (2007), arXiv:astro-ph/0701752 .
  * Riess  _et al._ (2016) A. G. Riess _et al._ , Astrophys. J. 826, 56 (2016), arXiv:1604.01424 [astro-ph.CO] .
  * Cao and Liang (2011) S. Cao and N. Liang, Research in Astronomy and Astrophysics 11, 1199 (2011), arXiv:1104.4942 [astro-ph.CO] .
  * Hees  _et al._ (2014) A. Hees, O. Minazzoli, and J. Larena, Phys. Rev. D 90, 124064 (2014), arXiv:1406.6187 [astro-ph.CO] .
  * Abbott  _et al._ (2017a) B. P. Abbott _et al._ (LIGO Scientific, Virgo), Phys. Rev. Lett. 119, 161101 (2017a), arXiv:1710.05832 [gr-qc] .
  * Abbott  _et al._ (2017b) B. P. Abbott _et al._ (LIGO Scientific, Virgo, Fermi GBM, INTEGRAL, IceCube, AstroSat Cadmium Zinc Telluride Imager Team, IPN, Insight-Hxmt, ANTARES, Swift, AGILE Team, 1M2H Team, Dark Energy Camera GW-EM, DES, DLT40, GRAWITA, Fermi-LAT, ATCA, ASKAP, Las Cumbres Observatory Group, OzGrav, DWF (Deeper Wider Faster Program), AST3, CAASTRO, VINROUGE, MASTER, J-GEM, GROWTH, JAGWAR, CaltechNRAO, TTU-NRAO, NuSTAR, Pan-STARRS, MAXI Team, TZAC Consortium, KU, Nordic Optical Telescope, ePESSTO, GROND, Texas Tech University, SALT Group, TOROS, BOOTES, MWA, CALET, IKI-GW Follow-up, H.E.S.S., LOFAR, LWA, HAWC, Pierre Auger, ALMA, Euro VLBI Team, Pi of Sky, Chandra Team at McGill University, DFN, ATLAS Telescopes, High Time Resolution Universe Survey, RIMAS, RATIR, SKA South Africa/MeerKAT), Astrophys. J. Lett. 848, L12 (2017b), arXiv:1710.05833 [astro-ph.HE] .
  * Fanizza  _et al._ (2020) G. Fanizza, G. Franchini, M. Gasperini, and L. Tedesco, Gen. Rel. Grav. 52, 111 (2020), arXiv:2010.06569 [gr-qc] .
  * Finke  _et al._ (2021) A. Finke, S. Foffa, F. Iacovelli, M. Maggiore, and M. Mancarella, (2021), arXiv:2101.12660 [astro-ph.CO] .
  * Mukherjee  _et al._ (2020) S. Mukherjee, B. D. Wandelt, and J. Silk,  (2020), 10.1093/mnras/stab001, arXiv:2012.15316 [astro-ph.CO] .
  * Minazzoli and Hees (2014) O. Minazzoli and A. Hees, Phys. Rev. D 90, 023017 (2014), arXiv:1404.4266 [gr-qc] .
  * Reitze  _et al._ (2019) D. Reitze _et al._ , Bull. Am. Astron. Soc. 51, 035 (2019), arXiv:1907.04833 [astro-ph.IM] .
  * Abbott  _et al._ (2018) B. P. Abbott _et al._ (KAGRA, LIGO Scientific, Virgo, VIRGO), Living Rev. Rel. 21, 3 (2018), arXiv:1304.0670 [gr-qc] .
  * Acernese  _et al._ (2015a) F. Acernese _et al._ (VIRGO), Class. Quant. Grav. 32, 024001 (2015a), arXiv:1408.3978 [gr-qc] .
  * Aso  _et al._ (2013) Y. Aso, Y. Michimura, K. Somiya, M. Ando, O. Miyakawa, T. Sekiguchi, D. Tatsumi, and H. Yamamoto (KAGRA), Phys. Rev. D 88, 043007 (2013), arXiv:1306.6747 [gr-qc] .
  * Somiya (2012) K. Somiya (KAGRA), Class. Quant. Grav. 29, 124007 (2012), arXiv:1111.7185 [gr-qc] .
  * Unnikrishnan (2013) C. S. Unnikrishnan, Int. J. Mod. Phys. D 22, 1341010 (2013), arXiv:1510.06059 [physics.ins-det] .
  * Saleem  _et al._ (2022) M. Saleem _et al._ , Class. Quant. Grav. 39, 025004 (2022), arXiv:2105.01716 [gr-qc] .
  * Adhikari  _et al._ (2020a) R. X. Adhikari _et al._ (LIGO), Class. Quant. Grav. 37, 165003 (2020a), arXiv:2001.11173 [astro-ph.IM] .
  * Evans  _et al._ (2021) M. Evans _et al._ , (2021), arXiv:2109.09882 [astro-ph.IM] .
  * Punturo  _et al._ (2010a) M. Punturo _et al._ , Class. Quant. Grav. 27, 194002 (2010a).
  * Punturo  _et al._ (2010b) M. Punturo _et al._ , Class. Quant. Grav. 27, 084007 (2010b).
  * Hild  _et al._ (2011) S. Hild _et al._ , Class. Quant. Grav. 28, 094013 (2011), arXiv:1012.0908 [gr-qc] .
  * Holanda  _et al._ (2012) R. F. L. Holanda, R. S. Gonçalves, and J. S. Alcaniz, JCAP 06, 022 (2012), arXiv:1201.2378 [astro-ph.CO] .
  * Yunes  _et al._ (2010) N. Yunes, F. Pretorius, and D. Spergel, Phys. Rev. D 81, 064018 (2010), arXiv:0912.2724 [gr-qc] .
  * Yunes  _et al._ (2016) N. Yunes, K. Yagi, and F. Pretorius, Phys. Rev. D 94, 084002 (2016), arXiv:1603.08955 [gr-qc] .
  * Vijaykumar  _et al._ (2021) A. Vijaykumar, S. J. Kapadia, and P. Ajith, Phys. Rev. Lett. 126, 141104 (2021), arXiv:2003.12832 [gr-qc] .
  * Aasi  _et al._ (2015) J. Aasi _et al._ (LIGO Scientific), Class. Quant. Grav. 32, 074001 (2015), arXiv:1411.4547 [gr-qc] .
  * Acernese  _et al._ (2015b) F. Acernese _et al._ (VIRGO), Class. Quant. Grav. 32, 024001 (2015b), arXiv:1408.3978 [gr-qc] .
  * Akutsu  _et al._ (2019) T. Akutsu _et al._ (KAGRA), Nature Astron. 3, 35 (2019), arXiv:1811.08079 [gr-qc] .
  * Adhikari  _et al._ (2020b) R. X. Adhikari _et al._ (LIGO), Class. Quant. Grav. 37, 165003 (2020b), arXiv:2001.11173 [astro-ph.IM] .
  * Abbott  _et al._ (2020) R. Abbott _et al._ (LIGO Scientific, Virgo), (2020), arXiv:2010.14533 [astro-ph.HE] .
  * Belgacem  _et al._ (2019) E. Belgacem, Y. Dirian, S. Foffa, E. J. Howell, M. Maggiore, and T. Regimbau, JCAP 08, 015 (2019), arXiv:1907.01487 [astro-ph.CO] .
  * Vangioni  _et al._ (2015) E. Vangioni, K. A. Olive, T. Prestegard, J. Silk, P. Petitjean, and V. Mandic, Mon. Not. Roy. Astron. Soc. 447, 2575 (2015), arXiv:1409.2462 [astro-ph.GA] .
  * Abbott  _et al._ (2021) R. Abbott _et al._ (LIGO Scientific, VIRGO, KAGRA), (2021), arXiv:2111.03634 [astro-ph.HE] .
  * Kalogera  _et al._ (2021) V. Kalogera _et al._ , (2021), arXiv:2111.06990 [gr-qc] .
  * Howell  _et al._ (2018) E. J. Howell, K. Ackley, A. Rowlinson, and D. Coward,  (2018), 10.1093/mnras/stz455, arXiv:1811.09168 [astro-ph.HE] .
  * Wanderman and Piran (2015) D. Wanderman and T. Piran, Mon. Not. Roy. Astron. Soc. 448, 3026 (2015), arXiv:1405.5878 [astro-ph.HE] .
  * Burns  _et al._ (2016) E. Burns, V. Connaughton, B.-B. Zhang, A. Lien, M. S. Briggs, A. Goldstein, V. Pelassa, and E. Troja, Astrophys. J. 818, 110 (2016), arXiv:1512.00923 [astro-ph.HE] .
  * Gupta  _et al._ (2019) A. Gupta, D. Fox, B. S. Sathyaprakash, and B. F. Schutz, The Astrophysical Journal 886, 71 (2019).
  * Li  _et al._ (2011) W. Li, R. Chornock, J. Leaman, A. V. Filippenko, D. Poznanski, X. Wang, M. Ganeshalingam, and F. Mannucci, MNRAS 412, 1473 (2011), arXiv:1006.4613 [astro-ph.SR] .
  * Abbott  _et al._ (2019) B. P. Abbott _et al._ (LIGO Scientific, Virgo), Phys. Rev. X 9, 031040 (2019), arXiv:1811.12907 [astro-ph.HE] .
  * Friedmann and Maoz (2018) M. Friedmann and D. Maoz, MNRAS 479, 3563 (2018), arXiv:1803.04421 [astro-ph.GA] .
  * Girardi  _et al._ (2002) M. Girardi, P. Manzato, M. Mezzetti, G. Giuricin, and F. Limboz, Astrophys. J. 569, 720 (2002), astro-ph/0112534 .
  * Aghanim  _et al._ (2020) N. Aghanim _et al._ (Planck), Astron. Astrophys. 641, A6 (2020), arXiv:1807.06209 [astro-ph.CO] .
  * Borhanian (2020) S. Borhanian, (2020), arXiv:2010.15202 [gr-qc] .
  * Cutler and Flanagan (1994) C. Cutler and E. E. Flanagan, Phys. Rev. D 49, 2658 (1994), arXiv:gr-qc/9402014 .
  * Poisson and Will (1995) E. Poisson and C. M. Will, Phys. Rev. D 52, 848 (1995), arXiv:gr-qc/9502040 .
  * Marković (1993) D. Marković, Phys. Rev. D 48, 4738 (1993).
  * Amanullah  _et al._ (2010) R. Amanullah, C. Lidman, D. Rubin, G. Aldering, P. Astier, K. Barbary, M. S. Burns, A. Conley, K. S. Dawson, S. E. Deustua, M. Doi, S. Fabbro, L. Faccioli, H. K. Fakhouri, G. Folatelli, A. S. Fruchter, H. Furusawa, G. Garavini, G. Goldhaber, A. Goobar, D. E. Groom, I. Hook, D. A. Howell, N. Kashikawa, A. G. Kim, R. A. Knop, M. Kowalski, E. Linder, J. Meyers, T. Morokuma, S. Nobili, J. Nordin, P. E. Nugent, L. Östman, R. Pain, N. Panagia, S. Perlmutter, J. Raux, P. Ruiz-Lapuente, A. L. Spadafora, M. Strovink, N. Suzuki, L. Wang, W. M. Wood-Vasey, N. Yasuda, and T. Supernova Cosmology Project, Astrophys. J. 716, 712 (2010), arXiv:1004.1711 [astro-ph.CO] .
  * (105) 11.5 Constraining the Dark Energy Equation of State, https://www.lsst.org/sites/default/files/docs/sciencebook/SB_11.pdf.
  * Riess  _et al._ (2021) A. G. Riess _et al._ , (2021), arXiv:2112.04510 [astro-ph.CO] .
  * Kashyap  _et al._ (2019) R. Kashyap, G. Raman, and P. Ajith, Astrophys. J. Lett. 886, L19 (2019), arXiv:1908.02168 [astro-ph.SR] .
  * Coughlin  _et al._ (2020) M. W. Coughlin, T. Dietrich, J. Heinzel, N. Khetan, S. Antier, M. Bulla, N. Christensen, D. A. Coulter, and R. J. Foley, Phys. Rev. Res. 2, 022006 (2020), arXiv:1908.00889 [astro-ph.HE] .
  * Abbott  _et al._ (2017c) B. P. Abbott _et al._ (LIGO Scientific, Virgo, 1M2H, Dark Energy Camera GW-E, DES, DLT40, Las Cumbres Observatory, VINROUGE, MASTER), Nature 551, 85 (2017c), arXiv:1710.05835 [astro-ph.CO] .
  * Kasen  _et al._ (2017) D. Kasen, B. Metzger, J. Barnes, E. Quataert, and E. Ramirez-Ruiz, Nature 551, 80 (2017), arXiv:1710.05463 [astro-ph.HE] .
  * Bulla (2019) M. Bulla, Mon. Not. Roy. Astron. Soc. 489, 5037 (2019), arXiv:1906.04205 [astro-ph.HE] .