Gravitational self-force on a particle in eccentric orbit around a Schwarzschild black hole

We present a numerical code for calculating the local gravitational self-force acting on a pointlike particle in a generic (bound) geodesic orbit around a Schwarzschild black hole. The calculation is carried out in the Lorenz gauge: For a given geodesic orbit, we decompose the Lorenz-gauge metric perturbation equations (sourced by the delta-function particle) into tensorial harmonics, and solve for each harmonic using numerical evolution in the time domain (in $1+1$ dimensions). The physical self-force along the orbit is then obtained via mode-sum regularization. The total self-force contains a dissipative piece as well as a conservative piece, and we describe a simple method for disentangling these two pieces in a time-domain framework. The dissipative component is responsible for the loss of orbital energy and angular momentum through gravitational radiation; as a test of our code we demonstrate that the work done by the dissipative component of the computed force is precisely balanced by the asymptotic fluxes of energy and angular momentum, which we extract independently from the wave-zone numerical solutions. The conservative piece of the self-force does not affect the time-averaged rate of energy and angular-momentum loss, but it influences the evolution of the orbital phases; this piece is calculated here for the first time in eccentric strong-field orbits. As a first concrete application of our code we recently reported the value of the shift in the location and frequency of the innermost stable circular orbit due to the conservative self-force [Phys. Rev. Lett. 102, 191101 (2009)]. Here we provide full details of this analysis, and discuss future applications.

Figure 1

Parameter space for bound geodesics in Schwarzschild spacetime. The (dimensionless) “semilatus rectum” $p$ and “eccentricity” $e$ are defined in Eq. (2.3). Bound geodesics have $e\ge 0$ and $p>6+2e$. Points along the separatrix $p=6+2e$ represent marginally unstable orbits. Stable circular orbits lie along the axis $e=0$ for $p\ge 6$. The point $(p,e)=(6,0)$ is the ISCO.

Figure 2

Numerical domain: a staggered $1+1$-dimensional mesh in null coordinates $v$, $u$. $r_{*}$ is the standard Schwarzschild “tortoise” radial coordinate. The dotted line represents the trajectory of a typical eccentric orbit. In actual implementation the mesh is, of course, much finer than it is depicted here.

Figure 3

Grid points involved in constructing our finite-difference scheme. The point C is the center of the cell over which we integrate the field equations, as described in the text. The dimensions of each grid cell are $\Delta v\times \Delta u=h\times h$.

Figure 4

Same as in Fig. 3, but now point C is located near the particle’s worldline, represented by the dashed line. The finite-difference scheme for this case is described in the text under “near-orbit cell.”

Figure 5

Illustration of the four cases considered in formulating the finite-difference equation for grid cells traversed by the particle’s worldline (here represented by dashed lines). The center of the cell is at $(u,v)=(u_c,v_c)$, and in each case we indicate the $u$, $v$ coordinates of the two points where the particle enters/leaves the cell.

Figure 6

The gravitational SF for a sample of eccentric orbits. We plot the temporal and radial components of the SF as functions of $t$ for the sample parameters $(p,e)=(10,0.2)$, (10, 0.5), and (10, 0.7). The solid and dashed lines show $F^t$ and $F^r$, respectively. The $\phi$ component of the SF can be trivially obtained using the orthogonality condition $F^{\alpha }{u}_{\alpha }=0$. In all graphs $t=0$ corresponds to a periapsis passage. The radial periods for the $e=0.2$, 0.5, 0.7 plots are, respectively, $T_r\simeq 328M$, $434M$, $693M$. We have cut out from these plots the early part of the numerical solution, where nonphysical initial spurious waves dominate; the plots display only the later, stationary part of the solutions, after the spurious waves have dissipated away. Note the small retardation manifest in the amplitude of the total SF (with reference to the orbital phase). A similar retardation is observed in the scalar and EM cases [12, 23].

Figure 7

The gravitational SF for a zoom-whirl orbit with parameters $(p,e)=(7.4001,0.7)$. In the upper panel we show the (contravariant) $t$, $r$, and $\phi$ components of the SF as functions of time $t$ along the orbit. The lower panel shows the radial motion of the particle, for reference. $t=0$ is periapsis, and the radial period is $\sim 789M$. While $F^t$ and $F^{\phi }$ quickly settle to a constant value during the quasicircular whirl episode (as one would expect), the radial component exhibits a peculiar linear behavior. This behavior is analyzed in the text (cf. Fig. 8).

Figure 8

Conservative and dissipative pieces of the radial SF for zoom-whirl orbits—left and right panels, respectively. Shown are results for four orbits with the same eccentricity ($e=0.7$) but with increasing proximity to the separatrix $p=6+2e$.

Figure 9

Derivation of the SF coefficients $F_{1\mathrm{is}}^r$ and $F_{\phi \mathrm{is}}^1$ entering the ISCO-shift formula (5.17), using method I. Plotted are the numerical values of ${\widehat{F}}_1^r(p,e)$ and ${\widehat{F}}_{\phi }^1(p,e)$ [see Eqs. (5.23, 5.24)] for a sequence of points in the $p\mathrm{\text{−}}e$ plane approaching the ISCO along the curve $p=6+\sqrt{e}$. Error bars indicate the estimated numerical error, which is evaluated as explained in the text. The solid curves are weighted polynomial fits used to extrapolate the values of ${\widehat{F}}_1^r$ and ${\widehat{F}}_{\phi }^1$ to $p\to 6$. The extrapolated values are the desired coefficients $F_{1\mathrm{is}}^r$ and $F_{\phi \mathrm{is}}^1$.

Figure 10

Derivation of $F_{1\mathrm{is}}^r$ and $F_{\phi \mathrm{is}}^1$, using method II. Plotted are the numerical values of $F_{\cos }^r(p)$ and $F_{\phi }^{\sin }(p)$ [see Eqs. (5.33, 5.35)] for a sequence of points approaching the ISCO along the $p$ axis in the $p\mathrm{\text{−}}e$ plane. Error bars indicate the estimated numerical error, which is evaluated as explained in the text. The solid curves are weighted polynomial fits used to extrapolate the values of $F_{\cos }^r$ and $F_{\phi }^{\sin }$ to $p\to 6$. The extrapolated values are the desired coefficients $F_{1\mathrm{is}}^r$ and $F_{\phi \mathrm{is}}^1$.