Transition from Inspiral to Plunge: A Complete Near-Extremal Trajectory and Associated Waveform
Ollie Burke Jonathan R. Gair Joan Simón


APS/123-QED





[]ollie.burke@aei.mpg.de

Max Planck Institute for Gravitational Physics (Albert Einstein Institute), Am Mühlenberg 1, Potsdam-Golm 14476, Germany

  School of Mathematics, University of Edinburgh, James Clerk Maxwell Building, Peter Guthrie Tait Road, Edinburgh EH9 3FD, UK



Max Planck Institute for Gravitational Physics (Albert Einstein Institute), Am Mühlenberg 1, Potsdam-Golm 14476, Germany

  School of Mathematics, University of Edinburgh, James Clerk Maxwell Building, Peter Guthrie Tait Road, Edinburgh EH9 3FD, UK



  School of Mathematics, University of Edinburgh, James Clerk Maxwell Building, Peter Guthrie Tait Road, Edinburgh EH9 3FD, UK




We extend the Ori and Thorne (OT) procedure to compute the transition from an adiabatic inspiral into a geodesic plunge for any spin, with emphasis on near-extremal ones. Our analysis revisits the validity of the approximations made in OT. In particular, we discuss possible effects coming from eccentricity and non-geodesic past-history of the orbital evolution. We find three different scaling regimes according to whether the mass ratio is much smaller, of the same order or much larger than the near extremal parameter describing how fast the primary black hole rotates. Eccentricity and non-geodesic past-history corrections are always sub-leading, indicating that the quasi-circular approximation applies throughout the transition regime. However, we show that the OT assumption that the energy and angular momentum evolve linearly with proper time must be modified in the near-extremal regime. Using our transition equations, we describe an algorithm to compute the full worldline in proper time for an extreme mass ratio inspiral (EMRI) and the resultant gravitational waveform in the high spin limit.












Valid PACS appear here
                             

                              
Transition from Inspiral to Plunge:
 A Complete Near-Extremal Trajectory and Associated Waveform
    Joan Simón
    
================================================================================================






§ INTRODUCTION

The LIGO observation of the transient gravitational wave (GW) signal from the collision of two stellar mass black holes <cit.> in September 2015 spectacularly opened the new field of gravitational wave astronomy. By the end of the O2 observing run in August 2017, the LIGO/Virgo detectors had observed ten binary black hole mergers and a single binary neutron star inspiral <cit.>. This handful of observations has already had a profound impact on our understanding of the astrophysics of compact objects and ruled out a number of modified theories of gravity <cit.>. During the ongoing O3 observing run new events are being reported at the rate of one per week, so these constraints are rapidly improving. However, the masses of the objects being observed are all in the range of 1–100M_⊙, which is determined by the frequency sensitivity of the instruments <cit.>. Black holes with much higher masses are expected to exist in the centres of most galaxies <cit.> and will be even stronger sources of GWs, but these waves will be at millihertz frequencies which are inaccessible to ground-based detectors due to the seismic noise background.

The launch of the Laser Interferometer Space Antennae (LISA) <cit.>, scheduled for 2034, will open the millihertz band from 10^-4–10^-1Hz for the first time. Expected sources in this frequency band include massive black hole binaries, cosmic strings and extreme mass ratio inspirals (EMRIs). Detection of these sources, and estimation of their parameters, will rely on the comparison of accurate theoretical models of the expected gravitational waveforms to the observed data. Building these models for LISA is extremely challenging, in particular for EMRIs, which are expected to have a very rich structure and to be observed for hundreds of thousands of waveform cycles prior to merger with the central object <cit.>. In this paper we focus on modelling of a particular class of EMRIs, in which the central black hole has very large angular momentum (spin). All of the LIGO observations to date are consistent with zero or small spin <cit.>, but the massive black holes that will be probed by LISA are a different population. These black holes are observed in high accretion states as quasars, and accretion tends to spin the black holes up. Semi-analytic models predict that the typical spins of these objects are a ≳ 0.95 <cit.>. 

The maximum spin of massive black holes is a quantity of fundamental interest for understanding the origin of black holes in the Universe. It was shown by Thorne <cit.> that the angular momentum of black holes being spun up through thin disc accretion saturates at a limit of a=0.998 where an equilibrium is reached between spin up by accreted material and spin down by captured retrograde photons. Black holes with higher spin could in principle be formed directly in the early Universe and for sufficiently high mass these black holes can retain spins above the Thorne limit for a Hubble time <cit.>. A direct observation of a system with spin above the Thorne limit would thus have profound implications for our understanding of the origin and growth of black holes. It is therefore important to understand how well observations of EMRIs can constrain the spin of near-extremal black holes and to determine this we first need to build accurate representations of the gravitational waves emitted by such systems.

The near extremal limit is also relevant for more theoretical considerations. Indeed, as the primary rotates faster, its Hawking's temperature decreases because the distance between the inner (_-) and outer (_+) horizons in Boyer-Lindquist (BL) coordinates reduces according to

    r̃_± = 1 + √(1 - a^2) = 1 + ϵ ,

where r̃ = r/M and a the dimensionless Kerr spin parameter. The existence of a double pole in the function determining the black hole horizons in this limit is responsible for an enhancement of symmetry in the near horizon geometry of the Kerr black hole <cit.>, a feature that remains true for any extremal black hole <cit.>. This enhancement of symmetry from time translations to the conformal group has allowed several groups to analytically solve the master Teukolsky equation in the presence of the in spiraling probe particle leading to an analytic expression for the energy fluxes carried by the gravitational waves generated by this source 
<cit.>. This provides a very exciting opportunity where analytic tools developed in the high energy theoretical physics community can provide accurate predictions to generate gravitational waveform templates. Future observations using such templates will be directly testing these theoretical predictions.



It has already been shown that gravitational waveforms emitted by these sources contain unique qualitative features that provide a smoking gun for the existence of near-extremal systems <cit.>. The amplitude of an EMRI waveform (averaged over a suitable amount of orbits) typically increases linearly in time for moderate spin a ≈ 0.9. It was shown in <cit.> that the amplitudes of these signals dampen in the high spin limit due to behaviour of the flux close to the horizon. There has been progress in modelling the inspiral from radial infinity to the innermost stable circular orbit (ISCO) <cit.>,
by integrating the geodesic equations in the near horizon geometry of the Kerr black hole <cit.> and exploiting the enhanced set of symmetries to compute the energy fluxes for more source trajectories <cit.>. However, no one has focused on providing a model which encapsulates the inspiral and plunge in the limit of high spins. This is precisely what this paper seeks to do.

In this work, we build such a model for an EMRI comprised of a small compact object of mass μ gravitationally bound to a supermassive Kerr black hole of mass M and study the transition from an adiabatic inspiral into a geodesic plunge for any spin of the primary black hole. This transition to plunge was originally discussed by Ori and Thorne (OT) <cit.> for moderate values of the spin in the limit of small mass ratios. A similar but independent analysis conducted by Buonanno and Damour in <cit.> solved the problem for Schwarzschild black holes with arbitrary (reduced) mass ratio. The technical reason why high spins require a separate discussion is because of the existence of a second independent small parameter competing with the mass ratio η=μ/M≪ 1. This new parameter is the near-extremal parameter ϵ = √(1-a^2) encoding the distance of the spin parameter a from its upper/lower bound, since Kerr black holes have spin parameters a ∈ [-1,+1]. Since the dynamical equations describing the transition depend on the spin, the near extremal limit, i.e. ϵ→ 0, modifies the original scaling discussed by OT. The transition to plunge for near-extremal EMRIs was previously considered in <cit.> and our work clarifies and extends those results in a number of ways. We point out the physical interpretation of the mathematical procedure used in that paper, identify a missing term in the near-extremal regime and incorporate recent analytic results for the near-extremal energy flux for the first time. 

In this paper, we will first review the treatment of the transition regime given by OT in <cit.>. We analyse their methodology and approximations and carefully estimate the scaling of terms that are being omitted. In each of <cit.> the notion of eccentricities and non-circular motion was ignored. We discuss the potential growth of eccentricities before and during the transition regime and find that corrections to our equations due to eccentric motion are sub-leading for any spin. We identify three separate transition regimes, each with a slightly different equation of motion: η≪ϵ, ϵ∼η and ϵ≪η.  
We then discuss a numerical algorithm to generate full inspiral trajectories in Boyer-Lindquist coordinates, alongside the corresponding evolution of the integrals of motion E(τ) and L(τ). Finally, we extend the waveform from the inspiral only results of <cit.> to include the plunge in the regime ϵ∼η.

This paper is organised as follows. In section <ref>, we review the properties of equatorial and circular orbits in the Kerr black hole and, in section  <ref>, we review and compare the results describing gravitational fluxes emitted by circular EMRIs as a function of the spin.  In section <ref> we set-up the master transition equation of motion in general and, in subsection <ref>, we estimate corrections due to eccentricity and non-geodesic past-history of the orbital evolution. The transition equations of motion in the three different sxcaling regimes are described in subsections <ref>, <ref> and <ref> respectively. The numerical scheme to integrate our transition equations of motion for the ϵ∼η regime is presented in section <ref>. We describe how to generate a near-extremal EMRI gravitational waveform encapsulating inspiral and plunge in subsections <ref> and <ref>. We finish with a summary of our main results in section <ref>. 

Notation: Any quantity carrying a tilde refers to a dimensionless quantity in units of the primary mass M, i.e.  r̃ = r/M, τ̃ = τ/M, t̃ = t/M, Ẽ = E/μ and L̃ = L/Mμ, but we keep a as the dimensionless Kerr spin parameter with no tilde.  Dotted quantities (eg Ẽ̇) denote coordinate time derivatives of that quantity. Finally, expressions A ∼𝒪(B) or, for brevity, A ∼ B stress that both A and B scale in the same way with the small parameters under consideration. We impose geometrized units by setting the constants G = c = 1. 



Note added. In the final stages of this work, we became aware of overlapping results that were independently obtained in <cit.>.



§ PRELIMINARIES




In Boyer-Lindquist (BL) coordinates (,ϕ,θ,t̃), the motion of a point particle with mass μ in a Kerr black hole on the equatorial plane  (θ = π/2) is given by <cit.>

    (dr̃/dτ̃)^2   = [Ẽ(r̃^2 + a^2) - aL̃]^2 - Δ [(L̃ - aẼ)^2 + r̃^2]/r̃^4
       = Ẽ^2 - V_eff(r̃,Ẽ,L̃,a) =  G(r̃,Ẽ,L̃,a) 
    dϕ/dτ̃   = -(aẼ - L̃) + a(Ẽ[r̃^2 + a^2] - aL̃)/r̃/r̃^2
       = Φ(r̃,Ẽ,L̃,a) 
    dt̃/dτ̃    = -Δ a(aẼ -L̃) + (r̃^2 + a^2)(Ẽ[r̃^2 + a^2] - aL̃)/Δr̃^2
        = T(r̃,Ẽ,L̃,a),

where the largest root of Δ = r̃^2 - 2r̃ + a^2 corresponds to the outer horizon _+

    r̃_+ = 1 + √(1-a^2) ,

=τ/M denotes proper time in units of the Kerr black hole mass M and a is the dimensionless spin parameter a ∈ [-1,1]. 

The particle is on a prograde (retrograde) orbit if it follows the same (opposite) direction as the rotation of the primary hole. Prograde (retrograde) orbits correspond to a > 0  (a < 0) while keeping the azimuthal component of the angular momentum >0. Since retrograde orbits do not reach the near horizon geometry of the primary black hole in the near-extremal limit (see appendix <ref>), while prograde orbits do, we only consider the latter from here on.

For an equatorial orbit to be circular,  the BL radial coordinate r̃ must be constant and to be stable, the latter must be at a minimum of the potential V_eff in (<ref>) so that 

    G = ∂G/∂r̃= 0,  and ∂^2G/∂r̃^2≥ 0 .

These conditions determine the energy Ẽ and angular momentum L̃ of these orbits to be <cit.> 

    Ẽ   = 1 - 2/r̃ + a/r̃^3/2/√(1 - 3/r̃ + 2a/r̃^3/2) , 
    L̃   = r̃^1/21 - 2a/r̃^3/2 + a^2/r̃^2/√(1 - 3/r̃ + 2a/r^3/2) .


Substituting (<ref>)-(<ref>) into (<ref>)-(<ref>) gives rise to

    dϕ/dτ̃   =   1/r̃^3/2√(1 - 3/r̃ + 2a/r̃^3/2) , 
    dt̃/dτ̃   =   1 + a/r̃^3/2/√(1 - 3/r̃ + 2a/r̃^3/2) ,

whose ratio defines the angular velocity Ω̃ of the particle

    dϕ/dt̃ =Ω̃ = (r̃^3/2 + a)^-1 .

  
Equatorial circular orbits are also known to satisfy the identity <cit.>

    ∂ G/∂() Ω̃() + ∂ G/∂() = 0 ,
 
where we stress the equality holds for any circular orbit labelled by (, , ). Differentiating (<ref>) with respect to , we can derive further equalities satisfied for any such orbits. The ones below

    ∂^2G/∂∂ Ω̃ + ∂^2G/∂∂   =    -∂Ω̃/∂∂ G/∂ , 
    
    
      -1/2(∂^3G/∂^2∂ Ω̃ + ∂^3G/∂^2∂)    =   ∂Ω̃/∂∂^2 G/∂∂
           + 1/2∂^2Ω̃/∂^2∂ G/∂

will play a role in our analysis later on.

The innermost stable circular orbit (ISCO) is the marginal circular stable orbit satisfying

    G|_ = ∂G/∂r̃|_= ∂^2G/∂r̃^2|_ = 0 .

The last equality, describing marginality, allows to solve for its radius as a function of the spin <cit.>

    r̃_   = 3 + Z_2 - [(3-Z_1)(3+Z_1+2Z_2)]^1/2
    
    Z_1   = 1 + (1-a^2)^1/3[(1+a)^1/3 + (1-a)^1/3]
    
    Z_2   = (3a^2 + Z_1^2)^1/2 .

This is the last radii before plunging into the horizon occurs. In appendix <ref>, we derive general formulas for (<ref>) and higher order derivatives of G(,,), which are valid for any spin a>0, when evaluated at ISCO that will be relevant in the rest of this work. 

For near-extremal Kerr black holes, it is natural to introduce the near extremal parameter

    ϵ = √(1-a^2)  for ϵ≪ 1 ,

to mathematically capture the large spin limit a→ 1. For these black holes, the ISCO location can be expanded in ϵ

    r̃_ = 1 + 2^1/3ϵ^2/3 + 𝒪(ϵ^4/3)

and the physical parameters of this marginal orbit reduce to

    Ẽ_   →1/√(3)(1 + 2^1/3ϵ^2/3)  + 𝒪(ϵ^4/3) , 
    L̃_   →2/√(3)(1 + 2^1/3ϵ^2/3)  + 𝒪(ϵ^4/3) , 
    Ω̃_   →1/2(1 -3/ 2^5/3ϵ ^2/3) + 𝒪(ϵ^4/3)  .

Notice |_ - r̃_+|∼𝒪(ϵ^2/3) for prograde orbits, whereas it is 𝒪(1) for retrograde ones, as mentioned below (<ref>). This further justifies our interest in prograde orbits in the near extremal limit.





 §.§ Gravitational Wave Flux


For equatorial orbits, the motion of a (point) particle on the Kerr spacetime background generates gravitational waves carrying energy and angular momentum, either escaping towards infinity or being absorbed by the horizon of the primary hole. 

Due to energy and angular momentum conservation, the orbit averaged rates of change ⟨⟩ and ⟨⟩ satisfy ⟨⟩ = -⟨_GW⟩ and ⟨⟩ = - ⟨_GW⟩, where the averaged gravitational wave dissipative fluxes are 

    ⟨Ẽ̇_GW⟩   = ⟨Ẽ̇_GW,H⟩ + ⟨Ẽ̇_GW,∞⟩ = ⟨ f^diss_t / u^t⟩ , 
    ⟨L̇̃̇_GW⟩   = ⟨L̇̃̇_GW,H⟩ +⟨L̇̃̇_GW,∞⟩ = -⟨ f^diss_ϕ / u^t⟩ ,

determined in terms of the time-averaged t and ϕ components of the dissipative self-force f̃^diss normalised by the t-component of the four velocity u^t. These quantities are discussed in more depth in the next subsection <ref> (see appendix <ref>, in particular the discussion around (<ref>), for a derivation of relations such as (<ref>)). Notice we split these fluxes into their horizon and asymptotic infinity contributions in the right hand side. We stress here this flux balance law only holds for adiabatically evolving binaries forcing the small mass ratio limit η→ 0. The (orbit-averaged and dissipative) gravitational wave fluxes are determined by solving the Teukolsky equation in the presence of the point particle source <cit.>.  From hereon, to avoid cumbersome notation, we shall drop the angled-brackets to denote time-averaging and simply write, for example, ⟨Ẽ̇⟩ = Ẽ̇. 



In <cit.>, Finn and Thorne (F&T) parameterise the energy flux as the (Peters and Mathews <cit.>) leading order post newtonian correction with an extra general relativistic correction ℰ̇ factor

    dẼ_GW/dt̃= 32/5η Ω̃^10/3ℰ̇(r̃) .

These fluxes are spin dependent and are typically computed through numerical means. See the tables in <cit.> for some of the values of these relativistic corrections.

As we increase the spin of the black hole, the two roots _± = 1 ±ϵ of the function Δ determining the outer and inner horizons of the rotating black hole coincide in the extremal limit ϵ=0. In this limit, the geometry close to the horizon of the black hole, which can be isolated using the change of coordinates

    - _+ = λ ρ ,     t̃ = T/λ ,     ϕ̃ = ϕ + t̃/2λ

with λ→ 0, has an enhancement of symmetry from ℝ×𝕌(1), i.e. time translations and rotational symmetry, to 𝕊𝕃(2,ℝ) ×𝕌(1). The resulting near horizon geometry is warped AdS_2 over a 2-sphere. The enhanced 𝕊𝕃(2,ℝ), the isometry group of AdS_2, includes the scaling symmetry ρ→ cρ and T→ T/c. This was already observed in the original work <cit.> and it is true for any extremal black hole <cit.>.

Larger symmetry in physics implies larger kinematic constraints which can provide further analytic control over the given problem, in this case, the calculation of the gravitational wave fluxes (<ref>). It is precisely the emergence of this conformal group (𝕊𝕃(2,ℝ)) and the use of asymptotic expansion matching methods that allowed to find analytic expressions for these energy and angular momentum fluxes for equatorial circular orbits close to the horizon <cit.>. This body of work led to the simple relationship for the flux given in <cit.>

    dẼ_GW/dt̃≈η(C̃_H + C̃_∞)r̃ - r̃_+/r̃_+ .

The quantities C̃_H and C̃_∞ are constants representing how much wave emission goes towards the horizon and infinity respectively. These constants are given analytically in equations (76) and (77) of <cit.>. Numerically evaluating them and summing the contribution of the first |m| ≤ l = 30 modes gives C̃_H = 0.987 and C̃_∞ = -0.133.

The flux (<ref>) is only reliable when working with near extremal black holes and interested in near horizon physics.  This fact can be checked by comparing the exact fluxes (<ref>), using exact results found in the black hole perturbation toolkit (http://bhptoolkit.orgBHPT), with the near extremal approximation (<ref>). This comparison is shown in figure <ref>. Fixing the radial coordinate to
r̃=r̃_ and varying the spin parameter a, we observe in Table (<ref>) that as a → 1, the NHEK flux (<ref>) converges towards the exact value computed using the http://bhptoolkit.orgBHPT. Furthermore, fixing the spin parameter to a = 1-10^-9, as in figure (<ref>), the NHEK flux (<ref>) provides a nearly-perfect agreement up to a coordinate radii r̃≈ 1.012. The reason for the (extremely small) discrepancy at the ISCO is because Eq.(<ref>) is only valid for ϵ→ 0 and we consider ϵ≈ 10^-5.  


Thus we can use (<ref>) to build a trajectory throughout the adiabatic inspiral regime. Then, as we near the ISCO, we can use the powerful analytic result given by Eq.(<ref>). Using Eq.(<ref>) allows for a more analytic treatment of the analysis of the transition regime. 









§ THE TRANSITION EQUATION OF MOTION



In this section we revisit the earlier work by OT <cit.> describing how a small body following an initial equatorial circular orbit around the large black hole inspirals and eventually transitions into a plunging trajectory falling into the black hole. 

Our discussion is organised as follows. First, we analyse in section <ref> the effects arising from the radial self-force in the vicinity of the ISCO on the dynamics of this small body, justifying the starting point in OT. Second, assuming the dissipative fluxes of energy and angular momentum for quasi-circular and equatorial orbits are still related as in circular orbits <cit.> 

    Ẽ̇ = Ω̃(r̃)L̇̃̇ ,

we derive in section <ref> the transition equation for arbitrary black hole spins without the OT assumption that both energy Ẽ and angular momenta L̃ evolve linearly in proper time τ̃. Third, given the quasi-circular nature of our assumed orbits, we argue in section <ref> there can be corrections to (<ref>) of the form

    Ẽ̇ - Ω̃(r̃)L̇̃̇∼ηṙ̃̇ + η e^2 ,

whose scaling behaviour on the trajectory of the small body is determined. Finally, in sections <ref>-<ref>, we discuss in great detail the existence of three different scaling regimes in our transition equation, depending on the black hole spin, paying special attention to the near-extremal ones which contain new physics. We show the corrections due to (<ref>) are subleading in all the regimes.



 §.§ The self-force




 This subsection shows both that quasi-circular and equatorial orbits have vanishing dissipative effects and the conservative piece of the radial self-force can be neglected close to the ISCO. 

Consider the radial geodesic equation, Eq. (<ref>). Differentiating it with respect to proper time, one obtains

    d^2r̃/dτ̃^2 - 1/2∂ G/∂r̃ = 1/2(d/dτ̃∂ G/∂ + d/dτ̃∂ G/∂L̃)(dr̃/dτ̃)^-1 .

It is shown in Appendix <ref> this is equivalent to

    d^2/dτ̃^2 + Γ^_ρσdx̃^ρ/dτ̃dx̃^σ/dτ̃ = f̃^_diss .

Hence, the terms on the left hand side of eq. (<ref>) correspond to the usual ones for geodesic motion, whereas the one on the right hand side can be understood as the perturbing force f̃^r̃_diss exerted on the particle driving energy (dẼ/dτ̃) and angular momentum (dL̃/dτ̃) loss due to gravitational wave emission. 











In <cit.>, Mino recognised that the forcing term for general geodesic motion can perturbatively be split into a radiative reactive dissipative and a conservative piece at first order in the mass ratio η

    f̃^r̃ = η(f̃^r̃_(1)diss + f̃^_(1)cons) + 𝒪(η^2) .

More details on this splitting can be found in <cit.>.

For circular orbits, as considered by OT, the dissipative fluxes of energy and angular momentum are related by
<cit.> 

    Ẽ̇ = Ω̃(r̃)L̇̃̇ .

Hence, the dissipative part of the self-force f̃^r̃_(1)diss vanishes, leading to

    d^2r̃/dτ̃^2 - 1/2∂ G/∂r̃ = ηf̃^r̃_(1)cons + 𝒪(η^2)

which is precisely Eq.(3.10) of OT in <cit.>. A gauge invariant way to quantify the leading order in η effect on the trajectory due to f̃^r̃_(1)cons is to study the orbital velocity Ω̃ (see <cit.> for a review). This generates a shifted orbital velocity Ω̃^shifted_ with respect to the Kerr orbital velocity Ω̃_ at ISCO given by <cit.>

    (1 + η)Ω̃^shifted_ = Ω̃_(1 + η C_Ω̃(a)) + 𝒪(η^2)

with the quantity C_Ω̃(a) discussed in depth and independently (numerically) calculated in both <cit.>. According to <cit.>, C_Ω̃(a)∈(1.24,1.39). Hence, C_Ω̃(a)∼𝒪(1) for all spins, and since Ω̃_∼𝒪(1), it follows that for η≪ 1

    |Ω̃^shifted_ - Ω̃_| ≈η (C_Ω̃(a) - 1)Ω̃_∼η.

It is further shown in <cit.> that  C_Ω̃(a) → 1 + 1/2√(3) as a → 1 in an averaged sense[C_Ω̃(a) is shown to actually oscillate around this limiting value as a→ 1. This phenomenon is non-trivial and still not well understood today. See <cit.> for more details.]. Using Eq.(<ref>) in the high spin limit and η≪ 1, equation (<ref>) becomes

    |Ω̃^shifted_ - Ω̃_| = η/4√(3) + 𝒪(η^2) + 𝒪(ηϵ^2/3).

This implies that the change in the orbital velocity at the ISCO due to conservative self-force effects is an 𝒪(η) quantity. 

Since Ω̃^shifted_ is related to the shifted Boyer-Lindquist radial coordinate at the ISCO by

    r̃^shifted_ = (1/Ω̃^shifted_ - a)^2/3 ,

it follows, using eq.(<ref>), that the “radial thickness"

    r̃ - r̃^shifted_∼r̃ - r̃_ + 𝒪(η) ,

differs by an 𝒪(η) quantity when including the conservative self-force effects.

It will be shown in this paper that there are three different transition regimes depending on the ratio of ϵ and η. The “radial thickness" of the transition in each regime scales according to

    
  * For η≪ϵ⇒r̃ - r̃_∼η^2/5,
    
  * For η∼ϵ⇒r̃ - r̃_∼ (η/ϵ)^2/5ϵ^2/3∼η^2/3,
    
  * For η≫ϵ⇒r̃ - r̃_∼η^2/3.

Thus, the effect of the conservative piece of the self-force is subleading in all regimes. For this reason, like in the original OT analysis, we shall ignore these effects.




 §.§ Transition Equation - Generalities
 


To discuss the evolution of the orbit, we pursue the following strategy : we assume the corrections in (<ref>) are subleading, and once the scaling behaviour of the different dynamical regimes is identified, we double check the consistency of our original assumption.

To evolve the orbit, OT used the circular flux relationship (<ref>) and additionally assumed that the energy Ẽ and angular momentum L̃ evolve linearly in proper time τ̃ throughout the transition regime

    Ẽ - _   =  Ω̃_dL̃_/dτ̃(τ̃ - τ̃_) , 
    L̃ - _   =  dL̃_/dτ̃(τ̃ - τ̃_).



In our analysis of the transition, we will not assume a strict equality in Eq.(<ref>). Instead, we will keep track of the evolution of Ω̃_ ^-1(Ẽ - _) - (L̃ - _), as also considered in <cit.>. 

OT proposed to analyse the transition to the plunging geodesic by expanding (<ref>) around the 
ISCO trajectory (_,_,_), since the latter provides the natural starting point for the plunging trajectory for equatorial and circular orbits. It is physically natural to introduce the new variables δ E, δ L and R

    δ E    = Ω̃_^-1(Ẽ - _)   
    δ L    = L̃ - _
    
    R    = r̃ - _


to study the inspiral evolution of the small body perturbatively around the primary. The presence of Ω̃_ is for technical convenience.

Instead of expanding (<ref>), we find it more convenient to expand (<ref>). Our conclusions do not depend on this choice. The latter is given by


    (d/dτ̃)^2 = G(r̃_,_,_)  + ∑_i=1^∞1/i!∂^i G/∂r̃^i|_(r̃ - _)^i
    
    + ∑_i=0^∞1/i!(∂^i+1 G/∂r̃^i∂|_ ( - _) + ∂^i+1 G/∂r̃^i∂L̃|_ (L̃ - _))  (r̃ - _)^i
    1/2∑_i=0^∞1/i!(∂^i+2 G/∂r̃^i∂^2|_ ( - _)^2 + 2∂^i+2 G/∂r̃^i∂L̃∂Ẽ|_ (Ẽ - _)(L̃ - _) + ∂^i+2 G/∂r̃^i∂^2|_ (L̃ - _)^2)  (r̃ - _)^i .

Since G(,,) is quadratic in  and , we have ignored the derivatives

    ∂^nG/∂Ẽ^n = ∂^nG/∂L̃^n-k∂Ẽ^k = ∂^nG/∂Ẽ^n-k∂L̃^k = ∂^nG/∂^n = 0   for n≥ 3  and k<n.

Plugging the perturbative variables (<ref>), using the definition of the coefficients (<ref>) and the results in (<ref>)-(<ref>),
one can rewrite the general transition equation as

    (dR/dτ̃)^2 =  ∑_n=3^∞1/n! A_n R^n + δ L 
    ∑_n=1^∞1/n! B_n R^n + δ L^2/2∑_n=0^∞1/n! C_n R^n  
    + Γ_⊙ ,

where Γ_⊙ is defined by

    Γ_⊙   = 1/2∑_n=0^∞1/n!Ω̃_(δ E - δ L)(2∂^n+1G/∂r̃^n∂|_ + 2(∂^n+2G/∂r̃^n∂Ẽ∂L̃|_+ Ω̃_∂^n+2G/∂r̃^n∂Ẽ^2|_) δ L . 
       .+  Ω̃_∂^n+2G/∂r̃^n∂Ẽ^2|_( δ E - δ L) )R^n .


Notice that Γ_⊙∝δ E - δ L at leading order in R. Hence, it encodes the deviations from the OT approximation (<ref>).

The time evolution of δ E - δ L near _ is controlled by the fluxes and the angular velocity. Throughout a quasi-circular inspiral far from ISCO, the compact object inspirals on a sequence of circular geodesics defined by the constants of motion Ẽ(r̃_circ) = Ẽ_circ and L̃(r̃_circ) = L̃_circ, as given in Eq (<ref>) and Eq.(<ref>) respectively. The evolution of the constants of motion is linked through Eq. (<ref>) above, which simply states that circular geodesics evolve into circular geodesics. It can be shown that solutions to the Teukolsky equation for circular orbits obey this condition <cit.>. For circular evolutions we therefore see that

    d/d(δ E - δ L)    = Ω̃_^-1dẼ/dτ̃ - dL̃/dτ̃
       = (Ω̃_^-1Ω̃(r̃) - 1)dL̃/dτ̃
       ≈ -∂logΩ̃/∂r̃|_ηκ R 
    ⟹ (δ E - δ L)    ∼η R ,

where we expanded Ω̃(r) to first order in R and approximated dL̃/dτ̃≈ (dL̃/dτ)_ = -ηκ for κ constant defined by

    κ = (Ω̃^-1dt̃/dd_GW/dt̃)_∼𝒪(1)   for a∈[0,1].

Thus we deduce that δ E - δ L ∼η Rτ̃ for circular inspirals close to _. We shall see that these corrections are indeed subleading in the regime considered by OT <cit.>. However, they will not be negligible for near-extremal black holes. 




 §.§ Corrections arising from deviations from adiabatic nearly-circular inspiral



Given our assumption that the orbit is nearly circular when it reaches the transition regime, one expects corrections to the relation (<ref>) between the fluxes of energy and angular momentum satisfied for an exactly equatorial circular adiabatic inspiral. We discuss below two possible physical effects giving rise to such corrections : eccentricity and the non-geodesic past-history of the orbital evolution. These will give rise to the corrections (<ref>).

Eccentricity can lead to corrections to the transition equation which we will discuss further below, but eccentricity corrections to the fluxes tend to be suppressed during the transition regime. This is because the transition, for an arbitrary eccentric inspiral, corresponds to the orbit passing over the maximum of the effective potential given by Eq.(<ref>). 

The radial velocity throughout the transition regime is therefore always small, while the angular velocity remains 𝒪(1). Hence the orbit looks very much like a circular orbit, even if it is technically eccentric or even plunging. For nearly-circular transitions, the orbit is passing over a point of inflection of the effective potential and corrections to this approximately-circular assumption are even smaller. 

Corrections from non-geodesic past-history enter because the self-force acting on the small object at a particular time is generated by the intersection of the particle world line with gravitational perturbations generated by the orbital motion in the immediate past <cit.>. The self-force acting on the orbit when it is at a particular radius will therefore have corrections that depend on how far, in radius, the orbit has moved over the relevant past-history. The latter is determined by the dominant, azimuthal, timescale, and is an 𝒪(1) quantity, when expressed in coordinate time [If we are more conservative, we could assume that the timescale for radial oscillations is the appropriate averaging timescale. This is not 𝒪(1), but 𝒪(T), the scaling of the time coordinate in the transition zone. While this condition is more restrictive we will see below that even this condition does not change the conclusion that past-history corrections can be ignored in the transition zone.]. The orbital radius therefore changes by an amount of 𝒪(ṙ̃̇) over the relevant past-history. This is the scaling of the fractional change in the fluxes, and since Ẽ̇∼𝒪(ηϵ^2/3), the non-geodesic past-history corrections to the coordinate-time fluxes thus scale like ηϵ^2/3ṙ̃̇. In the regime η≪ṙ̃̇, considered by OT, and discussed in Section <ref>, ϵ can be considered 𝒪(1) and so the scaling of this correction is ηṙ̃̇. This is the first type of correction in Eq.(<ref>). In the adiabatic inspiral phase, these corrections are 𝒪(η^2) and form part of the second-order component of the self-force. However, in the transition phase these corrections can be larger.

We have argued above that eccentricity corrections to the fluxes should be suppressed in the transition regime. We now make this more concrete. Eccentricity corrections to the fluxes enter as fractional corrections of 𝒪(e^2), since corrections to the orbit at linear order in eccentricity are oscillatory and average to zero over a complete orbit <cit.>. The corrections to the coordinate time fluxes thus scale like ηϵ^2/3 e^2 (which is η e^2 in the OT regime discussed in Section <ref>). This is the second type of correction in Eq.(<ref>). If these corrections are to be small relative to the non-geodesic past-history corrections, we need e^2 < ṙ̃̇. In the transition zone we will see that the proper time scales like R^-1/2, where R =  - _ is the distance from the ISCO, regardless of the spin of the primary. For non near-extremal black holes, i.e., those with η≪ϵ, proper time and coordinate time scale in the same way and the scaling of ṙ̃̇ is therefore the same as that of R^3/2. The constraint we obtain on eccentricity is therefore e < R^3/4. However, there is also a geometric constraint, which is that the variation in the orbital radius due to eccentricity should be small compared to the variation due to radiation reaction through the transition zone. The latter is the scaling of R, while the former is a quantity of 𝒪(e), so we deduce an additional constraint e < R < R^3/4, the latter inequality following from the fact that R is a small quantity throughout the transition. We deduce that the geometrical constraint is stronger than the flux-correction constraint in the regime η≪ϵ. In the near-extremal regime, η≲ϵ, dt̃/dτ̃∼ϵ^-2/3 and so the constraint on the eccentricity changes to e < ϵ^1/3R^3/4 if these corrections are to be subleading. This is then more stringent than the geometric constraint. However, in this regime we will see below that eccentricity cannot grow until deep inside the transition zone, so even the more stringent constraint is easily satisfied.


Eccentricity during the transition can arise either from the presence of residual eccentricity prior to the start of the transition zone, or due to the excitation of eccentricity during the transition. The latter manifests itself as additional terms in the transition equation, the existence of which we will check for carefully in our analysis. To understand the former, we need to analyse the growth of eccentricity during the adiabatic inspiral. We will assume that at the beginning of the inspiral the orbit is nearly circular. It was shown in <cit.> that, for small eccentricity, the evolution of eccentricity under radiation reaction takes the form ė = f(r̃_0) e, where r̃_0 is the mean orbital radius and e is an eccentricity defined such that the orbital apoapsis is at r̃=r̃_0(1+e). For large r̃_0, f(r̃_0) < 0 and so the eccentricity decreases. In this regime any small eccentricity that is excited by small perturbations arising due to inspiral evolution or other effects is damped away and does not grow. However, for all spins a < 1, as the innermost stable circular orbit (or separatrix) is approached the sign of f(r̃_0) changes and is greater than zero in the vicinity of the ISCO. This means that orbits near to the separatrix are unstable to eccentricity growth. We would therefore expect any eccentricity that is excited to begin to grow.

Denoting ṽ^2=1/r̃_0, Kennefick <cit.> showed that the evolution of the orbital parameters, for small eccentricity, was governed by equations of the form

    ṙ̃̇_0/r̃_0   =-2(1-3ṽ^2+2aṽ^3)^3/2/ṽ^2(1-6ṽ^2+8aṽ^3-3a^2ṽ^4)_0
    ė/e   = 1/e^2(_0 - Ω(ṽ) _0) - j(ṽ) [Γ - h(ṽ) _0)] 
       
    
    j(ṽ)    = (1+aṽ^3)(1-2ṽ^2+a^2 ṽ^4) (1-3ṽ^2+2aṽ^3)^1/2/ṽ^2 (1-6ṽ^2+8aṽ^3-3a^2 ṽ^4)
    
    h(ṽ)    =  H(ṽ)(1+aṽ^3)^-1(1-2ṽ^2+a^2ṽ^4)^-2/2(1-6ṽ^2+8aṽ^3-3a^2ṽ^4)
     H(ṽ)   =1-12ṽ^2+66ṽ^4-108ṽ^6 + aṽ^3+8a^2ṽ^4
       -72aṽ^5-20a^2ṽ^6+204aṽ^7+38a^3ṽ^7-42a^2ṽ^8
       -9a^4ṽ^8-144a^3ṽ^9+116a^4ṽ^10-27a^5ṽ^11.

Both Γ and _0 are components of the self-force, which can be evaluated by solving the Teukolsky equation. The quantity Γ is in fact a linear combination of quantities that are time derivatives and so the above equation takes the same form for any choice of time coordinate with respect to which to evaluate the fluxes. Kennefick's analysis used coordinate time and so we make the same choice in the following discussion. An explicit expression for Γ is given in <cit.> and the quantity Ẽ̇_0 is the energy flux given in (<ref>). Numerical calculations show that these are finite quantities of 𝒪(1) throughout parameter space. The first term in the eccentricity evolution equation vanishes for evolution driven by gravitational radiation reaction, while the quantity h(ṽ) is singular at the ISCO. Therefore, close to ISCO the eccentricity evolution takes the form

    ė/e   ≈ j(ṽ) h(ṽ) _0 
    ⇒   r_0  dln e/ d r_0   ≈ - ṽ^2(1-6ṽ^2+8aṽ^3-3a^2ṽ^4) j(ṽ) h(ṽ)/2 (1-3ṽ^2+2aṽ^3)^3/2.

Notice that the expression is entirely geodesic and independent of the energy flux _0. For non-extremal spin, both j(ṽ) and h(ṽ) have simple poles at r̃=r̃_ and there is a simple zero in the term (1-6ṽ^2+8aṽ^3-3a^2ṽ^4) in the numerator. Therefore as ISCO is approached the eccentricity evolves as

    dln e/ dR≈ - k(a)/R  ⇒   e = e_0 ( R_0/R)^k(a)

with  R = r̃-r̃_ as before, and e_0 denotes the eccentricity when R=R_0 and _0≫_. The exponent k(a) is given by

    k(a)    =  H(ṽ_)/ D(ṽ_) 
     D(ṽ)    = 2 ṽ^2 (1 -2ṽ^2 + a ṽ^4) 
       × (12 ṽ - 24 a ṽ^2 + 12 a^2 ṽ^3) 
       ×(1-3ṽ^2+2aṽ^3)

and ṽ^2_=1/_. We find that k(a) = 1/4 for all a < 1. The behaviour for near-extremal black holes is slightly different, which we will discuss further below.

For extremal black holes the various factors in the expression for dln e/ dr̃_0 have repeated roots at the ISCO. To understand the behaviour for near-extremal black holes we therefore need to do an expansion in both R and ϵ. This takes the form

    dln e/ dR = a_0 ϵ^4 + a_1  ϵ^4 R+ ∑_i=2^5 a_i ϵ^2(6-i)/3 R^i + a_6 R^6 + ⋯/∑_i=1^6 b_i ϵ^2(7-i)/3 R^i + b_7 R^7 + ⋯

The terms omitted from both the numerator and denominator above are 𝒪(1) in ϵ. The ratio a_0/b_1 = -1/4, agreeing with the result for k(a) found above. However, for ϵ≪ R, the behaviour is not dominated by this term, but by the terms from a_6 in the numerator and from b_7 in the denominator. The leading order behaviour in this regime is therefore 




    dln e/ dR = a_6/b_71/R.

This is also exponential, but we find the ratio a_6/b_7 = 3/2, i.e., it is greater than zero and therefore the eccentricity decreases exponentially until we reach the regime R ∼ϵ. This is the statement that the critical curve, where the sign of the eccentricity evolution changes, is in the near-horizon region, which is consistent with results in <cit.>. We conclude that for near-extremal black holes, eccentricity can only grow once the inspiraling object is already very close to the ISCO, which is typically already inside the transition zone.

To complete this discussion we need to determine the scaling of the initial eccentricity e_0. If the orbit is truly circular then the eccentricity remains zero, so there must be some mechanism to excite an initial eccentricity which can then grow. Eccentricity can be excited by other physical processes, such as the presence of perturbing material, e.g., dust, or gravitational interactions with third bodies. Those processes are important, but in the pure-vacuum case eccentricity could still in principle be excited by the evolution under radiation reaction. We argued earlier that corrections to the fluxes far from the horizon scale like ηṙ̃̇ which is η^2 during the adiabatic inspiral. These corrections mean that the first term in Eq. (<ref>) is no longer exactly zero. Setting that term to η^2 we find an evolution equation of the form de^2/ dt̃∼η^2. After a few orbits the eccentricity is then 𝒪(η). This eccentricity induced by second order corrections to the evolution is damped by the process described above, until we reach the critical curve where it grows, eventually exponentially near the ISCO. This suggests appropriate initial conditions are e_0 ∼η and R_0 ∼𝒪(1).[A natural continuation of this argument would be to say that the second-order self-force induced corrections continue to drive eccentricity growth, over the whole of the inspiral, lasting a coordinate time ∼η^-1, leading to a final eccentricity of 𝒪(η^1/2), which can be larger than the eccentricity grown through the mechanism discussed here. However, this assumes that the eccentricity grows coherently and monotonically. In practice, once the eccentricity is 𝒪(η), the radial motion due to eccentricity becomes larger than the amount the radius evolves over the relevant past-history that determines the self-force and so the argument that the latter is the dominant contribution to corrections no longer applies. Knowledge of the second-order self-force would be required to fully explore the further evolution of the eccentricity and this is not currently available. However, we expect that the growth of initial eccentricities of O(η) through the instability mechanism will be the dominant contributor to the residual eccentricity in the transition zone.] We note that this mechanism could also excite eccentricity during the transition zone itself, but this would be of order e^2 ∼ηṙ̃̇ϵ^2/3 and hence no larger than the non-geodesic past-history corrections described above. If eccentricity grew coherently throughout the transition zone, the eccentricity induced by this process would be no larger than e^2 ∼ηṙ̃̇ϵ^2/3 T, where T is the coordinate time elapsed through the transition zone, which is typically smaller than the eccentricity grown during adiabatic inspiral prior to the start of the transition zone.

To summarise, we expect corrections to the evolution equations that arise from higher-order terms in the flux to scale like ηϵ^2/3 (which is η in the OT regime discussed in Section <ref>), and we expect residual eccentricity in the transition zone to be no more than e ∼η R^-k(a). In the non-near extreme case, these eccentricity corrections will be important when e > R, which implies R < η^1/(1+k(a)). In the near-extreme case, the corrections only become important when R∼ϵ, so we simply need to check that this is well inside the transition zone. In the analysis that follows we will evaluate the scaling of these terms and show that they are sub-dominant for inspirals into near-extremal black holes.



 §.§ Ori and Thorne regime



Consider non-extremal black holes, i.e. rotating black holes where the extremality parameter ϵ is not close to zero so that η≪ϵ. In this regime of spins and according to the discussion below (<ref>)-(<ref>), 
 all the coefficients controlling the general transition equation (<ref>) and (<ref>) are 𝒪(1). This is the regime originally discussed in <cit.>. 

Omitting coefficients of order one, the dominant contributions to the transition equation are

    (dR/dτ̃)^2   ∼ R^3 + R δ L + Γ_⊙
    Γ_⊙   ∼δ E- δ L ,

where we also omitted any further terms from (<ref>) and (<ref>) since they are subleading. Looking for a scaling solution R∼η^p and ∼η^q, it follows, using equation (<ref>) that δ L∼η^1+q. Requiring all dominant terms to have the same scaling fixes p=2/5 and q=-1/5, so that

    R = η^2/5ℛ ,   τ = η^-1/5𝒯 ,   δ L =  η^4/5δℒ .


Notice the overall scaling of the transition equation is (dr̃/dτ̃)^2∼η^6/5. The remaining question is whether the dominant term in Γ_⊙∼δ E- δ L is subleading or not. 
From (<ref>), it follows δ E - δ L ∼η^6/5 in this regime, suggesting the change of variables

    Γ_⊙ = η^6/5𝒴 .

This allows to write the schematic transition equation as

    (dℛ/d𝒯)^2∼ℛ^3 + ℛδℒ + 𝒴 .

Terms in Eqs. (<ref>) and 
(<ref>) that have been dropped can be seen to scale like the above terms multiplied by additional powers of ℛ or δℒ. Since both ℛ and δℒ are small quantities in the transition zone, these terms are sub-leading and we can ignore them.



The above scaling analysis proves the dominant terms in (<ref>) in the regime η≪ϵ are captured by


    (dR/dτ̃)^2≃ -2/3α R^3 + 2β δ L R + γ δ L R^2 + Γ_⊙ + …



    (dR/dτ̃)^2≃ -2/3α R^3 + 2β δ L R + Γ_⊙ + …

where we neglected all subleading corrections, kept the same original notation as in OT <cit.> for the coefficients 


    α   = -1/4∂^3G/∂r̃^3|_
    β   = 1/2(∂^2G/∂r̃∂ẼΩ̃ + ∂^2G/∂r̃∂L̃)_
    γ   = 1/2(∂^3G/∂r̃^2∂ẼΩ̃ + ∂^3G/∂r̃^2∂L̃)_ ,



    α   = -1/4∂^3G/∂r̃^3|_
    β   = 1/2(∂^2G/∂r̃∂ẼΩ̃ + ∂^2G/∂r̃∂L̃)_

and the dominant contribution to (<ref>) reduces to

    Γ_⊙≃Ω̃_(δ E - δ L) ∂ G/∂Ẽ|_ + …


Keeping all coefficients of order one, the natural scaled variables to introduce are

    R    = η^2/5α^-3/5(βκ)^2/5X 
    τ̃ - τ̃_   = η^-1/5(αβκ)^-1/5T 
    δ E - δ L    = η^6/5 Y  
    δ L    = -η^4/5(αβ)^-1/5κ^4/5 T

where

    κ = (Ω̃^-1dt̃/ddẼ_GW/dt̃)_ .

Plugging this into (<ref>), one obtains

    (dX/dT)^2 = -2/3X^3 - 2XT + C_0(Ω̃∂ G/∂ E)_Y + 𝒪(η^2/5)

where we defined C_0 = α^4/5(κβ)^-6/5. From now on, we ignore the subleading corrections.

The analogue of the acceleration equation (<ref>) reduces to 

    d^2X/dT^2   = -X^2 - T 
        - 1/2(dX/dT)(2X - C_0[Ω̃∂ G/∂Ẽ]_dY/dT)

This depends on the time evolution of the circularity deviation parameter Y, whose dominant contribution is derived in (<ref>). Inserting the re-scaled variables (<ref>) into Eq.(<ref>)

    dY/dT = -∂logΩ̃/∂r̃|_(β C_0)^-1X ,

leads to a transition equation 

    d^2X/dT^2 = -X^2 - T - 1/2(dX/dT)(2X + β^-1[∂Ω̃/∂r̃∂ G/∂Ẽ]_X) .

Evaluating (<ref>) at ISCO, we find the term in square brackets equals

    (∂Ω̃/∂r̃∂ G/∂Ẽ)_ = -2β

and so the last term vanishes. This was inevitable, since this term is precisely the term that arises from the dissipative part of f^r̃ from Eq.(<ref>), as identified earlier. The leading order evolution of Y is driven by maintaining the circularity of the orbit and so with this condition we expect the radial self-force corrections to be subleading.


The resulting transition equation of motion in the regime of low spins η≪ϵ is

    ^2X/ T^2 = -X^2 - T

and Y is evolved through the ODE (<ref>). We note that the transition equation does not depend on Y in this regime. Corrections to this equation arising from evolution of Y enter at an order η^2/5 higher than leading and so are subdominant. As discussed earlier the evolution of Y is related to deviations from the linear-in-proper-time evolution of energy and angular momentum and so the fact that these corrections do not enter the transition equation for η≪ϵ demonstrate that the linear evolution assumed by OT is appropriate in this regime.

Let us check the self-consistency of the transition equation (<ref>) by verifying that all neglected corrections to it are indeed smaller when evaluated on the scaling regime (<ref>). First, as discussed in section <ref>, the corrections to the orbit due to the conservative piece in the self-force are order 𝒪(η), see (<ref>). This is indeed smaller than the "radial thickness" R∼η^2/5 in (<ref>). Second, corrections due to η, appearing in (<ref>), are 𝒪(η^8/5). Hence, these corrections are 𝒪(η^3/5) smaller than the dominant δ L and δ E scaling in (<ref>) [Using the more conservative assumption that the averaging timescale is determined by the period of radial oscillations, which scales with T∼η^-1/5, the corrections are still suppressed by a factor of η^2/5. Third, 
corrections to dΓ_⊙/dτ̃ arising from non-geodesic past history corrections to the fluxes scale like η^8/5 and those arising from additional terms in the expansion of the azimuthal frequency as a function of radius scale as η^9/5, which are both subdominant to the leading η^7/5 scaling, albeit only by a factor of η^1/5.]

Finally, corrections arising from eccentricity are subleading provided e < -_, as discussed in Section <ref>. In the non-extremal case we therefore need e < η^2/5, due to (<ref>). This yields the constraint

    η^1-2k/5 < η^2/5  ⇒   3 - 2 k > 0   ⇒   k < 3/2.

We saw previously that k=1/4 for all spins a<1, which satisfies this bound. We deduce that eccentricity corrections are subdominant in the non-near-extremal regime. 
















 §.§ General Transition Equation of Motion - Near-Extremal


Let us consider rapidly rotating black holes with spin parametrized by a = √(1-ϵ^2) for ϵ≪ 1, as in (<ref>). The discussion below equations (<ref>)-(<ref>) allows to identify the a priori dominant contributions to the transition equation 
(<ref>) as 

    (dR/dτ̃)^2   ∼ R^3 + R δ L ϵ^2/3 + R^2δ L + δ L^2 ϵ^4/3 + Γ_⊙
    Γ_⊙   ∼ (δ E- δ L) (ϵ^2/3 + R + δ L ϵ^2/3) .


Since the functional dependence of the above equation does not depend on η, we learn the η scaling should be the same as before if we keep the R^3 and R δ L terms. Hence, we are left to determine any possible ϵ scaling. Proceeding as before, we look for scalings of the form R ∼η^2/5ϵ^p and ∼η^-1/5ϵ^q. We learn from equation (<ref>) that δ L ∼η^4/5ϵ^q. Requiring these dominant terms to scale in the same way determines p=4/15 and q=-2/15, so that

    R    = η^2/5ϵ^4/15ℛ ,   τ̃ = η^-1/5ϵ^-2/15𝒯 , 
    δ L    =  η^4/5ϵ^-2/15δℒ .

Hence, if η∼ϵ, the term R^2δ L scales like the velocity squared (d/dτ̃)^2 ∼η^6/5ϵ^4/5∼ϵ^2 and must be kept in the transition equation, whereas the term δ L^2 ϵ^4/3 is 𝒪(ϵ^2/3) smaller and, consequently, subdominant.

The only remaining question is whether Γ_⊙ is relevant in this regime or not. Using (<ref>) and the scalings (<ref>), we infer (δ E - δ L)∼η^6/5ϵ^2/15.  Since in the regime η∼ϵ, R∼δ L ∼ϵ^2/3 we conclude
Γ_⊙∼η^6/5ϵ^4/5∼ (d/dτ̃)^2 and must be kept in the transition equation. Introducing the finite variable 𝒴

    Γ_⊙ = η^6/5ϵ^4/5𝒴 ,

the general transition equation in the η∼ϵ regime reduces to

    (dℛ/d𝒯)^2∼ℛ^3 + ℛδℒ + ℛ^2δℒ + 𝒴 .

 Notice the radial velocity throughout the transition regime scales like d/dτ̃∼η^3/5ϵ^2/5∼η in the regime ϵ∼η. This is as in the adiabatic regime, but smaller than in the OT regime where d/dτ̃∼η^3/5.



As a self-consistency check, we can write the radial geodesic equation using the change of variables (<ref>) and (<ref>)

    (dℛ/d𝒯)^2∼∑_i=3^∞η^2(i-3)/5ϵ^4(i-3)/15ℛ^i + δℒℛ + 
    ∑_m=2^∞(η/ϵ)^2(m-1)/5ϵ^2(m-2)/3ℛ^mδℒ + η^2/5ϵ^4/15δℒ^2
     + ∑_n=1^∞(η/ϵ)^2(n+1)/5ϵ^2(5n - 1)/15δℒ^2ℛ^n + 𝒴.

It is apparent that the dominant terms are the i=3 and m=2 terms, all others being subleading. 

The above scaling analysis proves the dominant terms in (<ref>) in the regime η∼ϵ are captured by

    (dR/dτ̃)^2≃ -2/3α R^3 + 2β δ L R + γ δ L R^2 + Γ_⊙ + …

where α and β are defined as in (<ref>)-(<ref>) and γ = B_2 in (<ref>). As shown in appendix <ref>, they are approximated by 

    α→ 1 ,  β→ 2^-2/3√(3)ϵ^2/3≡β̂ ϵ^2/3 ,   γ→√(3) .

Furthermore, the dominant contributions to Γ_⊙ are 

    Γ_⊙ = Ω̃_(δ E - δ L) (∂ G/∂Ẽ|_ + ∂^2 G/∂r̃∂Ẽ|_ R + …) .



Keeping all coefficients of order one, the natural scaled variables to introduce are

    R    = η^2/5ϵ^4/15α^-3/5(β̂κ)^2/5X , 
    τ̃ - τ̃_   = η^-1/5ϵ^-2/15(αβ̂κ)^-1/5T  , 
    δ E - δ L    = η^6/5ϵ^2/15 Y 
    δ L    = -η^4/5ϵ^-2/15(αβ̂)^-1/5κ^4/5 T .

Since η∼ϵ, it follows R∼ϵ^2/3. Hence, the near ISCO expansion corresponds to the near horizon geometry of the primary black hole since, in Boyer-Lindquist coordinates, |_ - r̃_+| ∼ϵ^2/3. As a result, we will be able to use the (leading order and analytic) expression for the energy flux due to gravitational radiation in (<ref>). This allows to compute κ in (<ref>) in this regime as 

    κ   = (Ω̃^-1dt̃/dτ̃dẼ/dt̃)_→8/√(3)(C̃_H + C̃_∞).

Notice κ∼𝒪(1) since C̃_H + C̃_∞∼𝒪(1). 

Ignoring subleading terms, the general transition equation (<ref>) reduces to

    (dX/dT)^2 = -2/3X^3 - 2XT -(η/ϵ)^2/5C_1TX^2 + Γ̃_⊙

with

    C_1   = γ (αβ̂κ)^-3/5κ
    Γ̃_⊙   = ϵ^-4/5η^-6/5α^4/5(β̂κ)^-6/5Γ_⊙.

Notice the appearance of the new term proportional to TX^2, compared to the OT regime, is due to the regime η∼ϵ.

Taking a further T derivative, we find the analogue of the acceleration equation (<ref>) in this regime

    d^2X/dT^2 = -X^2 - T - (η/ϵ)^2/5C_1XT + 
    1/2(dX/dT)(-2X - (η/ϵ)^2/5C_1X^2 + dΓ̃_⊙/dT)

The time evolution of Γ_⊙ in (<ref>) has two contributions : one proportional to dY/dT, which can be computed using (<ref>) and a second one proportional to Y(dX/dT). Altogether yields

    dΓ̃_⊙/dT =     2X - (η/ϵ)^2/5(αβ̂κ)^-3/5κ(∂Ω̃/∂r̃∂^2G/∂r̃∂Ẽ)_ X^2
        + (η/ϵ)^2/5α^1/5(β̂κ)^-4/5(Ω̃∂^2G/∂r̃∂Ẽ)_ Y dX/dT

where we used (<ref>) to simplify the first term. The latter cancels the -2X term in (<ref>). Using the dominant contribution to the identity (<ref>) evaluated at ISCO, the second term cancels the C_1X^2 term in (<ref>). Finally, the third term gives a non-trivial contribution to the acceleration equation

    d^2X/dT^2 = -X^2 - T - (η/ϵ)^2/5(C_1XT - C_2Y)

with constant defined by

    C_2 = 1/2α^1/5(β̂κ)^-4/5(Ω̃∂^2G/∂r̃∂Ẽ)_.

and evolution equation for Y such that

    dY/dT = -Λ ∂logΩ̃/∂r̃|_ X,  with Λ = α^-4/5κ^6/5β̂^1/5.


In our treatment of the OT regime (non near-extremal spins), the terms in Eq. (<ref>) were neglected since they scaled with η^2/5 and were subdominant. In the near-extremal case, one can clearly see that the XT and Y term are comparable to the (rescaled) radial acceleration provided η∼ϵ. As such, they must be included in the analysis. Our final transition equation of motion differs from Eq.(43) in <cit.>, which correctly included the Y term but missed the cross term XT, which is the same order as the terms being retained. Our analysis improves on <cit.> in two additional ways. Firstly, Y was introduced in <cit.> as a mathematical construct to ensure conservation of the four-velocity norm. The evolution equation for Y was derived by forcing the equation of motion obtained from differentiation of the kinetic energy equation, Eq. (<ref>), to agree with that obtained by expansion of the left-hand-side of the acceleration equation, Eq. (<ref>). This is equivalent to setting the radial self-force term to zero, which is equivalent to imposing the circular-to-circular condition. This physical interpretation of the procedure was not made clear in <cit.>, nor the interpretation of Y as representing departures from the linear-in-proper-time evolution. Secondly, the scaling of the flux given in Eq. (<ref>) was not known at that time and this was left as an unspecified power of ϵ. Now that we know this scaling we can do a more complete analysis of the near-extremal regime.

The quantities above can be computed in the near-extremal limit, ϵ→ 0,

    Λ    =  2^52/15(C̃_H + C̃_∞)^6/5/√(3) + 𝒪(ϵ^2/3)
     
    C_1    = 2^8/5(C̃_H + C̃_∞)^2/5+ 𝒪(ϵ^2/3)
    
    C_2    = 2^-13/15· 3^-1/2(C̃_H + C̃ _∞)^-4/5+ 𝒪(ϵ^2/3).

Equations (<ref>) and (<ref>) are a coupled set of ODEs which will link the adiabatic inspiral to a plunging geodesic.

As in the previous section we now consider the size of corrections to the transition equation. Corrections to the circular flux-balance law in the geodesic part of the transition equation scale like η  ϵ^2/3 according to (<ref>). These are 𝒪(ϵ) smaller than the terms kept in this regime [Using the conservative assumption about the averaging timescale, these corrections are sub-leading by a factor of η^2/3.]. Similarly, corrections to the linear-in-proper-time angular momentum evolution enter through corrections to δ E and δ L and scale like d/dτ̃ times terms that are being retained. These are therefore subdominant since d/dτ̃∼η^3/5ϵ^2/5≪ 1. These corrections also contribute additional terms through corrections to the radial self-force part of the transition equation. These are of order η·∂ G/∂Ẽ and η·∂ G/∂L̃, which scale like ηϵ^2/3 and so are a factor of (η/ϵ)^1/5ϵ^1/3 smaller than the leading order terms in the transition equation and are therefore sub-dominant.

Eccentricity corrections enter like fractional e^2 corrections to the fluxes, and are only more important than the corrections described above if e > R ∼η^2/5 or e^2 >. In the near-extremal regime ∼ϵ^2/3  dr/ dτ∼η^5/3 and so eccentricity corrections become important when e > η^5/6. However, as shown in Section <ref>, for near-extremal inspirals eccentricity can only grow once r̃-r̃_∼𝒪(ϵ). In the transition zone r̃-_∼ (η/ϵ)^2/5ϵ^2/3≫ϵ and so eccentricity has not started to grow when the transition zone is reached. Residual eccentricity from the adiabatic inspiral would be 𝒪(η) and eccentricity excited during the transition would be 𝒪(η^4/5ϵ^8/15) (or 𝒪(η^7/10ϵ^7/15) if it was coherently excited throughout the transition). These are smaller than the threshold η^5/6 at which the eccentricity corrections become more important than the non-geodesic past history corrections, which we have already shown to be sub-leading.










 §.§ General Transition Equation - Very Near-Extremal



The final regime concerns very rapidly rotating black holes, where ϵ≪η. Using the results in appendix <ref>, one can identify the a priori dominant contributions to the transition equation (<ref>) and (<ref>) to be (ignoring coefficients of 𝒪(1))

    (dR/dτ̃)^2   ∼ R^3 + R δ L ϵ^2/3 + R^2δ L + δ L^2 ϵ^4/3 + Γ_⊙
    Γ_⊙   ∼ (δ E- δ L) (ϵ^2/3 + R +  ϵ^2/3δ L) .


It is natural to expect that terms involving some explicit factors of ϵ should be sub-leading in this regime.
Assuming a scaling solution of the form R∼η^α and τ̃∼η^β, we learn using (<ref>) that δ L ∼η^β + 1. Imposing the dominant terms R^3 and R^2δ L scale like (dR/dτ̃)^2 yields the scaling solutions α = 2/3 and β = -1/3, so that

    R = η^2/3ℛ ,    = η^-1/3𝒯 ,   δ L =  η^2/3δℒ .

As a consistency check, notice the terms ϵ^2/3Rδ L ∼η^2(ϵ/η)^2/3 and ϵ^4/3δ L^2∼η^8/3(ϵ/η)^4/3 are sub-dominant compared to the leading scaling (d/dτ̃)^2∼η^2.

The remaining question is whether Γ_⊙ is negligible in this regime or not. Using the scalings (<ref>) together with (<ref>), we infer that δ E - δ L∼η^4/3. It follows Γ_⊙∼η^2 from the term linear in R in the second equation in  (<ref>). Introducing the finite variable 𝒴

    Γ_⊙ = η^2𝒴

leads to the transition equation of motion

    (dℛ/d𝒯)^2∼ℛ^3 + ℛ^2δℒ + 𝒴 .

Notice the radial velocity throughout the transition regime scales as d/dτ̃∼η, as it were throughout the adiabatic inspiral regime and in the near-extremal case [see sec.(<ref>)]. Thus the radial motion is fastest throughout the transition regime when the primary is of moderate spin: η≪ϵ.

 As a consistency check, we can substitute the scalings (<ref>) and (<ref>) into the general transition equation (<ref>)

    (dℛ/d𝒯)^2∼∑_i=3^∞η^2(i-3)/3ℛ^i + (ϵ/η)^2/3δℒℛ + 
    ∑_m=2^∞η^2(m-2)/3ℛ^mδℒ + ϵ^2/3(ϵ/η)^2/3δℒ^2 + 
    ϵ^4/3δℒ^2ℛ + ∑_n=2^∞η^2(n-1)/3δℒ^2ℛ^n + 𝒴.

Clearly the dominant terms occur when both i = 3 and m = 2 with the rest being subleading. 

The above scaling analysis proves the dominant terms in (<ref>) in the regime ϵ≪η are captured by

    (dR/dτ̃)^2≃ -2/3α R^3 + γ δ L R^2 + Γ_⊙ + …

where α and γ are given in Eq.(<ref>) with

    Γ_⊙≃Ω̃_(δ E - δ L) ∂^2 G/∂r̃∂Ẽ|_ R + … .


Keeping all coefficients of order one, the natural rescaled variables in this regime are

    R    = η^2/3α^-3/5κ^2/5X  
     - τ̃_   = η^-1/3(ακ)^-1/5T  
    δ E - δ L    = η^4/3 Y 
    δ L    = -η^2/3α^-1/5κ^4/5 T .

In these variables, the radial velocity equation (<ref>) can be expressed as

    (dX/dT)^2 = -2/3X^3 - K_1X^2T + Γ̃_⊙

with 

    K_1   = γα^-3/5κ^2/5 , 
    Γ̃_⊙   = η^-2α^4/5κ^-6/5Γ_⊙ ,

and κ as in (<ref>).

Taking a further derivative with respect to T yields the acceleration equation

    d^2X/dT^2 = -X^2 - K_1XT + 1/2(dX/dT)(dΓ̃_⊙/dT-K_1X^2).

Using (<ref>) together with (<ref>), one finds that 

    dΓ̃_⊙/dT   = α^1/5κ^-4/5Ω̃_∂^2G/∂r̃∂Ẽ|_dX/dTY  
        - α^-3/5κ^2/5(∂^2G/∂r̃∂Ẽ∂Ω̃/∂r̃)_X^2 .

Plugging this back in (<ref>) and using the dominant contribution to the identity (<ref>), the K_1X^2 term cancels and one is left with 

    d^2X/dT^2 = - X^2 - K_1XT + K_2Y

together with the evolution equation for Y(T) given by

    dY/dT = -α^-4/5κ^6/5∂logΩ̃/∂r̃|_X .

where 

    K_2 = 1/2α^1/5κ^-4/5Ω̃_∂^2G/∂r̃∂Ẽ|_ .

In the limit ϵ→ 0, the constants K_1 and K_2 approach the values

    K_1   → 2^6/53^3/10(C̃_H + C̃_∞)^2/5
     
       K_2   → 2^-7/53^-1/10(C̃_H + C̃_∞)^-4/5.

As argued in previous sections, corrections to the circular flux-balance law contribute terms to the transition equation which scale like d/dτ̃∼𝒪(η) times terms that are being retained and like ηϵ^2/3. Corrections to the linear-in-time angular momentum evolution enter with the same scaling as the former. The retained terms in the transition equation scale like η^4/3 in the very near extremal regime and so these corrections are both sub-leading. Eccentricity corrections enter like fractional e^2 corrections to the fluxes, but, as in the near-extremal case, eccentricity cannot grow until the transition zone has already been reached, and so these corrections are no larger than 𝒪(η^5/3 ϵ^2/3) and are also sub-leading. 







We conclude this subsection by noting that the transition equation of motion (<ref>) is perfectly well behaved in the limit ϵ→ 0 and can therefore be used to compute an inspiral into a maximally spinning black hole with a=1. In this case, the horizon coincides with the ISCO in BL coordinates. However, the proper distance is Δℓ = -M/3lnϵ (see Fig.2 in <cit.> together with explanations in <cit.> and more recently in appendix A of <cit.>). Hence, we terminate the integration of the ODE (<ref>) at r̃ = r̃_+, since our numerics are specific to BL coordinates. The presence of the horizon manifests itself in the transformation from proper time to coordinate time, which will be discussed, for non-extremal inspirals, in the next sub-section.













§ RESULTS




 §.§ Numerical Integration



We now seek to compute a full worldline r̃(τ̃) for ∞ > r̃≥r̃_+. Out of the three regimes just discussed, we restrict ourselves to the ϵ∼η one. This is because the η≪ϵ regime has already been considered in the literature <cit.> and the ϵ≪η regime has been argued to be inaccessible throughout the transition regime in <cit.>. The latter conclusion follows from the observation that the waves emanating from the secondary produce a spin down effect on the primary leading to a maximum attainable spin with ϵ∼η. For ϵ∼η, we try to find the solution to 

    d^2X/dT^2   = -X^2 - T - (η/ϵ)^2/5(C_1XT - C_2Y) 
    dY/dT   = -3/4Λ X

which deviates off the past adiabatic inspiral and evolves into a geodesic plunge. The constants in Eq.(<ref>) are given by Eq.(<ref>). We can derive an equation for an adiabatic inspiral in proper time by using the quasi-circular approximation. Using our far-horizon expression for the energy flux defined by Eq.(<ref>) with both equations (<ref>) and (<ref>), one derives

    dr̃/dτ̃ = -η64/5Ω̃^7/3(2a - 3r̃^1/2 + r̃^3/2)r̃/r̃^2 - 6r̃ + 8ar̃^1/2 - 3a^2ℰ̇(r̃).

This equation diverges at the ISCO which is a break down of the quasi-circular approximation. We shall use Eq.(<ref>) to smoothly transition from the adiabatic inspiral Eq.(<ref>) into a geodesic plunge to the horizon. We used a cubic spline to interpolate values for the relativistic correction ℰ̇() using exact flux data found in the http://bhptoolkit.orgBHPT. We then numerically integrate Eq.(<ref>) by stepping forwards in proper time until (τ̃) - _∼η^4/5ϵ^-2/15. We feel this criteria is suitable for turning on the transition equation of motion since our model for the flux is well represented during the transition regime. When this criteria is met we can be sure that our model for flux evolution throughout the transition regime is correct to leading order. Once this is satisfied, we stop integrating our adiabatic inspiral solution and begin integrating our transition equation of motion (<ref>).

Since we do not terminate our adiabatic inspiral solution at the ISCO, we do not know the precise proper time where the particle crosses the ISCO. As such, the variable T is not a good choice of variable to integrate on the right hand side of (<ref>). Instead, we substitute T for δ L from Eq.(<ref>) into our transition equation of motion, then

    d^2X/dT^2   = -X^2 + B_0δ L + (η/ϵ)^2/5(C_1B_0δ L + C_2Y) 
    dY/dT   = -3/4Λ X 
    dδ L/dT   = B_0^-1,   B_0 = -η^-4/5ϵ^2/15(αβ̂)^1/5κ^-4/5.

We use initial conditions determined by the end of the adiabatic inspiral Eq.(<ref>) at some time _init.

    X(T_init)    = η^-2/5ϵ^-4/15α^-3/5(β̂κ)^-2/5(r̃-_)
    dX/dT|_T_init   = η^-3/5ϵ^-2/5α^2/5(β̂κ)^-3/5dr̃/dτ̃|_τ̃_init
    
    Y(T_init)    = η^-6/5ϵ^-2/15(Ω̃_^-1
    δ E_init - δ L_init)
    δ L(T_init)    = L̃_circ(r̃_init) - L̃_.

Where L̃_circ(_init) corresponds to the circular angular momenta evaluated at the end of the inspiral, _init.  Using this prescription, we are able to integrate the coupled ODEs Eq.(<ref>) with initial conditions (<ref>) to obtain Fig.(<ref>).

The transition solution smoothly deviates away from the adiabatic inspiral (blue curve), passes through the ISCO and reaches the horizon where the solution terminates. The plot on the right shows the full worldline in proper time r̃(τ̃) where the inspiral starts at r̃ = 1.006 and terminates at the horizon. This method ensures that r̃(τ̃) is both continuous and once differentiable everywhere. 

Also, by our choice of integrating (<ref>) using the variable δ L, we ensure continuity but not differentiability in L̃  throughout the full inspiral. We note here that Apte and Hughes in <cit.> also found discontinuities in their evolution of both  and  and added corrections to ensure both (first order) differentiability and continuity at τ̃_init. We consider a correction of the form  

    _trans =  Δ_cor + _ + dL̃_/dτ̃(τ̃ - τ̃_).

We have discussed previously that the leading order term in L̃() - _ scales proportionally to η^4/5ϵ^-2/15. So we choose to add a constant offset Δ_cor∼η^6/5ϵ^2/15 to the angular momenta evolution to ensure continuity in the _trans evolution. 

To calculate the evolution in Ẽ, one computes Ẽ_circ given by Eq.(<ref>) during the adiabatic inspiral regime. Then, during the transition regime, one integrates

    Ẽ = Δ_cor + _ + ∫_τ̃_isco^τ̃_+Ω̃(r̃)L̇̃̇_trans dτ̃

from the flux balance law Ẽ̇ = Ω̃(r̃)L̇̃̇. The correction to Δ_cor is chosen to ensure continuity with the end of the inspiral energy given by Eq.(<ref>) as previously discussed after equation (<ref>).

Notice here that this ensures that the energy flux obeys Ẽ̇ = Ω̃(r̃)L̇̃̇_ and is thus not constant. This ensures that we are still granted a full cancellation of the dissipative part of the forcing term f̃^r̃ in Eq.(<ref>). This will yield a continuous evolution Ẽ at the matching point with a discontinuous first derivative. At this point we will have a full trajectory r̃(τ̃) with (continuous) integrals of motion in proper time Ẽ(τ̃) and L̃(τ̃). In each of  <cit.>, the authors compute three separate worldlines in proper time; Adiabatic inspiral, transition, geodesic plunge. Apte et al in <cit.>, provide an algorithm in which they freeze the constants of motion  and  when the extra terms in Eq.(<ref>) exceed the leading order terms X^2 and T by 5%. As one would expect, as one ventures farther from the ISCO, the Taylor expansion used to derive these transition equations of motion will break down. As such, it is very natural for each of the aforementioned authors to compute a geodesic plunge to complete their worldlines in proper time r̃(τ̃). Simply because, for moderate spins (non near-extremal), |_+ - _| ∼𝒪(1) η^2/5. For near-extreme black holes the ISCO is close to the horizon in Boyer Lindquist coordinates |r̃_ - r̃_+| ∼ϵ^2/3. The scaling of the near-extremal transition zone is also ϵ^2/3 and so the horizon is reached while the object is still in the transition zone. We therefore do not expect to need to add a geodesic plunge to compute full near-extremal inspirals. To verify this we numerically calculate  the extra terms in (<ref>), which are 

    C_3 X^3   ⇒ C_3 =1/12(η/ϵ)^2/5ϵ^2/3∂^4G/∂r̃^4|_α^-8/5(β̂κ)^2/5
    
    C_4XY    ⇒ C_4 = 1/2(η/ϵ)^4/5ϵ^2/3(Ω̃∂^3G/∂r̃^2∂Ẽ)_(αβ̂κ)^-2/5.

We compare the solution to (<ref>) when these terms are omitted or included in Figure <ref>. The difference is at most 1% even the horizon _+. We conclude that we can use the solution from (<ref>) throughout the plunging regime, for r̃(τ̃) ∈ [r̃_+,r̃_]. It would be useful in the future to compare our results with the analytic geodesic plunges found in <cit.>.



 §.§ Worldline in Boyer-Lindquist Coordinates


In the previous section, we computed the full worldline comprised of inspiral, transition and plunge parametrized as r̃(τ̃). We now intend to do the same but in coordinate time so that our worldline is in Boyer-Lindquist coordinates (t̃,r̃(t̃),θ = π/2,ϕ(t̃)). Loosely speaking, this is the time measured from Earth (at radial infinity) so is useful for observable purposes.

For the quasi-circular inspiral solution, we simply integrate the circular relation relating coordinate time to proper time via Eq.(<ref>)

    t̃ = ∫_0^τ̃_init1 + a/r̃^3/2/√(1 - 3/r̃ + 2a/r̃^3/2) dτ̃

where r̃(τ̃) is the worldline constructed by integrating Eq.(<ref>) up to some suitable point to begin the transition solution, in our case, r̃(τ̃_init) = r̃_init. To compute the trajectory in coordinate time r̃(t̃) throughout the transition regime, we must integrate 

    t̃ = t̃(_init) + ∫_τ̃_init^τ̃_+ T(r̃,Ẽ,L̃,a)dτ̃.

where T(r̃,Ẽ,L̃,a) is given by <ref> and t̃_insp is defined through t̃(τ̃_init). Throughout the transition regime, we use the model for both Ẽ(τ̃) and L̃(τ̃) given by Eq.(<ref>) and Eq.(<ref>).  This will yield the  r̃(t̃) throughout the transition regime. Combining these results yield a full trajectory from radial infinity to the horizon in coordinate time r̃(t̃).

To then calculate the orbital velocity dϕ/dt̃ = Ω̃ in coordinate time we substitute r̃(t̃) found previously into Eq.(<ref>). This now gives Ω̃(t̃) valid throughout the adiabatic inspiral regime. Using our solutions for Ẽ(τ̃) and L̃(τ̃) defined through Eq.(<ref>) and Eq.(<ref>) and r̃(t̃) throughout the transition regime, we calculate

    Ω̃ = dϕ/dt̃ = 2aẼr̃ - a^2 L̃ + Δ(r̃) L̃/Ẽ(r̃^2 + a^2)^2 - 2aL̃r̃ - Δ(r̃) a^2 Ẽ.

where Δ(r̃) = r̃^2 - 2r̃ + a^2. This algorithm will provide a worldline in coordinate time r̃(t̃) which will be used for our waveforms. We stress here that r̃(t̃) is continuous and (once) differentiable.






 §.§ Near-Extremal Waveform


Following <cit.>, the root mean square (rms) amplitude of gravitational waves emitted towards infinity at harmonic m is given by h_o,m = √(⟨ h_+,m^2 + h_×, m^2⟩). The plus and cross each represent individual transverse-traceless polarisations of the gravitational wave strain h. The amplitudes are averaged ⟨·⟩ over the direction and over the period of the waves.  Furthermore, the rms amplitude 
generated by a particle on an equatorial circular orbit in the limit η→ 0 is related to the outgoing radiation flux in harmonic m by 

    h_o,m = 2M√(ηẼ̇_∞,m)/mΩ̃D

for distance D and outgoing fluxes defined by

    Ẽ̇_∞,m = η𝒜_mΩ̃^2 + 2m/3ℰ̇_∞,m

where the amplitude 𝒜_m equals

    𝒜_m = 8(m+1)(m+2)(2m)! m^2m-1/(m-1)[2^mm!(2m+1)!!]^2 ,    m≥ 2

and ℰ̇_∞,m is the relativistic correction to Ẽ̇_∞,m at each harmonic m. 

An EMRI signal is a superposition of infinitely many harmonics of the fundamental frequency Ω̃

    h = ∑_m=2^∞h_o,msin(2πf̃_mt̃ + ϕ),

with 2π f̃_m= m·Ω̃. Recall that the total emission of radiation through gravitational waves is related to the outgoing and ingoing flux by 

    _GW   = _∞ + _H
       = ∑_m=2^∞(Ẽ̇_∞,m + Ẽ̇_H,m)

where Ẽ̇_H,m is the ingoing flux (towards the horizon) including the contribution from all l for each harmonic m. Using the exact results from the http://bhptoolkit.orgBHPT for a spin parameter of a = 1-10^-9, we constructed a cubic spline for each outgoing flux _∞,m. Our results are plotted in figure <ref>. It is clear that including the higher order modes become increasingly important as the spin parameter increases towards unity. This has already been observed in <cit.>. Hence, for near-extremal systems, only using the m=2 harmonic is not an accurate representation of the EMRI signal in general.


 

Figure <ref> suggests that truncating the sum at the eleventh harmonic in the outgoing flux (<ref>) is a good approximation to model a near-extremal waveform encapsulating quasi-circular inspiral, transition, and then plunge with suitable accuracy. The remaining difference from modes with m > 11 contributes a small difference in the amplitude of the waveform, but gravitational wave detectors are much less sensitive to amplitude corrections than corrections to the phase. The phase evolution is determined by the total flux, in which we are including all modes. We therefore believe that the approximate waveform with 11 modes is a sufficiently accurate model for parameter estimation studies and will use this henceforth

    h ≈∑_m=2^11h_o,msin(2πf̃_mt̃ + ϕ).


Once the ISCO is reached, we smoothly extrapolate each of the fluxes Ẽ̇_∞,m→ 0, as  r̃→r̃_+.  
This is a similar approach to that found in Taracchini et al in <cit.>. Using (<ref>) and the results obtained in this paper, we plot a near-extremal waveform, including the transition from inspiral to plunge, in Fig.(<ref>).
 
We notice that the waveform in Figure <ref> exhibits the usual dampening before the ISCO is reached as seen by Gralla et al in <cit.>. This, qualitatively, is a unique feature to near-extreme EMRIs as a gravitational wave source. 





§ CONCLUSION



This paper has presented a solution to the problem of the transition from inspiral to plunge, for any primary spin, for EMRIs on circular and equatorial orbits. This work has extended the treatment of Ori & Thorne <cit.> which was the first analysis of this problem but did not apply to systems with near-extremal spins. This work also extended the analysis of  <cit.> which did consider near extremal spins, by providing a better physical interpretation of the procedure, identifying a missing term in the analysis and updating the treatment to use recent calculations of the near-extremal energy flux. We have also carefully identified the scaling of the various higher order terms arising from effects such as eccentricity and non-geodesic past-history to carefully demonstrate that these are all sub-dominant. Previous treatments have assumed that the quasi-circular assumption holds throughout the inspiral, but without rigorous justification. We have demonstrated that initial eccentricities excited during the adiabatic inspiral regime grow by the time the transition regime is reached, but are still sufficiently small to be sub-dominant. We have shown that corrections to the flux balance law (<ref>) arising from eccentricity and from the non-geodesic past-history of the orbital evolution are also sub-dominant, if only marginally, but there are non-trivial deviations from the linear-in-proper-time evolution of energy and angular momentum in (<ref>) that was assumed in OT. These deviations are encoded in the evolution of the parameter Γ̃_⊙ through the transition regime. 

Based on these arguments, we have derived a transition equation for each of the three scaling regimes: η≪ϵ, η∼ϵ and η≫ϵ and described a numerical scheme to generate a full inspiral trajectory in coordinate time, from radial infinity to the horizon. For near-extremal black holes, we found that there was no need to attach a geodesic plunge onto the transition solution as the inspiraling object reaches the horizon while still within the transition regime. Finally, we used these inspiral trajectories to construct a near-extremal waveform exhibiting the transition and plunging dynamics using results from the http://bhptoolkit.orgBHPT <cit.>. 

The OT procedure is straightforward, but with surprisingly rich phenomenology. Through semi-analytic means, one is able to derive an equation which describes the dynamics within the vicinity of the ISCO. However, in practice, the OT theory has several shortcomings. The point at which the transition solution is taken to start has a significant influence on the time it takes the particle to reach the horizon and so the OT procedure does  not define a unique worldline given a particular set of parameters for the source. This is clearly not physical behaviour. We argued in section <ref> that if the switch from the adiabatic inspiral to the transition equation is made when the constraint δ L ∼η^4/5ϵ^-2/15 is satisfied, the solution will be almost unique. This was verified numerically and we found it leads to plunge times consistent within ± 0.5M. This very same problem was found in <cit.> but they saw no effect in their waveform analysis. For the η≪ϵ case, there is a further degree of freedom as to when to attach the geodesic plunge. To do so, one must “freeze" the integrals of motion Ẽ,L̃ at the end of the transition regime and integrate the Kerr geodesic equations forward in coordinate time. Attaching the geodesic plunge is discussed, at length, in <cit.> but does not have a unique solution. Care must be taken as to when the transition solution and the geodesic plunge is attached or comparatively different radial trajectories will be produced. Fortunately for the ϵ∼η case, there is no need to attach a geodesic plunge as shown earlier in section <ref>.

Another issue with the OT method is that it can lead to discontinuities in the constants of motion Ẽ(τ) and L̃(τ) if the OT equations are integrated backwards from the ISCO rather than forwards from the point of the switch from the adiabatic inspiral to transition regime. Discontinuities in the constants of motion lead to discontinuities in the coordinate time trajectories and in the waveforms which must be avoided if these waveforms are to give physically reasonable results in parameter estimation studies. Our solution, which was to integrate forward not backwards, yields continuous, but not first order differentiable, trajectories. The procedure described in <cit.> provides both. For parameter estimation studies we only require continuity of Ẽ and L̃ and first order differentiability of r̃(t̃) and so our procedure should be sufficient, although this should be examined more carefully.


There are natural extensions of this work. First, the waveforms constructed in this paper can be used to carry out a parameter estimation study to understand how well the parameters of near-extremal EMRIs can be measured with observations by LISA. Of particular interest is how well the spin can be determined, since the identification of an object that definitely has spin above the Thorne limit would be of profound significance. Second, our waveforms are missing the quasi-normal mode ringdown contribution. Hence, it would be very interesting to generate a full waveform taking these into account, together with the plunging dynamics discussed here. Details on how to construct the waveform including this effect were discussed in <cit.>.  Finally, it would also be of interest to extend the current analysis to inspirals that are not circular and equatorial. The extension of OT was first performed in  <cit.> who attempted to solve the problem for generic orbits; both eccentric and inclined orbits. Sundararajan's treatment was corrected by <cit.> in the case of arbitrary inclination. Hence, no one, as of yet, has considered the transition from inspiral to plunge in the case of eccentric orbits and inclined orbits. These orbits are expected for EMRIs formed through standard astrophysical channels <cit.>. The extension to eccentric orbits will require more careful modelling of the self-force and the use of the (eccentricity-dependent) separatrix in place of the ISCO among other complications. A model of the transition for inspirals on generic orbits into black holes arbitrary spin will be invaluable  for the analysis of future LISA EMRI observations and is an important future topic of study. 


This work made use of data hosted as part of the Black Hole Perturbation Toolkit at http://bhptoolkit.org/bhptoolkit.org. We wish to thank Maarten van de Meent for his comments on the final stages of this manuscript. We also give our thanks to Maarten for pointing out useful literature concerning the conservative part of the self-force and giving direction on why we could neglect it. We also wish to thank Niels Warburton for his guidance using the toolkit and Peter Zimmerman for his useful discussions regarding the behaviour of near-extremal Kerr Black holes. Finally, we thank Geoffrey Compère for comments on an earlier version of this manuscript.






§ THE INNERMOST STABLE CIRCULAR ORBIT



In this appendix we review the main properties of the function G(r̃,,) determining the radial geodesic (<ref>)

    G(r̃,Ẽ,L̃)    = Ẽ^2 - 1 + a^2(Ẽ^2-1)-L̃^2/r̃^2
       + 2 (a Ẽ-L̃)^2/r̃^3 +2/r̃.

together with its derivatives when evaluated at the ISCO orbit _. The spin dependence of these quantities will play a critical role in the identification of the different transition regimes discussed in section <ref>.

Remember the ISCO radial coordinate _ is characterised by marginal stability

    G(r̃,_,_)= . ∂ G/∂r̃|_ = . ∂^2 G/∂r̃^2|_   = 0 .

Labelling the energy and angular momentum of the ISCO orbit by _ and _, we can solve the second and third constraint equations by 

    _   = r̃_^2 - 3a^2 + 6r̃_/2√(3) r̃_ , 
    _   = 6r̃_ - 3a^2 - r̃_^2/2√(3) a r̃_ .

Plugging these into G(r̃_,_,_)= 0, one derives the relation 

    2/3r̃_ = 1-_^2 ,

which combined with (<ref>) yields

    r̃_^2 - 6r̃_ + 8a√(r̃_) - 3a^2 =0 ,

whose solution r_0(a) reproduces (<ref>) <cit.>. This equality allows to simplify the energy (<ref>) and angular momentum (<ref>)
of equatorial circular orbits when evaluated at ISCO to

    _   = 1-2/_ + a/_^3/2/√(1-3/_ + 2a/_^3/2) = 4√(_) - 3a/√(3)_
    _   = r̃^1/21 - 2a/r̃^3/2 + a^2/r̃^2/√(1 - 3/r̃ + 2a/r^3/2) = 2√(3) - 4a/√(3_) .


Armed with these identities, we move towards the evaluation of the derivatives controlling the expansions (<ref>) relevant to the transition regime. First, we introduce some notation

    A_n    = ∂^nG/∂r̃^n|_ ,
    
      B_n    = (∂^n+1 G/∂r̃^n∂Ω̃ + ∂^n+1G/∂r̃^n∂L̃)_ , 
    
       C_n    = (∂^n+2 G/∂r̃^n∂^2Ω̃^2 + 2∂^n+2 G/∂r̃^n∂L̃∂ẼΩ̃ + ∂^n+2 G/∂r̃^n∂^2)_

with Ω̃ as in(<ref>). Either by explicit calculation or by induction, one can prove for any integer n≥ 0


    ∂^nG/∂^n   =     (-1)^n((n+2)!(aẼ - L̃)^2/r̃^n+3 + (n+1)!(a^2(Ẽ^2 - 1) - L̃^2)/r̃^n+2 + 2n!/r̃^n+1) - δ_n0(1 - Ẽ^2) , 
    ∂^n+1G/∂r̃^n∂Ẽ   =    (-1)^n(2(n+2)!(a^2Ẽ - aL̃)/r̃^n+3 + 2(n+1)!a^2Ẽ/r̃^n+2) + 2δ_n0Ẽ , 
    ∂^n+1G/∂r̃^n∂L̃   =    -(-1)^n(2(n+2)!(aẼ - L̃)/r̃^n+3 + 2(n+1)!L̃/r̃^n+2) , 
    ∂^n+2G/∂r̃^n∂L̃^2   =    (-1)^n(2(n+2)!/r̃^n+3 - 2(n+1)!/r̃^n+2) , 
    ∂^n+2G/∂r̃^n∂L̃∂Ẽ   =    -(-1)^n(2a(n+2)!/r̃^n+3) , 
    ∂^n+2G/∂r̃^n∂Ẽ^2   =    (-1)^n(2a^2(n+2)!/r̃^n+3 - 2a^2(n+1)!/r̃^n+2) + 2δ_n0

where δ_n0 = 1 for n=0 and zero otherwise.  Finally, evaluating these derivatives at (_,_,_) and using the properties (<ref>)-(<ref>), we can derive the exact results

    A_n    =     (1 - δ_n0)(-1)^n(n-1)(n-2)n!/3_^1+n , 
    
      B_n    =    2(1 - δ_n0) (-1)^n (n+1)! n(a-√(_)) + a - 2√(_) +_^3/2/_^n√(3_)(a-√(_)) (a+_^3/2) ,
    
      C_n    =    2·δ_0n -(-1)^n(2a + √(r̃_)[r̃_ - 2 - n])(n+1)! /r̃_^(2n+1)/2(a + r̃_^3/2)^2 ,


Notice equations (<ref>)-(<ref>) recover the familiar identities for circular orbits 

    (∂ G/∂ẼΩ̃ + ∂ G/∂L̃)_ = 0, 
            G_ = ∂ G/∂r̃|_ = ∂^2G/∂r̃^2|_ = 0.


Let us study the behaviour of these derivatives for near extremal black holes, i.e. in the limit ϵ→ 0 as introduced in section <ref>. Remember _ is given by 

    _→ 1 + 2^1/3ϵ^2/3 + 7/4· 2^1/3ϵ^4/3 + 𝒪(ϵ^2) 
     .

Using this expansion together with a = √(1-ϵ^2), we can evaluate the leading terms of all previous derivatives to be 


    A_n    → (1 - δ_n0)(n-2) (n-1)(1/3 (-1)^n Γ(n+1) + 𝒪(ϵ^2/3)) ,
    
      B_n    → (1 - δ_n0)(-1)^nΓ(n+2)/√(3)(n-1 - 4n^2 + n + 1/2^5/3ϵ^2/3) + 𝒪(ϵ^4/3) ,
     
    C_n   →  -1/4 (-1)^n (n-1) (-2 + 2^1/3ϵ ^2/3[2 n +3]) (n+1)! + (-1)^n (4n^2 -3n -3) (n+2)!/2^10/3ϵ ^4/3
       + p_n0
    ∂^n+1G/∂r̃^n∂Ẽ|_   →2/√(3)(-1)^n+1(n+1)![(n+1) - 2^1/3(n^2 + 3n + 3)ϵ^2/3] + 2/√(3)(1 + 2^1/3ϵ^2/3)δ_n0 + 𝒪(ϵ^4/3)
    (∂^n+2G/∂r̃^n∂Ẽ + Ω̃∂^n+2G/∂r̃^n∂Ẽ^2)_   → (δ_0n - (-1)^n(n+1)^2n!) + (-1)^n(7 + 13n + 4n^2)(n+1)! - 3δ_0n/2^5/3ϵ^2/3 + 𝒪(ϵ^4/3)
    ∂^n+2G/∂r̃^n∂Ẽ^2|_   →  2((-1)^n(3+n)(n+1)! + δ_0,n)   -2^4/3(-1)^n(4+n)(n+2)!ϵ^2/3 + 𝒪(ϵ^4/3)


where we defined 

    p_n0 = 2-3· 2^1/3ϵ ^2/3/4δ _0n + 𝒪(ϵ^2).


What we learn is that A_n∼𝒪(1) for all n≥3, B_1∼ϵ^2/3, B_n∼𝒪(1) for n≥ 2, C_0∼ C_1∼ϵ^4/3 and C_n∼𝒪(1) for n≥ 2. Furthermore,
(<ref>) and (<ref>) are 𝒪(ϵ^2/3) for n=0 and 𝒪(1) for n≥ 1, whereas (<ref>) is always 𝒪(1).







§ RETROGRADE ORBITS



In this section, we will restrict our attention to retrograde orbits. That is, orbits opposing the direction with the primaries angular momenta. These orbits are of interest because the ISCO is much further away from the horizon, which implies that the radial distance travelled during plunge time is much longer. We plot the location of the ISCO as a function of spin a in figure (<ref>).
 
Due to frame-dragging, we expect the ISCO to be farther from the hole since the space is dragged in the opposite direction to the compact objects orbital direction. In our conventions, retrograde orbits correspond to a<0 and >0. Hence, near-extremal ones are characterised by a → -1, or equivalently, by


    a → -√(1-ϵ^2),  where ϵ≪ 1.


notice that the horizon takes the same form as in the case of prograde orbits

    r̃_+ = 1 + √(1-a^2) =  1 + ϵ

as to be expected. Using a spin parameter of negative parity, the expressions for Ẽ,L̃, Ω̃ and r̃_ remain the same. However, each quantity will be different at the ISCO of a retrograde orbit. By substituting Eq. (<ref>) into Eqs.(<ref>), (<ref>), (<ref>) and (<ref>) for small ϵ≪ 1, one finds

    r̃_   = 9 - 45/32ϵ^2 +𝒪(ϵ^4) 
    Ẽ_   = 5/3 √(3)-1/96 √(3)ϵ ^2 +𝒪(ϵ^4)
    L̃_   = 22/3 √(3)-3 √(3)/16ϵ ^2 +𝒪(ϵ^4)
    Ω̃_   = 1/26 + 373/43264ϵ ^2 +𝒪(ϵ^4).

Notice here that the expansion in ϵ is no longer increasing in powers of ϵ^2/3 and now in ϵ^2. Also notice that |_ - r̃_+|∼𝒪(1) rather than of order ϵ^2/3 like in the case of near-extremal prograde orbits. Like we have done previously, we consider the Kerr radial velocity expanded around the ISCO

    (dR/dτ̃)^2≃ -2/3α R^3 + 2βδ LR + γδ LR^2 + Γ_⊙ + ….

with small variables 

    Ẽ - _   = Ω̃_δ E 
    L̃ - _   = δ L 
    r̃ - _   = R.

The coeffcients in (<ref>) can be approximated for ϵ→ 0 under the retrograde condition Eq.(<ref>)

    α   = -1/4∂^3G/∂r̃^3|_→1/6561
    β   =  1/2(∂^2G/∂r̃∂L̃ + Ω̃∂^2G/∂r̃∂Ẽ)_→4/351√(3)
    γ   = 1/2(∂^3G/∂r̃^2∂L̃+ Ω̃∂^3G/∂r̃^2∂Ẽ)_→ -1/351√(3).

and Γ_⊙ in (<ref>) defined through equation (<ref>). Notice that none of the coefficients in our transition equation of motion depend on the extremality parameter ϵ. This gives us no reason to introduce any scalings on r̃,τ̃ and δ L like we did for prograde orbits around rapdily rotating black holes. As such, let us introduce similar scalings to OT

    R    = η^2/5α^-3/5(βκ)^2/5X 
     - τ̃_   = η^-1/5(αβκ)^-1/5T 
    δ E - δ L    = η^6/5 Y. 
    δ L    = -η^4/5(αβ)^-1/5κ^4/5 T

Substituting these results into Eq.(<ref>) we find that 

    (dX/dT)^2 = -2/3X^3 - 2XT + α^4/5(ηβκ)^-6/5Γ_⊙.

Since R∼η^2/5, we only need the first term of Γ_⊙ 

    Γ_⊙ = η^6/5(Ω̃∂ G/∂Ẽ)_Y.

and taking derivatives of Eq.(<ref>) and following an identical procedure to subsection (<ref>),

    d^2X/dT^2 = -X^2 - T

with evolution equation for Y

    dY/dT = -∂logΩ̃/∂r̃|_(C_1K_0)^-1 X

for K_0 = α^4/5(βκ)^-6/5. Which is precisely the equation of motion for the transition regime derived by OT in <cit.>. Although the quantities α,β and κ present in the change of coordinates are different, the physics and ultimate end goal are the same. As a result, we stop our analysis of retrograde orbits here since we feel that this problem has already been solved by the community for smaller spin values a ≥ -0.999. We conclude that, for near-extremal retrograde orbits, there is nothing new to learn about the transition regime. It can be solved in the matter of OT in <cit.>. We do remark that the quantity κ can no longer be computed using the near-extremal formula defined by Ẽ̇_GW = (C̃_H + C̃_∞)(r̃ - r̃_+)/r̃_+. This is because the transition region is far from the horizon of the primary hole [see figure (<ref>)]. Instead we have to use the numerical quantity

    κ   = (Ω̃^-1dt̃/dτ̃dẼ/dt̃)_
       = (-32/5Ω̃^7/31 + a/r̃^3/2/√(1 - 3/r̃ + 2a/r̃^3/2)ℰ̇(r̃))_.

Various results are tabulated (including retrograde orbits) in <cit.>. The downside of this equation is that it can only be evaluated numerically. 



§ OSCULATING ELEMENTS EQUATIONS



The proper time derivative of the radial geodesic equation (<ref>) yields (<ref>)

    d^2/dτ̃^2 - 1/2∂ G/∂r̃ = 1/2(d/dτ̃∂ G/∂ + d/dτ̃∂ G/∂L̃)(dr̃/dτ̃)^-1 .

The purpose of this appendix is to review how this equation is equivalent to the radial component of a forced geodesic equation

    u^ν∇_νu^r̃= d^2x̃^r̃/dτ̃^2 + Γ^r̃_ρσdx̃^ρ/dτ̃dx̃^σ/dτ̃= f̃^r̃ ,

where x̃^μ = (r̃,t̃,θ,ϕ), ∇_ν = ∇_x̃^ν, u^μ = dx̃^μ/dτ̃ is the four-velocity of the particle and f̃^μ a forcing term driving deviations from geodesic motion.

To show the equivalence between (<ref>) and (<ref>), we use the osculating elements formulation <cit.>. Since this method does not take into account conservative effects arising from the self-force 
<cit.>, the component f̃^r̃ in this appendix will only account for the dissipative piece in (<ref>). Its conservative piece is treated in more detail in the main text (See section <ref>).

Since the four velocity u^α is normalised, it follows f̃^α is normal to it by proper time differentiation 

    u^α u_α = -1   ⟹  f̃^α u_α = 0 .

Evaluating (<ref>) along the radial direction, solving (<ref>) for f̃^ and plugging it into (<ref>) yields

    d^2r̃/dτ̃^2 + Γ^r_ρσdx̃^ρ/dτ̃dx̃^σ/dτ̃ = -f̃^ϕũ_ϕ + f̃^t̃ũ_t̃/ũ_r̃ .

The left hand side

    Γ^r̃_ρσdx̃^ρ/dτ̃dx̃^σ/dτ̃ =  Γ^r̃_ r̃(dr̃/dτ̃)^2 + Γ^r̃_ ϕϕ(dϕ/dτ̃)^2 + 
    Γ^_ t̃t̃(dt̃/dτ̃)^2 + 2Γ^_ t̃ϕdt̃/dτ̃dϕ/dτ̃

is computed using the Kerr Christoffel symbols and the geodesic equations (<ref>)-(<ref>)

    Γ^r̃_ρσdx^ρ/dτ̃dx̃^σ/dτ̃    = 3(aẼ - L̃)^2/r̃^4 - a^2(Ẽ^2 - 1) - L̃^2/r̃^3 + 1/r^2
        = -1/2∂ G/∂r̃

To evaluate the right hand side, f̃^r̃, we first notice the existence of two Killing vectors : ξ^μ = ∂/∂t̃ and ψ^μ = ∂/∂ϕ, associated with time and angular translational invariance, respectively. There exists a conserved charge associated with each :

    Ẽ = -ξ^μu_μ,   L̃ = ψ^μu_μ.

It follows from Eq.(<ref>) that u_ϕ = L̃ and u_t̃= -Ẽ. Finally, we relate the proper time derivatives of these charges with the forcing terms in (<ref>). For example, consider  the proper time derivative of Ẽ 

    -dẼ/dτ̃   = u^β∇_β(ξ^αu_α) 
       = ξ^α(u^β∇_βu_α) + u^αu^β(∇_βξ_α)
        = f̃_t̃ ,

where Killing's equation was used in the last step. A similar calculation leads to dL̃/dτ̃ = f̃_ϕ. Solving the two equations dL̃/dτ̃ and dẼ/dτ̃ for f̃^ϕ and f̃^t̃ gives

    f̃^ϕ   = -1/Δ(g_t̃t̃dL̃/dτ̃ + g_ϕt̃dẼ/dτ̃) , 
    f̃^t̃   = 1/Δ(g_t̃ϕdL̃/dτ̃ + g_ϕϕdẼ/dτ̃).

where we used the identity (g_ϕt̃)^2 - g_ϕϕg_t̃t̃ = Δ for Δ = r̃^2 - 2r̃ + a^2. Since u_r̃ = g_r̃r̃(dr̃/dτ̃), it follows the right hand side of (<ref>) is

    f̃^r̃ = 1/Δ g_r̃r̃(dẼ/dτ̃[g_ϕt̃L̃ + g_ϕϕẼ] +  dL̃/dτ̃[g_t̃t̃L̃ + g_t̃ϕẼ])(dr̃/dτ̃)^-1

Noticing that

    1/Δ g_r̃r̃(g_ϕt̃L̃ + g_ϕϕẼ)    = 1
    /2∂ G/∂Ẽ
    1/Δ g_r̃r̃( g_t̃t̃L̃ + g_t̃ϕẼ)    = 1/2∂ G/∂L̃

we reach the desired conclusion

    f̃^ = 1/2(d/dτ̃∂ G/∂ + d/dτ̃∂ G/∂L̃)(dr̃/dτ̃)^-1.













§ THE SCALING OF ECCENTRICITY

We are interested in determining how the eccentricity e scales near the separatrix (analogue of ISCO). Kennefick derives an equation incorporating the rate of change of eccentricity with various orbital parameters

    ė/e   = -j(ṽ)[Γ - h(ṽ)Ẽ̇]
    
        j(ṽ)    = (1 + aṽ^3)(1 - 2ṽ^2 + a^2ṽ^4)^1/2/ṽ^2(1 - 6ṽ^2 + 8aṽ^3 - 3a^2ṽ^4
    
        h(ṽ)    = ℋ(ṽ)/𝒢(ṽ)

with

    ℋ = 1 - 12ṽ^2 + 66ṽ^4 - 108ṽ^6 + aṽ^3 + 8a^2ṽ^4
     - 72aṽ - 20a^2ṽ^6 + 204aṽ^7 + 38a^3ṽ^7
    
        -42a^2ṽ^8 - 9a^4ṽ^8 - 144a^3ṽ^9
     + 116a^4ṽ^10 - 27a^5ṽ^11


    𝒢 = 2(1 + aṽ^3)(1 - 2ṽ^2 + a^2ṽ^4)^2(1 - 6ṽ^2 + 8aṽ^3 - 3a^2ṽ^4).

with ṽ = 1/√(r̃) and Ẽ̇ is the (dimensionless) circular energy flux. What is important to appreciate is that both Γ and Ẽ̇ depend on data extracted from the numerical resolution of the Teukolsky problem. Kennefick later shows that close to ISCO the dominant term sourcing ė is Ẽ̇. Hence, following this observation, we approximate the rate of change in the eccentricity by

    ė/e∼ j(v)h(v) Ė_0 .     (close to ISCO)

Consider the following expansion 

    a    ∼ 1-1/2ϵ^2 , 
    r̃   ∼ 1 + R ϵ^2/3,    R≥2^1/3 .

Using these expansions, we can analyse the behaviour of <ref> in the limit as ϵ→ 0. Consider h(ṽ) for ϵ→ 0

    h(ṽ)   →9/164-R^3/R(R^3-2) ϵ^-2/3 , 
    
      j(ṽ)   → 2√(3)(R^3/R^3 - 2),

which gives the result valid close to the ISCO

    ė/e∼9√(3)/8 R^2(4-R^3)/(R^3-2)^2 Ẽ̇/ϵ^2/3 + 𝒪(1/[R^3 - 2]^3).

Using the circular relationship for the (dimensionless) energy of the particle at some r̃, one can compute

    ṙ/r = - 2(1-3ṽ^2+2aṽ^3)^3/2/ṽ^2(1-6ṽ^2+8aṽ^3-3a^2ṽ^4)Ẽ̇

and performing the same expansion, we find

    ṙ̃̇/r̃∼ -3√(3)/2 R^3/R^3-2Ė_0 + 𝒪(ϵ^2/3/[R^3 - 2]^3).

Finally, by taking the ratio of both quantities in (<ref>) and (<ref>) we find a relation which only depends on geodesic properties (absent of radiation reactive effects)

    ė/e/ṙ̃̇/r̃ = -3/44 - R^3/R(R^3 - 2)ϵ^-2/3 + 𝒪(1/[R^3 - 2]^3).

By using the change of variables r̃ = 1 + Rϵ^2/3, the quantity ϵ^-2/3 cancels from either side leaving an integral which can be done analyically. The result is

    e = e_0(R^6/R^3-2)^1/4

where e_0 is some initial eccentricity from the adiabatic inspiral regime (far from transition regime). Near R=2^1/3 one finds that

    e = ê(R - 2^1/3)^-1/4


which diverges for R = 2^1/3. There will be extra terms of greater negative powers of (R - 2^1/3) with higher (positive) powers of ϵ in Eq.(<ref>). We remark on this in a moment. First of all, we know that during the transition regime for ϵ∼η, we have that  (r̃-_)∼η^2/5ϵ^4/15 to show that (R-2^1/3)∼(η/ϵ)^2/5 and similarly for the regime ϵ≪η one finds that (R - 2^1/3)∼(η/ϵ)^2/3. Finally, since since ê is an eccentricity sourced from the radial evolution of the earlier inspiral, we must have that ê∼ṙ̃̇∼η. Combining these results provide a scaling relationship on the eccentricity

    e ∼η·(η/ϵ)^-1/10   ϵ∼η
    η·(η/ϵ)^-1/6   ϵ≪η.

This implies that once we enter the transition regime, the eccentricties induced are low enough in order so that they contribute sub-dominantly. To see this, we know that for approximately a circular orbit, one has correction to the circular flux balance law of the form

    Ẽ̇ - Ω̃L̇̃̇∼ηṙ + η e^2 + higher order corrections

The first term on the right hand side takes into account the past history of the orbit. The change in the fluxes Ẽ̇ and L̇̃̇ are η quantities and change by fractional amounts of the radial velocity ṙ̃̇. If the eccentricity scales as Eq.(<ref>) suggests, then we see that corrections to the circular relation <ref> are indeed subdominant and will not contribute throughout the transition regime. An important question to ask is when the eccentricity corrections are of the same order as the other terms in the master equation. Let us consider the ϵ∼η regime, we question when 

    d^2r̃/dτ̃^2∼η^4/3∼η· e^2.

and, given our solutions before, consider e∼η·η^-γ. What this suggests that if γ > 5/6 then the eccentricity induced during the transition regime will significantly contribute to the motion and must be taken into account. Notice that in the two near-extremal regimes we have considered, γ < 5/6 so the eccentricity is not significant throughout the transition regime. Let us now consider an inspiral then transition to plunge for a Schwarzschild black hole (a = 0). For a = 0 and around r̃ = 6, a short calculation shows that

    1/ede/dr = -9/r - 6 + 𝒪(r̃ - 6),

Which has solution

    e = ê_0(r - 6)^-9.

Applying similar arguments as before, the initial eccentricity is sourced through the adiabatic inspiral giving ê∼ṙ∼η. We also know that r̃ - _ = r̃ - 6 ∼η^2/5 yielding e ∼η·η^-18/5. Obviously γ = 5/6 < 18/5 which implies that the eccentricity correction are significant throughout the transition regime. A natural question; at what spin can we ignore eccentricity corrections during the transition?















 















IEEEtran


