Constraint on cosmological parameters by Hubble parameter from gravitational wave standard sirens of neutron star binary system

TL;DR

This paper proposes a novel method to measure the Hubble parameter using gravitational wave standard sirens from neutron star binaries, demonstrating its potential for precise cosmological constraints with future detectors.

Contribution

The paper introduces a new approach utilizing GW anisotropy of luminosity distance to measure H(z), validated through simulations for upcoming GW detectors like DECIGO.

Findings

DECIGO provides accurate H(z) measurements.

Mock data significantly tightens cosmological parameter constraints.

Future GW observations could vastly improve cosmological parameter estimation.

Abstract

In this paper, we present a new method of measuring Hubble parameter($H(z)$), making use of the anisotropy of luminosity distance($d_{L}$), and the analysis of gravitational wave(GW) of neutron star(NS) binary system. The method has never been put into practice before due to the lack of the ability of detecting GW. LIGO's success in detecting GW of black hole(BH) binary system merger announced the possibility of this new method. We apply this method to several GW detecting projects, including Advanced LIGO(Adv-LIGO), Einstein Telescope(ET) and DECIGO, finding that the $H(z)$ by Adv-LIGO and ET is of bad accuracy, while the $H(z)$ by DECIGO shows a good accuracy. We use the error information of $H(z)$ by DECIGO to simulate $H(z)$ data at every 0.1 redshift span, and put the mock data into the forecasting of cosmological parameters. Compared with the available 38 observed $H(z)$ data(OHD), mock data shows an obviously tighter constraint on cosmological parameters, and a concomitantly higher value of Figure of Merit(FoM). For a 3-year-observation by standard sirens of DECIGO, the FoM value is as high as 834.9. If a 10-year-observation is launched, the FoM could reach 2783.1. For comparison, the FoM of 38 actual observed $H(z)$ data is 9.3. These improvement indicates that the new method has great potential in further cosmological constraints.

Figures

Equations

$$\{\begin{array}{rl}& d_L^{(1)}=\frac{(1+z{)}^2}{H(z)}|v_0|\\ & \frac{\Delta H(z)}{H(z)}=\sqrt{3}{\left[\frac{d_L^{(1)}}{d_L^{(0)}}\right]}^{-1}\left[\frac{\Delta d_L^{(0)}}{d_L^{(0)}}\right]\end{array}$$

$${\Gamma }_{ab}=4Re{\int }_{f_{min}}^{f_{max}}\frac{{\partial }_a{\tilde{h}}_i^{\ast }(f){\partial }_b{\tilde{h}}_i(f)}{P(f)}df$$

$$\begin{array}{rl}{P}_1(f)=& 6.53\times 10^{-49}[1+(\frac{f}{7.36Hz}{)}^2\\ & +4.45\times 10^{-51}\times (\frac{f}{1Hz}{)}^{-4}\times \frac{1}{1+(\frac{f}{7.36Hz}{)}^2}\\ & +4.94\times 10^{-52}\times (\frac{f}{1Hz}{)}^{-4}]H{z}^{-1}\end{array}$$

$$\begin{array}{rl}{P}_2(f)=& 10^{-50}[2.39\times 10^{-27}(\frac{f}{100Hz}{)}^{-15.64}+0.349\times (\frac{f}{100Hz}{)}^{-2.145}\\ & +1.76\times (\frac{f}{100Hz}{)}^{-0.12}+0.409\times (\frac{f}{100Hz}{)}^{1.1}{]}^{2}H{z}^{-1}\end{array}$$

$$\begin{array}{rl}{P}_3(f)=& 10^{-49}[(\frac{f}{215Hz}{)}^{-4.14}-5\times (\frac{f}{215Hz}{)}^{-2}\\ & +111\times (\frac{1-(\frac{f}{215Hz}{)}^2+(\frac{f}{215Hz}{)}^4/2}{1+(\frac{f}{215Hz}{)}^2/2})]H{z}^{-1}\end{array}$$

$$f_{min}=0.233(\frac{1M_{⨀}}{M_z}{)}^{5/8}(\frac{1yr}{T_{obs}}{)}^{3/8}Hz$$

$$\Delta N(z)=T_{obs}{\int }_{z-\frac{\Delta z}{2}}^{z+\frac{\Delta z}{2}}\frac{\dot{n}(z^{\prime })}{1+z^{\prime }}\times 4\pi [\frac{d_L(z^{\prime })}{1+z^{\prime }}{]}^2\times \frac{1}{H(z^{\prime })}d{z}^{\prime }$$

$$P({\Omega }_m,{\Omega }_{\Lambda }|\{H_i\})=\int P({\Omega }_m,{\Omega }_{\Lambda },H_0|\{H_i\})d{H}_0$$

$$=\int \ell (\{H_i\}|{\Omega }_m,{\Omega }_{\Lambda },H_0)P(H_0)d{H}_0$$

$$\ell (\{H_i\}|{\Omega }_m,{\Omega }_{\Lambda },H_0)=\left(\prod _i\frac{1}{\sqrt{2\pi {\sigma }_i^2}}\right)exp\left(-\frac{{\chi }^2}{2}\right)$$

$${\chi }^2=\sum _i\frac{[H_{0}E(z)-H_i{]}^2}{{\sigma }_i^2}$$

Full text

dipole of luminosity distance — gravitational wave — Neutron star binary system — Hubble parameter — cosmological parameters

1 Introduction

In the twenty-first century, we witnessed the bloom of accurate cosmology. Accurate cosmology even ranked second on a list named ”Insights of the decade” from Science magazine in 2010. The key of accurate cosmology is to accurately constrain cosmological parameters and their state equations, which can lead us to a better understanding of the evolution of our universe. We mainly developed four observations to constrain cosmological parameters: Supernova(SN), Baryon Acoustic Oscillation(BAO), Galaxy Cluster(CL), Weak Lens(WL) so far(Albrecht et al. 2006). Actually, a relatively new tool, Hubble parameter( $H(z)$ ), is becoming increasingly popular these years because of its effective constraint on cosmological parameters. $H(z)$’s high efficiency lies on the fact that it is the only observation that can directly represent the expanding history of our universe. Compared with the Luminosity distance($d_{L}$) of SN, $H(z)$ contains no integral terms and directly connects with cosmic parameters, which makes it a powerful tool in constraining cosmological parameters, because the integral term can conceal many details and hide important information. As Ma and and Zhang reported, $H(z)$ constrains cosmological parameter much tighter than the same-number SN does. To achieve the same constraint effect of $H(z)$, ones need four times as many SNs as $H(z)$(Ma & Zhang 2011). Getting a more accurate measurement of $H(z)$ could encourage the development of accurate cosmology a lot.

There are various ways to detect $H(z)$, which can be generally classified into three types: 1, differential age method(Stern et al. 2010); 2, radial BAO method(Gaztañaga et al. 2009a); 3, standard sirens method(Gaztañaga et al. 2009b). The first two techniques have been employed in the past detection of $H(z)$, but the data of observed hubble parameter data(OHD) are still insufficient. We get only 38 OHDs so far, whose accuracies are far from desirable. Now with the development of GW detecting technology, it is time to look forwards to the third method: GW standard sirens. In 2015, even the second generation GW detector Adv-LIGO operated not at its design sensitivity, it still detected the first GW signal at its first run(Abbott et al. 2016). According to theoretical understanding, formula of GW of binary system encodes the information of luminosity $d_{L}$, providing an access to the direct measurement of $d_{L}$. Several frequency windows of GW are targeted by different detectors. The second generation detector are mainly aimed at frequency window $10-1000Hz$, such as LIGO and VIRGO. The next generation detector plan to reach lower frequency region. The underground project DECIGO was designed most sensitive at $0.1-10Hz$, while the Einstein Telescope(ET) may also reach the frequency of the order $1Hz$. The space-based eLISA can even detect GW of $10^{-4}-10^{-1}Hz$. In this paper, we make use of GW sirens to measure $H(z)$ by estimating the error of $d_{L}$. Because NS binary system are used as the source of GW in this paper, the frequency of the GW signal of whom mostly concentrate on $10-1000Hz$, we ignore the projects whose optimal sensitivity are far away from $10-1000Hz$, such eLISA, and choose the ones whose optimal sensitivity locate around $10-1000Hz$.

In 2006, a new way to narrow the relatively error of $H(z)$ by studying the dipole of $d_{L}$ has been proposed(Bonvin et al. 2006b). But the problem is that the new method needs plenty of SNs if we want to get a relatively accurate $H(z)$, which can not be met in reality. And this method has problem in detecting high-$z$ $H(z)$. Still, it is an instructive idea. Atsushi Nishizawa and Atsushi Tamga et al. gave us an alternative by pointing it out that we can get the error information of $d_{L}$ through the gravitational wave function of NS binary system, instead of SN (Nishizawa et al. 2011). As the technology developing, now Adv-LIGO already detected gravitational wave, which was really a big strick. It not only proved the general relativity, but also implied that the new method to measure $H(z)$ was feasible. We sense the possibility and potential from detecting gravitational wave. It is meaningful and cutting-edge to study GW detecting projects. In this paper, we focus mainly on two aspect: 1, How will it work out if we apply the new method to some other projects? 2, with the new $H(z)$ information we get, to what degree could we constrain cosmological parameters?

This paper is organized as follows. In sec 2, we are going to sketch the idea of GW standard sirens method by Atsushi Nishizawa et al.(Nishizawa et al. 2011), and apply it to some GW detecting projects. In sec 3, we simulate the $H(z)$ data, and analyze the constraining ability of the mock data. In sec 4, we discuss the result. All through this paper, we adopt the natural unit, $c=G=1$.

2 Method

2.1 dipole of luminosity distance

If the universe is completely homogeneous and isotropic on large scale, and the observer is relatively rest with the cosmic microwave background(CMB), the luminosity distance, $d_{L}$, would be just the same form and expression as in standard cosmology. But in fact, there are perturbations around ideal condition leading into the appearance of correction term, multiple of $d_{L}$ (Sasaki 1987).Then $d_{L}$ can be written as follow:

$d_{L}=d_{L}^{(0)}+d_{L}^{(1)}+$higher order terms

$d_{L}^{(0)}$represents the traditional meaning of luminosity distance in unperturbed Friedmann universe , also the average of $d_{L}$ on all direction. $d_{L}^{(1)}$ means the dipole of $d_{L}$. As to the ”higher order terms”, it is self-explantory. The contribution to ”higher order terms” coming from the weak gravitational lens effect is so small when compared with dipole that we are going to ignore them here(Bonvin et al. 2006a). The dipole is dominated by the peculiar velocity of observers. If you want to check it or feel intrigued by the theory, you can look up reference[Bonvin et al. 2006b] for the details. Here is the final result:

[TABLE]

where $|v_{0}|$, z, $H(z)$ respectively denotes the projection of observer peculiar velocity on the direction of sight, the redshift of the observed celestial body, the expanding rate at the redshift z, and $\Delta d_{L}^{(0)}$,$\Delta d_{L}^{(1)}$ means the error of $d_{L}^{(0)}$,$d_{L}^{(1)}$. The mean error of $H(z)$ will reduce to $\Delta H(z)/\sqrt{N}$ if we observe N independent sources at the given redshift. Thus, we can improve the accuracy of $H(z)$ by the observation of a large number of sources. From the equations above, given the value of $d_{L}^{(1)}/d_{L}^{(0)}$ and $\Delta d_{L}^{(0)}/d_{L}^{(0)}$, $\Delta H(z)/H(z)$ can be easily calculated by multiplication. We already know the meaning and expression of $d_{L}^{(0)}$ and $d_{L}^{(1)}$. The result of the term $d_{L}^{(1)}/d_{L}^{(0)}$ is shown in Fig. 1. To get $\Delta H(z)/H(z)$, the only remaining problem is to find out $\Delta d_{L}^{(0)}/d_{L}^{(0)}$, which can be solved by analyzing observed GW function in following subsection.

2.2 GW standard siren

One can use SN to illustrate the method of reducing the error of $H(z)$. But due to the number and distribution of SN, it works not that well, especially at high-z region. Considering the advantage of the larger number of observed sources, which can dramatically narrow down the error of $H(z)$, we choose NS binary system as an alternative of SN. What’s more, the distribution and property of NS binary system is easier to estimate than that of black hole binary system. We are more familiar with NS. And for black hole binary system, there is a serious problem: black hole seldom radiates electromagnetic wave, rendering it impossible to measure its corresponding redshift up to now. This is an important factor to choose NS binary system

In GW experiments, one can extract the property of the source and cosmological information by comparing detected waveform with theoretical template. That is exactly what LIGO team did when the first detected the GW of two back holes merged(Abbott et al. 2016). The typical Fourier transform of GW waveform can be expressed by the following formula of frequency:

$\widetilde{h}(f)=\dfrac{A}{d_{L}(z)}M_{z}^{5/6}f^{-7/6}e^{i\Psi(f)}$

$A=(\sqrt{6}\pi^{2/3})^{-1}$ is a constant which is already geometrically averaged over the inclination angle of a binary system. $d_{L}(z)$ is the luminosity distance at redshift z, and we can set it as $d_{L}^{(0)}$ cause we are going to observe plenty of source at the given redshift. $M_{z}=(1+z)\eta^{3/5}M_{t}$ with the definition of total mass $M_{t}=m_{1}+m_{2}$ and symmetric mass ratio $\eta=m_{1}m_{2}/M_{t}^{2}$. The last unknown function $\Psi(f)$ is a little intricate. It is the frequency-dependent phase caused by orbital evolution. Usually we deal with it by post-Newtanion(PN) approximation. Here we are not going to give too much explain. Because this term will be eliminated when we do the following math, so its concrete expression will not affect the final result. Here we just need to know that it is a function of the coalescence time $t_{c}$, the phase $\phi_{c}$ when emitted, $M_{z}$, $f$, $\eta$.

There are five unknown parameters, namely: $M_{z}$, $\eta$, $t_{c}$, $\phi_{c}$, $d_{L}$. $d_{L}$ is the only parameter that has nothing do with the own property of binary system. For the convenience of calculating, we just take account of equal mass NS binary system with $1.4M_{\bigodot}$ ,and set $t_{c}=0,\phi_{c}=0$. Then $M_{z}=1.22(1+z)M_{\bigodot},\eta=0.25$. The observation of GW can tell us no information about reshift, which means we should still resort to electromagnetic observation to find out corresponding redshift. Cutler and Holz has demonstrated its technological viability(Cutler & Holz 2009).

The estimate of error of $d_{L}$ is based on Fisher matrix. And the Fisher matrix here is given by

[TABLE]

$\partial_{a}$ means derivative with respect to parameter $\theta_{a}$. For DECIGO, who has eight interferometric signals, $\Gamma_{ab}$ should multiplied by 8. We have set values to parameters expect for $d_{L}$. So the only parameter in $\Gamma_{ab}$ is $d_{L}$. P(f) is the noise power spectrum. And the $P(f)$ for DECIGO, ET and Adv-LIGO is shown in Fig. 2. Here we give the expression of each detector’s noise curve,$P_{1}(f),P_{2}(f),P_{3}(f)$ respective for DECIGO, ET and Adv-LIGO.

DECIGO: DECIGO is the acronym for Deci-Hertz Interferometer Gravitational Wave Observatory. DECIGO is a planed space-based GW observation aimed at $0.1~{}10Hz$ frequency region. Its configuration is still to be decided. Here we adopt the following parameters in its configuration: the arm length $1000km$, the output laser power $10W$ with wavelength $\lambda=532nm$, the mirror diameter $1m$ with its mass $100kg$, and the finesse of FP cavity 10. Its noise curve is(Kawamura et al. 2006)

[TABLE]

ET: ET is a third generation GW detector, whose design is not finished. Here we just consider the simplest case with 10 km arms. We adopt the expression given by Keppel and Ajith’s fitting used(Keppel & Ajith 2010)

[TABLE]

Adv-LIGO: Adv-LIGO is an available second generation detector whose optimal sensitivity band match with the frequency window of GW from NS binary system . The first run of Adv-LIGO did not reach its design sensitivity. Here we use the noise curve fitted by [Arun et al. 2005]. It is not an accurate expression, but an approximation of the original curve given by [Cutler & Thorne 2002].

[TABLE]

In the expression of $\Gamma_{ab}$, the lower cutoff of frequency, $f_{min}$, is a function of observation time $T_{obs}$.

[TABLE]

In the case of our paper, for a given $T_{obs}$, $f_{min}$ changes little with $M_{z}$, which is always in the high strain noise region. It makes no big difference to the result of the integral. A reasonable appointment of the value of $f_{min}$ will work. But for prudence, we just take the original expression of $f_{min}$ when calculate the integral. And the higher cutoff, $f_{max}$, can be decided by the property of the integrand. When the value of integrand goes comparatively small, its contribution to $\Gamma_{ab}$ can be neglected. Setting a upper limit of the integrating region is all right for the calculation. For the reason of their property of integrand, we set the $f_{max}$ of DECIGO, ET, Adv-LIGO respectively as $100Hz$, $10000Hz$ and $10000Hz$

The 1-sigma error of parameter is

$\dfrac{\Delta d_{L}^{(0)}}{d_{L}^{(0)}}=\Delta\theta_{a}=\sqrt{\{\Gamma^{-1}\}_{aa}}$

For different observation time, ie: 1 year, 3 years, 10 years, the 1-sigma error estimate of $d_{L}$, which arising from instrumental noise, is showed in Fig. 3. we use $\sigma_{instr}$ to denote it. For a given device, no matter it is DECIGO, ET or Adv-LIGO, the accuracy of $d_{L}$ is all the same even for different observation time. It makes no difference for the error no matter how long the observation continues. It is mainly because that the error is due to the property of device, having nothing to do with observation time. Then we can calculate the $\Delta H(z)/H(z)$ by a simple multiply. The value of $\Delta H(z)/H(z)$ by analyzing the GW function of a NS binary system at a given redshift is shown is Fig. 4. As we can see, the accuracy is far from desirable. We need to take measure to narrow down the error.

2.3 H(z) error

In last subsection, we already calculate the $H(z)$ relative error from a given NS binary system. The mean error will statistically abate if we have many independent sources. We can observe many NS binary system at the same redshift, which will lead us to a remarkably reduced error of $H(z)$. But to what degree can we reduce the error? First we need to figure out the number distribution of NS binary system, namely N, at different redshift. The distribution of NS binary system can be described and estimated. We are going to make use of it to estimate the accuracy of $H(z)$ in this subsection. According to Ref[Cutler & Harms 2006], the fitting of NS-NS merger rate presents us the following mathematic describe:

$\dot{n}(z)=\dot{n}_{0}s(z)$, $s(z)=\begin{cases}1+2z,&z\leq 1\\ 0.75(5-z),&1<z<5\\ 0,&z\geq 5\end{cases}$

where the s(z) is estimated from star formation history inferred from UV luminosity, and the $\dot{n}_{0}$ represents the merger rate at present time. Then $\Delta N$, the number of NS-NS merger at redshift bin $\Delta z$, is expressed by:

[TABLE]

Current study doesn’t provide solid evidence of the exact value of $s(z)$ and $\dot{n}_{0}$. But we can take a reasonable estimate . After all, we are aimed at evaluating the method, not launching an actual observation here. It is rational to set that $\dot{n}_{0}$ equals to the most recent estimate, $10^{-6}Mpc^{-3}yr^{-1}$, and $\Delta z$ equals to 0.1. Thus we get the estimation of 10-year observed number of NS binary system merger at different redshift, which is shown in Fig. 5.

Given the fact that the total number of SN is just of hundred-magnitude by now, the observed number of NS-NS merger event is much larger than that of SN, showing a tremendous potential in reducing the mean error of $H(z)$. And from above equation, the number of NS-NS merger at fixed redshift increases with $T_{obs}^{1/2}$. The elongation of observation time can remedy the drawback of the device sensibility. This is also an advantage of the replacement of SN by NS-NS merger event. Combining with the information we get in last subsection, we are able to calculate the $H(z)$ error now. To make it intuitive, we convert relative error to absolute error of $H(z)$ under $\Lambda CDM$. We choose $\Lambda CDM$ as fiducial model is out of the consideration that $\Lambda CDM$ is mostly consistent with cosmological observation. The relative error of $H(z)$ by DECIGO ,ET and Adv-LIGO is shown in Fig. 6, Fig. 7, and Fig. 8, each for 1-year, 3-year, 10-year observation. The error by Adv-LIGO is a totally disaster, which basically has little application value in constraining cosmological parameters. The error by ET is a little better, especially at low redshift region, because ET is more sensitive than Adv-LIGO. DECIGO plays quite well in this method. For 10-year observation, the relative error is less than $1\%$ at the range $z\leq 0.5$. When redshift reach 3, Due to the decreasing of the number of observed NS-NS merger event with redshift, the relative error of $H(z)$ was magnified, but still quite small. And the elongation of $T_{obs}^{1/2}$ shows a great ability in narrowing down the error.

3 evaluation

3.1 simulate data

The new method for detecting $H(z)$ has been proposed and a little error analysis has been done. The problem is how accurate could $H(z)$ observed by this way constrain cosmological parameters? Now it is just a method, which could not operate so far. We get no actual OHD by this way. But it doesn’t necessarily mean we can do nothing about it. A reasonable and rational simulation would help our forecasting and evaluating a lot.

Since the $H(z)$ data by Adv-LIGO have a bad accuracy, we are just going to carry on no simulation and forecast for $H(z)$ data by Adv-LIGO here. ET can do some simulation and forecast. The problem is that the effect is a little bad, even worse than 38 OHDs. We do not plan to show it here, too. In following sections, DECIGO is the only one we discuss. we are going to follow the way Yuan Shuo once took to generate mock data(Yuan & Zhang 2015),

$H_{sim}=H_{\Lambda CDM}+H_{drift}$

We treat $H_{sim}$ as a drift, $H_{drift}$, based on the theoretical $H(z)$ value under $\Lambda$CDM, $H_{\Lambda CDM}$, caused by various errors. $H_{drift}$ is a random value under gaussian distribution, $N(0,\Delta H)$. $\Delta H$ is the calculated by relative error we get in last section. Using a piece of python code, we generate our mock $H_{sim}$ data of 3-year observation at very 0.1 z, shown in Fig. 9. We have got 38 OHDs up to now. The data were obtained by different ways from different groups(Jimenez et al. 2003;Simon et al. 2005;Stern et al. 2010;Moresco et al. 2012; Moresco et al. 2016;Zhang et al. 2014;Moresco 2015;Gaztañaga et al. 2009a;Blake et al. 2012; Samushia et al. 2013;Xu et al. 2013;Meng et al. 2015). And their value has a distribution. Fig. 10 shows the $H_{\Lambda CDM}$ at every redshift and the 38 OHDs so far. As we can see, the OHD value goes up and down around the $H_{\Lambda CDM}$ at the same redshift, which justifies the validity of our simulation.

3.2 forecasting of mock data

Now that we have got the 3-year-observation mock data, we can use them to constrain cosmological parameters. Before that, we need a criteria to evaluate the constraining ability of the dataset-Figure of Merit(FoM). We can define FoM in different ways, as long as its value can reflect how tightly or loosely the data constrain parameters. Here for the convenience of our analysis, we adopt the definition in Ref(Albrecht et al. 2006), the reciprocal of the area enclosed by the contour, coinciding with a specially appointed confidence region under gaussian distribution.

We choose the $\Lambda CDM$ as our prior model. In a standard $\Lambda CDM$ universe with a curvature term $\Omega_{k}=1-\Omega_{m}-\Omega_{\Lambda}$,the Hubble parameter is given by

$H(z)=H_{0}E(z);E(z)=\sqrt{\Omega_{m}(1+z)^{3}+\Omega_{\Lambda}+\Omega_{k}(1+z)^{2}}$

The determination of $H_{0}$ has been carried on in different projects. Its value varies slightly from one project to another, but always around $70kms^{-1}Mpc^{-1}$. Freedman suggested $H_{0}=72\pm 8kms^{-1}Mpc^{-1}$. For 7-year WMAP observation, $H_{0}=73\pm 3kms^{-1}Mpc^{-1}$(Spergel et al. 2007). In this paper, we take the most recent value $H_{0}=74.2\pm 3.6kms^{-1}Mpc^{-1}$(Riess et al. 2009). And the best value of $\Omega_{m}$,$\Omega_{\Lambda}$ we adopt is 0.27, 0.73 respectively, due to the coherence that they are consistent with the observations and the fact that we use these value to generate our simulation data. All the three parameters are assumed under gaussian distribution. By Bayes’ theorem, the posterior probability density function of parameters given the data set {$H_{i}$} is :

[TABLE]

[TABLE]

where $\ell$ is the likelihood and $P(H_{0})$ is the prior probability density function of $H_{0}$. And the expression of $\ell$ is given by

[TABLE]

;

[TABLE]

and $\sigma_{i}$ is the uncertainty of the data $H_{i}$. Then the integral can be worked out for a given $P(H_{0})$. There would be a point in the parameter space maximizing the probability density, $P_{max}$. Because of what we have described in last paragraph, such a point in this forecasting would be {0.27, 0.73, 74.2}. The formula

$P=P_{max}exp\left(-\dfrac{\Delta\chi^{2}}{2}\right)$

means the contour of a given confidence region, which corresponds to the value of $\Delta\chi^{2}$. We have three parameters, $\Omega_{m},\Omega_{\Lambda},H_{0}$, here. $\Delta\chi^{2}$ is statistically set to 2.3, 6.17, 11.8 respectively for $1\sigma,2\sigma,3\sigma$ confidence region. For a direct comparing and understanding, here we choose $2\sigma$ confidence region, namely $\Delta\chi^{2}=6.17$, when calculate FoM.

To estimate the FoM, we take the Fisher Matrix forecast technique in [Dodelson 2003], and

$F_{ij}=\dfrac{1}{2}\dfrac{\partial^{2}\chi^{2}}{\partial\theta_{i}\partial\theta_{j}}$

the value of matrix elements is taken at the most-likely value of parameters. It is easy to calculate the Fisher matrix, a $3\times 3$ matrix in parameter space { $\Omega_{m},\Omega_{\Lambda},H_{0}$}. But we need to marginalize it to obtain the FoM in subspace {$\Omega_{m},\Omega_{\Lambda}$}. The marginalization is not intricate: inverse the matrix and remove the row and column to be marginalized, then inverse the reduced matrix again. Let’s denote the marginalized Fisher matrix by $\tilde{F}$. Then the contour in subspace can be given by

$(\Delta\theta)^{T}\tilde{F}\Delta\theta=\Delta\chi^{2};\Delta\theta=\theta-\theta_{best-value}$

$\Delta\theta$ is the deviation from the beat value of the parameters. The contour is shown in Fig. 11. As we can see, the contour is an ellipse, which is consistent with the equation of $\tilde{F}$. For a more direct and concrete comparison, we cite the 38 OHDs here. Their constraint on $\Omega_{m}$ and $\Omega_{\Lambda}$ is shown in Fig. 12. Apparently, the constraint of the mock data on parameters is much tighter, compared with that of available OHDs, which has an significant improvement on ”accurate cosmology”. The simulation and forecasting of 10-year-observation just carries out in the same way. We are going to skip the elaborate explanation here. As Fig. 13 shows, its constraint on cosmological parameters is even much tighter, implying a consequent higher FoM value.

When calculating FoM, we take $\Delta\chi^{2}$ as 6.17. The enclosed area is $\pi/\sqrt{det(\tilde{F}/\Delta\chi^{2})}$. So FoM, the reciprocal of the area, is

$FoM=\dfrac{\sqrt{det(\tilde{F}/\Delta\chi^{2})}}{\pi}$

By this way, we have the FoM value of mock data, which is about 834.9, while the FoM of 38 OHDs is just about 9.3. It is a remarkable improvement. For 10-year-observation mock data, the FoM has a farther improvement, reaching 2783.1. We have enough reason to look forward to the excellent application of $H(z)$ data by this method.

4 Discussion and Conclusions

In this paper, we mainly evaluate the quality of $H(z)$ data by GW standard siren method of several GW detection plans, whose optimal frequency locate around the frequency window of GW from typical BN system. We calculate the relative error of $H(z)$ on three devices, DECIGO and ET and Adv-LIGO. Though the sensitivity of the three devices is almost of the same order of magnitude, the $H(z)$ error of DECIGO is quite optimistic while that of other two is far from satisfying. It is because of the term in the integrand of $\Gamma$, $f^{-7/3}/P(f)$. The most sensitive frequency of Adv-LIGO and ET is of  kHz, while DECIGO sensitive on  10Hz, which is the exact reason leading to a difference on the integral. In low frequency region, the $P(f)$ of DECIGO is smaller, its contribution to integral is bigger. When it goes to high frequency region. Though $P(f)$ of ET and Adv-LIGO is small, $f^{-7/3}$ diminishes the value contributed by high frequency region to $\Gamma$. So the integral of ET and Adv-LIGO is smaller, leading a comparatively bigger error. But it does not mean that $H(z)$ data by this method is a dead end or of no meaning, which is justified by the forecasting of DECIGO-based $H(z)$ data. If the sensitivity of Adv-LIGO or ET is sightly improved, or just move the most sensitive frequency to a lower region, the error of $H(z)$ will be comparable that by DECIGO.

Considering the absence of real $H(z)$ data by DECIGO, we simulate $H(z)$ and the data show an alluring constraining ability on cosmological parameters. After all, we are aimed at evaluating the viability and quality of $H(z)$ data by GW standard siren method, not putting the method into actual operation. We find that the FoM of mock data shows a huge improvement when compared with that of 38 actual OHDs. For contrast, the FoM is 9.3 for 38 OHDs, 834.9 for 3-year-observation, 2783.1 for 10-year-observation. The tight constraint of mock data and the FoM of the contour indicates a bright further of detecting $H(z)$ data by this method. It will bring ”accurate cosmology” to next stage.

It is well worth to point out that at current stage, the range of Adv-LIGO’s detecting ability for NS binary system is just 70 or 80Mpc. This range is much smaller than what we assumed here. In the further, if we want to detect $H(z)$ by this way, a farther detecting range is necessary, implying a lower strain noise. A lower strain noise will eventually lead to a more accurate $H(z)$ data. So if one day we can put GW standard sirens into practice, the $H(z)$ data we get should be more accurate than we simulate in this paper, and the constraining on cosmological parameters will be tighter.

Acknowledgements

We thank to the many helpful suggestion from teachers and classmates. Especially Duan Xiaowei and Yu Hai help a lot. The discuss with them improve this work. This work was supported by the National Science Foundation of China (Grants No. 11573006, 11528306), the Ministry of Science and Technology National Basic Science program (project 973) under grant No. 2012CB821804. And Y. K. Tang acknowledges financial support from the National Natural Science Foundation of China (Grant Nos. 11303003 and 11673005) and from the Shandong Provincial Natural Science Foundation of China (Grant No. ZR2013AM002)