Measuring the delay time distribution of binary neutron stars. II. Using the redshift distribution from third-generation gravitational wave detectors network

TL;DR

This study assesses how current and future gravitational wave detectors can measure the delay time distribution of binary neutron star mergers by analyzing their redshift distribution, highlighting the potential for third-generation detectors to constrain DTD parameters accurately.

Contribution

It demonstrates that third-generation GW detectors can precisely determine the BNS delay time distribution and formation efficiency, overcoming degeneracies present in current detector measurements.

Findings

Current detectors cannot uniquely determine DTD due to degeneracies.

Third-generation detectors can constrain DTD and efficiency to better than 10%.

Redshift uncertainties at high z are manageable for DTD inference.

Abstract

We investigate the ability of current and third-generation gravitational wave (GW) detectors to determine the delay time distribution (DTD) of binary neutron stars (BNS) through a direct measurement of the BNS merger rate as a function of redshift. We assume that the DTD follows a power law distribution with a slope $\Gamma$ and a minimum merger time $t_{\rm min}$, and also allow the overall BNS formation efficiency per unit stellar mass to vary. By convolving the DTD and mass efficiency with the cosmic star formation history, and then with the GW detector capabilities, we explore two relevant regimes. First, for the current generation of GW detectors, which are only sensitive to the local universe, but can lead to precise redshift determinations via the identification of electromagnetic counterparts and host galaxies, we show that the DTD parameters are strongly degenerate with the unknown mass efficiency and therefore cannot be determined uniquely. Second, for third-generation detectors such as Einstein Telescope (ET) and Cosmic Explorer (CE), which will detect BNS mergers at cosmological distances, but with a redshift uncertainty inherent to GW-only detections ($\delta(z)/z\approx 0.1z$), we show that the DTD and mass efficiency can be well-constrained to better than 10\% with a year of observations. This long-term approach to determining the DTD through a direct mapping of the BNS merger redshift distribution will be supplemented by more near term studies of the DTD through the properties of BNS merger host galaxies at $z\approx 0$ (Safarzadeh & Berger 2019).

Figures

Tables

Table 1. Table 1: Definition of w 𝑤 w and fitting parameters of Equation  A6 for CE, ET and CE+ET.

Network	$w\left(z\right)$	$A$	$B$
CE	${\rho }_T/{\rho }_{\text{opt}}\left(z\right)$	1.63	3.93
ET	${\rho }_T/{\rho }_{\text{opt}}\left(z\right)$	1.05	3.33
CE-ET	$z/z_{\text{horizon}}$	-0.58	-5.05

Table 2. Table 2: Definition of w 𝑤 w and fitting parameters of Equation  A7 for CE, ET and CE+ET.

Network	$w\left(z\right)$	$a_1$	$a_2$	$a_3$	$a_4$	$a_5$	$a_6$	$a_7$	$a_8$	$a_9$
CE	${\rho }_T/{\rho }_{\text{opt}}\left(z\right)$	0.038	-0.57	10.76	-49.76	143.1	-293.1	470.5	-512.4	307.5
ET	${\rho }_T/{\rho }_{\text{opt}}\left(z\right)$	0.22	-8.51	132.8	-966.0	3995.7	-9923.5	15017.1	-13487.3	6581.3
CE-ET	$z/z_{\text{horizon}}$	-0.15	5.48	-71.4	455.0	-1613.3	3363.5	-4119.8	2794.4	-877.2

Equations

$$\dot{n}(z)=$$

$$\psi (z)=0.015\frac{(1+z{)}^{2.7}}{1+[(1+z)/2.9{]}^{5.6}}{\mathrm{M}}_{\odot }\mathrm{y}{\mathrm{r}}^{-1}\mathrm{Mp}{\mathrm{c}}^{-3}.$$

$$\rho (z)=8\Theta \frac{r_0}{D_L}{\left(\frac{{\mathcal{M}}_z}{1.2M_{\odot }}\right)}^{5/6}\sqrt{\zeta (f_{\max })},$$

$$r_0^2$$

$$x_{7/3}$$

$$\zeta (f_{\max })$$

$$\mathcal{M}={\left(\frac{m_1{m}_2}{{(m_1+m_2)}^2}\right)}^{3/5}(m_1+m_2).$$

$$f_{\max }=\frac{785\text{ Hz}}{1+z}\left(\frac{2.8M_{\odot }}{M}\right),$$

$$R_D(z)=\frac{d{V}_c}{dz}\frac{\dot{n}(z)}{1+z}{P}_{\text{det}}(z),$$

$$\tilde{h}(f)=F_{+}{\tilde{h}}_{+}(f)+F_{\times }{\tilde{h}}_{\times }(f),$$

$$F_{+}$$

$$F_{\times }$$

$$P_{\text{det}}\left({\theta }_{\text{int}},z\right)=\int H\left(\frac{{\rho }_T}{{\rho }_{\text{net}}({\theta }_{\text{int}},z,\theta ,\varphi ,\psi ,\iota )}\right)d\cos \theta d\varphi d\psi d\cos \iota ,$$

$${\rho }^2=64{\Theta }^2{\left(\frac{r_0}{d_L}\right)}^2{\left(\frac{{\mathcal{M}}_z}{1.2{\mathrm{M}}_{\odot }}\right)}^{5/3}{\zeta }^2(f_{\text{max}}),$$

$$P_{\text{det}}(w;A,B)=\frac{\text{erf}\left(A-Bw\right)-\text{erf}\left(A-B\right)}{\text{erf}\left(A\right)-\text{erf}\left(A-B\right)},$$

$$P_{\text{det}}(w;a_k)=\sum _{k=1}^9{a}_k(1-w{)}^k+\left(1-\sum _{k=1}^9{a}_k\right){\left(1-w\right)}^{10},$$

Full text

1 Introduction

The joint gravitational wave (GW) and electromagnetic (EM) detections of the binary neutron star (BNS) merger, GW170817 (Abbott et al., 2017b), marked the dawn of multi-messenger astronomy. As the current generation of GW detectors increases in sensitivity and number, the local rate of BNS mergers will soon be determined accurately for the first time, providing initial insight into the formation channels of these binaries. Currently, the BNS merger rate is weakly constrained by the single detection of GW170817 ($1540^{+3200}_{-1220}$ Gpc*-3* yr*-1*; Abbott et al . 2017), by the small known sample of Galactic BNS systems ($21^{+28}_{-14}$ Myr*-1*; Kim et al. 2015), and by the beaming-corrected rate of short gamma-ray bursts ($270^{+1580}_{-180}$ Gpc*-3* yr*-1*; Fong et al. 2015). These rates are in broad agreement, but the uncertainties from all methods span at least two orders of magnitude.

Still, even when the local BNS merger rate is well determined, the more fundamental distribution of merger delay times may not be. The delay time distribution (DTD) encodes the time span between the formation of the BNS system (or alternatively the time since the formation of the parent stars) until the two neutron stars merge through the emission of gravitational waves. The DTD therefore provides fundamental insight into the evolutionary processes that govern the initial separation of the binaries, including poorly-understood effects such as common envelope evolution.

The DTD is usually parametrized as a power law distribution above some minimum merger timescale, $t_{\rm min}$, based on the following arguments: after the BNS formation, the binary’s orbit decays due to the emission of gravitational waves on a timescale that depends on the initial semi-major axis ($a$) as $t\propto a^{4}$. Therefore the resulting distribution of the merger times depends on the distribution of initial semi-major axes, $dN/da\propto a^{-\beta}$. The initial semi-major axis distribution of the O/B stellar progenitors is assumed to follow a power law $dN/da\propto a^{-1}$. If the binary experiences a common envelope phase, then the distribution becomes steeper. Therefore, the expected merger times follow $dN/dt_{\rm merge}\propto t^{-\beta/4-3/4}$, where we define $\Gamma\equiv-\beta/4-3/4$ (Belczynski et al., 2018).

Insight on the form of the DTD has been gained from studies of the small population of Galactic BNS systems (Vigna-Gómez et al., 2018), from the properties of SGRB host galaxies (Zheng & Ramirez-Ruiz, 2007; O’Shaughnessy et al., 2008; Leibler & Berger, 2010; Fong et al., 2013; Behroozi et al., 2014; Berger, 2014), and from arguments related to $r$-process enrichment (Matteucci et al., 2014; Komiya et al., 2014; van de Voort et al., 2015; Shen et al., 2015; Côté et al., 2018; Hotokezaka et al., 2018; Safarzadeh et al., 2018, 2019). These results point to the need for a fast merging channel if BNSs are assumed to be the primary source of $r$-process enrichment in the universe, which suggest that $t_{\rm min}$ may be rather small, $\lesssim 0.1$ Gyr or the slope of the DTD is steep. Population synthesis models have also made various predictions for the values of $\Gamma$ and $t_{\rm min}$, but those are dependent on uncertain binary evolution processes (Dominik et al., 2012). We stress that the local BNS merger rate in itself cannot fully characterize the DTD since it also depends on an additional unknown parameter, the efficiency of BNS formation per unit stellar mass, $\lambda$.

In a recent paper, Safarzadeh & Berger (2019) (hereafter, Paper I) showed that the mass distribution of BNS merger host galaxies at $z\approx 0$ can provide insight on $\Gamma$ and $t_{\rm min}$ with a sample size of $\mathcal{O}(10^{2}-10^{3})$. This is based on the fact that, on average, galaxy star formation histories (SFH) depend on their mass, and hence the convolution of the DTD and SFH leads to a specific prediction about the mass function of BNS merger host galaxies. Such an observational approach to determining the DTD is only feasible in the local universe due to the required detection of EM counterparts that will in turn lead to the identification of the host galaxies. It is anticipated that Advanced LIGO/Virgo, joined by KAGRA111https://gwcenter.icrr.u-tokyo.ac.jp/en and IndIGO222http://www.gw-indigo.org/tiki-index.php, can produce the required sample size within the next two decades.

Here, we instead explore how the DTD can be determined by directly observing the redshift distribution of BNS mergers well beyond the local universe. Mapping the rate of BNS mergers as a function of redshift can break the degeneracy between the shape of the DTD ($\Gamma$ and $t_{\rm min}$) and the BNS mass efficiency ($\lambda$) when comparing to the cosmic star formation history. This approach requires an order of magnitude increase in GW detector sensitivity to detect BNS mergers at cosmological distances. Such an improvement is expected for third-generation ground-based observatories such as Einstein Telescope333http://www.et-gw.eu (ET; Punturo et al. 2010) and Cosmic Explorer444http://www.cosmicexplorer.org (CE; Abbott et al. 2017a). However, at these distances, it is unlikely that EM counterparts will be detected for the majority of events, and therefore the distance (redshift) information will rely directly on the GW signal itself. We explore how the inherent distance-inclination degeneracy affects the ability to determine the DTD.

The structure of the paper is as follows: In §2 we delineate the method of estimating the observed redshift distribution of BNS mergers as a function of DTD and mass efficiency, for different GW interferometer networks; in §3 we show the results of DTD determination for existing GW detectors, which are only sensitive to the local universe; in §4 we expand our analysis to a future network of ET and CE, including a determination of the expected redshift uncertainties, and show the resulting constraints on the DTD and mass efficiency. We discuss some caveats and summarize the key results in §5. We adopt the Planck 2015 cosmological parameters (Planck Collaboration et al., 2016) where $\Omega_{M}=0.308$, $\Omega_{\Lambda}=0.692$, $\Omega_{b}=0.048$ are total matter, vacuum, and baryonic densities, in units of the critical density, $\rho_{c}$, $H_{0}=67.8$ km s*-1* Mpc*-1* is the Hubble constant, and $\sigma_{8}=0.82$ is the variance of linear fluctuations on the 8 $h^{-1}$ Mpc scale.

2 method

The BNS merger rate as a function of redshift is a convolution of the DTD with the cosmic star formation rate density:

[TABLE]

where $dt/dz=-[(1+z)E(z)H_{0}]^{-1}$, and $E(z)=\sqrt{{\Omega}_{m,0}(1+z)^{3}+{\Omega}_{k,0}(1+z)^{2}+{\Omega}_{\Lambda}(z)}$. Here, $\lambda$ is the currently unknown BNS mass efficiency (assumed not to evolve555Although the DTD for binary black holes is likely highly dependent on the metallicity, the DTD for BNS systems has been argued to be at most weakly dependent on metallicity (Dominik et al., 2012). with redshift) used as a free parameter that we try to recover alongside the parameters governing the DTD; $t_{b}$ is the time corresponding to the redshift $z_{b}$; $dP_{m}/dt$ is the DTD, parametrized to follow a power law distribution ($\propto t^{\Gamma}$) with a minimum delay time, $t_{\rm min}$ that refers to the time since birth of the ZAMS stars and not when the BNS system formed. Therefore, $t_{\rm min}$ corresponds to the sum of the nuclear lifetime of the lowest mass component of the binary system and the minimal gravitational delay that is induced by the existence of a minimal separation between the two newly born neutron stars. We also impose a maximum delay time of 10 Gyr for our fiducial case, although this does not affect our results, although we note that more than half of the observed BNS systems in the MW half merger times more than 10 Gyr (Pol et al., 2019). We adopt the cosmic star formation rate density666We neglect the uncertainties in the cosmic SFRD since the GW source redshift uncertainty (see Appendix B) dominates the overall error budget at cosmological redshifts. from Madau & Dickinson (2014):

[TABLE]

To determine the observed BNS merger rate as a function of redshift we need to consider the matched filtering signal-to-noise ratio as a function of GW detector sensitivity (Finn, 1996):

[TABLE]

where $\mathcal{M}_{z}=(1+z)\mathcal{M}$ is the redshifted chirp mass, $D_{L}$ is the luminosity distance, $\Theta$ is the orientation function, and

[TABLE]

where $2f_{\rm{max}}$ is the wave frequency at which the inspiral detection template ends, $r_{0}$ denotes the characteristic distance sensitivity, and $S_{h}(f)$ is the detector’s noise power spectral density. The intrinsic chirp mass, $\mathcal{M}$, is given in terms of the component masses by:

[TABLE]

Here we assume that both neutron stars have mass of $m_{1}=m_{2}=1.4$ ${\rm\,M_{\odot}}$. The frequency at the end of the inspiral (taken to correspond to the innermost stable circular orbit) is:

[TABLE]

where $M$ is the total mass of the binary. In Figure 1 we show the sensitivity curves for Advanced LIGO, ET, and CE. The substantial reduction in noise amplitude for the third-generation detectors with respect to Advanced LIGO leads to an increase in the typical values of $r_{0}$ from $\approx 0.1$ to $\approx 1.5$ Gpc. Finally, the observed BNS merger rate as a function of redshift is given by:

[TABLE]

where $P_{\text{det}}(z)$ is defined in Appendix A, and the redshift derivative of the comoving volume is given by $dV_{c}/dz=(4\pi c/H_{0})[D_{L}^{2}/(1+z)^{2}E(z)]$.

In Figure 2 we show the intrinsic merger rate density, $\dot{n}(z)$, for nine different choices of the DTD, with $\Gamma=[-1.5,-1,-0.5]$ and $t_{\rm min}=[10,100,1000]$ Myr, and a fixed mass efficiency value of $\lambda=10^{-5}$ M${}_{\odot}^{-1}$. For comparison, we also show the curve corresponding to no delay (i.e., the cosmic star formation rate density). Clearly, DTDs that prefer longer delays result in a merger distribution that is skewed to lower redshifts, with a higher merger rate at $z\approx 0$, but with some degeneracy between $\Gamma$ and $t_{\rm min}$. However, since the value of $\lambda$ is not presently known, all of the DTDs can reproduce the same local rate by simply scaling $\lambda$ appropriately. This is essentially why a local measurement of the merger rate cannot by itself constrain the DTD.

To explore how well current and third-generation GW detectors can determine the DTD, we inject a specific DTD model ($\Gamma$, $t_{\rm min}$, $\lambda$), generate the resulting redshift distribution with associated uncertainties, and then fit this distribution using an interpolation table that is based on the nine input DTDs. We fit for the input parameters using Markov Chain Monte Carlo (MCMC) sampling with emcee, a python based affine invariant sampler (Foreman-Mackey et al., 2013). The likelihood function is $ln(L)=-\chi^{2}/2$, with $\chi^{2}=\L\sum_{i=0}^{i=N}(R_{D,i}-\hat{R_{D,i}})/\sigma_{t,i}^{2}$. Here the summation is over all of the redshift bins; $R_{D,i}$ and $\hat{R_{D,i}}$ are the constructed and simulated detection rates at redshift bin $i$, respectively; $\sigma_{t,i}$ is the total error on the detection rate at redshift bin $i$, which is a combination of the Poisson error and the error due to the distance-inclination degeneracy from GW data, $\sigma_{t}^{2}=\sigma_{p}^{2}+\sigma_{z}^{2}$; $\sigma_{p}=\sqrt{N}$, where $N$ is the expected number of events at a given redshift during the integrated observation time of length $T_{\rm obs}$; and the redshift uncertainty ($\sigma_{z}$) for each redshift bin is estimated based on the vertical distance from the mean expected detection rate to the upper envelope corresponding to when the detections’ redshift are all biased high. We model the redshift uncertainty as $\delta z/z=0.1z$ based on re-scaled simulations of binary black hole redshift uncertainty estimates as detailed in Appendix B. We adopt a flat prior distribution for all of our parameters in the $\log\lambda\in[-7,-3]$, $\log t_{\rm min}\in[1,3]$, and $\Gamma\in[-1.5,-0.5]$.

3 results for Current GW Detectors

The predicted observed redshift distribution for the current generation of detectors at design sensitivity is shown in Figure 3. As expected, because the detection distance is limited to only a few hundred Mpc, all of the DTDs predict the same shape of observed distribution, with a simple change in scaling that can be accommodated by varying the unknown value of $\lambda$.

For the purpose of assessing the resulting constraints on the DTD and $\lambda$ we assume that BNS mergers from the current GW network will have precisely determined redshifts through associated EM counterparts and host galaxies. Therefore, the error budget is dominated by the Poisson error based on the detection rate. For the input model we assume $\Gamma=-0.6$, $t_{\rm min}=700$ Myr, and $\lambda=10^{-5}$ M${}_{\odot}^{-1}$. Using our MCMC approach we show the resulting constraints on the DTD parameters for a year of Advanced LIGO/Virgo operations at design sensitivity in Figure 4. The results indicate that the DTD remains largely unconstrained, with the posterior distributions strongly influenced by the flat priors. In particular, $t_{\rm min}$ is unconstrained, while $\Gamma$ and $\lambda$ show a strong degeneracy, with median values that are biased away from the injected model. We find the same result for a decade of Advanced LIGO/Virgo operations.

Our results for Advanced LIGO/Virgo suggest that even the proposed upgrades to the current facilities, such as A+ (Miller et al., 2015) and Voyager (Lantz et al., 2018) will not have a significant impact on the DTD since these facilities will still only detect BNS mergers in the local universe (see e.g., Figure 1, right panel of Reitze et al. 2019). As argued in Paper I, a more robust constraint on the DTD from the current generation of GW detectors may be achieved through the mass distribution of BNS merger host galaxies. However, even this approach leaves a lingering degeneracy between $\Gamma$ and $t_{\rm min}$.

4 results for Third-Generation Detectors

The situation is drastically different for the third-generation detectors, ET and CE. In Figure 5 we plot the expected detection rate as a function of redshift for ET, CE, and a network of ET+CE. Two improvements are readily apparent. First, the expected detection rate is about three orders of magnitude larger than for Advanced LIGO/Virgo. Second, the redshift range for BNS merger detections increases to $z\sim 5$ in the case of ET+CE. The latter improvement results in a clear difference between the redshift distributions of the various DTDs, while the former improvement provides the detection statistics needed to distinguish between the DTD models. The differences between the various DTDs can no longer be scaled away with a change in $\lambda$ (as is the case for Advanced LIGO/Virgo). In what follows we focus on the case of ET+CE as a realistic version of a third-generation detector network.

In Figure 6 we show the result of MCMC fitting for the same input model used in the previous section. We show the results for 1 day, 1 week, 1 month, and 1 year of observations. Unlike in the case of the current generation of GW detectors, we assume that the BNS mergers at cosmological distances will generally not have detectable EM counterparts (see Appendix C.). Instead we rely on distance information from the GW signal itself. We model the resulting redshift uncertainty as $\delta z/z=0.1z$; a detailed motivation for this parametrization is provided in Appendix B.

Our results show that an ET+CE network is able to constrain the DTD parameters and overcome the intrinsic degeneracy between $\Gamma$ and $t_{\rm min}$ within a year of observations. The values of $\Gamma$ and $t_{\rm min}$ can be determined to better than $10\%$ accuracy. We note that these numbers depend on the overall event rate, which is determined by $\lambda$; here we use an injected value of $10^{-5}{\rm\,M_{\odot}}^{-1}$, but the results can be rescaled for higher or lower values.

To assess the impact of our input model on the results, in Figure 7 we repeat the same exercise, but for two different DTDs that favor short merger timescales: $\Gamma=-1.2$ with $t_{\rm min}=30$ Myr, and $\Gamma=-1$ with $t_{\rm min}=100$ Myr. Although the power law index is recovered with the same accuracy as before, we find that $t_{\rm min}$ becomes more challenging to determine when its value is small. This is because the relative shift in the observed BNS merger redshift distribution becomes progressively smaller for small values of $t_{\rm min}$, which is challenging to detect in the presence of realistic GW redshift uncertainties.

5 Summary and Discussion

We investigated how well the DTD and mass efficiency of BNS systems can be determined through the redshift distribution of BNS mergers detected by current and future GW networks. We model the DTD as a power law with a minimum merger timescale, and leave the mass efficiency as a free parameter. While other DTDs have been proposed (e.g., Simonetti et al. 2019), our primary conclusions should not be affected by the exact form of the DTD.

We find that current GW detectors, which can only detect BNS mergers in the local universe, cannot directly constrain the DTD due to their limited sensitivity. In effect, the various DTDs, with an appropriate scaling of $\lambda$, predict the same BNS merger detection rate at $z\approx 0$. However, the situation is dramatically different for the anticipated third-generation detectors, which will be able to detect BNS mergers to $z\approx 5$. For this cosmological population, even in the presence of redshift uncertainties of $\delta z/z\approx 0.1z$ from the GW data, the large detection rate and broad redshift range will precisely determine the DTD within about a year of operations. It has been previously argued that the cosmological merger population uncovered by third-generation detectors will be able to constrain cosmological parameters (Sathyaprakash et al., 2009; Taylor et al., 2011; Taylor & Gair, 2012; Vitale & Farr, 2018); our results for the DTD further bolster the science case for third-generation GW detectors.

This work was supported by the National Science Foundation under grant AST14-07835 and by NASA under theory grant NNX15AK82G. The Berger Time-Domain Group at Harvard is supported in part by NSF under grant AST-1714498 and by NASA under grant NNX15AE50G. MTS is thankful to Harvard-Smithsonian Center for Astrophysics for hospitality which made this work possible. SV, KKYN and CW acknowledge support of the National Science Foundation, the LIGO Laboratory and the LIGO Data Grid clusters. LIGO was constructed by the California Institute of Technology and Massachusetts Institute of Technology with funding from the National Science Foundation and operates under cooperative agreement PHY-0757058.

Chapter \thechapter

Appendix A Detection probability of a network of detectors

The strain measured by a GW interferometer in frequency domain is given by

[TABLE]

where $\tilde{h}_{+,\times}$ are the $+,\times$-polarization bases and $F_{+,\times}$ are the corresponding beam pattern functions,

[TABLE]

where $\theta,\phi$ and $\psi$ are the zenith, azimuth and polarization angles respectively, and $g$ is a dimensionless coefficient determined by the geometry of an interferometer (Sathyaprakash & Schutz, 2009; Chen et al., 2017; Schutz, 2011). CE is a single interferometer with the angle between two arms equal to $90\deg$, hence $g_{\text{CE}}=1$ (Chen et al., 2017; Schutz, 2011). ET consists of 3 identical interferometers with the angle between two arms equal to $60\deg$, forming an equilateral triangle, hence $g_{\text{ET}}=\sqrt{3}/2$ for each interferometer in ET (Sathyaprakash et al., 2012; Punturo et al., 2010).

The detection probability, $P_{\text{det}}$, is defined as the probability of a detection with $\rho_{\text{net}}\geq\rho_{T}$, where $\rho_{T}$ is the SNR threshold of detection and $\rho^{2}_{\text{net}}=\sum_{i}{\rho^{2}_{i}}$ is the network SNR as a geometric sum of SNR of each inteferometer. Assuming isotropic sky locations, orbital orientation and polarization, (i.e., uniform distribution of $(\cos\theta,\phi,\psi,\cos\iota)$), $P_{\text{det}}$ is given analytically as

[TABLE]

where $\theta_{\text{int}}$ is the set of intrinsic parameters, which are fixed at $1.4-1.4{\rm\,M_{\odot}}$ and zero-spin for BNS systems, and $H(w(z,\theta,\phi,\psi,\iota))$ is the unitary step function defined in $w\in(0,1]$ for $w(z,\theta,\phi,\psi,\iota)=\rho_{T}/\rho_{\text{net}}(z,\theta,\phi,\psi,\iota)$.

For single CE or ET, we follow the inspiral approximation in Finn (1996). The SNR of a single interferometer is approximately

[TABLE]

where $\Theta^{2}=4\left[F_{+}^{2}\left(1+\cos^{2}\iota\right)^{2}+4F_{\times}^{2}\cos^{2}\iota\right]$, and proportional to $\Theta^{2}_{\text{net}}=\sum_{i}{\Theta^{2}_{i}}$ for the same signal strain observed by homogeneous detectors such as single CE or ET. Hence $P_{\text{det}}(w(z))$ is equivalently the survival function of $\Theta_{\text{net}}/\Theta^{\text{max}}_{\text{net}}$, where $\Theta^{\text{max}}_{\text{net}}$ is the maximum response of a particular network of detectors (Sathyaprakash & Schutz, 2009; Chen et al., 2017; Schutz, 2011; Finn, 1996). Here $w(z)=\rho_{T}/\rho_{\text{opt}}(z)$ where $\rho_{\text{opt}}$ is the optimal SNR at maximum $\Theta^{2}$. CE has $\Theta^{\text{max}}_{\text{CE}}=4$ and ET has $\Theta^{\text{max}}_{\text{ET}}=4\times\sqrt{3\times 3/4}=6$. We have verified that the inspiral approximation in 3G detectors only results in few-percent difference in SNR for $z\lesssim 6$, which is the region we are interested in.

Approximated forms of $P_{\text{det}}$ assuming uniformly distributed $(\cos\theta,\phi,\psi,\cos\iota)$ provided in Finn (1996) or Dominik et al. (2015) are only suitable for a single CE or a second-generation network. Therefore, we generate the survival function of $\Theta$ by drawing $10^{6}$ points of uniformly distributed $(\cos\theta,\phi,\psi,\cos\iota)$ in single CE or ET and fit the survival function with the following parametric form:

[TABLE]

where $\text{erf}(x)=\frac{2}{\sqrt{\pi}}\int_{0}^{x}e^{-t^{2}}dt$ is the error function. We also employ a 10th order polynomial,

[TABLE]

where $\{a_{k}\}$ are the polynomial coefficients.

The above simplification using distribution of $\Theta$ breaks down for a CE+ET network due to the different sensitivity and heterogeneous geometry of each interferometer. Instead, we calculate the integral in Equation A4 by simulating waveform and network SNR for each redshift. Then $w(z)=z/z_{\text{horizon}}$ and $z_{\text{horizon}}$ is the redshift of the detector horizon, which is $\sim 12.5$ in CE+ET for a BNS merger. Again we fit the simulation result for CE+ET using Equations A6 and A7.

Tables 1 and 2 show the definitions of $w$ and fitting parameters of Equations A6 and A7 in each network. Figure 8 shows the comparison of actual simulation and the parametric fits.

Appendix B redshift uncertainty from GW detections

In the case of third-generation detectors it is unlikely that most BNS merger detections will have EM counterparts. Instead, the distance information will need to be gleaned from the GW signal itself. There is an inherent degeneracy between a binary’s inclination angle in the sky with respect to a GW detector, $\theta_{\rm JN}$, and the luminosity distance (Schutz, 2011; Abbott et al . et al., 2016; Chen et al., 2018; Usman et al., 2018). Messenger & Read (2012) show that for a range of representative neutron star equations of state the redshift of such systems can be determined to an accuracy of $\sim 8-40\%$ for $z<1$ and $\sim 9-65\%$ for $1<z<4$. However, for a binary to be detectable at high redshifts, we expect the inclination to be close to face-on ($\theta_{\rm JN}=0^{o}$), or face-off ($\theta_{\rm JN}=180^{o}$), as most of the energy in gravitational waves is released along the angular momentum vector of the binary. Schutz (2011) derived an analytic formulation for the distribution of inclination angle of sources detectable by advanced detectors; $\lesssim 7\%$ ($\lesssim 3\%$) of detectable events will have viewing angles of $>70^{o}$ ($>80^{o}$) (Chen et al., 2018).

In principle, we can simulate BNS signals at different redshifts in third-generation detectors and estimate the evolution of uncertainty in parameter estimation using MCMC. Since BNS are a low-mass system with a long coalescence time, the parameter estimation is computationally expensive. To approximate the redshift uncertainty of a BNS at fixed redshift, we may extrapolate a binary black hole (BBH) signal by lowering the starting frequency from 10 Hz to 5 Hz, as to mimic the long inspiral phase in a BNS merger. We obtain the redshift uncertainties of BBHs in CE+ET from the simulations of Vitale & Whittle (2018), and fit the mean redshift uncertainties as a function of true redshifts. Then we calibrate our fit using the above extrapolation scheme to approximate the redshift uncertainties as $\delta z/z\approx 0.1z$.

Alternatively, Chen et al. (2018) developed a rapid algorithm that provides a luminosity distance uncertainty estimate for a large population of BNS merger detections. We use this algorithm to simulate 2000 BNS detections in a third-generation network, and compare the results to the extrapolation procedure above. We find that the two approaches yield consistent results. We therefore use the extrapolated distance uncertainty, and convert it to the redshift uncertainty through the adopted cosmology in this work.

Appendix C redshift uncertainty from joint SGRB detections

The distance-inclination degeneracy can be broken through the detections of an associated SGRB. The relative fraction of on-axis mergers is only a few percent. A $\gamma$-ray detection alone will thereby reduce the overall redshift uncertainty of at most a few percent of BNS merger detections, and likely over a restricted redshift range (perhaps to $z\sim 2$). An afterglow detection can further lead to a precise redshift determination through an associated host galaxy, but to date such detections have mainly been limited to $z\lesssim 1$ (Berger, 2014), which is generally not high enough to make an impact on the DTD determination (see Figure 5). Thus, it seems unlikely that associated SGRBs and their afterglows will substantially improve the constraints on the DTD.