Jianwei Zhang(张见微) Chengmin Zhang(张承民) Di Li(李菂)Xianghan Cui(崔翔翰) Wuming Yang(杨伍明) Dehua Wang(王德华)Yiyan Yang(杨佚沿) Shaolan Bi(毕少兰) and Xianfei Zhang(张先飞)

1CAS Key Laboratory of FAST,National Astronomical Observatories,Chinese Academy of Sciences,Beijing 100101,China

2University of Chinese Academy of Sciences,Beijing 100049,China

3Department of Astronomy,Beijing Normal University,Beijing 100875,China

4School of Physics and Electronic Science,Guizhou Normal University,Guiyang 550001,China

5School of Physics and Electronic Science,Guizhou Education University,Guiyang 550018,China

6NAOC-UKZN Computational Astrophysics Centre,University of KwaZulu-Natal,Durban 4000,South Africa

Keywords: gravitational waves,statistical methods,neutron stars,black holes

1. Introduction

Since the first direct detection of gravitational waves(GWs) from a binary black hole (BBH) merger on 14 September 2015, i.e., GW150914,[1]the GW detector network, i.e., advanced LIGO and Virgo (LIGO-Virgo),[2,3]has opened up a new window of multi-messenger astronomy. Later, LIGO-Virgo hunted a new type of astrophysical source of GWs from the collision of a binary neutron star (BNS) on 17 August 2017, i.e., GW170817.[4]During the first and second observing runs (O1/O2) of LIGOVirgo, a total of 11 confident detections (i.e., 10 BBHs,and 1 BNS) of GWs from the compact binary mergers were realized (see Ref. [5] and references therein). With the new round GW hunting (O3), started on 1 April 2019(https://www.ligo.caltech.edu/news/ligo20190320), scientists are expecting to observe more GW merger events.

To date, LIGO-Virgo has not yet convincingly hunted the GWs from the coalescence of a NS-BH binary, however,they might be hopefully observed in the near future. The exhibition of NS-BH mergers is still an open scientific question, and several possible astrophysical channels have been proposed. Channels include isolated binary evolution (e.g.,Refs. [6–9]), mergers in star clusters (e.g., Refs. [10–13]), in triple or quadruple systems(e.g.,Refs.[14–18]),and in active galactic nuclei(AGN)accretion discs(e.g.,Refs.[19,20]),and so on.

Due to our lacking knowledge of the actual population and mass distribution of NS-BH merger events, Ref. [21]simulated their chirp mass (M) distribution with a synthetic model,M=(m1m2)3/5(m1+m2)−1/5withm1andm2being the component masses of the binary system, where the BHs and NSs were inferred by LIGO-Virgo (O1/O2),[5]and suggested that the underlying NS-BH candidates may distribute in the range of 2.1M⊙

However,to the best of our knowledge,none of the existing works simulated the distribution of GW frequency(fGW)of NS-BH mergers using the above synthetic model of Ref.[21].Therefore, to fill this gap, in this paper we examined thefGWdistribution of NS-BH mergers within the model by the Monte Carlo method, as a follow-up work of Ref. [21]. This article is organized as follows: in Section 2,the procedure of the model construction is described. Then,the results of our random sampling process are shown in Section 3, as well as a brief summary.

2. Model construction

Usually,the coalescence of a compact binary system,e.g.,NS-BH,is characterized by three stages: inspiral,merger,and ringdown(see Refs.[22,23]and references therein for a pedagogical review). During the inspiral phase,the two compact objects lose their orbital energy due to GW emission,while the orbital separation shrinks, and the GW signal amplitude and GW frequencyfGW[as well as the orbital frequency (forb)]increase, then the amplitude reaches a maximum when the system fully approaches the merger phase. The evolution of this stage is usually described by the post-Newtonian (PN)methods(e.g.,Refs.[24–27])approximately,which will break down once the system starts the merger, where the numerical relativity is required(e.g.,Refs.[28–32]).

For brevity,here we adopted a basic binary system model in a Boyer–Lindquist coordinate system to describe the phenomenon of the inspiral stage of the NS-BH coalescence, in which the spin of the rotating BH(i.e.,Kerr BH[33])is perpendicular to the orbital plane. We noted that the NS was treated as a structureless point particle in our approximation,implying that the effect of the tidal interaction between the two objects was not taken into account in our model. Therefore, thefGWof the NS-BH binary when just before merger can be approximately calculated as twice the orbital frequency,[22,23,34–37]

whereM=mNS+mBHis the total mass, andmNSandmBHare the masses of NS and BH,Gis the gravitational constant,andrISCOis the radius of the innermost stable circular orbit(ISCO),aandrg=GmBH/c2(cis the speed of light)are the dimensionless spin parameter and the gravitational radius(one half of the Schwarzschild radius)of the BH,respectively.

In order to simulate the distribution offGWthrough Eq.(1)using the Monte Carlo method,[38,39]one needs to prepare the corresponding samples of four related quantities, i.e.,mNS,mBH,aandrISCO.

On the selection of NS and BH massesmNSandmBH,we adopted the samples of the synthetic model (i.e., model LIGO)used by Ref.[21].That is to say,themNSwas randomly sampled from a Gaussian distribution with the mean value ofµ=1.4M⊙and a standard deviation ofσ=0.2M⊙(e.g.,Refs.[4,40–46]),while themBHwas randomly sampled from the masses of 20 BHs inferred by LIGO-Virgo(O1/O2),[5]and the corresponding spin parameterafor each BH was also sampled along withmBH,as shown in Table 1.

Table 1. The component masses mi,[5] and dimensionless spin magnitude ai[47] of 10 BBHs observed by the LIGO-Virgo GW detections(O1/O2),with 90%credible intervals.

ForrISCO,we can calculate it from the associated BH spin which is usually characterized by the dimensionless spin parameter, i.e.,a=cJ/(Gm2BH), whereJandcare the spin angular momentum and speed of light, respectively, as the following equations:[48]

in whichrISCO=6rgfor a non-rotating BH(i.e.,Schwarzschild BH[49,50]) and it increases (decreases) asaincreases for a counter-rotating(co-rotating)orbit,torISCO=9rg(rISCO=rg)for a maximally rotating Kerr BH with spin parametera=1,as shown in Fig.1.

Fig. 1. Radius of the ISCO of a Kerr BH as a function of the spin parameter a, the scatters with error bars that encompass 90% credible intervals correspond to 20 BHs observed by the LIGO-Virgo GW detections(O1/O2). For guides of the eyes,the upper and lower dashed lines refer to the counter- and co- rotating orbits, respectively, represented from Ref.[48].

As remarked, with the current detector sensitivity of LIGO-Virgo, it is difficult to measure the spins of the individual BHs (see Refs. [51–54]), and Refs. [5,55] introduced the dimensionless effective aligned spin parameter,[56–58]i.e.,χeff=(m1a1+m2a2)·ˆLN/M,where ˆLNis the Newtonian angular momentum of the orbital system that is normal to the orbital plane,to describe the dominant spin effects,which show that most BBHs are consistent withχeffbeing near zero,suggesting that the progenitor BHs are mostly non-spinning or their spins lie in the orbital planes of the binary systems.

Regardless of the uncertainties of the BH spins, here we employed the spin estimations of the individual BHs by Ref.[47],as listed in Table 1,to calculate the associatedrISCOthrough Eqs. (2)–(4), as shown in Fig. 1. It is obvious that the 90%credible regions of the spin estimates are comparable with their astrophysical meaningful ranges (i.e., 0≤a ≤1),as well as therISCO, indicating that they are somewhat poorly constrained based on the current observations. For completeness,the case of the non-rotating BH(Schwarzschild BH)withrISCO= 6rgwas also considered in our simulation. As discussed above,we have prepared the samples of four quantities,i.e.,mNS,mBH,aandrISCO,which are required to simulate GW frequencyfGWdistribution of the NS-BH merges via Eq.(1).

3. Results

As described in Section 2, we simulated the GW frequencyfGWdistribution for the NS-BH merges in three different spin cases: (i) non-rotating, (ii) counter-rotating, (iii)co-rotating. In each case, we generated 1 000 000 synthetic NS-BH systems using Monte Carlo random sampling,[38,39]then analyzed the correspondingfGWdistribution, the results are shown in Fig.2 and Table 2. As is expected,the values of the median(50%)and the 90%upper and lower limits onfGWfor non-rotating case are in-between those of the counter-and co-rotating cases,and the magnitudes of the medians for three cases are consistent with each other,i.e.,nearly a few hundreds of Hz. In addition, we note that the 90% upper limit for corotating case is much larger than those of the other two cases,roughly as 1065 Hz,because the lower value ofrISCOleads to the higherfGWaccording to Eq.(1),thus there are more events reaching the high frequency range,as shown in Fig.2,but our results cannot be applied to investigate the question of probability of extremes due to the approximations of our model,which may raise the over-fitting problem.

Overall, from the conservative point of views, we prefer using non-rotating case to predict the GW frequencyfGWdistribution of NS-BH events when they just before merger, in which the median is 165 Hz, and the 90% credible interval is 101≤fGW≤640 Hz, then the GW frequency is expected to sweep upwards from the initialfGWto a few times of it in the merger stage. This simulated GW frequency range lies in the frequency band of ground-based interferometric detectors (e.g., LIGO-Virgo[2,3]), i.e., about 15 Hz to a few kHz,therefore the GW signals of NS-BH mergers may be hopefully hunted by LIGO-Virgo in the ongoing observation run O3. Moreover, it is out of the frequency range of the spacebased laser interferometers(e.g.,LISA[59,60]or DECIGO[61]),which worked in the band of 0.1–100 mHz for the mergers of the binary systems with BHs of∼106M⊙,and it is also out of that of the pulsar timing array experiments,which is working in the band of 1–100 nHz for mergers of binaries with huge BHs of∼109M⊙.[62]

Fig.2. The GW frequency fGW distribution of our simulation for three different spin cases,and the vertical dashed lines represent the 90%upper and lower limits,respectively.

Table 2. The simulation results of GW frequency fGW for three different spin cases,and the reported values are the medians(50%)with 90%confidence intervals.

It was noted that although both BHs and NSs are compact objects, their space-time geometries are different, which will affect the characteristics of the frequency and the structure of the spectrum for the GWs;for instance,if we consider the internal structure of the NS in NS-BH mergers, the tidal disruption event(TDE)may occur before the NS arrives at the ISCO of the BH. Therefore, we plan to investigate the frequency,amplitude,and the structure of the spectrum for GWs,i.e., the three key parameters in the GW detection, by taking into account the detailed space-time geometries of BHs and NSs in our future research.

4. Conclusion

As a short summary,our simulation results would play an important role in predicting the GW frequency(fGW)distribution of the NS-BH mergers in the upcoming GW observations,and it is shown that the mean value with 90% credible interval is 16Hz for GW signals at the moment of systems just before merger,and it increases several times in the merger phase. The simulated GW frequency range lies in the GW frequency band of LIGO-Virgo,thus it provides an useful reference for the future GW detections of NS-BH mergers.