Analytical effective one-body formalism for extreme-mass-ratio inspirals with eccentric orbits

The EOB formalism was originally introduced in [23, 24] to describe the evolution of binary systems. We start by considering an EMRI system with a central Kerr black hole m₁ and an inspiraling object m₂ (we assume it is nonspinning, for simplicity) which is restricted to the equatorial plane of m₁ (m₂ ≪ m₁). For the moment, we neglect the radiation reaction effects and focus on purely geodesic motion. The conservative orbital dynamics is derived via Hamilton's equations using the EOB Hamiltonian ${H}_{\mathrm{EOB}}\,=M\sqrt{1+2\nu ({\hat{H}}_{\mathrm{eff}}-1)}$, where M = m₁ + m₂, ν = m₁ m₂/M², μ = ν M, and ${\hat{H}}_{\mathrm{eff}}={H}_{\mathrm{eff}}/\mu$. The deformed Kerr metric is given by [42]

$$\begin{eqnarray}&&{g}^{{tt}}=-\displaystyle \frac{{{\rm{\Lambda }}}_{t}}{{{\rm{\Delta }}}_{t}{\rm{\Sigma }}},\end{eqnarray} \tag{ 1 }$$

$$\begin{eqnarray}&&{g}^{{rr}}=\displaystyle \frac{{{\rm{\Delta }}}_{r}}{{\rm{\Sigma }}},\end{eqnarray} \tag{ 2 }$$

$$\begin{eqnarray}&&{g}^{\theta \theta }=\displaystyle \frac{1}{{\rm{\Sigma }}},\end{eqnarray} \tag{ 3 }$$

$$\begin{eqnarray}&&{g}^{\phi \phi }=\displaystyle \frac{1}{{{\rm{\Lambda }}}_{t}}\left(-\displaystyle \frac{{\widetilde{\omega }}_{{fd}}^{2}}{{{\rm{\Delta }}}_{t}{\rm{\Sigma }}}+{\rm{\Sigma }}\right),\end{eqnarray} \tag{ 4 }$$

$$\begin{eqnarray}&&{g}^{t\phi }=-\displaystyle \frac{{\widetilde{\omega }}_{{fd}}}{{{\rm{\Delta }}}_{t}{\rm{\Sigma }}}.\end{eqnarray} \tag{ 5 }$$

The quantities Σ, Δ_t , Δ_r , Λ_t , and ${\widetilde{\omega }}_{\mathrm{fd}}$ in equations (1)–(5) are given by

$$\begin{eqnarray}&&{\rm{\Sigma }}={r}^{2},\end{eqnarray} \tag{ 6 }$$

$$\begin{eqnarray}&&{{\rm{\Delta }}}_{t}={r}^{2}\left[A(u)+\displaystyle \frac{{a}^{2}}{{M}^{2}}\,{u}^{2}\right],\end{eqnarray} \tag{ 7 }$$

$$\begin{eqnarray}&&{{\rm{\Delta }}}_{r}={{\rm{\Delta }}}_{t}\,{D}^{-1}(u),\end{eqnarray} \tag{ 8 }$$

$$\begin{eqnarray}&&{{\rm{\Lambda }}}_{t}={\left({r}^{2}+{a}^{2}\right)}^{2}-{a}^{2}\,{{\rm{\Delta }}}_{t},\end{eqnarray} \tag{ 9 }$$

$$\begin{eqnarray}&&{\widetilde{\omega }}_{\mathrm{fd}}=2a\,M\,r+{\omega }_{1}^{\mathrm{fd}}\,\nu \,\displaystyle \frac{{{aM}}^{3}}{r}+{\omega }_{2}^{\mathrm{fd}}\,\nu \,\displaystyle \frac{{{Ma}}^{3}}{r},\end{eqnarray} \tag{ 10 }$$

where a = ∣ S _Kerr∣/M is the effective Kerr parameter and u = M/r. The values of ${\omega }_{1}^{\mathrm{fd}}$ and ${\omega }_{2}^{\mathrm{fd}}$ are about −10 and 20, respectively, as given by a preliminary comparison of the EOB model with numerical relativity results [43, 44]. The metric potentials A and D for the EOB model are given in appendix A. Though the spin of the small body is zero, its effective spin is nonzero. The effective spins of the particle and the deformed Kerr black hole are [45]

$$\begin{eqnarray}&&{S}^{* }={S}_{1}\displaystyle \frac{{m}_{2}}{{m}_{1}}+{S}_{2}\displaystyle \frac{{m}_{1}}{{m}_{2}}+{{\rm{\Delta }}}_{{\sigma }^{* }}^{(1)}+{{\rm{\Delta }}}_{{\sigma }^{* }}^{(2)}+{{\rm{\Delta }}}_{{\sigma }^{* }}^{(3)},\end{eqnarray} \tag{ 11a }$$

$$\begin{eqnarray}&&{S}_{\mathrm{Kerr}}={S}_{1}+{S}_{2}.\end{eqnarray} \tag{ 11b }$$

For S₂ = 0, it is clear that S_Kerr = S₁ and ${S}^{* }\sim { \mathcal O }(\nu )$. The last three terms in the equation for S^* are at least ${ \mathcal O }({\nu }^{2})$ [42, 45, 46]. In this paper, considering the small mass ratio, we can omit these three terms for simplicity, and retain sufficient accuracy. The effective spins can then be rewritten as

$$\begin{eqnarray}&&{S}^{* }={aM}\displaystyle \frac{{m}_{2}}{{m}_{1}},\end{eqnarray} \tag{ 12a }$$

$$\begin{eqnarray}&&{S}_{\mathrm{Kerr}}={aM}.\end{eqnarray} \tag{ 12b }$$

The effective Hamiltonian should then describe a spinning test particle in the deformed Kerr metric, and reads

$$\begin{eqnarray}&&{H}_{\mathrm{eff}}={H}_{\mathrm{NS}}+{H}_{{\rm{S}}},\end{eqnarray} \tag{ 13 }$$

where H_S is the Hamiltonian caused by the spin of the effective particle. In this work, to develop the model step by step, we ignore the spin of the small compact object, though some researchers have argued that the spin is non-negligible for EMRIs [36, 47]. In this situation, considering that the effective spin parameter of the small object s ∼ μ a/M should be very small for EMRIs, we plan to omit this term in this work, for simplification. In table 1, the ratios of H_S to H_NS for different spin and orbital radii for small mass-ratio (equal or less than 10⁻¹ ) circular orbits are listed. We find that the magnitude of this spin Hamiltonian term is a fraction of only about 10⁻² ν of the nonspinning term; even for extreme Kerr parameters and close to the horizon, it is just 10⁻¹ ν, and the relative error in the orbital frequency due to omitting H_S is also ≲ 10⁻² ν (for smaller spin and larger orbital radii, this error will be reduced). Therefore, for the EMRIs described in this paper, we omit this nonlinear term, and will include it in the next work on the topic of comparable binaries.

Table 1. Ratios of H_S to H_NS with different spin and orbital radii (left); frequency errors compared to the full EOB Hamiltonian (right).

$\tfrac{{H}_{{\rm{S}}}}{{H}_{\mathrm{NS}}}\left(/\nu \right)$	${r}_{{}_{\mathrm{ISCO}}}$	r = 10	r = 20
a = 0.9	1.3 × 10−1	3.8 × 10−3	6.9 × 10−4
a = 0.5	2.1 × 10−2	2.1 × 10−3	4.1 × 10−4
a = 0.1	2.4 × 10−3	5.3 × 10−4	8.8 × 10−5

$\tfrac{{\rm{\Delta }}{\omega }_{\phi }}{{\omega }_{\phi }}\left(/\nu \right)$	${r}_{{}_{\mathrm{ISCO}}}$	r = 10	r = 20
a = 0.9	1.3 × 10−1	4.5 × 10−2	1.8 × 10−2
a = 0.5	6.8 × 10−2	2.8 × 10−2	1.1 × 10−2
a = 0.1	1.2 × 10−2	6.0 × 10−3	2.3 × 10−3

The most dominant part, H_NS, is the Hamiltonian for a nonspinning test particle of mass μ in the deformed Kerr metric. Using equation (2.26) of [48], we have

$$\begin{eqnarray}\begin{array}{rcl}{H}_{\mathrm{NS}} & = & \displaystyle \frac{{g}^{t\phi }}{{g}^{{tt}}}{P}_{\phi }+\displaystyle \frac{1}{\sqrt{-{g}^{{tt}}}}\\ & & \times \sqrt{{\mu }^{2}+\left[{g}^{\phi \phi }-\displaystyle \frac{{\left({g}^{t\phi }\right)}^{2}}{{g}^{{tt}}}\right]{P}_{\phi }^{2}+{g}^{{rr}}{P}_{r}^{2}+\displaystyle \frac{{Q}_{4}{P}_{r}^{4}}{{\mu }^{2}}},\end{array}\end{eqnarray} \tag{ 14 }$$

where P_r and P_ϕ are the radial and azimuthal angular momenta. The functions ${Q}_{4}=\tfrac{2(4-3\nu )\nu {M}^{2}}{{r}^{2}}$, and ${{ \mathcal Q }}_{4}\,={Q}_{4}{P}_{r}^{4}/{\mu }^{2}$ [49] represent a non-geodesic term that appears at the third PN order.

The energy of the system is given by

$$\begin{eqnarray}&&E={H}_{\mathrm{EOB}},\end{eqnarray} \tag{ 15 }$$

which implies the relation

$$\begin{eqnarray}&&{\hat{H}}_{\mathrm{eff}}(E)=1+\displaystyle \frac{1}{2\nu }\left(\displaystyle \frac{{E}^{2}}{{M}^{2}}-1\right).\end{eqnarray} \tag{ 16 }$$

The canonical EOB dynamics without the radiation reaction are

$$\begin{eqnarray}\begin{array}{rcl}\dot{r} & = & \displaystyle \frac{\partial {H}_{\mathrm{EOB}}}{\partial {P}_{r}},{\dot{P}}_{r}=-\displaystyle \frac{\partial {H}_{\mathrm{EOB}}}{\partial r},\\ \dot{\phi } & = & \displaystyle \frac{\partial {H}_{\mathrm{EOB}}}{\partial {P}_{\phi }},{\dot{P}}_{\phi }=-\displaystyle \frac{\partial {H}_{\mathrm{EOB}}}{\partial \phi }=0.\end{array}\end{eqnarray} \tag{ 17 }$$

The effective Hamiltonian associated with the deformed Kerr metric has the form

$$\begin{eqnarray}&&{\hat{H}}_{\mathrm{eff}}=\displaystyle \frac{\sqrt{{{\rm{\Delta }}}_{t}\left({{\rm{\Lambda }}}_{t}\left({r}^{2}+{{\rm{\Delta }}}_{r}{\hat{P}}_{r}^{2}+{r}^{2}{Q}_{4}{\hat{P}}_{r}^{4}\right)+{r}^{4}{\hat{P}}_{\phi }^{2}\right)}+{\widetilde{\omega }}_{\mathrm{fd}}{\hat{P}}_{\phi }}{{{\rm{\Lambda }}}_{t}},\end{eqnarray} \tag{ 18 }$$

where we define the reduced momenta ${\hat{P}}_{r}={P}_{r}/\mu$ and ${\hat{P}}_{\phi }={P}_{\phi }/\mu$. Solving equation (18) for P_r in terms of (E, P_ϕ , r) leads to

$$\begin{eqnarray}&&{\hat{P}}_{r}^{2}=\displaystyle \frac{1}{2{r}^{2}{Q}_{4}}\left[\sqrt{{{\rm{\Delta }}}_{r}^{2}-\displaystyle \frac{4{r}^{2}{Q}_{4}\left({{\rm{\Lambda }}}_{t}\left(2{\widetilde{\omega }}_{\mathrm{fd}}{\hat{H}}_{\mathrm{eff}}{\hat{P}}_{\phi }-{{\rm{\Lambda }}}_{t}{\hat{H}}_{\mathrm{eff}}^{2}+{r}^{2}{{\rm{\Delta }}}_{t}\right)+{\hat{P}}_{\phi }^{2}\left({r}^{4}{{\rm{\Delta }}}_{t}-{\widetilde{\omega }}_{\mathrm{fd}}^{2}\right)\right)}{{{\rm{\Delta }}}_{t}{{\rm{\Lambda }}}_{t}}}-{{\rm{\Delta }}}_{r}\right],\end{eqnarray} \tag{ 19 }$$

The conservative EOB equations of equatorial motion without the radiation reaction can then be written as

$$\begin{eqnarray}&&\dot{r}=\displaystyle \frac{{{\rm{\Delta }}}_{t}\left({{\rm{\Delta }}}_{r}{\hat{P}}_{r}+2{r}^{2}{Q}_{4}{\hat{P}}_{r}^{3}\right)}{E\sqrt{{{\rm{\Delta }}}_{t}\left({{\rm{\Lambda }}}_{t}\left({r}^{2}+{{\rm{\Delta }}}_{r}{\hat{P}}_{r}^{2}+{r}^{2}{Q}_{4}{\hat{P}}_{r}^{4}\right)+{r}^{4}{\hat{P}}_{\phi }^{2}\right)}},\end{eqnarray} \tag{ 20a }$$

$$\begin{eqnarray}&&\dot{\phi }=\displaystyle \frac{1}{E\ {{\rm{\Lambda }}}_{t}}\left[{\widetilde{\omega }}_{\mathrm{fd}}+\displaystyle \frac{{r}^{4}{{\rm{\Delta }}}_{t}{\hat{P}}_{\phi }}{\sqrt{{{\rm{\Delta }}}_{t}\left({{\rm{\Lambda }}}_{t}\left({r}^{2}+{{\rm{\Delta }}}_{r}{\hat{P}}_{r}^{2}+{r}^{2}{Q}_{4}{\hat{P}}_{r}^{4}\right)+{r}^{4}{\hat{P}}_{\phi }^{2}\right)}}\right],\end{eqnarray} \tag{ 20b }$$

Equations (20) are more convenient than the original canonical EOB equations (17), because the dependence on ${\hat{P}}_{r}$ (19) has been eliminated by the energy E, which is a constant in a conservative system and only changes due to the radiation reaction if GW fluxes are considered.

The constants of motion and the dynamic equations in the last subsection can be written in terms of the geometrized orbital elements, i.e., the semilatus rectum, p, and the eccentricity, e. This will make the description of the system more intuitive. For an eccentric orbit, periastron and apastron points exist, which can be expressed as

$$\begin{eqnarray}&&{r}_{1}=\displaystyle \frac{{pM}}{1-e},\ \ \ \ {r}_{2}=\displaystyle \frac{{pM}}{1+e},\end{eqnarray} \tag{ 21 }$$

where r_1,2 are the turning points of the radial motion, i.e., the apastron and periastron, respectively. By setting the radial equation of motion (20a ) equal to zero, i.e., ${\hat{P}}_{r}=0$, we can solve for the two points. Furthermore, by taking r₁, r₂ into equation (19) to make ${\hat{P}}_{r}=0$ again, we finally obtain the constants of motion $(E,{\hat{P}}_{\phi })$ in terms of (p,e)

$$\begin{eqnarray}&&{\hat{P}}_{\phi }^{2}=\displaystyle \frac{{\left({a}_{1}-{a}_{2}\right)}^{2}\left({b}_{1}^{2}+{b}_{2}^{2}\right)-({b}_{1}^{2}-{b}_{2}^{2})({b}_{1}^{2}{c}_{1}-{b}_{2}^{2}{c}_{2})-2({a}_{1}-{a}_{2}){b}_{1}{b}_{2}\sqrt{{\left({a}_{1}-{a}_{2}\right)}^{2}-\left({b}_{1}^{2}-{b}_{2}^{2}\right)\left({c}_{1}-{c}_{2}\right)}}{{\left({\left({a}_{1}-{a}_{2}\right)}^{2}-\left({b}_{1}^{2}{c}_{1}-{b}_{2}^{2}{c}_{2}\right)\right)}^{2}},\end{eqnarray} \tag{ 22a }$$

$$\begin{eqnarray}&&\displaystyle \frac{{E}^{2}}{{M}^{2}}=1+2\nu \left({a}_{1}{\hat{P}}_{\phi }+\sqrt{\displaystyle \frac{{c}_{1}{\hat{P}}_{\phi }^{2}+1}{{b}_{1}}}-1\right),\end{eqnarray} \tag{ 22b }$$

where the coefficients are

$$\begin{eqnarray}&&{a}_{1}=\displaystyle \frac{a{\left(1-e\right)}^{3}\left({\left(1-e\right)}^{2}\nu \left({\omega }_{1}^{\mathrm{fd}}+{a}^{2}{\omega }_{2}^{\mathrm{fd}}\right)+2{p}^{2}\right)}{{p}^{5}-{a}^{2}{p}^{3}(A({r}_{1})-2){\left(1-e\right)}^{2}},\end{eqnarray} \tag{ 23a }$$

$$\begin{eqnarray}&&{a}_{2}=\displaystyle \frac{a{\left(1+e\right)}^{3}\left({\left(1+e\right)}^{2}\nu \left({\omega }_{1}^{\mathrm{fd}}+{a}^{2}{\omega }_{2}^{\mathrm{fd}}\right)+2{p}^{2}\right)}{{p}^{5}-{a}^{2}{p}^{3}(A({r}_{2})-2){\left(1+e\right)}^{2}},\end{eqnarray} \tag{ 23b }$$

$$\begin{eqnarray}&&{b}_{1}=\sqrt{\displaystyle \frac{{a}^{2}{\left(1-e\right)}^{2}+A({r}_{1}){p}^{2}}{{p}^{2}-{a}^{2}(A({r}_{1})-2){\left(1-e\right)}^{2}}},\end{eqnarray} \tag{ 23c }$$

$$\begin{eqnarray}&&{b}_{2}=\sqrt{\displaystyle \frac{{a}^{2}{\left(1+e\right)}^{2}+A({r}_{2}){p}^{2}}{{p}^{2}-{a}^{2}(A({r}_{2})-2){\left(1+e\right)}^{2}}},\end{eqnarray} \tag{ 23d }$$

$$\begin{eqnarray}&&{c}_{1}=\displaystyle \frac{{\left(1-e\right)}^{2}}{{p}^{2}-{a}^{2}(A({r}_{1})-2){\left(1-e\right)}^{2}},\end{eqnarray} \tag{ 23e }$$

$$\begin{eqnarray}&&{c}_{2}=\displaystyle \frac{{\left(1+e\right)}^{2}}{{p}^{2}-{a}^{2}(A({r}_{2})-2){\left(1+e\right)}^{2}}.\end{eqnarray} \tag{ 23f }$$

The above formalism for a Kerr black hole is much more complicated than the Schwarzschild ones described in [33]. Obviously, for the test-particle limit ν → 0, the above results will return the geodesic motion of a test particle in Kerr space-time.

The description of geodesic motion around BHs is based on the semilatus rectum, p, and the eccentricity, e, together with two phase variables associated with the spatial geometries of the radial and azimuthal motions, denoted by ξ and ϕ. These variables are defined by expressing the radial motion as

$$\begin{eqnarray}&&r=\displaystyle \frac{{pM}}{1+e\cos \xi },\ \end{eqnarray} \tag{ 24 }$$

so that the periastron and apastron correspond to $\xi =(0,\pi )\mathrm{mod}\,2\pi$, respectively. Taking a derivation of equation (24), we can obtain the evolution equation for the phase variable ξ

$$\begin{eqnarray}&&\dot{\xi }=\frac{{\left(1+e\cos \xi \right)}^{2}}{{epM}\sin \xi }\dot{r}+\frac{\cot \xi }{e}\dot{e}-\frac{1+e\cos \xi }{{ep}\sin \xi }\dot{p}.\end{eqnarray} \tag{ 25 }$$

For a conservative system, $\dot{p}=\dot{e}=0$. However, if we take the radiation reaction of the GWs into account, the rate of change of the energy and the reduced angular momentum are given by

$$\begin{eqnarray}&&\displaystyle \frac{{\rm{d}}E}{{\rm{d}}t}=\displaystyle \frac{\partial E}{\partial p}\dot{p}+\displaystyle \frac{\partial E}{\partial e}\dot{e},\end{eqnarray} \tag{ 26a }$$

$$\begin{eqnarray}&&\displaystyle \frac{{\rm{d}}{\hat{P}}_{\phi }}{{\rm{d}}t}=\displaystyle \frac{\partial {\hat{P}}_{\phi }}{\partial p}\dot{p}+\displaystyle \frac{\partial {\hat{P}}_{\phi }}{\partial e}\dot{e}.\end{eqnarray} \tag{ 26b }$$

Using equations (26), we can express the evolution of (p, e) in terms of the energy and angular momentum fluxes of gravitational radiation

$$\begin{eqnarray}&&\dot{e}={c}_{{Ep}}\displaystyle \frac{{\rm{d}}{\hat{P}}_{\phi }}{{\rm{d}}t}-{c}_{{Lp}}\displaystyle \frac{{\rm{d}}E}{{\rm{d}}t},\qquad \dot{p}={c}_{{Le}}\displaystyle \frac{{\rm{d}}E}{{\rm{d}}t}-{c}_{{Ee}}\displaystyle \frac{{\rm{d}}{\hat{P}}_{\phi }}{{\rm{d}}t},\qquad \end{eqnarray} \tag{ 27a }$$

where the coefficients are given by

$$\begin{eqnarray}&&{c}_{{Cb}}=\displaystyle \frac{\partial C/\partial b}{(\partial E/\partial p)(\partial {\hat{P}}_{\phi }/\partial e)-(\partial E/\partial e)(\partial {\hat{P}}_{\phi }/\partial p)},\end{eqnarray} \tag{ 27b }$$

where $C=\left\{E,{\hat{P}}_{\phi }\right\}$, $b=\left\{p,e\right\}$, and the derivatives can be computed from the expressions in equation (22)

The final set of EOB equations of motion with the radiation reaction are equations (27a ) together with the evolution of the phases described by equations (25) and (20b ), and the radius of motion at an arbitrary moment, as given by equation (24). All the equations of motion are now expressed only in terms of the geometric parameters (p, e, ξ) and the effective Kerr parameter a.

Finally, we can get the coordinates of the effective test particle in terms of ξ only:

$$\begin{eqnarray}&&t(\xi )={\int }_{0}^{\xi }\displaystyle \frac{1}{{ \mathcal P }}{\rm{d}}\xi ,\end{eqnarray} \tag{ 28a }$$

$$\begin{eqnarray}&&r(\xi )=\displaystyle \frac{{pM}}{1+e\cos \xi },\end{eqnarray} \tag{ 28b }$$

$$\begin{eqnarray}&&\theta =\displaystyle \frac{\pi }{2},\end{eqnarray} \tag{ 28c }$$

$$\begin{eqnarray}&&\phi (\xi )={\int }_{0}^{\xi }\displaystyle \frac{\dot{\phi }}{{ \mathcal P }}{\rm{d}}\xi ,\end{eqnarray} \tag{ 28d }$$

where ${ \mathcal P }(e,p,\xi )$ is the first term on the right-hand side of equation (25)

$$\begin{eqnarray}\begin{array}{rcl}{ \mathcal P }(e,p,\xi ) & = & \displaystyle \frac{{\left(1+e\cos \xi \right)}^{2}}{{ep}\sin \xi }\\ & & \times \displaystyle \frac{{{\rm{\Delta }}}_{t}\left({{\rm{\Delta }}}_{r}{\hat{P}}_{r}+2{r}^{2}{Q}_{4}{\hat{P}}_{r}^{3}\right)}{E\sqrt{{{\rm{\Delta }}}_{t}\left({{\rm{\Lambda }}}_{t}\left({r}^{2}+{{\rm{\Delta }}}_{r}{\hat{P}}_{r}^{2}+{r}^{2}{Q}_{4}{\hat{P}}_{r}^{4}\right)+{r}^{4}{\hat{P}}_{\phi }^{2}\right)}}.\end{array}\end{eqnarray} \tag{ 28e }$$

We now calculate the fundamental orbital frequencies in terms of the geometric parameters. This has been done in the test-particle (ν = 0) limit due to the analytical integrals of the geodesic in Kerr space-time [50, 51]. However, in the EOB formalism with the mass-ratio correction, the situation becomes complicated. First, we express the radial frequency ω_r which reflects the period of the radial motion from the periastron to the apastron and back to the periastron again, and the orbital period T_r can be calculated by taking ξ from 0 to 2π, then

$$\begin{eqnarray}&&{\omega }_{r}=\displaystyle \frac{2\pi }{{\int }_{0}^{2\pi }\tfrac{1}{{ \mathcal P }}{\rm{d}}\xi }.\end{eqnarray} \tag{ 29 }$$

In equation (20b ), we have already written the variation of ϕ with coordinate time t. Strictly speaking, the radial motion is the real periodic motion. When the particle passes through the periastron twice, Δϕ will be larger than 2π because of the periastron precession caused by the relativistic effect. We can then define the frequency of the azimuthal motion to be Δϕ/T_r . Therefore, we can compute the azimuthal frequency from the orbital average of the ϕ motion as

$$\begin{eqnarray}&&{\omega }_{\phi }=\langle \dot{\phi }\rangle =\displaystyle \frac{{\int }_{0}^{2\pi }\tfrac{\dot{\phi }}{{ \mathcal P }}{\rm{d}}\xi }{{\int }_{0}^{2\pi }\tfrac{1}{{ \mathcal P }}{\rm{d}}\xi },\end{eqnarray} \tag{ 30 }$$

where $\dot{r}$ is given in equation (20a ), and P_r , E, and P_ϕ are expressed in equations (19) and (22). By just replacing r in these equations with ξ, we cause the above two integrals to only contain the argument ξ, meaning that they can easily be integrated. We do not write down the fully expanded expressions for these two integrals here, as they is direct and trivial.

With the expressions of the two fundamental frequencies at hand, we now investigate how the effective Kerr parameter a and the symmetric mass ratio ν affect the features of the radial and azimuthal frequencies and compare them with the case of the test-particle limit in conservative dynamics. The effect of the spin parameter on these frequencies for various values of the semilatus rectum, p, and the eccentricity, e, is shown in figure 1. Obviously, as the semilatus rectum increases, we can say that the mean orbital radius correspondingly increases, leading to a decrease in the periastron precession ω_ϕ /ω_r − 1, which approaches the Newtonian limit ω_r = ω_ϕ . Interestingly, as the spin of the central SMBH increases, the periastron precession decreases. The effect of eccentricity on the precession becomes obvious only when the separation, p, and spin, a, are both small (see the top three curves in figure 1).

Zoom In Zoom Out Reset image size

We now consider the radial and azimuthal frequencies' shifts due to the mass ratio of the binary in the absence of the radiation reaction. It is known that the motion of the test particle has a precise analytical solution. For a binary system, where the mass ratio cannot be omitted, we have given the analytical EOB orbital solution equations (28). In order to observe the impact of the mass ratio on the radial and azimuthal frequencies under different conditions with various orbital parameters, we compare the test-particle results and the EOB results in table 2. The rightmost column of this table shows the relative frequency shift divided by the mass ratio: Δω/(ω ν), where Δω/ω is the relative difference of the radial/azimuthal frequency between the EOB frequencies (ω_r0, ω_ϕ0) and the test-particle frequencies, i.e., (ω_r − ω_r0)/ω_r0, (ω_ϕ − ω_ϕ0)/ω_ϕ0.

Table 2. Relative differences in the orbital frequencies ω_r , ω_ϕ between the EOB model with small mass ratios and the test-particle approximation.

a/M	p/M	e		Test-particle	ν = 10−2	ν = 10−3	ν = 10−4	ν = 10−5	ν = 10−6	$\tfrac{{\rm{\Delta }}\omega }{\omega }(/\nu )$
0.99	5	0.1	ωr	0.050841032	0.052967383	0.051054792	0.050862411	0.050843170	0.050841250	4.2
			ωϕ	0.081375480	0.080796016	0.081318545	0.081369800	0.081374912	0.081375431	0.70
0.99	10	0.1	ωr	0.023863900	0.023996289	0.023876964	0.023865208	0.023864071	0.023863903	0.55
			ωϕ	0.030284712	0.030221400	0.030278593	0.030284103	0.030284660	0.030284683	0.20
0.99	5	0.6	ωr	0.031648541	0.032582060	0.031743453	0.031658044	0.031649491	0.031648640	3.0
			ωϕ	0.051669258	0.050254323	0.051523452	0.051654638	0.051667794	0.051669119	2.8
0.99	20	0.6	ωr	0.0052501837	0.0052537928	0.0052505388	0.0052502166	0.0052501873	0.0052501841	0.07
			ωϕ	0.0059449013	0.0059371561	0.0059441512	0.0059448231	0.0059448938	0.0059449005	0.13
0.5	5	0.1	ωr	0.030270602	0.032859397	0.030538233	0.030297445	0.030273290	0.030270870	8.9
			ωϕ	0.085402159	0.084958240	0.085355747	0.085397498	0.085401693	0.085402113	0.55

The frequency shift Δω/(ω ν) due to the mass ratio is almost independent of the mass ratio itself, and the relative shift Δω/ω is in the range of about [0.1ν, 10ν], based on table 2. This indicates that we must consider the influence of the mass ratio in the conservative dynamics of EMRIs, because the cycles of typical EMRI waves are ∼ 1/ν in the LISA waveband, and if the relative frequency error reaches ∼ ν, the dephasing will accumulate to a few radians at the end of the evolution. This may cause a failure to detect EMRIs using the test-particle (ν = 0) model.

Figure 2 illustrates the effect of eccentricity on the radial frequency shift (left panel) and the azimuthal frequency shift (right panel) for various values of the semilatus rectum, p, and the spin parameter, a. It shows that the shifts of ω_r and ω_ϕ due to the mass ratio have different effects versus eccentricity. The results for different mass ratios are also plotted (triangles, solid lines, and points represent ν = 10⁻⁶, 10⁻⁴, and 10⁻², respectively), clearly showing that Δω/(ω ν) is not sensitive to the mass ratio, except for the radial frequency when the trajectory approaches the vicinity of the innermost stable bound orbit (ISBO) and the mass ratio becomes 0.01 (blue lines in the left panel).

Zoom In Zoom Out Reset image size

Figure 3 demonstrates the effect of the semilatus rectum, p, on the radial frequency shift (left panel) and the azimuthal frequency shift (right panel) for various values of the spin a and the eccentricity e. When p becomes smaller, the frequency shifts become larger. The sudden growth of the frequency shift is due to the orbit approaching the ISBO. Figure 4 illustrates the effect of the BH's spin on the radial frequency shift (left panel) and the azimuthal frequency shift (right panel) for various values of the semilatus rectum, p, and the eccentricity, e. For the cases of p = 5, when a becomes small, the orbits are very close to the ISBO, and the frequency shifts grow very rapidly.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The semilatus rectum, p, of the ISBO is the separatrix of the bound orbits, and in the test-particle limit, it is given by the analytic expression from equation (24) of [52]

$$\begin{eqnarray}&&{p}_{s}=(6+2a)M\mp 8a\sqrt{\displaystyle \frac{1+e}{2e+6}}+{ \mathcal O }({a}^{2}).\end{eqnarray} \tag{ 31 }$$

However, the above equation is approximate for the ISBO, even for test particles. Considering the influence of the small mass on the background, the ISBOs of EMRIs should deviate from the test-particle ISBO. Figure 5 shows the effect of the mass ratio on the separatrix value p_s (e). We can see that the EOB's ISBO can deviate by as much as 10% from the test-particle limit. In addition, the error of approximate expression (31) becomes large for a rapidly spinning black hole.

Zoom In Zoom Out Reset image size

The analytic 4PN ${ \mathcal O }({e}^{6})$ formula for energy and angular momentum fluxes in Boyer-Lindquist coordinates is given by [50] in terms of the parameter $\upsilon \equiv \sqrt{1/p}$. For convenience, we use p and new expressions, as follows

$$\begin{eqnarray*}\begin{array}{rcl}\langle { \mathcal F }{\rangle }_{4\mathrm{PN}} & = & \displaystyle \frac{32{\nu }^{2}{\left(1-{e}^{2}\right)}^{3/2}}{5{p}^{5}}\left\{1+\displaystyle \frac{73}{24}{e}^{2}+\displaystyle \frac{37}{96}{e}^{4}\right.\\ & & -\displaystyle \frac{1}{p}\left[\displaystyle \frac{1247}{336}+\displaystyle \frac{9181{e}^{2}}{672}-\displaystyle \frac{809{e}^{4}}{128}-\displaystyle \frac{8609{e}^{6}}{5376}\right]\\ & & +\displaystyle \frac{1}{{p}^{3/2}}\left[\pi \left(4+\displaystyle \frac{1375{e}^{2}}{48}+\displaystyle \frac{3935{e}^{4}}{192}+\displaystyle \frac{10007{e}^{6}}{9216}\right)\right.\\ & & -\left.q\left(\displaystyle \frac{73}{12}+\displaystyle \frac{823{e}^{2}}{24}+\displaystyle \frac{949{e}^{4}}{32}+\displaystyle \frac{491{e}^{6}}{192}\right)\right]\\ & & -\displaystyle \frac{1}{{p}^{2}}\left[\displaystyle \frac{44711}{9072}+\displaystyle \frac{172157{e}^{2}}{2592}\right.\\ & & +\displaystyle \frac{2764345{e}^{4}}{24192}-\displaystyle \frac{3743{e}^{6}}{2304}\\ & & -\left.{q}^{2}\left(\displaystyle \frac{33}{16}+\displaystyle \frac{359{e}^{2}}{32}+\displaystyle \frac{1465{e}^{4}}{128}+\displaystyle \frac{883{e}^{6}}{768}\right)\right]\\ & & -\displaystyle \frac{1}{{p}^{5/2}}\left[\pi \left(\displaystyle \frac{8191}{672}+\displaystyle \frac{44531{e}^{2}}{336}\right.\right.\end{array}\end{eqnarray*}$$

$$\begin{eqnarray*}\begin{array}{rcl} & & +\left.\displaystyle \frac{4311389{e}^{4}}{43008}-\displaystyle \frac{15670391{e}^{6}}{387072}\right)\\ & & -\left.q\left(\displaystyle \frac{3749}{336}+\displaystyle \frac{1759{e}^{2}}{56}-\displaystyle \frac{111203{e}^{4}}{1344}-\displaystyle \frac{49685{e}^{6}}{448}\right)\right]\\ & & +\displaystyle \frac{1}{{p}^{3}}\left[\displaystyle \frac{6643739519}{69854400}+\displaystyle \frac{43072561991{e}^{2}}{27941760}\right.\\ & & +\displaystyle \frac{919773569303{e}^{4}}{279417600}+\displaystyle \frac{308822406727{e}^{6}}{186278400}\\ & & +\mathrm{log}(p)\left(\displaystyle \frac{856}{105}+\displaystyle \frac{7276{e}^{2}}{63}+\displaystyle \frac{553297{e}^{4}}{2520}+\displaystyle \frac{187357{e}^{6}}{2520}\right)\\ & & -\gamma \left(\displaystyle \frac{1712}{105}+\displaystyle \frac{14552{e}^{2}}{63}+\displaystyle \frac{553297{e}^{4}}{1260}+\displaystyle \frac{187357{e}^{6}}{1260}\right)\\ & & -\mathrm{log}\left(2\right)\left(\displaystyle \frac{3424}{105}+\displaystyle \frac{13696{e}^{2}}{315}\right.\\ & & +\left.\displaystyle \frac{12295049{e}^{4}}{1260}-\displaystyle \frac{24908851{e}^{6}}{252}\right)\end{array}\end{eqnarray*}$$

$$\begin{eqnarray*}\begin{array}{rcl} & & -\mathrm{log}\left(3\right)\left(\displaystyle \frac{234009{e}^{2}}{560}-\displaystyle \frac{2106081{e}^{4}}{448}+\displaystyle \frac{864819261{e}^{6}}{35840}\right)\\ & & +{\pi }^{2}\left(\displaystyle \frac{16}{3}+\displaystyle \frac{680{e}^{2}}{9}+\displaystyle \frac{5171{e}^{4}}{36}+\displaystyle \frac{1751{e}^{6}}{36}\right)\\ & & -q\pi \left(\displaystyle \frac{169}{6}+\displaystyle \frac{4339{e}^{2}}{16}+\displaystyle \frac{42271{e}^{4}}{96}+\displaystyle \frac{4867907{e}^{6}}{27648}\right)\\ & & +{q}^{2}\left(\displaystyle \frac{3419}{168}+\displaystyle \frac{50271{e}^{2}}{224}+\displaystyle \frac{340141{e}^{4}}{896}+\displaystyle \frac{1013347{e}^{6}}{5376}\right)\\ & & -\left.\mathrm{log}(5)\displaystyle \frac{5224609375{e}^{6}}{193536}\right]\\ & & -\displaystyle \frac{1}{{p}^{7/2}}\left[\pi \left(\displaystyle \frac{16285}{504}+\displaystyle \frac{22798583{e}^{2}}{48384}+\displaystyle \frac{65448785{e}^{4}}{48384}\right.\right.\\ & & +\left.\displaystyle \frac{41758768871{e}^{6}}{41803776}\right)\\ & & -{q}^{2}\pi \left(\displaystyle \frac{65}{8}+\displaystyle \frac{2277{e}^{2}}{32}+\displaystyle \frac{103229{e}^{4}}{768}+\displaystyle \frac{1307875{e}^{6}}{18432}\right)\\ & & -q\left(\displaystyle \frac{83819}{1296}+\displaystyle \frac{12203083{e}^{2}}{18144}\right.\\ & & +\left.\displaystyle \frac{3918011{e}^{4}}{2268}+\displaystyle \frac{1325729{e}^{6}}{1728}\right)\\ & & +\left.{q}^{3}\left(\displaystyle \frac{151}{12}+\displaystyle \frac{6497{e}^{2}}{48}+\displaystyle \frac{14041{e}^{4}}{64}+\displaystyle \frac{12789{e}^{6}}{128}\right)\right]\end{array}\end{eqnarray*}$$

$$\begin{eqnarray*}\begin{array}{rcl} & & -\displaystyle \frac{1}{{p}^{4}}\left[\displaystyle \frac{323105549467}{3178375200}+\displaystyle \frac{13084171033763{e}^{2}}{3178375200}\right.\\ & & +\displaystyle \frac{454079900391707{e}^{4}}{25427001600}+\displaystyle \frac{12759433101997{e}^{6}}{847566720}\\ & & +\mathrm{log}(p)\left(\displaystyle \frac{232597}{8820}+\displaystyle \frac{831307{e}^{2}}{1176}+\displaystyle \frac{25695073{e}^{4}}{11760}\right.\end{array}\end{eqnarray*}$$

$$\begin{eqnarray*}\begin{array}{rcl} & & +\left.\displaystyle \frac{14109647{e}^{6}}{14112}\right)\\ & & -\mathrm{log}(2)\left(\displaystyle \frac{39931}{294}-\displaystyle \frac{2373293{e}^{2}}{1260}+\displaystyle \frac{2273961523{e}^{4}}{17640}\right.\\ & & -\left.\displaystyle \frac{503388564173{e}^{6}}{317520}\right)\\ & & +\mathrm{log}(3)\left(\displaystyle \frac{47385}{1568}-\displaystyle \frac{55105839{e}^{2}}{15680}+\displaystyle \frac{3074548023{e}^{4}}{71680}\right.\\ & & -\left.\displaystyle \frac{107130980133{e}^{6}}{1003520}\right)\\ & & +\mathrm{log}(5)\left(\displaystyle \frac{15869140625{e}^{4}}{903168}-\displaystyle \frac{10089048828125{e}^{6}}{16257024}\right)\\ & & -\gamma \left(\displaystyle \frac{232597}{4410}+\displaystyle \frac{831307{e}^{2}}{588}\right.\\ & & +\left.\displaystyle \frac{25695073{e}^{4}}{5880}+\displaystyle \frac{14109647{e}^{6}}{7056}\right)\end{array}\end{eqnarray*}$$

$$\begin{eqnarray}\begin{array}{rcl} & & +{\pi }^{2}\left(\displaystyle \frac{1369}{126}+\displaystyle \frac{83045{e}^{2}}{252}+\displaystyle \frac{56023{e}^{4}}{56}+\displaystyle \frac{327895{e}^{6}}{1008}\right)\\ & & -q\pi \left(\displaystyle \frac{3883}{168}+\displaystyle \frac{57405{e}^{2}}{224}\right.\\ & & -\left.\displaystyle \frac{15198859{e}^{4}}{43008}-\displaystyle \frac{59177891{e}^{6}}{64512}\right)\\ & & +{q}^{2}\left(\displaystyle \frac{124091}{9072}-\displaystyle \frac{6211679{e}^{2}}{9072}\right.\\ & & -\left.\displaystyle \frac{32107339{e}^{4}}{18144}-\displaystyle \frac{16038995{e}^{6}}{10368}\right)\\ & & -\left.\left.{q}^{4}\left(\displaystyle \frac{17}{16}+\displaystyle \frac{279{e}^{2}}{16}+\displaystyle \frac{5399{e}^{4}}{192}+\displaystyle \frac{2513{e}^{6}}{192}\right)\right]\right\},\end{array}\end{eqnarray} \tag{ 32 }$$

$$\begin{eqnarray*}\begin{array}{rcl}\langle {{ \mathcal G }}^{z}{\rangle }_{4\mathrm{PN}} & = & \displaystyle \frac{32{\nu }^{2}{\left(1-{e}^{2}\right)}^{3/2}}{5{p}^{7/2}}\left\{1+\displaystyle \frac{7}{8}{e}^{2}\right.\\ & & -\displaystyle \frac{1}{p}\left[\displaystyle \frac{1247}{336}+\displaystyle \frac{425{e}^{2}}{336}-\displaystyle \frac{10751{e}^{4}}{2688}\right]\\ & & +\displaystyle \frac{1}{{p}^{3/2}}\left[\pi \left(4+\displaystyle \frac{97{e}^{2}}{8}+\displaystyle \frac{49{e}^{4}}{32}-\displaystyle \frac{49{e}^{6}}{4608}\right)\right.\\ & & -\left.q\left(\displaystyle \frac{61}{12}+\displaystyle \frac{119{e}^{2}}{8}+\displaystyle \frac{183{e}^{4}}{32}\right)\right]\\ & & -\displaystyle \frac{1}{{p}^{2}}\left[\displaystyle \frac{44711}{9072}+\displaystyle \frac{302893{e}^{2}}{6048}\right.\\ & & +\displaystyle \frac{701675{e}^{4}}{24192}-\displaystyle \frac{162661{e}^{6}}{16128}\\ & & -\left.{q}^{2}\left(\displaystyle \frac{33}{16}+\displaystyle \frac{95{e}^{2}}{16}+\displaystyle \frac{311{e}^{4}}{128}\right)\right]\\ & & -\displaystyle \frac{1}{{p}^{5/2}}\left[\pi \left(\displaystyle \frac{8191}{672}+\displaystyle \frac{48361{e}^{2}}{1344}\right.\right.\end{array}\end{eqnarray*}$$

$$\begin{eqnarray*}\begin{array}{rcl} & & -\left.\displaystyle \frac{1657493{e}^{4}}{43008}-\displaystyle \frac{5458969{e}^{6}}{774144}\right)\\ & & -q\left(\displaystyle \frac{417}{56}-\displaystyle \frac{5441{e}^{2}}{672}\right.\\ & & -\left.\left.\displaystyle \frac{1097{e}^{4}}{24}-\displaystyle \frac{153605{e}^{6}}{5376}\right)\right]\\ & & +\displaystyle \frac{1}{{p}^{3}}\left[\displaystyle \frac{6643739519}{69854400}+\displaystyle \frac{6769212511{e}^{2}}{8731800}\right.\\ & & +\displaystyle \frac{4795392143{e}^{4}}{7761600}+\displaystyle \frac{31707715321{e}^{6}}{186278400}\\ & & +\mathrm{log}(p)\left(\displaystyle \frac{856}{105}+\displaystyle \frac{24503{e}^{2}}{420}\right.\\ & & +\left.\displaystyle \frac{11663{e}^{4}}{280}+\displaystyle \frac{2461{e}^{6}}{1120}\right)\\ & & -\gamma \left(\displaystyle \frac{1712}{105}+\displaystyle \frac{24503{e}^{2}}{210}\right.\\ & & +\left.\displaystyle \frac{11663{e}^{4}}{140}+\displaystyle \frac{2461{e}^{6}}{560}\right)\end{array}\end{eqnarray*}$$

$$\begin{eqnarray*}\begin{array}{rcl} & & -\mathrm{log}\left(2\right)\left(\displaystyle \frac{3424}{105}-\displaystyle \frac{1391{e}^{2}}{30}+\displaystyle \frac{418049{e}^{4}}{84}-\displaystyle \frac{94138279{e}^{6}}{2160}\right)\\ & & -\mathrm{log}\left(3\right)\left(\displaystyle \frac{78003{e}^{2}}{280}-\displaystyle \frac{3042117{e}^{4}}{1120}+\displaystyle \frac{42667641{e}^{6}}{3584}\right)\\ & & +{\pi }^{2}\left(\displaystyle \frac{16}{3}+\displaystyle \frac{229{e}^{2}}{6}+\displaystyle \frac{109{e}^{4}}{4}+\displaystyle \frac{23{e}^{6}}{16}\right)\\ & & -q\pi \left(\displaystyle \frac{145}{6}+\displaystyle \frac{409{e}^{2}}{3}+\displaystyle \frac{22631{e}^{4}}{192}+\displaystyle \frac{69887{e}^{6}}{6912}\right)\\ & & +{q}^{2}\left(\displaystyle \frac{799}{56}+\displaystyle \frac{6213{e}^{2}}{56}+\displaystyle \frac{143159{e}^{4}}{1344}+\displaystyle \frac{14447{e}^{6}}{448}\right)\\ & & -\left.\mathrm{log}(5)\displaystyle \frac{1044921875{e}^{6}}{96768}\right]\\ & & -\displaystyle \frac{1}{{p}^{7/2}}\left[\pi \left(\displaystyle \frac{16285}{504}+\displaystyle \frac{396733{e}^{2}}{1008}+\displaystyle \frac{27639271{e}^{4}}{48384}\right.\right.\\ & & -\left.\displaystyle \frac{658044449{e}^{6}}{20901888}\right)\\ & & -{q}^{2}\pi \left(\displaystyle \frac{65}{8}+\displaystyle \frac{705{e}^{2}}{16}+\displaystyle \frac{11531{e}^{4}}{256}+\displaystyle \frac{10525{e}^{6}}{2304}\right)\end{array}\end{eqnarray*}$$

$$\begin{eqnarray*}\begin{array}{rcl} & & -q\left(\displaystyle \frac{90337}{1512}+\displaystyle \frac{6143581{e}^{2}}{18144}+\displaystyle \frac{10846363{e}^{4}}{24192}+\displaystyle \frac{409361{e}^{6}}{48384}\right)\\ & & +\left.{q}^{3}\left(\displaystyle \frac{505}{48}+\displaystyle \frac{1225{e}^{2}}{16}+\displaystyle \frac{8065{e}^{4}}{128}+\displaystyle \frac{881{e}^{6}}{64}\right)\right]\\ & & -\displaystyle \frac{1}{{p}^{4}}\left[\displaystyle \frac{323105549467}{3178375200}+\displaystyle \frac{4079638426121{e}^{2}}{4237833600}\right.\\ & & +\displaystyle \frac{9439621323313{e}^{4}}{12713500800}+\displaystyle \frac{5778193989793{e}^{6}}{50854003200}\end{array}\end{eqnarray*}$$

$$\begin{eqnarray*}\begin{array}{rcl} & & +\mathrm{log}(p)\left(\displaystyle \frac{232597}{8820}+\displaystyle \frac{2548127{e}^{2}}{8820}\right.\\ & & +\left.\displaystyle \frac{8214023{e}^{4}}{35280}-\displaystyle \frac{5314333{e}^{6}}{70560}\right)\\ & & -\mathrm{log}(2)\left(\displaystyle \frac{39931}{294}-\displaystyle \frac{571519{e}^{2}}{294}\right.\\ & & +\left.\displaystyle \frac{228770411{e}^{4}}{3528}-\displaystyle \frac{31677796739{e}^{6}}{45360}\right)\\ & & +\mathrm{log}(3)\left(\displaystyle \frac{47385}{1568}-\displaystyle \frac{8684577{e}^{2}}{3920}\right.\\ & & +\left.\displaystyle \frac{2303677179{e}^{4}}{100352}-\displaystyle \frac{22448097729{e}^{6}}{501760}\right)\\ & & +\mathrm{log}(5)\left(\displaystyle \frac{3173828125{e}^{4}}{301056}-\displaystyle \frac{314655859375{e}^{6}}{1161216}\right)\end{array}\end{eqnarray*}$$

$$\begin{eqnarray}\begin{array}{rcl} & & -\gamma \left(\displaystyle \frac{232597}{4410}+\displaystyle \frac{2548127{e}^{2}}{4410}\right.\\ & & +\left.\displaystyle \frac{8214023{e}^{4}}{17640}-\displaystyle \frac{5314333{e}^{6}}{35280}\right)\\ & & +{\pi }^{2}\left(\displaystyle \frac{1369}{126}+\displaystyle \frac{13759{e}^{2}}{126}+\displaystyle \frac{9901{e}^{4}}{504}-\displaystyle \frac{29635{e}^{6}}{336}\right)\\ & & -q\pi \left(\displaystyle \frac{3883}{168}-\displaystyle \frac{21013{e}^{2}}{336}-\displaystyle \frac{269503{e}^{4}}{448}-\displaystyle \frac{110189735{e}^{6}}{193536}\right)\\ & & +{q}^{2}\left(\displaystyle \frac{5717}{2268}-\displaystyle \frac{10073201{e}^{2}}{18144}\right.\\ & & -\left.\displaystyle \frac{11932475{e}^{4}}{16128}-\displaystyle \frac{40782425{e}^{6}}{96768}\right)\\ & & -\left.\left.{q}^{4}\left(\displaystyle \frac{17}{16}+\displaystyle \frac{197{e}^{2}}{16}+\displaystyle \frac{1209{e}^{4}}{128}+\displaystyle \frac{211{e}^{6}}{128}\right)\right]\right\}\,,\end{array}\end{eqnarray} \tag{ 33 }$$

where q = a/M is the dimensionless spin. The averages of the fluxes of the 1PN order, given in terms of the quantities (e, p) and the mass ratio ν in the EOB gauge without spin, are [33]:

$$\begin{eqnarray}\begin{array}{rcl}\langle { \mathcal F }\rangle & = & \displaystyle \frac{32{\nu }^{2}{\left(1-{e}^{2}\right)}^{3/2}}{5{p}^{5}}\left\{1+\displaystyle \frac{73}{24}{e}^{2}+\displaystyle \frac{37}{96}{e}^{4}\right.\\ & & -\displaystyle \frac{1}{p}\left[\displaystyle \frac{1247}{336}+\displaystyle \frac{5\nu }{4}+{e}^{2}\left(\displaystyle \frac{9181}{672}+\displaystyle \frac{325\nu }{24}\right)\right.\\ & & -\left.\left.{e}^{4}\left(\displaystyle \frac{809}{128}-\displaystyle \frac{435\nu }{32}\right)-{e}^{6}\left(\displaystyle \frac{8609}{5376}-\displaystyle \frac{185\nu }{192}\right)\right]\right\},\end{array}\end{eqnarray} \tag{ 34 }$$

$$\begin{eqnarray}\begin{array}{rcl}\langle {{ \mathcal G }}^{z}\rangle & = & \displaystyle \frac{32{\nu }^{2}{\left(1-{e}^{2}\right)}^{3/2}}{5{p}^{7/2}}\left\{1+\displaystyle \frac{7}{8}{e}^{2}\right.\\ & & -\displaystyle \frac{1}{p}\left[\displaystyle \frac{1247}{336}+\displaystyle \frac{7\nu }{4}+{e}^{2}\left(\displaystyle \frac{425}{336}+\displaystyle \frac{401\nu }{48}\right)\right.\\ & & -\left.\left.{e}^{4}\left(\displaystyle \frac{10751}{2688}-\displaystyle \frac{205\nu }{96}\right)\right]\right\}\,.\end{array}\end{eqnarray} \tag{ 35 }$$

We find that the coefficients of equations (34)–(35) coincide with equations (32)–(33), which are in the harmonic coordinates at the 1PN order. This may hint that the PN fluxes can be used in the EOB equations for the extreme-mass-ratio limit. Therefore, by combining the PN fluxes of a test particle orbiting a Kerr black hole and the 1PN fluxes of a nonspinning binary with a mass ratio, we derive expressions for energy and angular momentum fluxes containing both spin and mass ratio:

$$\begin{eqnarray}\begin{array}{rcl}\langle { \mathcal F }{\rangle }_{2\mathrm{PN}} & = & \displaystyle \frac{32{\nu }^{2}{\left(1-{e}^{2}\right)}^{3/2}}{5{p}^{5}}\left\{1+\displaystyle \frac{73}{24}{e}^{2}+\displaystyle \frac{37}{96}{e}^{4}\right.\\ & & -\displaystyle \frac{1}{p}\left[\displaystyle \frac{1247}{336}+\displaystyle \frac{5\nu }{4}+{e}^{2}\left(\displaystyle \frac{9181}{672}+\displaystyle \frac{325\nu }{24}\right)\right.\\ & & -\left.{e}^{4}\left(\displaystyle \frac{809}{128}-\displaystyle \frac{435\nu }{32}\right)-{e}^{6}\left(\displaystyle \frac{8609}{5376}-\displaystyle \frac{185\nu }{192}\right)\right]\\ & & +\displaystyle \frac{1}{{p}^{3/2}}\left[\pi \left(4+\displaystyle \frac{1375{e}^{2}}{48}\right.\right.\\ & & +\left.\displaystyle \frac{3935{e}^{4}}{192}+\displaystyle \frac{10007{e}^{6}}{9216}\right)\\ & & -\left.q\left(\displaystyle \frac{73}{12}+\displaystyle \frac{823{e}^{2}}{24}+\displaystyle \frac{949{e}^{4}}{32}+\displaystyle \frac{491{e}^{6}}{192}\right)\right]\\ & & -\displaystyle \frac{1}{{p}^{2}}\left[\displaystyle \frac{44711}{9072}+\displaystyle \frac{172157{e}^{2}}{2592}\right.\\ & & +\displaystyle \frac{2764345{e}^{4}}{24192}-\displaystyle \frac{3743{e}^{6}}{2304}\\ & & -\left.\left.{q}^{2}\left(\displaystyle \frac{33}{16}+\displaystyle \frac{359{e}^{2}}{32}+\displaystyle \frac{1465{e}^{4}}{128}+\displaystyle \frac{883{e}^{6}}{768}\right)\right]\right\},\end{array}\end{eqnarray} \tag{ 36 }$$

$$\begin{eqnarray}\begin{array}{rcl}\langle {{ \mathcal G }}^{z}{\rangle }_{2\mathrm{PN}} & = & \displaystyle \frac{32{\nu }^{2}{\left(1-{e}^{2}\right)}^{3/2}}{5{p}^{7/2}}\left\{1+\displaystyle \frac{7}{8}{e}^{2}\right.\\ & & -\displaystyle \frac{1}{p}\left[\displaystyle \frac{1247}{336}+\displaystyle \frac{7\nu }{4}+{e}^{2}\left(\displaystyle \frac{425}{336}+\displaystyle \frac{401\nu }{48}\right)\right.\\ & & -\left.{e}^{4}\left(\displaystyle \frac{10751}{2688}-\displaystyle \frac{205\nu }{96}\right)\right]\\ & & +\displaystyle \frac{1}{{p}^{3/2}}\left[\pi \left(4+\displaystyle \frac{97{e}^{2}}{8}+\displaystyle \frac{49{e}^{4}}{32}-\displaystyle \frac{49{e}^{6}}{4608}\right)\right.\\ & & -\left.q\left(\displaystyle \frac{61}{12}+\displaystyle \frac{119{e}^{2}}{8}+\displaystyle \frac{183{e}^{4}}{32}\right)\right]\\ & & -\displaystyle \frac{1}{{p}^{2}}\left[\displaystyle \frac{44711}{9072}+\displaystyle \frac{302893{e}^{2}}{6048}\right.\\ & & -\displaystyle \frac{701675{e}^{4}}{24192}+\displaystyle \frac{162661{e}^{6}}{16128}\\ & & -\left.\left.{q}^{2}\left(\displaystyle \frac{33}{16}+\displaystyle \frac{95{e}^{2}}{16}+\displaystyle \frac{311{e}^{4}}{128}\right)\right]\right\}.\end{array}\end{eqnarray} \tag{ 37 }$$

Here, for simplicity, we just write the 2PN formalism. Note that we only add the mass-ratio terms into the 1PN. Indeed, there are ν, ν², and ν q terms in the second PN that may also be important for the flux calculation. The impact of the second-order self-forces based on the PN fluxes [53] was comprehensively studied by Isoyama et al [54] for circular orbits with exponential resummation. The next-to-leading tail-induced spin-orbit effects on the fluxes were given in [55]. Compared with these PN fluxes, the above formalism may miss the ν, ν², and ν q terms. In the next step, we will include these missing terms in our fluxes and improve the orbital evolution. In the adiabatic limit, the evolution of energy and angular momentum is driven by the orbit-averaged radiation reaction forces, so that

$$\begin{eqnarray}&&\dot{E}=- \langle { \mathcal F } \rangle ,\quad \dot{{P}_{\phi }}=- \langle {{ \mathcal G }}^{z} \rangle .\end{eqnarray} \tag{ 38 }$$

Using equations (27a ) and (16), we obtain the evolution of e and p due to gravitational radiation

$$\begin{eqnarray}&&\dot{e}=\displaystyle \frac{(E/M)(\partial {\hat{P}}_{\phi }/\partial p)(\langle {{ \mathcal F }}^{z}\rangle /\mu )-(\partial {\hat{H}}_{\mathrm{eff}}/\partial p)(\langle {{ \mathcal G }}^{z}\rangle /\mu )}{(\partial {\hat{H}}_{\mathrm{eff}}/\partial p)(\partial {\hat{P}}_{\phi }/\partial e)-(\partial {\hat{H}}_{\mathrm{eff}}/\partial e)(\partial {\hat{P}}_{\phi }/\partial p)},\end{eqnarray} \tag{ 39a }$$

$$\begin{eqnarray}&&\dot{p}=\displaystyle \frac{(\partial {\hat{H}}_{\mathrm{eff}}/\partial e)(\langle {{ \mathcal G }}^{z}\rangle /\mu )-(E/M)(\partial {\hat{P}}_{\phi }/\partial e)(\langle {{ \mathcal F }}^{z}\rangle /\mu )}{(\partial {\hat{H}}_{\mathrm{eff}}/\partial p)(\partial {\hat{P}}_{\phi }/\partial e)-(\partial {\hat{H}}_{\mathrm{eff}}/\partial e)(\partial {\hat{P}}_{\phi }/\partial p)}.\end{eqnarray} \tag{ 39b }$$

The derivatives are computed from the expressions in equation (22). The evolution of the auxiliary phase of radial motion can now be calculated by

$$\begin{eqnarray}&&\dot{\xi }=\displaystyle \frac{{\left(1+e\cos \xi \right)}^{2}}{{epM}\sin \xi }\dot{r}+\displaystyle \frac{\cot \xi }{e}\dot{e}-\displaystyle \frac{1+e\cos \xi }{{ep}\sin \xi }\dot{p}.\end{eqnarray} \tag{ 40 }$$

Figure 6 demonstrates the evolution of eccentricity and the semilatus rectum with different orbital parameters and mass ratios. We find that the evolution of ν = 10⁻³ deviates from the other two evolutions with smaller mass ratios when the evolution is close to its end. Figure 7 illustrates the evolution of orbits with 2PN (left panels) and 3PN fluxes (right panels). Figure 8 shows the performances of various PN fluxes in orbital evolutions with very extreme spin. Similarly to [51], in the following section, we will use 2PN fluxes to calculate the orbital evolution.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In this subsection, we calculate the waveforms by solving the Teukolsky equations [35] to demonstrate the advantages of our analytical EOB solution. Our method is based on frequency-domain decomposition, and has been developed in previous works [21, 36–38], in which the gravitational waveform from an eccentric EMRI with a total mass M at a distance R , a latitudinal angle Θ, and an azimuthal angle Φ could be written as

$$\begin{eqnarray}&&{h}_{+}-{\mathrm{ih}}_{\times }=\displaystyle \frac{2}{{\rm{R}}}\sum _{\mathrm{lmk}}\displaystyle \frac{{{\rm{Z}}}_{\mathrm{lmk}}^{{\rm{H}}}}{{\omega }_{\mathrm{mk}}^{2}}{{}_{-2}{\rm{S}}}_{\mathrm{lmk}}^{{\rm{a}}\omega }\left({\rm{\Theta }}\right){{\rm{e}}}^{-{\rm{i}}{\phi }_{\mathrm{mk}}+\mathrm{im}{\rm{\Phi }}},\end{eqnarray} \tag{ 41 }$$

where l, m, and k are the harmonic numbers and ϕ_mk ≡ ∫ω_mk (t)dt, ${}_{-2}{S}_{{lmk}}^{a\omega }\left({\rm{\Theta }}\right)$ denotes spin-weighted spheroidal harmonics that depend on the polar angles Θ of the observer's direction of sight and the direction of orbital angular momentum of the source. ${Z}_{{lmk}}^{{\rm{H}}}$ describes the amplitude of each mode, which can be calculated by the radial component of the Teukolsky equation (see appendix C for details). In this article, we set Θ = 0 ('face on'), Φ = 0, and ω_mk is

$$\begin{eqnarray}&&{\omega }_{{mk}}=m{\omega }_{\phi }+k{\omega }_{r},\end{eqnarray} \tag{ 42 }$$

where ω_r and ω_ϕ denote the orbital frequencies of the radial and azimuthal directions, respectively, which are given in equations (29), (30). Due to our analytical solution for orbits and frequencies given in the previous section, the calculation of the Teukolsky-base waveform becomes very convenient and accurate.

As an example, figure 9 illustrates the numerical waveforms of four evolutionary stages of an EMRI with initial parameters of p₀ = 10 M, e₀ = 0.6, and the mass ratio ν = 10⁻³. We find that, as time passes, since the semilatus rectum, p, becomes smaller, the GW frequency becomes larger. At the same time, because of the decrease in the eccentricity e, the GW strains no longer vary strongly when the small object passes through the periastron and the apastron.

Zoom In Zoom Out Reset image size

Matched filtering [56] is widely used in GW detection in the LIGO and Virgo data analyses and will also be used in future spaceborne detectors. We employ this technology to quantitatively analyze the influence of the mass ratio on the EMRI waveforms. Given time series a(t) and b(t), the overlap of the two series is

$$\begin{eqnarray}&&\mathrm{Overlap}=\displaystyle \frac{(a| b)}{\sqrt{(a| a)(b| b)}},\end{eqnarray} \tag{ 43 }$$

where the inner product between a time series signal a(t) and a template b(t) is

$$\begin{eqnarray}&&\left(a| b\right)=2{\int }_{0}^{\infty }\displaystyle \frac{{\tilde{a}}^{* }(f)\tilde{b}(f)+\tilde{a}(f){\tilde{b}}^{* }(f)}{{S}_{n}(f)}{\rm{d}}f.\end{eqnarray} \tag{ 44 }$$

Here, $\tilde{a}(f)$ is the Fourier transform of the time series signal a(t), ${\tilde{a}}^{* }(f)$ is the complexity conjugate of $\tilde{a}(f)$, and S_n (f) is the power spectral density of the gravitational wave detectors' noise. Throughout this paper, the power spectral density is taken to be the LISA noise.

The standard fitting factor should be maximized over both time (t_s ) and the phase shift (ϕ_s ),

$$\begin{eqnarray}&&\mathrm{FF}(a,b)=\mathop{{\bf{\max }}}\limits_{{t}_{s},{\phi }_{s}}\displaystyle \frac{(a(t)| b(t+{t}_{s}){{\rm{e}}}^{{\rm{i}}{\phi }_{s}})}{\sqrt{(a| a)(b| b)}}.\end{eqnarray} \tag{ 45 }$$

We now use equation (45) to calculate the 'match' between the two waveforms. Figure 10 illustrates the fitting factors of the waveforms with the mass-ratio correction, in both conservative dynamics and fluxes, and the waveforms generated by the test-particle (ν = 0) limit with the same initial orbital parameters. It should be noted that, here, the test-particle waveform does not mean that we did not consider the radiation reaction of GWs. It just means that in the motion equation, ν → 0, but the radiation reaction is still included. We employ six EMRIs with mass ratios ranging from 10⁻⁶ to 10⁻³. The results show that for a mass ratio as low as 10⁻⁶, the test-particle model with adiabatic orbital evolution may be still accurate in waveform template calculations. However, for a mass ratio ∼ 10⁻⁵, if the small body is in an extremely relativistic orbit around the central BH, the waveform templates of the test-particle approximation will at least lead to relative errors of the order of the mass ratio in the parameter estimations after a few months of evolution. For a mass ratio of ∼ 10⁻³, even for a large orbital separation, the match of the two waveforms rapidly degrades after just two months. We may conclude that the mass-ratio correction in the EMRI waveform model should be important in the detection of this kind of system by LISA. It is worth pointing out that variations of the initial p, e, and other parameters may still match the true waveform.

Zoom In Zoom Out Reset image size

The method described in this paper could be useful in the development of an efficient waveform model for mock data analysis for spaceborne projects. Though we have omitted the nonlinear terms in the EOB Hamiltonian here, in extreme mass-ratio cases, we still get the orbital frequency and then the frequencies of the GWs with sufficient accuracy. The frequency is a very sensitive detector of EMRIs. Once we include the effective spin terms, we can expect to extend this model to comparable binaries.

However, the evolution of the orbital parameters in this work may still not be enough. This is due to our use of the PN formalism for GW fluxes. From figure 8, one can see that higher PN orders do not give more convergent results. Of course, one could numerically calculate the energy fluxes accurately using the Teukolsky equation, but this would be computationally expensive. Some researchers have developed analytical or fitting formalisms for the Teukolsky-based fluxes with very high accuracies, but only for circular orbits in the Schwarzschild [57] and Kerr [58] cases. A recent study developed a model that evolved the EMRI orbits precisely and rapidly using the gravitational self-force method [59]. In our next work, we may try to employ this fast orbital evolution.

In addition, a waveform model for the mock data analysis should be efficient, i.e., it should generate the waveforms quickly. For a typical EMRI with a 10⁶M_⊙ SMBH and a mass ratio of 10⁻⁵, our model takes about 522 seconds to generate one-year waveforms using a single CPU. The majority of the CPU time is used to solve the Teukolsky equation and to sum the harmonic number k from −1 to 8. The well-known AK and AAK templates take 229 and 209 seconds, respectively, but the NK template consumes much more CPU time. The XSPEG [17], an EMRI template that was also developed by some of us, needs 37.5 minutes to finish the evolution.

There is still some room to speed up our numerical codes. The orbital evolution is very fast, because we just need to integrate the two ordinary differential equations for p and e. However, the numerical calculation of the Teukolsky equation (even though we use a semi-analytical algorithm [60]) requires more CPU time. Fortunately, although the frequency needs to be updated at every step, the amplitude of the GW does not. If we had a fitting formula for accurate Teukolsky-base waveforms, the computation time would be greatly reduced. As a conclusion, our model may be a potential candidate for the mock data analysis.