# Forecast for cosmological parameter estimation with gravitational-wave
standard sirens from the LISA-Taiji network

Ling-Feng Wang  Department of Physics, College of Sciences, Northeastern
University, Shenyang 110819, China  Shang-Jie Jin  Department of Physics,
College of Sciences, Northeastern University, Shenyang 110819, China  Jing-Fei
Zhang  Department of Physics, College of Sciences, Northeastern University,
Shenyang 110819, China  Xin Zhang111Corresponding author
<EMAIL_ADDRESS>Department of Physics, College of Sciences,
Northeastern University, Shenyang 110819, China  Key Laboratory of Data
Analytics and Optimization for Smart Industry (Northeastern University),
Ministry of Education, Shenyang 110819, China

###### Abstract

LISA and Taiji are expected to form a space-based gravitational-wave (GW)
detection network in the future. In this work, we make a forecast for the
cosmological parameter estimation with the standard siren observation from the
LISA-Taiji network. We simulate the standard siren data based on a scenario
with configuration angle of $40^{\circ}$ between LISA and Taiji. Three models
for the population of massive black hole binary (MBHB), i.e., pop III, Q3d,
and Q3nod, are considered to predict the events of MBHB mergers. We find that,
based on the LISA-Taiji network, the number of electromagnetic (EM)
counterparts detected is almost doubled compared with the case of single Taiji
mission. Therefore, the LISA-Taiji network’s standard siren observation could
provide much tighter constraints on cosmological parameters. For example,
solely using the standard sirens from the LISA-Taiji network, the constraint
precision of $H_{0}$ could reach $1.3\%$. Moreover, combined with the CMB
data, the GW-EM observation based on the LISA-Taiji network could also tightly
constrain the equation of state of dark energy, e.g., the constraint precision
of $w$ reaches about $4\%$, which is comparable with the result of CMB+BAO+SN.
It is concluded that the GW standard sirens from the LISA-Taiji network will
become a useful cosmological probe in understanding the nature of dark energy
in the future.

gravitational wave, space-based, network

LISA-Taiji network, gravitational-wave standard siren, cosmological parameter
estimation, massive black bole binary, electromagnetic counterpart

###### pacs:

95.36.+x, 98.80.Es, 98.80.-k, 04.80.Nn, 95.55.Ym

##  I Introduction

The precise measurement of the cosmic microwave background (CMB) anisotropies
initiated the era of precision cosmology Spergel  _et al._ (2003); Bennett
_et al._ (2003). Constraining the standard cosmological model, i.e., the
$\Lambda$CDM model, with the high-precision CMB observations enables
cosmologists to have a comprehensive understanding of the evolution history of
the universe. However, the accurate measurements also led to some puzzling
issues. For example, there is a 4.4$\sigma$ tension between the $H_{0}$ values
inferred from the CMB observation Aghanim  _et al._ (2020) and the distance
ladder measurement Riess  _et al._ (2019). Essentially, the $H_{0}$ tension
reflects an inconsistency of measurements between the early universe and the
late universe Verde  _et al._ (2019); Riess (2019). In addition to facing the
challenge of the Hubble tension, the $\Lambda$CDM model also has some
theoretical problems, such as the “fine-tuning” and “cosmic coincidence”
problems Weinberg (1989); Sahni and Starobinsky (2000); Bean  _et al._ (2005),
which implies that the $\Lambda$CDM model needs to be further adjusted.
Therefore, the current development of cosmology can be divided into two main
aspects: (i) further extending the standard $\Lambda$CDM model Guo  _et al._
(2019a); Zhao  _et al._ (2020a); Di Valentino  _et al._ (2021); Li  _et al._
(2019); Vagnozzi (2020); Guo  _et al._ (2019b); Li  _et al._ (2020a); Feng
_et al._ (2020a); Zhang  _et al._ (2020a); Feng  _et al._ (2020b); Lin  _et
al._ (2020); Guo  _et al._ (2020); Hryczuk and Jodłowski (2020); Zhang  _et
al._ (2020b); Gao  _et al._ (2021), and (ii) developing more low-redshift
observation projects aimed at precisely measuring the late universe Zhang
(2019); Xu and Zhang (2020); Li and Zhang (2020); Zhang  _et al._ (2021); Wang
_et al._ (2021). For the second aspect, the gravitational-wave (GW) standard
siren method Schutz (1986); Holz and Hughes (2005) is one of the most
promising options and has been widely discussed Holz and Hughes (2005); Cai
_et al._ (2018a); Di Valentino and Melchiorri (2018); Wei (2018); Zhao  _et
al._ (2018); Di Valentino  _et al._ (2018); Du  _et al._ (2019); Mifsud and
van de Bruck (2019); Wei (2019); Gray  _et al._ (2020); Howlett and Davis
(2020); Chassande-Mottin  _et al._ (2019); Doctor (2020); Fu  _et al._ (2019);
Mukherjee  _et al._ (2021); Abbott  _et al._ (2021); Palmese  _et al._ (2020);
Chen (2020); Wang  _et al._ (2020a); Chen  _et al._ (2021); Zhao  _et al._
(2011); Cai and Yang (2017); Cai  _et al._ (2018b); Yang  _et al._ (2019a);
Cai and Yang (2018); Wang  _et al._ (2018); Zhang  _et al._ (2019a); Li  _et
al._ (2020b); Zhang  _et al._ (2019b, 2020c); Yan  _et al._ (2019); Yang  _et
al._ (2019b); Jin  _et al._ (2020); Qi  _et al._ (2021).

The amplitude of GW generated by the merger of compact binary encodes
luminosity distance and chirp mass of the source. The absolute luminosity
distance could be obtained if the amplitude is measured precisely, and the
chirp mass could be inferred from the variation of GW frequency. The relation
between luminosity distance $d_{\rm L}$ and redshift $z$ can be established,
once the electromagnetic (EM) counterpart of a GW source is detected by
optical observatories. Since the $d_{\rm L}$–$z$ relation is determined by the
expansion history of the universe, cosmological models could be constrained by
this relation. This method is usually referred to as the “standard siren”
method Schutz (1986); Holz and Hughes (2005). GW170817 Abbott  _et al._
(2017a, b, c), the first detected binary neutron star (BNS) merger event with
an EM counterpart (GRB 170817A), has provided an independent measurement of
the Hubble constant, giving the result of $H_{0}=70^{+12}_{-8}~{}{\rm
km}~{}{\rm s}^{-1}~{}{\rm Mpc}^{-1}$ Abbott  _et al._ (2017d). Recently, the
event ZTF19abanrhr Graham  _et al._ (2020) reported by the Zwicky Transient
Facility is regarded as a candidate of the first plausible optical EM
counterpart to the binary black hole (BBH) merger event GW190521 Abbott  _et
al._ (2020). Chen _et al._ Chen  _et al._ (2020) and Mukherjee _et al._
Mukherjee  _et al._ (2020) have given the constraints on the $\Lambda$CDM and
$w$CDM models, assuming ZTF19abanrhr is the actual EM counterpart to GW190521.
Even if BBH merger events are expected to have no EM counterparts, these “dark
sirens” could also be used in cosmological fits using the statistical method
discussed in Refs. Del Pozzo (2012); Chen  _et al._ (2018); Fishbach  _et al._
(2019); Feeney  _et al._ (2019); Ding  _et al._ (2019); Soares-Santos  _et
al._ (2019); Lagos  _et al._ (2019); Yu  _et al._ (2020).

The potential of standard siren method in constraining cosmological parameters
has been forecasted in Refs. Zhao  _et al._ (2011); Cai and Yang (2017); Cai
_et al._ (2018b); Yang  _et al._ (2019a); Cai and Yang (2018); Wang  _et al._
(2018); Zhang  _et al._ (2019a); Li  _et al._ (2020b); Zhang  _et al._ (2019b,
2020c); Yan  _et al._ (2019); Yang  _et al._ (2019b); Jin  _et al._ (2020); Qi
_et al._ (2021), based on future ground-based GW detectors, e.g., Einstein
Telescope Punturo  _et al._ (2010); ET-  and Cosmic Explorer Abbott  _et al._
(2017e); CE- . Several mechanisms for producing fast radio bursts (FRBs) by
the mergers of binaries, such as charged black holes or neutron stars, are
proposed in Refs. Totani (2013); Mingarelli  _et al._ (2015); Wang  _et al._
(2016); Liu  _et al._ (2016); Zhang (2016); Yamasaki  _et al._ (2018), and
subsequently, GW/FRB association systems as a complementary cosmological probe
are discussed in Refs. Wei  _et al._ (2018); Cai  _et al._ (2019) (see also
Ref. Zhao  _et al._ (2020b)). The GW sources detected by future ground-based
GW detectors are BNSs or stellar-mass BBHs, which are mainly distributed at
$z<3$. In addition, the GWs produced by the massive black hole binaries
(MBHBs) with EM counterparts are also expected to serve as standard sirens
Holz and Hughes (2005). The MBHB standard sirens may be detected at the
redshift up to $z\simeq 10$ in the future Caldwell  _et al._ (2019), providing
a promising method of measuring the expansion history of the universe back to
a much earlier time. Space-based GW detection missions have been proposed and
implemented to detect the GWs produced by MBHBs.

The Laser Interferometer Space Antenna (LISA) LIS ; Armano  _et al._ (2016);
Amaro-Seoane  _et al._ (2017); Armano  _et al._ (2018); Abich  _et al._
(2019); Speri  _et al._ (2021) is a European space-based GW observatory, with
three identical drag-free spacecraft that form an equilateral triangle with
arm length of $2.5\times 10^{6}$ km. Taiji Wu (2018); Ruan  _et al._ (2020a);
Hu and Wu (2017) and TianQin Luo  _et al._ (2020); Wang  _et al._ (2020b); Liu
_et al._ (2020); Milyukov (2020); Mei  _et al._ (2020); Fan  _et al._ (2020)
are two space-based GW observatories proposed by Chinese researchers. Taiji is
a LISA-like space-based GW observatory proposed by the Chinese Academy of
Sciences, also with a triangle of three satellites but with arm length of
$3\times 10^{6}$ km. Some forecasts for the capability of LISA and Taiji in
cosmological parameter estimation have been discussed in Refs. Tamanini  _et
al._ (2016); Belgacem  _et al._ (2019); Zhao  _et al._ (2020c). Multiple GW
observatories can constitute a network to improve the measurement precision of
source parameters Abbott  _et al._ (2017f); Fan  _et al._ (2019), by measuring
phase differences and amplitude ratios of GWs in different detectors. LISA’s
orbit is proposed to be at the ecliptic plane behind the Earth with a
$20^{\circ}$ trailing angle, and Taiji is planed to be localized in front of
the Earth with a $20^{\circ}$ leading angle, so that LISA and Taiji could form
a space-based network aimed at detecting GW signals within an mHz range (i.e.
$10^{-4}$ Hz to $10^{-1}$ Hz) Ruan  _et al._ (2020b); Wang  _et al._ (2020c);
Hu  _et al._ (2021); Wang  _et al._ (2020d); Omiya and Seto (2020); Orlando
_et al._ (2021); Wang and Han (2021).

Recently, the LISA-Taiji network’s capability of localizing GW sources was
detailedly discussed in Refs. Ruan  _et al._ (2020b); Wang  _et al._ (2020c).
Wang _et al._ Wang  _et al._ (2020d) showed that within 5-year operation time,
the LISA-Taiji network is able to constrain the Hubble constant within $1\%$
accuracy via dark sirens. Omiya and Seto Omiya and Seto (2020) explored the
detectability of vector and scalar polarization modes in a stochastic
gravitational wave background (SGWB) around 1 mHz with the LISA-Taiji network.
Orlando _et al._ Orlando  _et al._ (2021) proposed that the chirality of an
isotropic SGWB can be detected by cross-correlating the data streams of LISA
and Taiji. Wang and Han Wang and Han (2021) showed that the LISA-Taiji network
could improve the observations on the anomalous polarization predicted by the
theories beyond general relativity.

In this work, we focus on the LISA-Taiji network’s capability of improving the
constraint accuracies of cosmological parameters. We first use the Fisher
information matrix method to estimate the uncertainty of luminosity distance
$d_{\rm L}$, and then simulate 5-year standard siren data based on the LISA-
Taiji network. Then, we constrain three typical dark-energy cosmological
models, i.e., the $\Lambda$CDM, $w$CDM, and CPL models, using the standard
siren (joint GW-EM detection) mock data. We mainly analyze the improvement on
the constraint accuracy of equation of state (EoS) of dark energy, and show
that the LISA-Taiji network will play an important role in exploring the
properties of dark energy in the future.

The rest of this paper is organized as follows. In Sec. II.1, we describe the
GW waveform and the detector response. In Sec. II.2, we describe the Fisher
matrix analysis for GW parameter estimation. In Sec. II.3, we discuss the
identifications of EM counterparts. In Sec. II.4, we introduce the methods of
simulating the standard siren catalog. In Sec. III, we display the constraint
results and make some discussions. The conclusion is given in Sec. IV. Unless
otherwise specified, we adopt the system of units in which $c=G=1$ throughout
this paper.

##  II Simulation of standard siren observation

###  II.1 GW waveform and detector response

The GW signal from the inspiral of a non-spinning MBHB can be modeled by the
restricted post-Newtonian (PN) waveform. The GW strain $h(t)$ can be described
by two independent polarizations ${h_{+,\times}}(t)$ in the transverse-
traceless gauge,

$\displaystyle
h(t)={F_{+}}(t;\theta,\phi,\psi){h_{+}}(t)+{F_{\times}}(t;\theta,\phi,\psi){h_{\times}}(t)\,,$
(1)

where $F_{+,\times}$ are antenna pattern functions, $(\theta,\phi)$ denote the
source’s polar angle and azimuthal angle in the ecliptic frame, and $\psi$ is
the polarization angle of GW. We can separate the antenna pattern function
into a polarization angle part and a $D_{+,\times}$ part that describes the
dependence of time,

$\displaystyle F_{+}(t)=$ $\displaystyle\frac{1}{2}\Big{(}{\rm
cos}(2\psi)D_{+}(t)-{\rm sin}(2\psi)D_{\times}(t)\Big{)},$ (2) $\displaystyle
F_{\times}(t)=$ $\displaystyle\frac{1}{2}\Big{(}{\rm sin}(2\psi)D_{+}(t)+{\rm
cos}(2\psi)D_{\times}(t)\Big{)}.$ (3)

For the inspiral process, the specific forms of $D_{+,\times}$ with the low-
frequency approximation are given in Ref. Ruan  _et al._ (2020b),

$\displaystyle D_{+}(t)=$ $\displaystyle\frac{\sqrt{3}}{64}\bigg{[}-36{\rm
sin}^{2}\theta\,{\rm sin}\big{(}2\alpha(t)-2\beta\big{)}+\big{(}3+{\rm
cos(2\theta)}\big{)}$ $\displaystyle\times\bigg{(}{\rm
cos}(2\phi)\Big{(}9\sin(2\beta)-{\rm
sin}\big{(}4\alpha(t)-2\beta\big{)}\Big{)}$ $\displaystyle+{\rm
sin}(2\phi)\Big{(}{\rm
cos}\big{(}4\alpha(t)-2\beta\big{)}-9\cos(2\beta)\Big{)}\bigg{)}$
$\displaystyle-4\sqrt{3}{\rm sin}(2\theta)\Big{(}{\rm
sin}\big{(}3\alpha(t)-2\beta-\phi\big{)}-3{\rm sin}\big{(}\alpha(t)$
$\displaystyle-2\beta+\phi\big{)}\Big{)}\bigg{]}\,,$ (4) $\displaystyle
D_{\times}(t)=$ $\displaystyle\frac{1}{16}\bigg{[}\sqrt{3}{\rm
cos}\theta\Big{(}9{\rm cos}(2\phi-2\beta)-{\rm cos}\big{(}4\alpha(t)-2\beta$
$\displaystyle-2\phi\big{)}\Big{)}-6{\rm sin}\theta\Big{(}{\rm
cos}\big{(}3\alpha(t)-2\beta-\phi\big{)}$ $\displaystyle+3{\rm
cos}\big{(}\alpha(t)-2\beta+\phi\big{)}\Big{)}\bigg{]}\,,$ (5)

where $\alpha=2\pi f_{m}t+\kappa$ is the orbital phase of the guiding center,
and $\beta$ is the initial orientation of the constellation. We simply set
$\beta=0$ in our simulation. The triangular GW detectors with three arms, such
as LISA and Taiji, can be equivalent to two independent
$90^{\circ}$-interferometers (i.e., “L-shaped” interferometers). The second
interferometer is equivalent to the first one rotated by $\pi/4$ radians. The
response functions of the two interferometers are
${F_{+,\times}}(t;\theta,\phi,\psi)$ and
${F_{+,\times}}(t;\theta,\phi-\pi/4,\psi)$ Cutler (1998). Here $\kappa$ is the
initial ecliptic longitude of the guiding center and $f_{m}=1/{\rm yr}$. We
assume that $\kappa=0$ for LISA and $\kappa=40^{\circ}$ for Taiji, so that the
separation angle between LISA and Taiji is $40^{\circ}$.

For the sake of describing GW signals in the Fourier space, the observation
time $t$ in Eqs. (4) and (5) is replaced by Krolak  _et al._ (1995); Buonanno
_et al._ (2009)

$t(f)=t_{\rm c}-\frac{5}{256}M_{\rm c}^{-5/3}(\pi f)^{-8/3},$ (6)

where $t_{\rm c}$ is the coalescence time of MBHB. The Fourier transformation
of the strain can be obtained, i.e.,

$\displaystyle\tilde{h}(f)=-\left(\frac{5\pi}{24}\right)^{1/2}M_{\rm
c}^{5/6}\left[\frac{(\pi f)^{-7/6}}{D_{{\rm eff}}}\right]e^{-i\Psi}.$ (7)

The effective luminosity distance, $D_{\rm eff}$, is defined as

$D_{{\rm eff}}=d_{\rm L}\left[F^{2}_{+}\left(\frac{1+{\rm
cos}^{2}\iota}{2}\right)^{2}+F^{2}_{\times}{\rm cos}^{2}\iota\right]^{-1/2},$
(8)

where $d_{\rm L}$ is the luminosity distance to a GW source and $\iota$ is the
inclination angle between the orbital angular momentum and the line of sight.
$\Psi$ can be written to the second PN order as

$\displaystyle\Psi(f;M_{c},\eta)=$ $\displaystyle 2\pi
ft_{0}-2\phi_{0}-\frac{\pi}{4}+\frac{3}{128\eta}\bigg{[}\nu^{-5}+\Big{(}\frac{3715}{756}$
$\displaystyle+\frac{55}{9}\eta\Big{)}\nu^{-3}-16\pi\nu^{-2}+\Big{(}\frac{15293365}{508032}$
$\displaystyle+\frac{27145}{504}\eta+\frac{3085}{72}{\eta}^{2}\Big{)}\nu^{-1}\bigg{]},$
(9) $\displaystyle\nu=$ $\displaystyle\Big{(}\frac{G\pi
M}{c^{3}}f\Big{)}^{1/3},$ (10)

where $t_{0}=t_{c}+\tau(t)$ is the coalescence time at a detector. According
to the forward modeling of LISA described in Ref. Rubbo  _et al._ (2004), to
linear order in eccentricity, the time delay $\tau(t)$ and the phase
$\phi_{0}$ take the forms

$\displaystyle\tau(t)=$ $\displaystyle-\frac{R}{c}{\rm sin}\theta\,{\rm
cos}(\alpha-\phi)-\frac{1}{2}e\frac{R}{c}{\rm
sin}\theta\big{[}\cos(2\alpha-\phi-\beta)$
$\displaystyle-3\cos(\phi-\beta)\big{]},$ (11) $\displaystyle
2\phi_{0}=2\phi_{c}-{\rm
arctan}\left(\frac{F_{\times}(\theta,\phi,\iota,\psi;t)}{F_{+}(\theta,\phi,\iota,\psi;t)}\frac{2{\rm
cos}\iota}{1+{\rm cos}^{2}\iota}\right),$ (12)

where $R=1\ {\rm AU}$, and $e$ is the eccentricity of detector’s orbit. In the
observer’s reference frame, $M_{\rm c}=(1+z)\eta^{3/5}M$ is the redshifted
chirp mass. $M={M_{1}}+{M_{2}}$ is the total mass of MBHB with
${M_{1}}>{M_{2}}$, and ${\eta=M_{1}M_{2}/M^{2}}$ is the symmetric mass ratio.

###  II.2 Fisher matrix analysis for GW parameter estimation

For a network including $N$ independent detectors, the Fisher information
matrix can be written as

$\bm{F}_{ij}=\left(\frac{\partial\bm{h}(f)}{\partial\theta_{i}}\bigg{|}\frac{\partial\bm{h}(f)}{\partial\theta_{j}}\right),$
(13)

with $\bm{h}$ being given by

$\bm{h}(f)=\left[\frac{\tilde{h}_{1}(f)}{\sqrt{S_{\rm
n}(f)}},\frac{\tilde{h}_{2}(f)}{\sqrt{S_{\rm
n}(f)}},\cdots,\frac{\tilde{h}_{N}(f)}{\sqrt{S_{\rm n}(f)}}\right]^{\rm T},$
(14)

where $\theta_{i}$ denotes nine parameters ($d_{L}$, $M_{c}$, $\eta$,
$\theta$, $\phi$, $\iota$, $t_{c}$, $\phi_{c}$, $\psi$) for a GW event. Here,
$S_{\rm n}(f)$ is the noise power spectral density. The specific forms of
$S_{\rm n}(f)$ for Taiji and LISA are obtained from Refs. Ruan  _et al._
(2020a); Klein  _et al._ (2016). The bracket in Eq. (13) for two functions
$a(t)$ and $b(t)$ is defined as

$(a|b)=4\int_{f_{\rm low}}^{f_{\rm
up}}\frac{\tilde{a}(f)\tilde{b}^{*}(f)+\tilde{a}^{*}(f)\tilde{b}(f)}{2}{\rm
d}f,$ (15)

where “$\sim$” above a function denotes the Fourier transform of the function.
The upper limit of the integral is set to the innermost stable circular orbit
(ISCO) frequency $f_{\rm ISCO}=c^{3}/(6\sqrt{6}\pi GM)$ Feng  _et al._ (2019),
and the lower frequency cutoff is set to $f_{\rm low}=10^{-4}\,\rm Hz$. The
signal-to-noise ratio (SNR) of a GW event is given by

$\rho^{2}=(\bm{h}|\bm{h}),$ (16)

and we consider the SNR threshold of 8 in our simulation.

The Fisher matrix of the LISA-Taiji network is the sum of the Fisher matrix of
LISA and the one of Taiji, which can be expressed as

$F_{{\rm network}}=F_{{\rm LISA}}+F_{{\rm Taiji}}.$ (17)

Then, the errors of GW parameters can be estimated by the Fisher information
matrix,

$\Delta\theta_{i}=\sqrt{(F^{-1})_{ii}}.$ (18)

In our analysis, we take into account nine parameters ($d_{L}$, $M_{c}$,
$\eta$, $\theta$, $\phi$, $\iota$, $t_{c}$, $\phi_{c}$, $\psi$) in the Fisher
matrix. The error of luminosity distance, $\Delta d_{\rm L}$, and the angular
resolution, $\Delta\Omega$, could be calculated by the Fisher matrix. Here,
$\Delta\Omega$ is given by
$\Delta\Omega=2\pi|\textrm{sin}\theta|\sqrt{\langle\Delta\theta^{2}\rangle\langle\Delta\phi^{2}\rangle-\langle\Delta\theta\Delta\phi\rangle^{2}}$,
with $\langle\Delta\theta^{2}\rangle$, $\langle\Delta\phi^{2}\rangle$, and
$\langle\Delta\theta\Delta\phi\rangle$ being given by the inverse of the
9-parameter Fisher matrix Zhao and Wen (2018). Improving the angular
resolutions of GW sources is helpful to identify EM counterparts, so
$\Delta\Omega$ is important for determining the number of GW-EM events. In
addition, $\Delta d_{\rm L}$ could directly affect the constraint accuracies
of cosmological parameters. Therefore, before making a cosmological analysis,
it is necessary to study the reductions of $\Delta\Omega$ and $\Delta d_{\rm
L}$ made by the LISA-Taiji network.

For the sake of clearly showing the reductions of $\Delta\Omega$ and $\Delta
d_{\rm L}$, we plot $\Delta\Omega$ and $\Delta d_{\rm L}$ as functions of
redshift in Figure 1. In this figure, we simulate 500 GW events to show
statistical distributions of $\Delta\Omega$ and $\Delta d_{\rm L}$. We choose
$\kappa=0$ for LISA and $\kappa=40^{\circ}$ for Taiji. The mass of MBHs, the
sky location ($\theta$, $\phi$), the binary inclination $\iota$, the
polarization angle $\psi$, and the coalescence phase $\phi_{c}$ are randomly
chosen in the ranges of $[10^{4},10^{7}]M_{\odot}$, $[0,\pi]$, $[0,2\pi]$,
$[0,\pi]$, $[0,2\pi]$, and $[0,2\pi]$, respectively. We can clearly see that
the LISA-Taiji network could reduce $\Delta\Omega$ by several orders of
magnitude compared with the single Taiji mission, which implies that the LISA-
Taiji network could greatly improve the capability of locating GW sources, and
thus could increase the detection number of GW-EM events. For the uncertainty
of luminosity distance, $\Delta d_{\rm L}$, it is not reduced as much as
$\Delta\Omega$, but is still reduced by a factor of a few. We have a more
specific analysis on these aspects in Sec. II.4.

Figure 1: The uncertainties of angular resolutions and luminosity distances as
functions of redshift. Here we choose $\kappa=0$ for LISA, and
$\kappa=40^{\circ}$ for Taiji.

Actually, $\Delta d_{\rm L}$ estimated by the Fisher matrix is the
instrumental error, $\sigma_{d_{\rm L}}^{\rm inst}$. The total measurement
error of luminosity distance $\sigma_{d_{\rm L}}$ also consists of the lensing
error, the peculiar velocity error, and the redshift measurement error, which
can be expressed as Wang  _et al._ (2020b)

$\displaystyle(\sigma_{d_{\rm L}})^{2}=(\sigma_{d_{\rm L}}^{\rm
inst})^{2}+(\sigma_{d_{\rm L}}^{\rm lens})^{2}+(\sigma_{d_{\rm L}}^{\rm
pv})^{2}+(\sigma_{d_{\rm L}}^{\rm reds})^{2}.$ (19)

The main systematic error caused by weak lensing is adopted from the fitting
formula Tamanini  _et al._ (2016),

$\displaystyle\sigma_{d_{\rm L}}^{\rm lens}(z)=d_{\rm L}(z)\times
0.066\bigg{[}\frac{1-(1+z)^{-0.25}}{0.25}\bigg{]}^{1.8},$ (20)

and the error caused by the peculiar velocity of a source should also be
included Kocsis  _et al._ (2006),

$\displaystyle\sigma_{d_{\rm L}}^{\rm pv}(z)=d_{\rm
L}(z)\times\bigg{[}1+\frac{c(1+z)^{2}}{H(z)d_{\rm
L}(z)}\bigg{]}\frac{\sqrt{\langle v^{2}\rangle}}{c},$ (21)

where the peculiar velocity $\sqrt{\langle v^{2}\rangle}$ of the source with
respect to the Hubble flow is roughly set to $500\,\mathrm{km\,s^{-1}}$.

The error from the redshift measurement of the EM counterpart could be ignored
if the redshift is measured spectroscopically. But when using photometric
redshift for a distant source, this factor should be taken into account. We
estimate the error on the redshift measurement as

$\displaystyle\sigma_{d_{\rm L}}^{\rm reds}=\frac{\partial d_{\rm L}}{\partial
z}(\Delta z)_{n},$ (22)

with $(\Delta z)_{n}\simeq 0.03(1+z_{n})$ Ilbert  _et al._ (2013).

###  II.3 Identifications of electromagnetic counterparts

In addition to the luminosity distance obtained from the GW waveform, the
realization of standard siren also requires the redshift information.
Moreover, the number of GW-EM events has a great impact on the estimations of
cosmological parameters, so we give an analysis on the EM counterpart below.

In the process of the merger of MBHB with external magnetic fields, it is
assumed that EM radiations could be emitted in both the radio and optical
bands Palenzuela  _et al._ (2010); O’Shaughnessy  _et al._ (2011); Moesta  _et
al._ (2012); Kaplan  _et al._ (2011); Shi  _et al._ (2012); Blandford and
Znajek (1977); Meier (2001); Dotti  _et al._ (2012). The radiation in the
radio frequency band can be detected using SKA SKA , to identify the host
galaxy of the GW source. Subsequently, the radiation in the optical band can
be measured spectroscopically or photometrically through optical/IR projects,
such as the Vera C. Rubin Observatory Ivezić  _et al._ (2019) and the European
Extremely Large Telescope (E-ELT) ELT , to obtain the redshift information.

In practice, before the EM projects make observations, it is necessary to
accurately locate the GW source. Because the fields of view of the EM
projects, such as SKA and ELT, are about 10 ${\rm deg}^{2}$, in our
simulation, we choose those GW events with $\Delta\Omega<10~{}{\rm deg}^{2}$.
From our analyses in Sec. II.2, one can see that the most obvious advantage of
a GW detection network over a single detector is that it can greatly improve
the capability of locating GW sources, which is also discussed in Ref. Ruan
_et al._ (2020b). It should be mentioned that we use a conservative scenario
to estimate the location parameter by using only the inspiral phase Tamanini
_et al._ (2016). A more optimistic scenario is to include the merger and
ringdown phases, which could lead to more GW-EM events. We leave this issue
for future research.

After locating the GW source within $10~{}{\rm deg}^{2}$ by the GW detectors,
we need to further uniquely identify the host galaxy by the EM counterpart. Of
course, if the host galaxy cannot be uniquely identified, statistical methods
can also be applied in estimating cosmological parameters. We leave the
relevant discussion in the future work. To simulate the EM counterpart, it is
necessary to understand the formation mechanism of MBHB and its external
environment. Regarding these aspects, different theoretical models have been
proposed Klein  _et al._ (2016). In this paper, we use an analysis method
similar to that used in Refs. Tamanini  _et al._ (2016); Yang (2021). Next, we
briefly introduce this method. More details can be found in Ref. Tamanini  _et
al._ (2016).

We first simulate the EM radiation in the radio band, which will be used to
uniquely identify the host galaxy. The total luminosity $L_{\rm radio}$ in the
radio band consists of two parts. One part is the dual jet $L_{\rm flare}$
emitted when a binary is close to merging, which is caused by the twisting of
external magnetic field lines by the rapidly inspiralling MBHB Palenzuela  _et
al._ (2010); O’Shaughnessy  _et al._ (2011); Moesta  _et al._ (2012); Kaplan
_et al._ (2011); Shi  _et al._ (2012), being given by

$\displaystyle L_{\rm flare}=\epsilon_{\rm edd}\epsilon_{\rm radio}(v/v_{\rm
max})^{2}q^{2}L_{\rm edd}.$ (23)

The factor $(v/v_{\rm max})^{2}$ describes the luminosity evolution as the
MBHB inspirals ($v_{\rm max}=c/\sqrt{3}$ is the circular speed at the
innermost stable circular orbit for a binary of BHs, and $v$ is the binary’s
coordinate relative circular velocity O’Shaughnessy  _et al._ (2011)).
$\epsilon_{\rm radio}$ is the fraction of EM radiations emitted in the radio
band (i.e., a radio-to-bolometric luminosity correction), and is set to a
fiducial value of 0.1. $\epsilon_{\rm edd}$ is the Eddington ratio, which is
calculated according to the formulas shown in Appendix A of Ref. Tamanini  _et
al._ (2016). $q=M_{2}/M_{1}\leq 1$ is the binary’s mass ratio. The other part
of $L_{\rm radio}$ is the standard radio jet due to the Blandford-Znajeck
effect Blandford and Znajek (1977); Meier (2001), with luminosity dependent on
the mass accretion rate. Following Ref. Tamanini  _et al._ (2016), we use the
jet luminosity

$\displaystyle L_{\rm jet}=\begin{cases}10^{42.7}{\rm erg}~{}{\rm
s}^{-1}(\frac{\alpha}{0.01})^{-0.1}m_{9}^{0.9}(\frac{\dot{m}}{0.1})^{6/5}(1+1.1a_{1}+0.29a_{1}^{2}),~{}{\rm
if}~{}10^{-2}\leq\epsilon_{\rm edd}\leq 0.3,\\\ 10^{45.7}{\rm erg}~{}{\rm
s}^{-1}(\frac{\alpha}{0.3})^{-0.1}m_{9}(\frac{\dot{m}}{0.1})g^{2}(0.55f^{2}+1.5fa_{1}+a_{1}^{2}),~{}{\rm
otherwise.}\end{cases}$ (24)

We assume the Shakura-Sunyev viscosity parameter $\alpha=0.1$; $m_{9}$ is
defined as $m_{9}=M_{1}/(10^{9}M_{\odot})$; $\dot{m}$ is the central accretion
rate, which is calculated from Appendix A of Ref. Tamanini  _et al._ (2016);
$a_{1}$ is the spin parameter of the BH with the mass of $M_{1}$; $f$ and $g$
are dimensionless quantities regulating the angular velocity and the azimuthal
magnetic field, respectively, and are set to $f=1$ and $g=2.3$ Meier (2001).
The total luminosity in the radio band is given by $L_{\rm radio}=L_{\rm
flare}+L_{\rm jet}$. For the GW event satisfying $L_{\rm radio}\geq 4\pi
d_{\rm L}^{2}F^{\rm SKA}_{\rm min}$ O’Shaughnessy  _et al._ (2011), its EM
radiation in the radio band is expected to be detected by SKA, thus its host
galaxy could be uniquely identified. Here, $F^{\rm SKA}_{\rm min}=\nu_{\rm
SKA}F^{\rm SKA}_{\nu,{\rm min}}$ is the detector’s flux limit, with $\nu_{\rm
SKA}\simeq 1.4~{}{\rm GHz}$ and $F^{\rm SKA}_{\nu,{\rm min}}\simeq 1~{}\mu{\rm
Jy}$. It should be noted that we assume that the radio radiation is isotropic,
according to Ref. Tamanini  _et al._ (2016). Actually, the synchrotron
emission within the jet is beamed along the jet, which makes it impossible to
detect the events with jets being not towards the earth, thus reducing the
number of GW-EM events. However, collimation also implies a larger flux for a
given source luminosity, which makes some intrinsically fainter sources being
observable, thus increasing the number of GW-EM events. Because these two
factors have opposite effects on the number of GW-EM events, and may
counteract each other, we assume that the radiation is isotropic for the
purpose of simplification.

Only the radio identification cannot complete the measurement of the redshift,
so the optical/IR facilities are needed to observe the spectral features to
obtain the redshift. The host galaxy’s luminosity in the $K$-band, $L_{k}$, is
computed by converting the host total stellar mass into luminosity Tamanini
_et al._ (2016). According to the results of Ref. Bruzual and Charlot (2003),
for young stellar populations at moderate redshift, the mass-to-light ratio
$M/L_{k}$ falls in the range $0.01-0.05$. We assume a fiducial $M/L_{k}=0.03$
in the simulation. By converting $L_{k}$ into apparent magnitude $m_{\rm
gal}$, we assume that the redshifts of MBHB merger events that satisfy the
following relationship can be measured by ELT,

$\displaystyle m_{\rm gal}=82.5-\frac{5}{2}{\rm
log}_{10}\left(\frac{L_{k}}{3.02}\frac{{\rm s}}{{\rm erg}}\right)+5{\rm
log}_{10}\left(\frac{d_{\rm L}}{{\rm pc}}\right)\leq m_{{\rm ELT}},$ (25)

with the detection threshold $m_{{\rm ELT}}$ being set to 31.3, which is the
photometric limiting magnitude of ELT corresponding to $J$-band and $H$-band
Davies  _et al._ (2010). In principle, the detection threshold should be set
to 30.2 that is the limiting magnitude of $K$-band, because the host galaxy’s
luminosity in $K$-band is used to calculate apparent magnitude. Actually,
MICADO (Multi-AO Imaging Camera for Deep Observations) on ELT will cover the
wavelength range of 1000–2400 nm ($J$-band to $K$-band), so we simply choose
the highest limiting magnitude, 31.3, as the detection threshold. This
simplification may lead to a very small overestimation of the detection
threshold, but will have no an obvious effect on the number of GW-EM events
and the cosmological analysis. Based on the above criterion, we can pick out
the GW events whose redshifts can be determined. For the redshift error
$\sigma_{d_{\rm L}}^{\rm reds}$, the authors of Ref. Tamanini  _et al._ (2016)
take into account it for the GW-EM events satisfying $27.2<m_{\rm gal}<31.3$,
with 27.2 being the spectroscopy limiting magnitude of ELT. In Ref. Speri  _et
al._ (2021), the authors use a more simple method, namely taking into account
the redshift error for all the GW-EM events of $z>2$. Namely, it is assumed
that the redshifts of GW events with $z<2$ are measured by spectroscopy, while
those with $z>2$ are measured by photometry, because the spectroscopic
redshift in the range of $z>2$ is usually unavailable Dahlen  _et al._ (2013);
Speri  _et al._ (2021). In this work, we adopt the simplified method in Ref.
Speri  _et al._ (2021). Actually, one can see from Figure 2 that most GW-EM
events are at $z>2$. Thus, in our analysis, the redshift error is actually
considered for most data points.

###  II.4 Standard siren catalog

Based on the methods discussed in the previous sections, we construct the
standard siren catalogs in preparation for cosmological parameter constraints.
In this work, we discuss three population models of MBHB Madau and Rees
(2001); Volonteri  _et al._ (2008), i.e., pop III, Q3d, and Q3nod, which are
proposed Klein  _et al._ (2016) for the simulation of standard sirens,
according to the birth mechanism of MBHs and whether there exists a delay
between the mergers of MBHB and their host galaxies. For the total numbers of
the MBHB merger events within 5 years, we estimate them based on the event
rates given in Table I of Ref. Klein  _et al._ (2016). Specifically, the
numbers are 877, 41, and 610 for pop III, Q3d, and Q3nod, respectively. For
the redshift distribution and the mass distribution of MBHBs, we give
numerical fitting formulas of the curves shown in Figure 3 of Ref. Klein  _et
al._ (2016). The predicted MBHB merger rates as functions of redshift are
given by

$\displaystyle R(z)_{{\rm pop~{}III}}=$
$\displaystyle\begin{cases}2.11z,&0\leq z\leq 9,\\\ -1.8z+35.2,&9<z\leq 19,\\\
\end{cases}$ (26) $\displaystyle R(z)_{\rm{Q3d}}=$
$\displaystyle\begin{cases}0.43z,&0\leq z\leq 3.5,\\\ -0.18z+2.12,&3.5<z\leq
12,\\\ \end{cases}$ (27) $\displaystyle R(z)_{\rm{Q3nod}}=$
$\displaystyle\begin{cases}1.67z,&0\leq z\leq 6,\\\ -0.69z+14.62,&6<z\leq
19.\end{cases}$ (28)

The predicted MBHB merger rates as functions of the total redshifted mass,
$M_{z}$, are given by

$\displaystyle R(M_{z})_{{\rm pop~{}III}}=$
$\displaystyle\begin{cases}10^{-10.6}\left(\frac{M_{z}}{M_{\odot}}\right)^{3.2},&10^{3}M_{\odot}\leq
M_{z}\leq 10^{4}M_{\odot},\\\
10^{5.32}\left(\frac{M_{z}}{M_{\odot}}\right)^{-0.78},&10^{4}M_{\odot}<M_{z}\leq
10^{8}M_{\odot},\\\ \end{cases}$ (29) $\displaystyle R(M_{z})_{\rm{Q3d}}=$
$\displaystyle\begin{cases}10^{-4.4}\left(\frac{M_{z}}{M_{\odot}}\right)^{{0.81}},&10^{4}M_{\odot}\leq
M_{z}\leq 10^{6.3}M_{\odot},\\\
10^{6.65}\left(\frac{M_{z}}{M_{\odot}}\right)^{-0.94},&10^{6.3}M_{\odot}<M_{z}\leq
10^{8}M_{\odot},\end{cases}$ (30) $\displaystyle R(M_{z})_{\rm{Q3nod}}=$
$\displaystyle\begin{cases}10^{-5.44}\left(\frac{M_{z}}{M_{\odot}}\right)^{1.2},&10^{4}M_{\odot}\leq
M_{z}\leq 10^{6.2}M_{\odot},\\\
10^{9.75}\left(\frac{M_{z}}{M_{\odot}}\right)^{-1.25},&10^{6.2}M_{\odot}<M_{z}\leq
10^{8}M_{\odot},\end{cases}$ (31)

with $M_{z}=M(1+z)$. $R(z)$ and $R(M_{z})$ are in units of ${\rm yr}^{-1}$.
The sky location ($\theta$, $\phi$), the binary inclination $\iota$, the
polarization angle $\psi$, and the coalescence phase $\phi_{c}$ are randomly
chosen in the ranges of $[0,\pi]$, $[0,2\pi]$, $[0,\pi]$, $[0,2\pi]$, and
$[0,2\pi]$, respectively.

Figure 2: The standard siren catalogs simulated from Taiji and the LISA-Taiji
network within 5-year operation time based on the pop III, Q3d, and Q3nod
models.

After selecting the GW events satisfying SNR $>$ 8 and $\Delta\Omega<10~{}{\rm
deg}^{2}$, and further selecting the GW events whose redshifts can be detected
by SKA and ELT, we obtain the useful standard sirens. For each standard siren,
we calculate its luminosity distance $d_{\rm L}$ and the error of luminosity
distance $\sigma_{d_{\rm L}}$. The fiducial values of cosmological parameters
are set to the best-fit values of the Planck 2018 results Aghanim  _et al._
(2020). Then, for each MBHB model, we can construct a standard siren catalog
including the redshift $z$, luminosity distance $d_{\rm L}$, and error of
luminosity distance $\sigma_{d_{\rm L}}$ of MBHBs.

We show the simulated standard sirens in Figure 2. Firstly, it is found that
the numbers of standard sirens detected by the LISA-Taiji network are much
more than those detected by the single Taiji mission, which also can be seen
in Table 1. For each MBHB model, the LISA-Taiji network approximately doubles
the detection number compared to the single Taiji mission. This is due to the
fact that the network could improve SNR and the location accuracy of the GW
event. Secondly, the network can detect the GW events at higher redshifts due
to the improvement of SNR. For example, from the middle panel of Figure 2, we
can see that for the Q3d model, the redshifts of the standard sirens detected
by the network could reach $z\sim 9$. Thirdly, the LISA-Taiji network can
improve the measurement accuracy of luminosity distance to some extent, which
can be seen from the error bars in Figure 2. These improvements make us expect
that the LISA-Taiji network can greatly improve the capability of constraining
cosmological parameters.

Table 1:  The numbers of the standard sirens simulated from Taiji and the LISA-Taiji network within 5-year operation time, based on the pop III, Q3d, and Q3nod models of MBHB population. Model | Number of standard sirens
---|---
| Taiji | network
pop III |  | 25 | 50
Q3d |  | 12 | 20
Q3nod |  | 24 | 44

##  III Cosmological parameter estimation

In this section, we shall report the constraint results of cosmological
parameters. In theory, the luminosity distance $d_{\rm L}$ of a GW source at
redshift $z$ is determined by a specific cosmological model. The $\Lambda$CDM
model [$w(z)=-1$], the $w$CDM model [$w(z)=\rm{constant}$], and the
Chevallier–Polarski–Linder (CPL) model [$w(z)=w_{\rm{0}}+w_{\rm{a}}z/(1+z)$]
Chevallier and Polarski (2001); Linder (2003) are considered in this paper.
For the CMB data, we employ the “Planck distance priors” from the Planck 2018
observation Chen  _et al._ (2019). We use $\sigma(\xi)$ and $\varepsilon(\xi)$
to represent the absolute error and the relative error of the parameter $\xi$,
respectively, with $\varepsilon(\xi)$ defined as
$\varepsilon(\xi)=\sigma(\xi)/\xi$.

Table 2:  The absolute errors ($1\sigma$) and the relative errors of the cosmological parameters in the $\Lambda$CDM, $w$CDM, and CPL models using the mock data of the LISA-Taiji network. Here, $\sigma(\xi)$ and $\varepsilon(\xi)$ represent the absolute and relative errors of the parameter $\xi$, respectively. Note also that $\varepsilon(\xi)$ is defined as $\varepsilon(\xi)=\sigma(\xi)/\xi$, and $H_{0}$ is in units of ${\rm km\ s^{-1}\ Mpc^{-1}}$. Error |  $\Lambda$CDM |  |  $w$CDM |  | CPL |
---|---|---|---|---|---|---
pop III | Q3d | Q3nod |  | pop III | Q3d | Q3nod |  | pop III | Q3d | Q3nod |
$\sigma(\Omega_{\rm m})$ | $0.025$ | $0.073$ | $0.027$ |  | $0.036$ | $0.078$ | $0.037$ |  | $0.054$ | $0.099$ | $0.056$ |
$\sigma(H_{0})$ | $0.86$ | $3.25$ | $0.94$ |  | $1.85$ | $10.10$ | $1.75$ |  | $2.40$ | $10.00$ | $2.10$ |
$\sigma(w)$ | $-$ | $-$ | $-$ |  | $0.245$ | $0.735$ | $0.230$ |  | $-$ | $-$ | $-$ |
$\sigma(w_{0})$ | $-$ | $-$ | $-$ |  | $-$ | $-$ | $-$ |  | $0.380$ | $0.965$ | $0.340$ |
$\sigma(w_{a})$ | $-$ | $-$ | $-$ |  | $-$ | $-$ | $-$ |  | $1.95$ | $-$ | $2.05$ |
$\varepsilon(\Omega_{\rm m})$ | $0.078$ | $0.212$ | $0.084$ |  | $0.113$ | $0.258$ | $0.114$ |  | $0.165$ | $0.304$ | $0.171$ |
$\varepsilon(H_{0})$ | $0.013$ | $0.049$ | $0.014$ |  | $0.027$ | $0.141$ | $0.026$ |  | $0.036$ | $0.141$ | $0.031$ |
$\varepsilon(w)$ | $-$ | $-$ | $-$ |  | $0.229$ | $0.525$ | $0.215$ |  | $-$ | $-$ | $-$ |
$\varepsilon(w_{0})$ | $-$ | $-$ | $-$ |  | $-$ | $-$ | $-$ |  | $0.409$ | $0.832$ | $0.378$ |
Table 3:  The absolute errors ($1\sigma$) and the relative errors of the
cosmological parameters in the $\Lambda$CDM, $w$CDM, and CPL models using the
CMB, Taiji(Q3nod), network(Q3nod), CMB+Taiji(Q3nod), and CMB+network(Q3nod)
data. Here, $\sigma(\xi)$ and $\varepsilon(\xi)$ represent the absolute and
relative errors of the parameter $\xi$, respectively. Note also that
$\varepsilon(\xi)$ is defined as $\varepsilon(\xi)=\sigma(\xi)/\xi$, and
$H_{0}$ is in units of ${\rm km\ s^{-1}\ Mpc^{-1}}$.

Error |  $\Lambda$CDM |  |  $w$CDM |  | CPL
---|---|---|---|---|---
CMB | Taiji | network | CMB+Taiji | CMB+network |  | CMB | Taiji | network | CMB+Taiji | CMB+network |  | CMB | Taiji | network | CMB+Taiji | CMB+network
$\sigma(\Omega_{\rm m})$ | $0.009$ | $0.068$ | $0.027$ | $0.008$ | $0.007$ |  | $0.057$ | $0.087$ | $0.037$ | $0.026$ | $0.010$ |  | $0.059$ | $0.086$ | $0.056$ | $0.032$ | $0.018$
$\sigma(H_{0})$ | $0.61$ | $2.60$ | $0.94$ | $0.58$ | $0.46$ |  | $6.15$ | $3.90$ | $1.75$ | $2.80$ | $1.00$ |  | $6.25$ | $4.05$ | $2.10$ | $3.25$ | $1.90$
$\sigma(w)$ | $-$ | $-$ | $-$ | $-$ | $-$ |  | $0.215$ | $0.695$ | $0.230$ | $0.097$ | $0.042$ |  | $-$ | $-$ | $-$ | $-$ | $-$
$\sigma(w_{0})$ | $-$ | $-$ | $-$ | $-$ | $-$ |  | $-$ | $-$ | $-$ | $-$ | $-$ |  | $0.575$ | $0.885$ | $0.340$ | $0.510$ | $0.230$
$\sigma(w_{a})$ | $-$ | $-$ | $-$ | $-$ | $-$ |  | $-$ | $-$ | $-$ | $-$ | $-$ |  | $-$ | $-$ | $2.05$ | $1.80$ | $0.67$
$\varepsilon(\Omega_{\rm m})$ | $0.027$ | $0.201$ | $0.084$ | $0.026$ | $0.020$ |  | $0.178$ | $0.263$ | $0.114$ | $0.082$ | $0.030$ |  | $0.183$ | $0.246$ | $0.171$ | $0.097$ | $0.056$
$\varepsilon(H_{0})$ | $0.009$ | $0.039$ | $0.014$ | $0.009$ | $0.007$ |  | $0.090$ | $0.057$ | $0.026$ | $0.041$ | $0.015$ |  | $0.092$ | $0.059$ | $0.031$ | $0.049$ | $0.028$
$\varepsilon(w)$ | $-$ | $-$ | $-$ | $-$ | $-$ |  | $0.211$ | $0.489$ | $0.215$ | $0.097$ | $0.042$ |  | $-$ | $-$ | $-$ | $-$ | $-$
$\varepsilon(w_{0})$ | $-$ | $-$ | $-$ | $-$ | $-$ |  | $-$ | $-$ | $-$ | $-$ | $-$ |  | $0.871$ | $0.651$ | $0.378$ | $0.823$ | $0.242$

Firstly, let us take the simplest $\Lambda$CDM model as an example to discuss
the constraint results for the three MBHB models. From Figure 3, we find that
the pop III model gives the best constraints because it has the maximal event
number among the three MBHB models. The Q3nod model is similar with the pop
III model, but the Q3d model gives the worst constraints because it has the
fewest event number. Quantitatively, the pop III model gives the relative
errors $\varepsilon(\Omega_{\rm m})=7.8\%$ and $\varepsilon(H_{0})=1.3\%$, the
Q3nod model gives $\varepsilon(\Omega_{\rm m})=8.4\%$ and
$\varepsilon(H_{0})=1.4\%$, and the Q3d model gives $\varepsilon(\Omega_{\rm
m})=21.2\%$ and $\varepsilon(H_{0})=4.9\%$, as shown in Table 2. Here, we
notice that solely using the standard sirens from the LISA-Taiji network could
provide a tight constraint on $H_{0}$, with the precision close to 1%, which
is the standard of precision cosmology. Compared with the single Taiji
mission, the LISA-Taiji network reduces the absolute errors of $\Omega_{\rm
m}$ and $H_{0}$ by 60.3% and 63.8%, respectively, which can be seen from Table
3. This indicates that the GW detection network could provide a more accurate
measurement of $H_{0}$ than a single detector, which is helpful for solving
the Hubble constant tension.

Now let us have a look at the improvements of the constraints on the EoS
parameter of dark energy. We first discuss the $w$CDM model that has only one
dark-energy EoS parameter. In Figure 4, we show the two-dimensional posterior
contours in the $\Omega_{\rm m}-w$ and $w-H_{0}$ planes using the data of
Taiji, network, CMB, CMB+Taiji, and CMB+network. From the purple contour and
the green contour, we can see that the network can give much better
constraints than the single Taiji mission. The gray contour represents the
constraint from the CMB data, and this contour is almost orthogonal to the
green contour, which indicates that the standard siren data could
significantly break the parameter degeneracies. By comparing the red contour
and the blue contour, we can see that, although both the single Taiji mission
and the network can break degeneracies, the network can improve the constraint
accuracies much better. As we mentioned above, this is because the LISA-Taiji
network can detect more GW-EM events than Taiji, and can detect the events at
higher redshifts. In addition, the network also improves SNR, thus reducing
the error of $d_{\rm L}$. Concretely, the LISA-Taiji network could reduce the
absolute error of $w$ by 66.9% compared with the single Taiji mission. The
data combination CMB+network could reduce the absolute error of $w$ by 56.7%
compared with CMB+Taiji. It is worth emphasizing that using the CMB+network
data, the constraint precision of $w$ could reach $4.2\%$, which is comparable
with the result of Planck 2018 TT,TE,EE+lowE+lensing+SNe+BAO Aghanim  _et al._
(2020). Note also that here we use only the Planck distance priors, but not
the Planck full data of CMB power spectra.

For the CPL model that has two dark-energy EoS parameters $w_{0}$ and $w_{a}$,
we also find that the LISA-Taiji network data give better constraints compared
with the single Taiji mission. The detailed results are shown in Table 3. For
the parameter $w_{0}$, Taiji gives the relative error
$\varepsilon(w_{0})=65.1\%$, and the LISA-Taiji network gives the relative
error $\varepsilon(w_{0})=37.8\%$. Compared with the single Taiji mission, the
LISA-Taiji network can reduce the absolute error of $w_{0}$ by 61.6%. For the
parameter $w_{a}$, Taiji cannot constrain $w_{a}$ well, but the LISA-Taiji
network gives the absolute error $\sigma(w_{a})=2.05$. When the CMB data are
combined, CMB+Taiji gives $\varepsilon(w_{0})=82.3\%$ and
$\sigma(w_{a})=1.80$, and CMB+network gives $\varepsilon(w_{0})=24.2\%$ and
$\sigma(w_{a})=0.67$. Compared with CMB+Taiji, CMB+network could reduce the
absolute errors of $w_{0}$ and $w_{a}$ by 41.5% and 54.9%, respectively. This
indicates that the future LISA-Taiji network combined with the CMB observation
will play an important role in constraining the EoS parameter of dark energy.

Figure 3: The two-dimensional marginalized contours (68.3% and 95.4%
confidence level) in the $\Omega_{m}$–$H_{0}$ plane considering three MBHB
models for the $\Lambda$CDM model.

Figure 4: The two-dimensional marginalized contours (68.3% and 95.4%
confidence level) in the $\Omega_{\rm m}-w$ and $w-H_{0}$ planes for the
$w$CDM. The mock data of Taiji and the LISA-Taiji network are simulated based
on the Q3nod model.

It should be noted that we make some assumptions and approximations in this
work. For the redshift distribution and the mass distribution of MBHBs, we use
the numerical fitting formulas. For the redshift measurement, we assume that
the redshifts of GW events satisfying $z>2$ are measured by photometry, so we
take into account the redshift errors for these events. In addition, we simply
assume that the radio radiation in the EM counterpart is isotropic. We have
discussed the rationalities of these assumptions in the previous sections. The
main purpose of this work is to make a preliminary forecast on the capability
of the LISA-Taiji network of improving the estimations of cosmological
parameters, compared with a single GW detector. For this purpose, these
assumptions do not affect our main conclusions. For a more specific
investigation on these issues, we leave it for future works.

##  IV Conclusion

In this work, we forecast the capability of the future space-based GW
detection network to constrain cosmological parameters. We consider a
detection network composed of the European LISA mission and the Chinese Taiji
mission. The configuration angle between the LISA and Taiji is considered to
be $40^{\circ}$. Three models for MBHB, i.e., pop III, Q3d, and Q3nod, are
used to simulate the EM counterpart detection and estimate the number of GW-EM
events. Three typical cosmological models, i.e., the $\Lambda$CDM, $w$CDM, and
CPL models, are chosen as representatives.

We find that the LISA-Taiji network could significantly improve the constraint
accuracies of cosmological parameters compared with the single Taiji mission.
This is mainly due to three aspects: (i) the LISA-Taiji network could increase
the number of the GW-EM events; (ii) the LISA-Taiji network could detect the
GW events at higher redshift; and (iii) the LISA-Taiji network could reduce
the error of luminosity distance. Taking the Q3nod model as an example, the
LISA-Taiji network increases the GW-EM number from 24 to 44, compared with the
singe Taiji mission. The redshifts of GW events detected by the LISA-Taiji
network can reach to $z\sim 8$. The errors of luminosity distances are also
reduced to some extent, which can be seen from the error bars in Figure 2.

For the simplest $\Lambda$CDM model, we find that the pop III model gives the
best constraints due to its maximal event number among these three MBHB
models. The Q3nod model is similar with pop III, but the Q3d model gives the
worst constraints due to its fewest event number. The pop III model gives the
relative errors $\varepsilon(\Omega_{\rm m})=7.8\%$ and
$\varepsilon(H_{0})=1.3\%$; the Q3nod model gives $\varepsilon(\Omega_{\rm
m})=8.4\%$ and $\varepsilon(H_{0})=1.4\%$; the Q3d model gives
$\varepsilon(\Omega_{\rm m})=21.2\%$ and $\varepsilon(H_{0})=4.9\%$. So, we
find that solely using the standard sirens from the LISA-Taiji network, the
constraint precision of $H_{0}$ could reach about 1%, which is the standard of
precision cosmology. In addition, compared with the single Taiji mission, the
LISA-Taiji network could reduce the relative error of $H_{0}$ by 63.8%. This
indicates that the LISA-Taiji network is helpful in addressing the $H_{0}$
tension.

For the dark-energy EoS parameter, we first discuss it in the $w$CDM model.
Taking the Q3nod model as an example, we find that the LISA-Taiji network
reduces the absolute error of $w$ by 66.9% compared with the single Taiji
mission. We also see that the LISA-Taiji network improves the capability of
breaking the parameter degeneracies inherent in the CMB data. Specifically,
the CMB+network data reduces the absolute error of $w$ by 56.7% compared with
CMB+Taiji. The constraint precision of $w$ could reach $4.2\%$ using the
CMB+network data, which is comparable with the result of $Planck$ 2018
TT,TE,EE+lowE+lensing+SNe+BAO. For the CPL model that has two dark-energy EoS
parameters $w_{0}$ and $w_{a}$, the LISA-Taiji network can reduce the absolute
error of $w_{0}$ by 61.6%, compared with the single Taiji mission. For the
parameter $w_{a}$, Taiji cannot constrain $w_{a}$ well, but the LISA-Taiji
network can give the absolute error $\sigma(w_{a})=2.05$. Compared with
CMB+Taiji, CMB+network could reduce the absolute errors of $w_{0}$ and $w_{a}$
by 41.5% and 54.9%, respectively. This indicates that the future LISA-Taiji
network combined with the CMB observation will play an important role in
constraining the EoS parameters of dark energy.

In the next few decades, GW detectors are expected to form powerful detection
networks. This allows researchers to use multiple detectors to jointly detect
a GW source, thereby improving the measurement accuracies of source
parameters. At the same time, the detections of different GW sources in
multiple GW frequency bands can help us test the theory of gravity, and can
also lead to a more comprehensive understanding of the properties of compact
objects and the evolution history of the universe. We expect that the LISA-
Taiji network could play an important role in these aspects, and we will
further discuss them in depth in our future work.

###### Acknowledgements.

We are very grateful to Antoine Klein, Alberto Mangiagli, Alberto Sesana, and
Nicola Tamanini for fruitful discussions on the identifications of EM
counterparts, and also grateful to Tao Yang, Wen-Hong Ruan, and Ze-Wei Zhao
for discussions on the detections of GW sources. This work was supported by
the National Natural Science Foundation of China (Grants Nos. 11975072,
11835009, 11875102, and 11690021), the Liaoning Revitalization Talents Program
(Grant No. XLYC1905011), the Fundamental Research Funds for the Central
Universities (Grant No. N2005030), and the National Program for Support of
Top-Notch Young Professionals (Grant No. W02070050).

## References

  * Spergel  _et al._ (2003) D. N. Spergel _et al._ (WMAP), Astrophys. J. Suppl. 148, 175 (2003), arXiv:astro-ph/0302209 .
  * Bennett  _et al._ (2003) C. L. Bennett _et al._ (WMAP), Astrophys. J. Suppl. 148, 1 (2003), arXiv:astro-ph/0302207 .
  * Aghanim  _et al._ (2020) N. Aghanim _et al._ (Planck), Astron. Astrophys. 641, A6 (2020), arXiv:1807.06209 [astro-ph.CO] .
  * Riess  _et al._ (2019) A. G. Riess, S. Casertano, W. Yuan, L. M. Macri, and D. Scolnic, Astrophys. J. 876, 85 (2019), arXiv:1903.07603 [astro-ph.CO] .
  * Verde  _et al._ (2019) L. Verde, T. Treu, and A. G. Riess, Nature Astron. 3, 891 (2019), arXiv:1907.10625 [astro-ph.CO] .
  * Riess (2019) A. G. Riess, Nature Rev. Phys. 2, 10 (2019), arXiv:2001.03624 [astro-ph.CO] .
  * Weinberg (1989) S. Weinberg, Rev. Mod. Phys. 61, 1 (1989).
  * Sahni and Starobinsky (2000) V. Sahni and A. A. Starobinsky, Int. J. Mod. Phys. D 9, 373 (2000), arXiv:astro-ph/9904398 .
  * Bean  _et al._ (2005) R. Bean, S. M. Carroll, and M. Trodden, (2005), arXiv:astro-ph/0510059 .
  * Guo  _et al._ (2019a) R.-Y. Guo, J.-F. Zhang, and X. Zhang, JCAP 02, 054 (2019a), arXiv:1809.02340 [astro-ph.CO] .
  * Zhao  _et al._ (2020a) M. Zhao, R. Guo, D. He, J. Zhang, and X. Zhang, Sci. China Phys. Mech. Astron. 63, 230412 (2020a), arXiv:1810.11658 [astro-ph.CO] .
  * Di Valentino  _et al._ (2021) E. Di Valentino, O. Mena, S. Pan, L. Visinelli, W. Yang, A. Melchiorri, D. F. Mota, A. G. Riess, and J. Silk,  (2021), 10.1088/1361-6382/ac086d, arXiv:2103.01183 [astro-ph.CO] .
  * Li  _et al._ (2019) H.-L. Li, L. Feng, J.-F. Zhang, and X. Zhang, Sci. China Phys. Mech. Astron. 62, 120411 (2019), arXiv:1812.00319 [astro-ph.CO] .
  * Vagnozzi (2020) S. Vagnozzi, Phys. Rev. D 102, 023518 (2020), arXiv:1907.07569 [astro-ph.CO] .
  * Guo  _et al._ (2019b) R.-Y. Guo, L. Zhang, J.-F. Zhang, and X. Zhang, Sci. China Phys. Mech. Astron. 62, 30411 (2019b), arXiv:1801.02187 [astro-ph.CO] .
  * Li  _et al._ (2020a) H.-L. Li, J.-F. Zhang, and X. Zhang, Commun. Theor. Phys. 72, 125401 (2020a), arXiv:2005.12041 [astro-ph.CO] .
  * Feng  _et al._ (2020a) L. Feng, D.-Z. He, H.-L. Li, J.-F. Zhang, and X. Zhang, Sci. China Phys. Mech. Astron. 63, 290404 (2020a), arXiv:1910.03872 [astro-ph.CO] .
  * Zhang  _et al._ (2020a) M. Zhang, J.-F. Zhang, and X. Zhang, Commun. Theor. Phys. 72, 125402 (2020a), arXiv:2005.04647 [astro-ph.CO] .
  * Feng  _et al._ (2020b) L. Feng, H.-L. Li, J.-F. Zhang, and X. Zhang, Sci. China Phys. Mech. Astron. 63, 220401 (2020b), arXiv:1903.08848 [astro-ph.CO] .
  * Lin  _et al._ (2020) M.-X. Lin, W. Hu, and M. Raveri, Phys. Rev. D 102, 123523 (2020), arXiv:2009.08974 [astro-ph.CO] .
  * Guo  _et al._ (2020) R.-Y. Guo, J.-F. Zhang, and X. Zhang, Sci. China Phys. Mech. Astron. 63, 290406 (2020), arXiv:1910.13944 [astro-ph.CO] .
  * Hryczuk and Jodłowski (2020) A. Hryczuk and K. Jodłowski, Phys. Rev. D 102, 043024 (2020), arXiv:2006.16139 [hep-ph] .
  * Zhang  _et al._ (2020b) J.-F. Zhang, B. Wang, and X. Zhang, Sci. China Phys. Mech. Astron. 63, 280411 (2020b), arXiv:1907.00179 [astro-ph.CO] .
  * Gao  _et al._ (2021) L.-Y. Gao, S.-S. Xue, and X. Zhang, (2021), arXiv:2101.10714 [astro-ph.CO] .
  * Zhang (2019) X. Zhang, Sci. China Phys. Mech. Astron. 62, 110431 (2019), arXiv:1905.11122 [astro-ph.CO] .
  * Xu and Zhang (2020) Y. Xu and X. Zhang, Sci. China Phys. Mech. Astron. 63, 270431 (2020), arXiv:2002.00572 [astro-ph.CO] .
  * Li and Zhang (2020) H. Li and X. Zhang, Sci. Bull. 65, 1419 (2020), arXiv:2005.10458 [astro-ph.CO] .
  * Zhang  _et al._ (2021) M. Zhang, B. Wang, J.-Z. Qi, Y. Xu, J.-F. Zhang, and X. Zhang, (2021), arXiv:2102.03979 [astro-ph.CO] .
  * Wang  _et al._ (2021) L.-F. Wang, D.-Z. He, J.-F. Zhang, and X. Zhang, (2021), arXiv:2102.09331 [astro-ph.CO] .
  * Schutz (1986) B. F. Schutz, Nature 323, 310 (1986).
  * Holz and Hughes (2005) D. E. Holz and S. A. Hughes, Astrophys. J. 629, 15 (2005), arXiv:astro-ph/0504616 .
  * Cai  _et al._ (2018a) R.-G. Cai, T.-B. Liu, and S.-J. Wang, Phys. Rev. D 97, 023027 (2018a), arXiv:1710.02425 [hep-ph] .
  * Di Valentino and Melchiorri (2018) E. Di Valentino and A. Melchiorri, Phys. Rev. D 97, 041301 (2018), arXiv:1710.06370 [astro-ph.CO] .
  * Wei (2018) J.-J. Wei, Astrophys. J. 868, 29 (2018), arXiv:1806.09781 [astro-ph.CO] .
  * Zhao  _et al._ (2018) W. Zhao, B. S. Wright, and B. Li, JCAP 10, 052 (2018), arXiv:1804.03066 [astro-ph.CO] .
  * Di Valentino  _et al._ (2018) E. Di Valentino, D. E. Holz, A. Melchiorri, and F. Renzi, Phys. Rev. D 98, 083523 (2018), arXiv:1806.07463 [astro-ph.CO] .
  * Du  _et al._ (2019) M. Du, W. Yang, L. Xu, S. Pan, and D. F. Mota, Phys. Rev. D 100, 043535 (2019), arXiv:1812.01440 [astro-ph.CO] .
  * Mifsud and van de Bruck (2019) J. Mifsud and C. van de Bruck, Mon. Not. Roy. Astron. Soc. 487, 900 (2019).
  * Wei (2019) J.-J. Wei, Astrophys. J. 876, 66 (2019), arXiv:1902.00223 [astro-ph.CO] .
  * Gray  _et al._ (2020) R. Gray _et al._ , Phys. Rev. D 101, 122001 (2020), arXiv:1908.06050 [gr-qc] .
  * Howlett and Davis (2020) C. Howlett and T. M. Davis, Mon. Not. Roy. Astron. Soc. 492, 3803 (2020), arXiv:1909.00587 [astro-ph.CO] .
  * Chassande-Mottin  _et al._ (2019) E. Chassande-Mottin, K. Leyde, S. Mastrogiovanni, and D. Steer, Phys. Rev. D 100, 083514 (2019), arXiv:1906.02670 [astro-ph.CO] .
  * Doctor (2020) Z. Doctor, Astrophys. J. Lett. 892, L16 (2020), arXiv:1912.12218 [astro-ph.HE] .
  * Fu  _et al._ (2019) X. Fu, L. Zhou, and J. Chen, Phys. Rev. D 99, 083523 (2019), arXiv:1903.09913 [gr-qc] .
  * Mukherjee  _et al._ (2021) S. Mukherjee, G. Lavaux, F. R. Bouchet, J. Jasche, B. D. Wandelt, S. M. Nissanke, F. Leclercq, and K. Hotokezaka, Astron. Astrophys. 646, A65 (2021), arXiv:1909.08627 [astro-ph.CO] .
  * Abbott  _et al._ (2021) B. P. Abbott _et al._ (LIGO Scientific, Virgo), Astrophys. J. 909, 218 (2021), arXiv:1908.06060 [astro-ph.CO] .
  * Palmese  _et al._ (2020) A. Palmese _et al._ (DES), Astrophys. J. Lett. 900, L33 (2020), arXiv:2006.14961 [astro-ph.CO] .
  * Chen (2020) H.-Y. Chen, Phys. Rev. Lett. 125, 201301 (2020), arXiv:2006.02779 [astro-ph.HE] .
  * Wang  _et al._ (2020a) B. Wang, Z. Zhu, A. Li, and W. Zhao, Astrophys. J. Suppl. 250, 6 (2020a), arXiv:2005.12875 [gr-qc] .
  * Chen  _et al._ (2021) H.-Y. Chen, P. S. Cowperthwaite, B. D. Metzger, and E. Berger, Astrophys. J. Lett. 908, L4 (2021), arXiv:2011.01211 [astro-ph.CO] .
  * Zhao  _et al._ (2011) W. Zhao, C. Van Den Broeck, D. Baskaran, and T. Li, Phys. Rev. D 83, 023005 (2011), arXiv:1009.0206 [astro-ph.CO] .
  * Cai and Yang (2017) R.-G. Cai and T. Yang, Phys. Rev. D 95, 044024 (2017), arXiv:1608.08008 [astro-ph.CO] .
  * Cai  _et al._ (2018b) R.-G. Cai, T.-B. Liu, X.-W. Liu, S.-J. Wang, and T. Yang, Phys. Rev. D 97, 103005 (2018b), arXiv:1712.00952 [astro-ph.CO] .
  * Yang  _et al._ (2019a) T. Yang, R. Holanda, and B. Hu, Astropart. Phys. 108, 57 (2019a), arXiv:1710.10929 [astro-ph.CO] .
  * Cai and Yang (2018) R.-G. Cai and T. Yang, EPJ Web Conf. 168, 01008 (2018), arXiv:1709.00837 [astro-ph.CO] .
  * Wang  _et al._ (2018) L.-F. Wang, X.-N. Zhang, J.-F. Zhang, and X. Zhang, Phys. Lett. B 782, 87 (2018), arXiv:1802.04720 [astro-ph.CO] .
  * Zhang  _et al._ (2019a) X.-N. Zhang, L.-F. Wang, J.-F. Zhang, and X. Zhang, Phys. Rev. D 99, 063510 (2019a), arXiv:1804.08379 [astro-ph.CO] .
  * Li  _et al._ (2020b) H.-L. Li, D.-Z. He, J.-F. Zhang, and X. Zhang, JCAP 06, 038 (2020b), arXiv:1908.03098 [astro-ph.CO] .
  * Zhang  _et al._ (2019b) J.-F. Zhang, M. Zhang, S.-J. Jin, J.-Z. Qi, and X. Zhang, JCAP 09, 068 (2019b), arXiv:1907.03238 [astro-ph.CO] .
  * Zhang  _et al._ (2020c) J.-F. Zhang, H.-Y. Dong, J.-Z. Qi, and X. Zhang, Eur. Phys. J. C 80, 217 (2020c), arXiv:1906.07504 [astro-ph.CO] .
  * Yan  _et al._ (2019) C. Yan, W. Zhao, and Y. Lu,  (2019), 10.3847/1538-4357/ab60a6, arXiv:1912.04103 [astro-ph.GA] .
  * Yang  _et al._ (2019b) W. Yang, S. Vagnozzi, E. Di Valentino, R. C. Nunes, S. Pan, and D. F. Mota, JCAP 07, 037 (2019b), arXiv:1905.08286 [astro-ph.CO] .
  * Jin  _et al._ (2020) S.-J. Jin, D.-Z. He, Y. Xu, J.-F. Zhang, and X. Zhang, JCAP 03, 051 (2020), arXiv:2001.05393 [astro-ph.CO] .
  * Qi  _et al._ (2021) J.-Z. Qi, S.-J. Jin, X.-L. Fan, J.-F. Zhang, and X. Zhang, (2021), arXiv:2102.01292 [astro-ph.CO] .
  * Abbott  _et al._ (2017a) B. Abbott _et al._ (LIGO Scientific, Virgo), Phys. Rev. Lett. 119, 161101 (2017a), arXiv:1710.05832 [gr-qc] .
  * Abbott  _et al._ (2017b) B. 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] .
  * Abbott  _et al._ (2017c) B. Abbott _et al._ (LIGO Scientific, Virgo, Fermi-GBM, INTEGRAL), Astrophys. J. Lett. 848, L13 (2017c), arXiv:1710.05834 [astro-ph.HE] .
  * Abbott  _et al._ (2017d) B. Abbott _et al._ (LIGO Scientific, Virgo, 1M2H, Dark Energy Camera GW-E, DES, DLT40, Las Cumbres Observatory, VINROUGE, MASTER), Nature 551, 85 (2017d), arXiv:1710.05835 [astro-ph.CO] .
  * Graham  _et al._ (2020) M. Graham _et al._ , Phys. Rev. Lett. 124, 251102 (2020), arXiv:2006.14122 [astro-ph.HE] .
  * Abbott  _et al._ (2020) R. Abbott _et al._ (LIGO Scientific, Virgo), Phys. Rev. Lett. 125, 101102 (2020), arXiv:2009.01075 [gr-qc] .
  * Chen  _et al._ (2020) H.-Y. Chen, C.-J. Haster, S. Vitale, W. M. Farr, and M. Isi, (2020), arXiv:2009.14057 [astro-ph.CO] .
  * Mukherjee  _et al._ (2020) S. Mukherjee, A. Ghosh, M. J. Graham, C. Karathanasis, M. M. Kasliwal, I. Magaña Hernandez, S. M. Nissanke, A. Silvestri, and B. D. Wandelt, (2020), arXiv:2009.14199 [astro-ph.CO] .
  * Del Pozzo (2012) W. Del Pozzo, Phys. Rev. D 86, 043011 (2012), arXiv:1108.1317 [astro-ph.CO] .
  * Chen  _et al._ (2018) H.-Y. Chen, M. Fishbach, and D. E. Holz, Nature 562, 545 (2018), arXiv:1712.06531 [astro-ph.CO] .
  * Fishbach  _et al._ (2019) M. Fishbach _et al._ (LIGO Scientific, Virgo), Astrophys. J. Lett. 871, L13 (2019), arXiv:1807.05667 [astro-ph.CO] .
  * Feeney  _et al._ (2019) S. M. Feeney, H. V. Peiris, A. R. Williamson, S. M. Nissanke, D. J. Mortlock, J. Alsing, and D. Scolnic, Phys. Rev. Lett. 122, 061105 (2019), arXiv:1802.03404 [astro-ph.CO] .
  * Ding  _et al._ (2019) X. Ding, M. Biesiada, X. Zheng, K. Liao, Z. Li, and Z.-H. Zhu, JCAP 04, 033 (2019), arXiv:1801.05073 [astro-ph.CO] .
  * Soares-Santos  _et al._ (2019) M. Soares-Santos _et al._ (DES, LIGO Scientific, Virgo), Astrophys. J. Lett. 876, L7 (2019), arXiv:1901.01540 [astro-ph.CO] .
  * Lagos  _et al._ (2019) M. Lagos, M. Fishbach, P. Landry, and D. E. Holz, Phys. Rev. D 99, 083504 (2019), arXiv:1901.03321 [astro-ph.CO] .
  * Yu  _et al._ (2020) J. Yu, Y. Wang, W. Zhao, and Y. Lu, Mon. Not. Roy. Astron. Soc. 498, 1786 (2020), arXiv:2003.06586 [astro-ph.CO] .
  * Punturo  _et al._ (2010) M. Punturo _et al._ , Class. Quant. Grav. 27, 194002 (2010).
  * (82) “ET,” www.et-gw.eu/.
  * Abbott  _et al._ (2017e) B. P. Abbott _et al._ (LIGO Scientific), Class. Quant. Grav. 34, 044001 (2017e), arXiv:1607.08697 [astro-ph.IM] .
  * (84) “CE,” https://cosmicexplorer.org/.
  * Totani (2013) T. Totani, Pub. Astron. Soc. Jpn. 65, L12 (2013), arXiv:1307.4985 [astro-ph.HE] .
  * Mingarelli  _et al._ (2015) C. M. F. Mingarelli, J. Levin, and T. J. W. Lazio, Astrophys. J. Lett. 814, L20 (2015), arXiv:1511.02870 [astro-ph.HE] .
  * Wang  _et al._ (2016) J.-S. Wang, Y.-P. Yang, X.-F. Wu, Z.-G. Dai, and F.-Y. Wang, Astrophys. J. Lett. 822, L7 (2016), arXiv:1603.02014 [astro-ph.HE] .
  * Liu  _et al._ (2016) T. Liu, G. E. Romero, M.-L. Liu, and A. Li, Astrophys. J. 826, 82 (2016), arXiv:1602.06907 [astro-ph.HE] .
  * Zhang (2016) B. Zhang, Astrophys. J. Lett. 827, L31 (2016), arXiv:1602.04542 [astro-ph.HE] .
  * Yamasaki  _et al._ (2018) S. Yamasaki, T. Totani, and K. Kiuchi, Publ. Astron. Soc. Jap. 70, Publications of the Astronomical Society of Japan, Volume 70, Issue 3, 1 June 2018, 39, https://doi.org/10.1093/pasj/psy029 (2018), arXiv:1710.02302 [astro-ph.HE] .
  * Wei  _et al._ (2018) J.-J. Wei, X.-F. Wu, and H. Gao, Astrophys. J. Lett. 860, L7 (2018), arXiv:1805.12265 [astro-ph.CO] .
  * Cai  _et al._ (2019) R.-G. Cai, T.-B. Liu, S.-J. Wang, and W.-T. Xu, JCAP 09, 016 (2019), arXiv:1905.01803 [astro-ph.CO] .
  * Zhao  _et al._ (2020b) Z.-W. Zhao, Z.-X. Li, J.-Z. Qi, H. Gao, J.-F. Zhang, and X. Zhang, Astrophys. J. 903, 83 (2020b), arXiv:2006.01450 [astro-ph.CO] .
  * Caldwell  _et al._ (2019) R. Caldwell _et al._ , (2019), arXiv:1903.04657 [astro-ph.CO] .
  * (95) “LISA,” lisa.nasa.gov/.
  * Armano  _et al._ (2016) M. Armano _et al._ , Phys. Rev. Lett. 116, 231101 (2016).
  * Amaro-Seoane  _et al._ (2017) P. Amaro-Seoane _et al._ (LISA), (2017), arXiv:1702.00786 [astro-ph.IM] .
  * Armano  _et al._ (2018) M. Armano _et al._ , Phys. Rev. Lett. 120, 061101 (2018).
  * Abich  _et al._ (2019) K. Abich _et al._ , Phys. Rev. Lett. 123, 031101 (2019), arXiv:1907.00104 [astro-ph.IM] .
  * Speri  _et al._ (2021) L. Speri, N. Tamanini, R. R. Caldwell, J. R. Gair, and B. Wang, Phys. Rev. D 103, 083526 (2021), arXiv:2010.09049 [astro-ph.CO] .
  * Wu (2018) Y.-L. Wu, Int. J. Mod. Phys. A 33, 1844014 (2018), arXiv:1805.10119 [physics.gen-ph] .
  * Ruan  _et al._ (2020a) W.-H. Ruan, Z.-K. Guo, R.-G. Cai, and Y.-Z. Zhang, Int. J. Mod. Phys. A 35, 2050075 (2020a), arXiv:1807.09495 [gr-qc] .
  * Hu and Wu (2017) W.-R. Hu and Y.-L. Wu, Natl. Sci. Rev. 4, 685 (2017).
  * Luo  _et al._ (2020) J. Luo _et al._ , Class. Quant. Grav. 37, 185013 (2020), arXiv:2008.09534 [physics.ins-det] .
  * Wang  _et al._ (2020b) L.-F. Wang, Z.-W. Zhao, J.-F. Zhang, and X. Zhang, JCAP 11, 012 (2020b), arXiv:1907.01838 [astro-ph.CO] .
  * Liu  _et al._ (2020) S. Liu, Y.-M. Hu, J.-d. Zhang, and J. Mei, Phys. Rev. D 101, 103027 (2020), arXiv:2004.14242 [astro-ph.HE] .
  * Milyukov (2020) V. Milyukov, Astron. Rep. 64, 1067 (2020).
  * Mei  _et al._ (2020) J. Mei _et al._ (TianQin),  (2020), 10.1093/ptep/ptaa114, arXiv:2008.10332 [gr-qc] .
  * Fan  _et al._ (2020) H.-M. Fan, Y.-M. Hu, E. Barausse, A. Sesana, J.-d. Zhang, X. Zhang, T.-G. Zi, and J. Mei, Phys. Rev. D 102, 063016 (2020), arXiv:2005.08212 [astro-ph.HE] .
  * Tamanini  _et al._ (2016) N. Tamanini, C. Caprini, E. Barausse, A. Sesana, A. Klein, and A. Petiteau, JCAP 04, 002 (2016), arXiv:1601.07112 [astro-ph.CO] .
  * Belgacem  _et al._ (2019) E. Belgacem _et al._ (LISA Cosmology Working Group), JCAP 07, 024 (2019), arXiv:1906.01593 [astro-ph.CO] .
  * Zhao  _et al._ (2020c) Z.-W. Zhao, L.-F. Wang, J.-F. Zhang, and X. Zhang, Sci. Bull. 65, 1340 (2020c), arXiv:1912.11629 [astro-ph.CO] .
  * Abbott  _et al._ (2017f) B. Abbott _et al._ (LIGO Scientific, Virgo), Phys. Rev. Lett. 119, 141101 (2017f), arXiv:1709.09660 [gr-qc] .
  * Fan  _et al._ (2019) X. Fan, J. Li, X. Li, Y. Zhong, and J. Cao, Sci. China Phys. Mech. Astron. 62, 969512 (2019), arXiv:1811.01380 [astro-ph.IM] .
  * Ruan  _et al._ (2020b) W.-H. Ruan, C. Liu, Z.-K. Guo, Y.-L. Wu, and R.-G. Cai, Nature Astron. 4, 108 (2020b), arXiv:2002.03603 [gr-qc] .
  * Wang  _et al._ (2020c) G. Wang, W.-T. Ni, W.-B. Han, S.-C. Yang, and X.-Y. Zhong, Phys. Rev. D 102, 024089 (2020c), arXiv:2002.12628 [gr-qc] .
  * Hu  _et al._ (2021) Q. Hu, M. Li, R. Niu, and W. Zhao, Phys. Rev. D 103, 064057 (2021), arXiv:2006.05670 [gr-qc] .
  * Wang  _et al._ (2020d) R. Wang, W.-H. Ruan, Q. Yang, Z.-K. Guo, R.-G. Cai, and B. Hu, (2020d), arXiv:2010.14732 [astro-ph.CO] .
  * Omiya and Seto (2020) H. Omiya and N. Seto, Phys. Rev. D 102, 084053 (2020), arXiv:2010.00771 [gr-qc] .
  * Orlando  _et al._ (2021) G. Orlando, M. Pieroni, and A. Ricciardone, JCAP 03, 069 (2021), arXiv:2011.07059 [astro-ph.CO] .
  * Wang and Han (2021) G. Wang and W.-B. Han, Phys. Rev. D 103, 064021 (2021), arXiv:2101.01991 [gr-qc] .
  * Cutler (1998) C. Cutler, Phys. Rev. D 57, 7089 (1998), arXiv:gr-qc/9703068 .
  * Krolak  _et al._ (1995) A. Krolak, K. D. Kokkotas, and G. Schaefer, Phys. Rev. D 52, 2089 (1995), arXiv:gr-qc/9503013 .
  * Buonanno  _et al._ (2009) A. Buonanno, B. Iyer, E. Ochsner, Y. Pan, and B. Sathyaprakash, Phys. Rev. D 80, 084043 (2009), arXiv:0907.0700 [gr-qc] .
  * Rubbo  _et al._ (2004) L. J. Rubbo, N. J. Cornish, and O. Poujade, Phys. Rev. D 69, 082003 (2004), arXiv:gr-qc/0311069 .
  * Klein  _et al._ (2016) A. Klein _et al._ , Phys. Rev. D 93, 024003 (2016), arXiv:1511.05581 [gr-qc] .
  * Feng  _et al._ (2019) W.-F. Feng, H.-T. Wang, X.-C. Hu, Y.-M. Hu, and Y. Wang, Phys. Rev. D 99, 123002 (2019), arXiv:1901.02159 [astro-ph.IM] .
  * Zhao and Wen (2018) W. Zhao and L. Wen, Phys. Rev. D 97, 064031 (2018), arXiv:1710.05325 [astro-ph.CO] .
  * Kocsis  _et al._ (2006) B. Kocsis, Z. Frei, Z. Haiman, and K. Menou, Astrophys. J. 637, 27 (2006), arXiv:astro-ph/0505394 .
  * Ilbert  _et al._ (2013) O. Ilbert _et al._ , Astron. Astrophys. 556, A55 (2013), arXiv:1301.3157 [astro-ph.CO] .
  * Palenzuela  _et al._ (2010) C. Palenzuela, L. Lehner, and S. L. Liebling, Science 329, 927 (2010), arXiv:1005.1067 [astro-ph.HE] .
  * O’Shaughnessy  _et al._ (2011) R. O’Shaughnessy, D. L. Kaplan, A. Sesana, and A. Kamble, Astrophys. J. 743, 136 (2011), arXiv:1109.1050 [astro-ph.CO] .
  * Moesta  _et al._ (2012) P. Moesta, D. Alic, L. Rezzolla, O. Zanotti, and C. Palenzuela, Astrophys. J. Lett. 749, L32 (2012), arXiv:1109.1177 [gr-qc] .
  * Kaplan  _et al._ (2011) D. L. Kaplan, R. O’Shaughnessy, A. Sesana, and M. Volonteri, Astrophys. J. Lett. 734, L37 (2011), arXiv:1105.3653 [astro-ph.HE] .
  * Shi  _et al._ (2012) J.-M. Shi, J. H. Krolik, S. H. Lubow, and J. F. Hawley, Astrophys. J. 749, 118 (2012), arXiv:1110.4866 [astro-ph.HE] .
  * Blandford and Znajek (1977) R. D. Blandford and R. L. Znajek, Mon. Not. Roy. Astron. Soc. 179, 433 (1977).
  * Meier (2001) D. L. Meier, Astrophys. J. Lett. 548, L9 (2001), arXiv:astro-ph/0010231 .
  * Dotti  _et al._ (2012) M. Dotti, A. Sesana, and R. Decarli, Adv. Astron. 2012, 940568 (2012), arXiv:1111.0664 [astro-ph.CO] .
  * (139) “SKA,” www.skatelescope.org.
  * Ivezić  _et al._ (2019) v. Ivezić _et al._ (LSST), Astrophys. J. 873, 111 (2019), arXiv:0805.2366 [astro-ph] .
  * (141) “E-ELT,” www.eso.org/sci/facilities/eelt/.
  * Yang (2021) T. Yang, JCAP 05, 044 (2021), arXiv:2103.01923 [astro-ph.CO] .
  * Bruzual and Charlot (2003) G. Bruzual and S. Charlot, Mon. Not. Roy. Astron. Soc. 344, 1000 (2003), arXiv:astro-ph/0309134 .
  * Davies  _et al._ (2010) R. Davies _et al._ (MICADO Team), Proc. SPIE Int. Soc. Opt. Eng. 7735, 77352A (2010), arXiv:1005.5009 [astro-ph.IM] .
  * Dahlen  _et al._ (2013) T. Dahlen _et al._ , Astrophys. J. 775, 93 (2013), arXiv:1308.5353 [astro-ph.CO] .
  * Madau and Rees (2001) P. Madau and M. J. Rees, Astrophys. J. Lett. 551, L27 (2001), arXiv:astro-ph/0101223 .
  * Volonteri  _et al._ (2008) M. Volonteri, G. Lodato, and P. Natarajan, Mon. Not. Roy. Astron. Soc. 383, 1079 (2008), arXiv:0709.0529 [astro-ph] .
  * Chevallier and Polarski (2001) M. Chevallier and D. Polarski, Int. J. Mod. Phys. D 10, 213 (2001), arXiv:gr-qc/0009008 .
  * Linder (2003) E. V. Linder, Phys. Rev. Lett. 90, 091301 (2003), arXiv:astro-ph/0208512 .
  * Chen  _et al._ (2019) L. Chen, Q.-G. Huang, and K. Wang, JCAP 02, 028 (2019), arXiv:1808.05724 [astro-ph.CO] .