Effects of the merger history on the merger rate of primordial black hole binaries

We develop a formalism to calculate the merger rate of primordial black hole binaries with a general mass function, by taking into account the merger history of primordial black holes. We apply the formalism to two specific mass functions, monochromatic and power-law cases. In the former case, the merger rate is dominated by the single-merger events, while in the latter case, the contribution of the multiple-merger events on the merger rate can not be ignored. The effects of the merger history on the merger rate depend on the mass function.


I. INTRODUCTION
Various astrophysical and cosmological observations provide substantial evidences firmly establishing the existence of dark matter (DM) in our Universe. However, the nature of DM remains one of the major unsolved problems in fundamental physics. Primordial black holes (PBHs) produced in the early Universe due to the collapse of large energy density fluctuations, as a promising candidate for DM, have recently attracted much attention [1][2][3][4][5][6][7][8][9][10][11].
Two neighboring PBHs can form a binary in the early Universe and coalesce within the age of the Universe. The merge rate of PBH binaries was first estimated through the three-body interaction for the case where all PBHs have the same mass [12,13]. In the PBH binary formation scenario, the gravitational wave event GW150914 detected by LIGO [14] and the merger rate estimated by the LIGO-Virgo Collaboration can be explained by the coalescence of PBH binaries if PBHs have the mass about 30M ⊙ and constitute a tiny fraction of DM [15]. The binary formation was extended to account for an arbitrary PBH mass function based on the three-body approximation [16] or to account for the torque from the surrounding PBHs as well as standard large-scale adiabatic perturbations assuming a monochromatic mass function [17]. The mechanism has recently been developed for a general mass function by taking into account the torque from the surrounding PBHs [18][19][20][21].
However, these studies ignore the possibility that a PBH binary merges into a new black hole which together with another PBH form a new PBH binary. Such a second-merge event can in principle be detected by LIGO-Virgo at the present time. In this paper, we develop an analytic formalism to work out the merger rate of PBH binaries with a general mass function, by taking into account the merger history of PBHs. * liulang@itp.ac.cn † guozk@itp.ac.cn ‡ cairg@itp.ac.cn The paper is organized as follows. In the next section, we summarize the basic equation for the primordial input parameters of PBHs and revisit the merger rate for a monochromatic mass function as the first-merger process. In Sec. III, we develop a formalism to calculate the merger rate of PBH binaries with a general mass function, by taking into account the merger history of PBHs. In Sec. IV, we consider two specific examples, monochromatic and power-law mass functions, to investigate the effects of the merger history on the merger rate of PBH binaries. The final section is devoted to conclusions.
In this paper, we use units of c = G = 1. Whenever relevant, we adopt the values of cosmological parameters from the Planck measurements [22]. The scale factor is normalized to unity at the present time.

II. SINGLE-MERGER EVENTS
Let us start with reviewing the basic equation of the merger rate of PBH binaries. It could be easily checked that the gravitational attraction between two approximately isolated PBHs dominates their dynamics if their average mass is bigger than the background mass contained in a comoving sphere whose radius equals to their conformal distance. Considering the different scaling with time of the two competing effects (their gravitational attraction versus the expansion of the Universe) in the equation of motion for their separation [17]. Following Ref. [15], in this section, we assume that all PBHs have the same mass, M , and PBH binaries decouple from the expansion of the Universe during radiation domination provided that their comoving separation, x, approximately satisfies where f pbh is the fraction of PBHs in DM, n pbh denotes the comoving average number density of PBHs and ρ dm denotes the present energy density of DM. The redshift z dec at which the binary decoupling occurs is given by where z eq ≃ 3400 is the redshift at matter-radiation equality, assuming negligible initial peculiar velocities here and throughout. Therefore, given PBH mass M and the initial comoving distance of PBHs x, the decoupling time is determined by PBHs. In this work, we assume that accretion and evaporation are negligible before the epoch of binary formation. When two PBHs come closer and closer, the surrounding PBHs, especially the nearest PBH, will exert torques on the PBH binary.
As the result, two PBHs avoid a head-on collision with each other. The tidal force will provide an angular momentum to prevent this system from direct coalescence. The major and minor axes are given by (denoted by a and b, respectively) where y is the comoving distance to the third PBH, A and B are numerical factors of O(1). A detailed investigation of the dynamics of the binary formation suggests A = 0.4 and B = 0.8 [13]. To be exact, in the following calculation, we adopt A = 0.4 and B = 0.8. The dimensionless angular momentum of PBH binaries is given by where e is the eccentricity of the binary at the formation time. Once two PBHs form a binary, they gradually shrink through the emission of gravitational radiation and eventually merge at the time τ after its formation, which can be estimated as [23] τ ≃ 3a 4 j 7 170M 3 .
To calculate the merger rate of PBH binaries, we have to know the spatial distribution of PBHs. Assuming that the spatial distribution of PBHs is random one, the probability that the comoving distances, x and y, are in the intervals (x, x + dx) and (y, y + dy) is given by To deal with this probability distribution, we can rewrite Eq. (7) as follows where y max = (4πn pbh /3) −1/3 , which is adopted in [12]. In Fig. 1 we show the merger rate estimated by using the initial distribution (7) and the simplified distribution (8), which indicates that the difference between the two cases is insignificant compared to the uncertainty of the merger rate estimated by the LIGO-Virgo Collaboration.
The fraction of PBHs which have merged before the time t is given by The fraction of PBHs which merge in the interval (t, t + dt) is G (t + dt) − G (t). Therefore, the merger rate of PBH binaries per unit volume per unit time (at the time t) can be easily obtained by where the factor 1/2 accounts for that each merger event involves two PBHs. From Eq. (10), the final result is given by which can be interpreted as the merger rate in Gpc −3 yr −1 and f c is given by We show the single-merger rate of PBH binaries as a function of the PBH abundance in Fig The merger rate R = 12 − 213 Gpc −3 yr −1 estimated by the LIGO-Virgo Collaboration is shown as the shaded region colored orange [24].

III. MULTIPLE-MERGER EVENTS
So far, several gravitational wave events from black hole binary mergers have been detected by the LIGO-Virgo collaboration, such as GW150914 (36 +5 [27]. These events detected by LIGO-Virgo suggest that the black holes should have an extended mass function. In this section, we calculate the merger rate distribution for PBH binaries with a general mass function by taking into account the effect of merger history on the merger rate of PBH binaries.
First of all, we consider the condition that two neighboring PBHs with the masses m i and m j decouple from the expansion of the Universe and form a bound system. Their comoving separation, x, approximately satisfies where m b = m i +m j is the total mass of the PBH binary. When two PBHs come closer, the nearest PBH with the mass m l , exert torque on the bound system. As a result, the two PBHs avoid a head-on collision and form a highly eccentric binary. The major axis a of the binary orbit and the dimensionless angular momentum are given by where y is the comoving distance to the third PBH with the mass m l . Once the PBHs decouple from the expansion of the Universe and form a binary, they gradually shrink by gravitational radiation and finally merge. The coalescence time is given by [23], The two neighboring PBHs with the masses m i and m j merge into a bigger black hole. The mass is given by where E GW is the energy of gravitational wave and γ is a factor of O(1). In the monochromatic case, γ = 0.95 is adopted in [28]. For simplicity, in this paper, we take γ = 1, which means we assume that the energy of gravitational wave is zero. The probability distribution function of PBHs P (m) is normalized to be dmP (m) = 1.
The abundance of PBHs in the mass interval (m, m+dm) is given by where f is a fraction of PBHs in non-relativistic matter including DM and baryons. The fraction of PBHs in DM is related to f by f pbh ≡ Ω pbh /Ω dm ≈ f /0.85. At the present time, the average number density of PBHs in the mass interval (m, m + dm) is given by where ρ m is the total energy density of matter and the present total average number density of PBHs, n T , is given by For simplicity, we define m pbh as We define F (m) as which is the fraction of the average number density of PBHs with the mass m in the total average number density of PBHs. The result in [17] indicates that in the case of f pbh < f c , the effects of the linear density perturbations on the merger rate of PBH binaries is significant. Here, we only consider the the case of f pbh > f c which is shown to be relevant to the LIGO observations [15]. In other words, we ignore the bound (13).
The only essential ingredient that we need is the spatial distribution of PBHs. We firstly consider the spatial distribution of two PBHs. The probability distribution of the comving separation x between two nearest PBHs with the masses m i and m j and without other PBHs in the comving volume of 4πx 3 /3 is given by Clearly, in the non-monochromatic case, to calculate the merger rate in the first-merger process, the differential probability distribution is given by where x is the comoving separation between two nearest PBHs with the masses m i and m j and y is the comoving distance to the third PBH with the mass m l which provides the angular momentum for the bound system. Similarly, to deal with this probability distribution, we can rewrite Eq. (25) as where y max = 4π 3 n T −1/3 . The fraction of PBHs that have merged before the time t is given by R 1 (t, m i , m j , m l ) is given by where the factor 1/2 accounts for that each merger event involves two PBHs. From Eq. (28), one has The merger rate density of PBH binaries with the masses m i and m j in the first-merger process is Let us estimate the merger rate in the second-merger process. In the first-merger process, two neighboring PBHs decouple from the expansion of the Universe and then merge into a new black hole with the mass m i + m j . In the second-merger process, the new black hole and the nearest PBH with mass m k form a new binary. The merge event of the new binary is detected by LIGO-Virgo at the time t. Statistically, the second coalescence time is larger than the first one, therefore, we can ignore the first coalescence time. The differential probability distribution is given by Similarly, to deal with this probability distribution, we can rewrite Eq. (31) as So, the fraction of PBHs that have merged in the secondmerger process is given by τ (y, z)).
(33) R 2 (t, m i , m j , m l ) is given by where the factor 1/3 accounts for that each merger event in second-merger process involves three PBHs. From Eq. (34), the final result is given by The merger rate density of PBH binaries with the masses m i and m j in the second-merger process is given by R 2 (t, m i , m j ) = 1 2 dm l dm e R 2 (t, m i − m e , m e , m j , m l ) Similarly, R 3 is given by pbh .
The merger rate density of PBH binaries with the masses m i and m j in the third-merger process is given by The total merger rate density of PBH binaries with the masses m i and m j detected by LIGO-Virgo is given by In the single-merger case, we have α = −(m i + m j ) 2 ∂ 2 ln R(t, m i , m j )/∂m i ∂m j = 36/37 which is independent of the PBH mass function. It is consistent with the result obtained in [18]. However, by taking account into the merger history of PBHs, α depends on the PBH mass function, which could help us reconstruct the mass function of PBHs.

IV. APPLICATIONS
Let us consider two typical PBH mass functions. One is a monochromatic function and the other is a power law function.

A. Monochromatic mass function
In this subsection, we consider the following monochromatic mass function [15,29,30] In this case, we can rewrite (22) and (23) as which is consistent with (10). Similarly, the merger rate of PBH binaries in the second-merger process is given by In Fig. 4, we show the merger rate of PBH binaries in the third-merger process as a function of f pbh , which scales as f Therefore, in the monochromatic case, the effect of the merger history on the merger rate of PBH binaries is negligible. FIG. 5. Event rate of first-merger (green), second-merger (blue) and third-merger (red) of PBH binaries with the mass 30M⊙ at the present time as a function of the PBH abundance. The merger rate R = 12 − 213 Gpc −3 yr −1 estimated by the LIGO-Virgo Collaboration is shown as the shaded region colored orange [24].

V. CONCLUSIONS
We have developed the formalism to calculate the merger rate of PBH binaries with a general mass function, by taking into account the merger history of PBHs. In the monochromatic case, we find that R 1 ≫ R 2 ≫ R 3 , which is independent f pbh . Therefore, the effect of the merger history on the merger rate of PBH binaries is negligible. However, the multiple-merger events may play an important role in the merger rate density of PBH binaries in the non-monochromatic case. For example, for the power-law mass function the effect of the merger history on the merger rate of PBH binaries could not be negligible. In the future, more and more coalescence events of black hole binaries will be detected by LIGO-Virgo. This will provide more rich information on the merger rate distribution of black hole binaries to test the PBH scenario.
The effects of the tidal field from the smooth halo, the encountering with other PBHs, the baryon accretion and present-day halos, are carefully investigated in [17]. It is found in [17] that these effects make no significant contributions to the overall merger rate. We therefore neglected these subdominant effects throughout our computation.