Extreme-mass-ratio inspirals produced by tidal capture of binary black holes

Extreme-mass-ratio inspirals (EMRIs) are important gravitational-wave (GW) sources for future space-based detectors. The standard model consists of one stellar-mass black hole spiraling into a supermassive one, and such a process emits low-frequency (~10−3 Hz) GWs, which contain rich information about the space–time geometry around the central massive body. Here we show that the small bodies in EMRIs, in fact, could be binary black holes, which are captured by the massive black holes during earlier close encounters. About 30% of the captured binaries coalesce due to the perturbation by the massive bodies, resulting in a merger rate of 0.03 Gpc3 yr−1 in the most optimistic scenario. The coalescence generates also high-frequency (~102 Hz) GWs detectable by ground-based observatories, making these binary-EMRIs ideal targets for future multi-band GW observations. An extreme-mass-ratio inspiral, generally consists of a stellar-mass black hole and a supermassive black hole. The authors propose an alternative scenario where the small black hole is replaced by a binary black hole, and show how likely their gravitation wave signal can be detected, simultaneously, by LISA and LIGO.

A n extreme-mass-ratio inspiral (EMRI) normally consists of a compact stellar object, such as a stellar-mass black hole (BH), and a supermassive black hole (SMBH). It is an important target for future space-borne, milli-Hz gravitationalwave (GW) detectors, such as the laser interferometer space antenna (LISA 1 ): an EMRI could dwell in the LISA band for years [2][3][4] , accumulating as many as 10 3 -10 4 GW cycles in the data stream. Such a long waveform contains rich information about the space-time, as well as the astrophysical environment at the immediate exterior of a SMBH [5][6][7] . To decode this information, our model of an EMRI has to be highly accurate, both mathematically and physically.
In the canonical model of an EMRI, the stellar object is captured by the SMBH in two possible ways 3 : (i) it is scattered by other stars, a process known as relaxation, to such a small distance to the SMBH that the stellar object loses a significant amount of its orbital energy through GW radiation and becomes bound to the big BH 8,9 . The event rate of this type of EMRIs is difficult to estimate because it depends on factors that are poorly constrained by observations. The current estimation lies in a broad range between 10 −9 and 10 −6 per galaxy per year 4,10-15 ; (ii) the small body could come from a binary which is also scattered to the vicinity of the SMBH. If the distance between the SMBH and the centre-of-mass of the binary becomes smaller than the tidal radius, R t ≡ a(M 3 /m 12 ) 1/3 , where a and m 12 are the semi-major axis and total mass of the binary and M 3 the mass of the SMBH, the interaction in general ejects the lighter member of the binary and leaves the other, more massive member on a bound orbit around the SMBH [16][17][18] . The corresponding event rate could be comparable to that produced by the aforementioned capture of individual BHs from eccentric orbits 18 .
In this article, we point to a third possibility. We show that a BH binary (BHB) could be tidally captured by a SMBH to a bound orbit. This could happen when the binary passes by the SMBH at such a close distance that the tidal interaction transforms a fraction of the kinetic energy into the internal potential energy of the binary 17,19 . We refer to the resulting triple system as the "binary EMRI" (b-EMRI) and investigate its long-term evolution using numerical simulations. We find that the captured binary, as a single unit, could survive around the SMBH for a long time. A significant fraction of the BHBs coalesce before their orbits around the SMBHs completely circularize. The result could be a burst of high-frequency (10 2 Hz) GWs happening at the same time and in the same sky position of a low-frequency (10 −3 Hz) EMRI.

Results
Formation of a b-EMRI. Our BHB starts most likely around the influence radius of the SMBH, R inf ≡ GM 3 /σ 2 , where σ is the onedimensional velocity dispersion of the stars surrounding the SMBH (or "nuclear star cluster", NSC 20 ). Deeper inside this radius the BHB is susceptible to ionization by energetic encounters with interlopers 21,22 and further outside the gravitational influence of the SMBH becomes negligible. Using the typical values M 3 = 10 6 M ⊙ and σ = 60 km s −1 for NSCs, we find R inf ≃ 1 pc.
Because the specific energy of the BHB initially is 3σ 2 /2 − GM 3 /R inf = σ 2 /2, the orbit is a hyperbola with an asymptotic velocity of v 0 = σ at infinity. Previous studies focus on parabolic orbits 19 and have shown that tidal capture happens when the pericentre distance of the parabola, R p , becomes comparable to the tidal radius R t , i.e., where 1 ≲ ξ ≲ 5. In the last equation, R g = GM 3 /c 2 is the gravitational radius of the SMBH, c is the speed of light, m 1 and m 2 are the masses of the two stellar BHs, q ≡ m 2 /m 1 is the mass ratio assuming that m 1 ≥ m 2 , and r g = Gm 1 /c 2 is the gravitational radius of the bigger stellar BH. In the following we assume q ≃ 1 because star clusters produce, most likely, equalmass binaries [23][24][25] . We also scale a with 10 5 r g and the reason will become clear later. The lifetime of these BHBs can been calculated in the Keplerian approximation 26 , where _ a is the decay rate of the semi-major axis due to GW radiation, e is the orbital eccentricity, and F(e) = (1 − e 2 ) 7/2 (1 + 73e 2 /24 + 37e 4 /96) −1 .
To prove that the above criterion for tidal capture also applies to hyperbolic orbits, we conduct scattering experiments using the numerical tool FEWBODY 27 . The default parameters are M 3 = 10 6 M ⊙ , m 1 = m 2 = 10M ⊙ , e = 0.1, and a 0 = 5 × 10 4 r g ≃ 0.005 AU. Initially we place the BHB at infinity with a velocity of v 0 = 60 km s −1 . The direction of the velocity is chosen in such a way that the hyperbolic orbit initially has a pericentre distance of ξR t . We then numerically integrate the triple system and, after the first pericentre passage, record the relative energy change of the BHB, η ≡ (E − E 0 )/|E 0 |, where E is the final energy of the BHB and E 0 = −Gm 1 m 2 /(2a 0 ) is the initial energy.
If η ≥ 1, the binary is tidally disrupted. Otherwise, the binary survives the pericentre passage. Since the BHB initially has a kinetic energy of K 0 ¼ m 12 v 2 0 =2, we know that the binary becomes bound to the SMBH if E − E 0 > K 0 , i.e., η ≳ a 0 v 2 0 =ðGμÞ ≃ 0.0041, where μ = m 1 m 2 /m 12 is the reduced mass of the binary. If η < 0.0041, the BHB is re-ejected to infinity. Figure 1 shows the cumulative distribution of η for different ξ, where we have randomized the initial inclination of the BHB and repeated the experiment 10 3 times for each value of ξ. We find that the captured fraction is f cap = (10-30)% when ξ varies between 0.7 and 2.5. Therefore, b-EMRIs can form via capturing the BHBs on hyperbolic orbits.
The orbital elements of the b-EMRIs can be derived from two conservation laws: (i) because of the conservation of energy, the binding energy (relative to the SMBH) of a captured BHB is , where the last approximation uses the previous result that η~0.1. From E 3 we derive a semi-major axis of R ≃ (a 0 /η)(M 3 /μ). (ii) The pericentre remains to be R p because of the conservation of angular momentum. Consequently, the eccentricity e 3 satisfies the condition Interestingly, it does not depend on our assumption of a 0 . The small value of 1 − e 3 indicates that a newly formed b-EMRI, in general, is very eccentric.
The high eccentricity will affect the stability of a b-EMRI in two ways. On the one hand, the orbital motion of the BHB around the SMBH is more susceptible to perturbations by the surrounding stars because the angular momentum, ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi , is small. To be more quantitative, suppose T rlx characterizes the typical timescale for stellar relaxation processes to completely alter the orbital elements of a circular orbit 3 , the timescale to spoil an orbit of an eccentricity of e 3 is only T rlx 1 À e 2 3 À Á . On the other hand, the orbit circularizes very fast due to GW radiation because the associated timescale T gw (during which 1 − e 3 increases) is proportional to R 4 1 À e 2 3 À Á 7=2 . For our purpose, we use the Based on the understanding of the above two effects, we find that the b-EMRIs with 1 À e 2 3 À Á T rlx >T gw are not any more affected by the relaxation of the background stars. Together with Eqs. (4) and (6), we find that these "clean" b-EMRIs satisfy This is the reason that we scaled a with 10 5 r g in the previous equations. In the last equation we have adopted a typical value of 10 9 years for T rlx . We choose this value because our b-EMRI resides in the phase space of (1 − e 3 )~10 −4 and R = R p /(1 − e 3 ) 0.02 pc (assuming a = a cri ), a region where two-body scattering dominates the relaxation process 28,29 .
Long-term evolution. Although the long-term interaction between a SMBH and a binary has been studied previously 21,[30][31][32][33][34][35][36][37] , these earlier works focus on a binary that is far away from the SMBH, so that there is no energy exchange between the "inner binary", i.e., the BHB, and the "outer binary", i.e. the pair formed by the BHB and the SMBH. Such triple systems are called "hierarchical triples". Moreover, because of the large distance, the GW radiation of the outer binary is also unimportant to the overall dynamics. Under these circumstances, the evolution reduces to a periodic oscillation of the inner binary in the eccentricity space, which is known as the "Lidov-Kozai" cycle 38,39 . In our problem, however, the BHB is much closer to the SMBH so that the energies of the inner and outer binaries are no longer conserved and the GW radiation from the outer binary is not negligible. These are the key differences between our problem and other studies.
The b-EMRIs, in general, do not satisfy the criterion of hierarchical triple. In our problem, since M=m 12 ) 1 and 1 + e 3 ≃ 2, the criterion for hierarchical triple 40 reduces to Our b-EMRIs do not satisfy this condition because tidal capture requires that 1 ≲ R p /R t ≲ 2.5 (see Fig. 1). For this reason, we use again Fewbody to evolve the triple system numerically. Besides adopting the default parameters, we also choose ξ = 2 and 1 − e 3 = 10 −4 as the initial conditions and start the BHB at the apocentre with a random inclination. Moreover, we notice that T gw is of the same order of t gw . If the eccentricity of the inner binary gets excited by the tidal force of the SMBH, the BHB may coalesce. This would terminate the b-EMRI. To include this effect in the simulation, we implement post-Newtonian (PN) corrections to the equations of motion of the inner binary, as has been done by one of the authors in a previous work 24 . Figure 2 shows an example of a coalescence. It also shows that the energy of the inner binary is not conserved, confirming our prediction that the system is not a hierarchical triple.
With the above setup, we run 10 4 simulations to get reasonable statistics. Figure 3 summarizes the outcomes. Interestingly, in about 30% of the cases, the BHBs coalesce well before the time limit t = T gw /2 (for circularization) is reached. Because coalescing BHBs generate high-frequency GWs (10 2 Hz) and their orbital motion around the SMBH produces low-frequency ones (10 −3 Hz), we conclude that b-EMRIs, in principle, could be detected by both ground-based and space-based detectors. Figure 4 illustrates this idea and we will elaborate this point in Discussion.
Furthermore, about 0.6% of our b-EMRIs survive till the time t = T gw /2. In principle they could have significantly circularized, i.e., both the semi-major axis R and the eccentricity e 3 of the outer binary could have significantly decreased. However, our code does not include PN corrections to the outer binary and this shortcoming prevents us from predicting the fate of these longlived b-EMRIs. One way of circumventing this issue is to repeat our b-EMRI simulations with smaller e, i.e., e = 0.999 or 0.99, to mimic the process of circularization, but the integration time exceeds our current computational capacity. We plan to resolve this issue in a future work. We note that PN corrections would not significantly affect the previous results about the capturing process because during the first periapsis passage the periastron shift, which is of the order of 2R g /R p~0 .027, is small and the energy lost via GW radiation is negligible relative to the dynamical energy exchange |E − E 0 |. Formation rate. The formation rate of b-EMRIs can be estimated with Γ = pf cap Γ BHB . Here Γ BHB denotes the formation rate of the BHBs with a ≃ a cri in NSCs, p is the probability that such a BHB is on a hyperbolic orbit with ξ ≲ 2, and f cap ≃ 0.3 is the probability for capture (from Fig. 1).
To quantify Γ BHB , we notice that the BHBs of our interest are short-lived (see Eq. (3)) relative to the age of a galaxy, which is about 10 10 years. Therefore, Γ BHB is equal to the event rate of BHB mergers in a galaxy. Recent calculations 22,41 suggest that Γ BHB~f ew × 10 −7 galaxy −1 year −1 . This rate dose not account for the enhancement due to the interaction between BHBs and SMBHs 21 . When this enhancement is considered, the rate could be as high as Γ BHB = (0.8 − 2) × 10 −6 galaxy −1 year −1 according to the most optimistic estimations 33,35 .
To quantify p, we first recall that our BHBs start from a rather large distance, R inf , from the SMBH, where the binaries have a maximum angular momentum of L c ' ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi GM 3 R inf p (i.e., the angular momentum for a circular orbit of the same energy). To reach a pericentre distance of R p , the BHBs would have to diffuse in the angular-momentum space (L-space) due to relaxation to reach a region where L ≲ L p $ ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi 2GM 3 R p p . If the relaxation process is incoherent 42 , like the aforementioned two-body relaxation, the probability that the BHB appears at L ≲ L p is proportional to L 2 . With this consideration and assuming a cri = 5 × 10 4 r g , we derive p ≃ (L p /L c ) 2 ≃ 4 × 10 −6 . Such a small probability renders the b-EMRI rate formidably small, i.e. Γ( 10 −12 -10 −11 ) galaxy −1 year −1 .
However, if the relaxation process is coherent, e.g., due to a massive perturber or an axisymmetric gravitational potential [43][44][45][46] , the probability p scales with L, so that p ≃ L p /L c ≃ 2 × 10 −3 . With this number, we find Γ ranges from 10 −10 galaxy −1 year −1 to as large as few × 10 −9 galaxy −1 year −1 in the most optimistic estimation. These numbers should be regarded as upper limits because the coherent relaxation processes normally affect only part of the angular-momentum space so that the real value of p could be significantly smaller than L p /L c .
Comparing these rates with the event rate of standard EMRIs 4,10-15 , which falls in the range of 10 −9 -10 −6 yr −1 galaxy −1 , we conclude that in the most optimistic case the b-EMRI rate could approach the lower boundary of the EMRI rate. Furthermore, since there are on average 2 × 10 7 galaxies per Gpc 3 (earlier studies 41 use the same galaxy density), we finally find that the above formation rate of b-EMRIs is equivalent to (10 −5 -10 −4 ) Gpc 3 year −1 in the pessimistic case and 0.1 Gpc 3 year −1 in the most optimistic one.
Although the formation rate is relatively low, each b-EMRI has a non-negligible lifetime. The extended lifetime leads to a substantial duty cycle, D bemri , for a galaxy to host a b-EMRI. To estimate D bemri , we use (i) the formation rate Γ and (ii) the probability of surviving b-EMRIs as a function of time, f bemri (t), which is shown in Fig. 3 as the black solid curve, to derive D bemri ≃ Γtf bemri (t). In the previous equation, the maximum value of tf bemri (t) appears at t ≃ 7600 years, where f bemri (t) ≃ 0.1. Therefore, we find that D bemri ≃ (10 −9 -10 −8 ) galaxy −1 in our pessimistic scenario and D bemri~1 0 −6 galaxy −1 in the most optimistic one. This result means that within a volume of 1 Gpc 3 (corresponding to a distance of about 600 Mpc), there are, on average, N bemri~0 .02-20 b-EMRIs.

Discussion
We have seen that the BHBs in b-EMRIs, in general, have a ≲ a cri . Interestingly, these BHBs are detectable by LISA. This is so because LISA is sensitive to the GWs with a frequency of f~10 −3 Hz and the strongest GW mode that a BHB emits 47,48 is of a frequency of f ' π À1 ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi Gm 12 =½að1 À eÞ 3 q . As a result, those BHBs with a semi-major axis of a LISA $ 3:4 10 4 r g f À2=3 À3 ð1 þ qÞ 1=3 1 À e are inside the LISA band, where f −3 := f/(10 −3 Hz). This critical semi-major axis is comparable to a cri . It is important to note that the BHBs in b-EMRIs are accelerating in the gravitational potential of the SMBHs. This condition could induce detectable shift to the phase of the GW inspiral waveform 49,50 . We have seen that within a distance of about 600 Mpc there are at most 20 accelerating BHBs for LISA to detect. In our default model, the BHB is at a typical distance of R ≃ 0.02 pc from the central SMBH. The corresponding acceleration is GM 3 / R 2 and it could induce a phase shift as large as 10 3 π per year according to the formula derived in the earlier studies 49,50 .
Furthermore, the orbital motion of the BHBs around the SMBHs also generates GWs. It is important to understand whether this radiation is detectable. Similar to the previous analysis for BHBs, we calculate the frequency of the strongest GW mode From Eqs. (1) and (7), we find that our b-EMRIs satisfy the condition Therefore, they are indeed inside the LISA band: Each pericentre passage will generate a burst of GWs detectable by LISA [51][52][53][54] . These b-EMRIs are particularly interesting from the observational point of view because they generate two types of GWs, i.e., the continuous waves from the BHBs and bursts from the orbital motion of the BHBs around the SMBHs. Moreover, the two types of GWs are emitted at the same time in the same band. LISA can detect a single burst of the second type out to a distance of about 200(m 12 /20 M ⊙ ) Mpc 52 . If the BHB has a total mass of about 60M ⊙ , like what the Laser Interferometer Gravitational-wave Observatory (LIGO) first detected 23 , the detection horizon would become as far as 600 Mpc. We already know that there are at most N bemri ≃ 20 b-EMRIs within this distance. Now we estimate the chance of catching a GW burst at the moment of the pericentre passage. Since the majority of the b-EMRIs cannot circularize (see Fig. 3), the successive pericentre passages are separated by a long orbital period of the outer binary, which is about P 3 ≃ 200 years in our fiducial model. Given that LISA is designed to have a mission duration of Δt = 5 years 1 , the chance of catching a b-EMRI at its pericentre passage is about (Δt/P)N bemri~5 0% in our most optimistic scenario.
Probably the most interesting feature about b-EMRIs is that the BHBs have a 30% chance to coalesce. The coalescence gives rise to a LIGO/Virgo event. The event rate is 0.3Γ, which is 0.03 Gpc 3 year −1 in the most optimistic scenario and two orders of magnitude smaller in the most pessimistic one. If detected, the LIGO/ Virgo event is likely separated from a LISA event, i.e., the GW burst generated during the pericentre passage, by about half of the orbital period of the outer binary, or 100 years in our fiducial model, because the coalescence happens most likely at the apocentre where the passing time is the longest. For this reason, it would be difficult to identify the LISA counterpart to the LIGO/ Virgo event.
However, exceptions might exist. If the BHB in a b-EMRI could circularize around the SMBH, the orbital period P 3 would be much shorter and hence the LISA and LIGO/Virgo events would appear much closer in time. According to our calculation, the event rate of such an exceptional case is at most Γ f bemri (t = T gw /2)~6 × 10 −4 Gpc 3 year −1 .
In summary, we have shown that the small bodies in EMRIs could, in fact, be BH binaries. Such b-EMRIs could form via tidal capture and a significant fraction of the captured BHBs could merge due to the perturbation by the SMBHs. Although rare, they are ideal targets for future multi-band GW observations 55,56 . First, the merger generates high-frequency GWs, i.e. a LIGO/ Virgo event that is in the same sky location of a LISA EMRI event. Second, the high-frequency GWs could be redshifted because they are generated very close to a SMBH 57 , providing a rare opportunity of studying the propagation of GWs in the regime of strong gravity. Third, the merger also induces a kick to the BH remnant 58 . This kick causes a glitch in the EMRI waveform, which, through a careful analysis, is discernible in the data stream 59 .

Low frequency GW
High frequency GW (i) (ii) (iii) Fig. 4 Three evolutionary stages of a binary extrme-mass-ratio inspiral (b-EMRI). First, in (i), a compact black hole binary (BHB) is captured to a bound orbit around a supermassive black hole (SMBH) because the pericentre distance becomes comparable to the tidal-disruption radius, R t . Next, in (ii), the outer binary circularizes due to gravitational wave (GW) radiation, and the GW frequency lies in the band of a space-borne GW detector. Finally, in (iii), the tidal force of the SMBH becomes strong enough to excite the eccentricity of the inner BHB and drive it to merge. The merger produces high-frequency GWs

Data availability
The data sets generated during the current study are available from the corresponding author on reasonable request.