The Orbital Lense-Thirring Precession in a Strong Field

TL;DR

This paper analytically demonstrates that the orbital angular momentum of particles around a Kerr black hole precesses in a way that closely resembles a circle, extending the weak-field Lense-Thirring precession concept to strong gravitational fields.

Contribution

It provides an exact analytical description of orbital angular momentum evolution in Kerr spacetime, confirming the near-circular precession pattern for a wide class of orbits and deriving parameters for unstable spherical orbits.

Findings

Orbital hodographs are nearly circular with deviations under 10%.

Nutation accurately approximates deviations, especially for larger deviations.

Lense-Thirring precession remains valid in strong gravitational fields.

Abstract

We study the exact evolution of the orbital angular momentum of a massive particle in the gravitational field of a Kerr black hole. We show analytically that, for a wide class of orbits, the angular momentum's hodograph is always close to a circle. This applies to both bounded and unbounded orbits that do not end up in the black hole. Deviations from the circular shape do not exceed $\approx10\%$ and $\approx7\%$ for bounded and unbounded orbits, respectively. We also find that nutation provides an accurate approximation for those deviations, which fits the exact curve within $\sim 0.01\%$ for the orbits of maximal deviation. Remarkably, the more the deviation, the better the nutation approximates it. Thus, we demonstrate that the orbital Lense-Thirring precession, originally obtained in the weak-field limit, is also a valid description in the general case of (almost) arbitrary exact orbits. As a by-product, we also derive the parameters of unstable spherical timelike orbits as a function of their radii and arbitrary rotation parameter $a$ and Carter's constant $Q$. We verify our results numerically for all the kinds of orbits studied.

Tables

Table 1. Table 1: Major notation.

Notation	Explanation
$(t,r,\theta ,\varphi )$	the Boyer–Lindquist coordinates of the Kerr metric
$M$	the mass of a Kerr black hole
$a$	the angular momentum of a Kerr black hole (in the units of $M$)
	frequently used combinations
$E$	the conserved energy of a particle (in the units of the particle’s rest mass)
$L_x,L_y,L_z$	component of the orbital angular momentum of a particle
$Q$	Carter’s constant
$\dot{\left(\dots \right)}\equiv \frac{d\left(\dots \right)}{d\tau }$	derivative with respect to proper time $\tau$ of a particle
$\gamma$: $d\gamma =Ed\tau /{\rho }^2$	auxiliary parameter along geodesics
$ℒ\equiv \sqrt{L_x^2+L_y^2}$	radius of an angular momentum’s hodograph in $(L_x,L_y)$-plane
$\Delta Q\equiv a^2\left(E^2-1\right){\cos }^2\theta$	variation of ${ℒ}^2$
$r_{\mathrm{isco}}^{\pm }$, $r_{\mathrm{mb}}^{\pm }$, $r_{\mathrm{ph}}^{\pm }$	radii of innermost stable circular orbit, marginally stable orbit, and photon orbit, respectively (plus and minus stand for co- and counterrotating orbits, correspondingly)
	auxiliary notation
$y_0$	inverse radius of a spherical orbit with $L_z-aE=x\sqrt{Q}=0$
$y_{\ast }$	inverse radius of a spherical orbit with $L_z=0$ in the limit $Q\to +\infty$

Table 2. Table 2: Initial conditions for orbits with noncircular hodographs.

Figure	$a$	$E$	$L_z$	$Q$	$r_0$	${\theta }_0$	${\varphi }_0$
2	1	0.59	1.18026	0.00437	1.03	$\pi /2$	0
3 (left)	0.999	100	207.8128	0.1728	10	$\pi /2$	0
3 (right)	0.999	19.5232	40.5399	50	10	$\pi /2$	0
5	1	0.59	1.18026	0.00437	1.03	$\pi /2$	0
6	0.999	19	40.5399	50	10	$\pi /2$	0

Equations

$$L_x$$

$$L_y$$

$$L_z$$

$$\mathrm{d}{s}^2$$

$$Q=p_{\theta }^2+{\cos }^2\theta \left[a^2(1-E^2)+\frac{L_z^2}{{\sin }^2\theta }\right].$$

$${\rho }^4{\dot{r}}^2$$

$${\rho }^4{\dot{\theta }}^2$$

$${\rho }^2\dot{\varphi }$$

$${\rho }^2\dot{t}$$

$$U_{\theta }(\theta ,L_z,E)={\cos }^2\theta \left[a^2(1-E^2)+\frac{L_z^2}{{\sin }^2\theta }\right],$$

$$R(r)=rS(E-U_{-})(E-U_{+}),$$

$$U_{\pm }(r,L_z,Q)\equiv \frac{2a{L}_z\pm \sqrt{\Delta \left[r^2{L}_z^2+(r+Q/r)S\right]}}{S},$$

$$S\equiv r^3+a^{2}r+2a^2,$$

$${\mathcal{L}}^2\equiv L_x^2+L_y^2=Q+\Delta Q,\Delta Q=a^2(E^2-1){\cos }^2\theta .$$

$${{\cos }^2\theta |}_{\max }=\frac{Q}{L_z^2+a^2(1-E^2)}+\mathcal{O}\left(Q^2\right),Q\to 0,$$

$$\frac{|\Delta Q|}{Q}=\frac{a^2(1-E^2)}{L_z^2+a^2(1-E^2)}+\mathcal{O}\left(Q\right)\lesssim 1,$$

$$E(y,a)$$

$$L_z(y,a)$$

$$r_{\mathrm{isco}}^{\pm }$$

$$Z_1$$

$$Z_2$$

$$r_{\mathrm{mb}}^{\pm }$$

$$r_{\mathrm{ph}}^{\pm }$$

$$\frac{a^2(1-E^2)}{L_z^2+a^2(1-E^2)}\equiv 1-\frac{1}{1+a^2(1-E^2)/L_z^2}$$

$$\frac{|\Delta Q{|}_{\mathrm{ma}{\mathrm{x}}^{+}}}{Q}=\frac{1}{3}\Rightarrow \frac{{\mathcal{L}}_{\max }-{\mathcal{L}}_{\min }}{{\mathcal{L}}_{\max }+{\mathcal{L}}_{\min }}={\left(\sqrt{3}-\sqrt{2}\right)}^2\approx 0.10.$$

$$\frac{|\Delta Q{|}_{\mathrm{ma}{\mathrm{x}}^{-}}}{Q}=\frac{1}{243}\Rightarrow \frac{{\mathcal{L}}_{\max }-{\mathcal{L}}_{\min }}{{\mathcal{L}}_{\max }+{\mathcal{L}}_{\min }}\approx 0.0010.$$

$$\frac{\Delta Q_{\mathrm{ma}{\mathrm{x}}^{+}}}{Q}$$

$$\frac{\Delta Q_{\mathrm{ma}{\mathrm{x}}^{-}}}{Q}$$

$$R_Q(y)$$

$$\frac{1}{2}\frac{\mathrm{d}{R}_Q}{\mathrm{d}y}$$

$$\frac{1}{2}\frac{{\mathrm{d}}^2{R}_Q}{\mathrm{d}{y}^2}$$

$$ϵ\equiv \frac{E}{\sqrt{Q}},x\equiv {\lambda }_z-aϵ,{\lambda }_z(y,a)\equiv \frac{L_z}{\sqrt{Q}}.$$

$$x=-\frac{2a^2{y}^2-3y+1}{\sqrt{4a^2{y}^3-9y^2+6y-1}},1/r_{\mathrm{ph}}^{-}<y<1/r_{\mathrm{ph}}^{+}.$$

$$ϵ(y,a)$$

$${\lambda }_z(y,a)$$

$$\frac{{\mathrm{d}}^2{R}_Q}{\mathrm{d}{y}^2}$$

$$\frac{\Delta Q}{Q}=\frac{a^2(E^2-1){{\cos }^2\theta |}_{\max }}{Q}=a^2({ϵ}^2-q){{\cos }^2\theta |}_{\max }=$$

$$=-\frac{1}{2}\left({\lambda }_z^2+1-a^2({ϵ}^2-q)-\sqrt{{\left[{\lambda }_z^2+1-a^2({ϵ}^2-q)\right]}^2+4a^2({ϵ}^2-q)}\right)=$$

Full text

K****eywords Black holes  $\cdot$ Geodesics  $\cdot$ Lense–Thirring precession

1 Introduction

Solutions to the two-body problem in Newtonian celestial mechanics are well known and exhaustively described by two properties: the orbits remain in one plane and fall in one of the three classes: ellipse, parabola, or hyperbola. As general relativity (GR) comes into play, those properties are violated (except for special cases). In the gravitational field of a Schwarzschild black hole, although the orbits are still confined to a plane, their shape becomes increasingly complicated as dimensionless parameter $GM/(r_{\rm m}c^{2})$ grows ($G$ is the gravitational constant, $c$ is the speed of light; $M$, the mass of the attracting center, and $r_{\rm m}$ is a minimal distance to that center). When the parameter is small, for bounded orbits the approximation of a precessing ellipse works well (e.g. [1] or see [2] for derivation in terms of the Laplace–Runge–Lenz vector). However, in the fully relativistic regime this approximation is not valid anymore.

In the field of a rotating Kerr black hole, nor remain the orbits in a single plane (one exception is equatorial orbits). We say that a particle’s orbit is confined to one plane if the particle’s total angular momentum $\mathbf{L}\equiv(L_{x},L_{y},L_{z})$ is constant, projections of the angular momentum being

[TABLE]

where $\phi$ and $\theta$ are the angles of either the conventional curvature coordinates of the Schwarzschild metric [3, 4], or Boyer–Lindquist coordinates of the Kerr metric [5], and $p_{\phi}$ , $p_{\theta}$ are the respective covariant components of the particle’s momentum. This is a convenient representation of an apparent shape of the orbit plotted against “Cartesian” axes $(x,y,z)$ which are related to the angles $\theta,\phi$ and to Schwarzschild or Boyer–Lindquist radius $r$ through the conventional formulae of a transformation to spherical coordinates 111From the viewpoint of intrinsic geometry, the plane is, of course, not flat and has a nonvanishing curvature. Also note that the “Cartesian” coordinates are different from the Kerr–Schild coordinates.. Although the components $L_{x}$, $L_{y}$, and $L_{z}$ do not come from any Killing vector in the Kerr space–time, they do quantify the measurable effect of the change of a particle’s orbital plane with respect to a distant observer. As we are to show, the evolution of this plane in terms of these components is quite simple (though the derivation of the result is rather lengthy).

Nonetheless, as long as the Kerr’s field is weak enough, that out-of-the-plane motion can be attributed to the Lense–Thirring precession of the angular momentum [6, 7], that is, of the orbital plane 222The precession of the axis of a gyroscope [8] which is also referred to as the Lense–Thirring effect is not pertinent to this Paper (see also a discussion of certain analogies between the two effects in Conclusion).. When the field is fully relativistic, the orbits can be quite intricate (see Figure 1 for an example) and, in particular, essentially spherical [9, 10, 11]. Also, a number of rigorous results on geodesics in Lorentz manifolds were obtained by means of geometric analysis, e.g. [12, 13, 14]. On the other hand, a visual reduction of these seemingly unarranged orbits to a superposition of simpler motions would contribute to intuitive understanding of the celestial mechanics in a Kerr black hole’s gravitational field.

In this Paper, we propose such a reduction and show that any orbit of a massive particle that does not end up in the black hole results from the superposition of a motion in a plane and the precession and nutation of that plane. Although, strictly speaking, this reduction is not exact, the superposition approximates the orbit within $\sim 0.01\%$. In order to show it, we study the exact evolution of orbital angular momentum and then confirm our results numerically. We find that for most orbits the precession alone approximates the exact motion reasonably well. The discrepancy can be further accounted for by introducing nutation, the “second-order precession” around a precessing axis. And as it turns out, the larger the discrepancy, the better the nutation accounts for it. The superposition of precession and nutation is reminiscent of Ptolemaic epicycles, but remarkable is the fact that there is no need to consider higher-order precessions. Therefore, unlike the orbits themselves, their orbital momenta demonstrate a relatively simple behavior.

The Paper is organized as follows. In Sec. 2 we review the geodesic equations of the Kerr space–time and provide a formula for the radius of a hodograph of the orbital angular momentum (recall that a hodograph is the locus of the end of a vector). Bounded and unbounded orbits are considered in Sec. 3 and Sec. 4, respectively. In Sec. 5 we introduce nutation to account for deviations of the hodographs from a circle and estimate its magnitude in Sec. 6. In Sec. 7 we discuss our results.

In what follows, the signature $(-\,+\,+\,+)$ and units $G=c=1$ are used. All lengths and times are in units of $M$, the mass parameter of a black hole. Major notations are summarized in Table 1 while minor notations are given along the text.

2 Geodesic equations

Consider the Kerr space–time in Boyer–Lindquist coordinates, e.g. [5]:

[TABLE]

As is known [15], a massive particle’s free fall in this space–time is characterized by a few constants of motion involving covariant components of 4-momentum $p_{\mu}\equiv(p_{t},p_{r},p_{\theta},p_{\phi})$. Those are energy at infinity $E=-p_{t}$ conventionally measured in units of the particle’s rest mass, a component of the orbital angular momentum $L_{z}=p_{\phi}$ that is parallel to the rotation axis of the black hole, and the so-called Carter’s constant [16] defined as

[TABLE]

In addition, conserved is the norm of the 4-momentum vector, $p_{\mu}p^{\mu}=-1$.

Given a set of parameters $(E,L_{z},Q,a)$, geodesic equations can be partially integrated to yield [15]:

[TABLE]

where the dot stands for a derivative with respective to proper time $\tau$.

In this form, there is a partial reduction of the particle’s motion into independent evolution along $r$ and $\theta$. In the $\theta$-direction, the particle is moving in angular potential

[TABLE]

with the Carter’s constant $Q$ playing the role of “angular energy”. Also, since function $R(r)$ can be factorized [9] as follows 333Another factorization proposed in [17] can be convenient if the sum $Q+(L_{z}-aE)^{2}$ is used instead of the Carter’s constant.:

[TABLE]

with

[TABLE]

the particle’s motion in the $r$-direction occurs in the radial potential $U_{+}(r,L_{z},Q)$.

Since $L_{z}$ is an integral of motion, the evolution of the angular momentum is completely given by the latter’s hodograph in the $(L_{x},L_{y})$-plane. Consider the squared radius of this hodograph, which, from (1)–(3) and (5), reads:

[TABLE]

Note that the first term in the right-hand side is constant. Therefore, $\mathcal{L}$ varies only due to the second term $\Delta Q$, and the ratio of the two terms determines whether the hodograph will deviate from a circle significantly.

We proceed by separately considering bounded and unbounded orbits and assume that neither of them enter the black hole.

3 Bounded orbits

In this case, $E<1$ and $Q\geq 0$ (see [9] and eq. (5)). Hence, $|\Delta Q|<1$ and, if $Q>>1$, a hodograph’s radius $\mathcal{L}$ is approximately constant, i.e. a circle is a good approximation for the hodograph.

If $Q<<1$, $Q$ and $|\Delta Q|$ may be comparable. Indeed, as long as $E<1$, $\theta=\pi/2$ is a stable equilibrium point [17] of angular potential $U_{\theta}$. Therefore, for a given $Q$, the maximal value of $\cos^{2}{\theta}$ is on the order of $Q$, see (7). Namely,

[TABLE]

such that

[TABLE]

where $L_{z}$ and $E$ are the parameters of equatorial orbits $Q=0$.

Seemingly, in the last case $|\Delta Q|$ could be equal to $Q$ when $L_{z}=0$ (even exactly, according to (7)). However, this is not realized as long as the orbit is bounded. Moreover, in order to estimate the r.h.s. of (15) from above, it is sufficient to consider orbits of constant $r$ (and equatorial ones in the case in question), because, for a given $L_{z}$, to maximize the r.h.s., $E$ must be as low as possible, i.e. equal to the minimum of $U_{+}(r)$.

Recall that the parameters of equatorial circular orbits are given by relations [18, 19]:

[TABLE]

where $y\equiv 1/r$, $r$ being the radius of a circular orbit, and the upper and lower signs stand for co- and counterrotation, correspondingly.

These orbits are divided into the following subclasses:

• •

$y\in(0,1/r_{\rm isco}^{\pm})$, $r_{\rm isco}^{\pm}$ is the radius of the Innermost Stable Circular Orbit (hereafter, “$\pm$” in the index also stand for co- and counterrotation, respectively); orbits are stable,

• •

$y\in(1/r_{\rm isco}^{\pm},1/r_{\rm mb}^{\pm})$, $r_{\rm mb}^{\pm}$ is the radius of the Marginally Bound orbit; orbit are unstable with $E<1$ (an inward perturbation sends the particle into the black hole and an outward perturbation accompanied by a decrease in energy sends the particle into a noncircular bounded orbit),

• •

$y\in(1/r_{\rm mb}^{\pm},1/r_{\rm ph}^{\pm})$, $r_{\rm ph}^{\pm}$ is the radius of the circular PHoton orbit; orbits are unstable with $E>1$ (an inward perturbation sends the particle into the black hole and, under an outward perturbation, the particle escapes to infinity).

The characteristic radii are ordered, $r_{\rm isco}^{\pm}\geq r_{\rm mb}^{\pm}\geq r_{\rm ph}^{\pm}$, the equalities holding for an extremely rotation black hole, $a=1$, and corotating orbits, and read:

[TABLE]

Note that both $E(y,a)$ and $L_{z}^{2}(y,a)$ as functions of $y$ have minima at $y=1/r_{\rm isco}^{\pm}$ (see Appendix B). Therefore, the leading term in the r.h.s. of (15) has a maximum at that point, because

[TABLE]

and $(1-E^{2})>0$ for bounded orbits as well as $L_{z}^{2}>0$. Also, the value of this maximum grows as $a\to 1$ and is the highest for corotating orbits as $r_{\rm isco}^{+}\to 1$. In that limit  [18], $E\to 1/\sqrt{3}$ and $L_{z}\to 2/\sqrt{3}$, and from (13) and (15) we obtain:

[TABLE]

For counterrotating orbits, $E\to 5\sqrt{3}/9$, $L_{z}\to-22\sqrt{3}/9$, and

[TABLE]

Thus, the hodographs of bounded orbits deviate from a circle by no more than $\approx 10\%$. Figure 2 shows such a close-to-extreme hodograph.

4 Unbounded orbits

A particle on an unbounded orbit that escapes to infinity after having been reflected by radial potential $U_{+}$ has $E\geq 1$. Therefore, there are two possibilities for $Q$ and $\Delta Q$ to be comparable. The first possibility is similar to the case $Q<<1$ for bounded orbits. The second possibility arises, because one might expect $\Delta Q>>Q$ if $E>>1$. Let us consider these possibilities in more detail.

$Q<<1$. In the range between $1/r_{\rm mb}^{\pm}$ and $1/r_{\rm ph}^{\pm}$, function $a^{2}(E^{2}(y,a)-1)$ is less than function $L_{z}^{2}(y,a)$, their ratio being maximal at $y=1/r_{\rm ph}^{\pm}$ (notice that either function tends to infinity at that point). Therefore, $\theta=\pi/2$ continues to be a stable center of the angular potential $U_{\theta}(\theta)$.

Also, that maximal value is the highest for an extremely rotating black hole. In the case of corotation [18], $(E/L_{z})^{2}\to 1/4$ at $y=1/r_{\rm ph}^{+}$ as $a\to 1$. Hence, similarly to (22),

[TABLE]

For counterrotating orbits, $(E/L_{z})^{2}\to 1/49$ at $y=1/r_{\rm ph}^{-}$ as $a\to 1$. Therefore,

[TABLE]

$Q>>1$. Consider also $E>>1$, because only in this case may we expect $\Delta Q\sim Q$. Hence, $\Delta Q/Q\lesssim a^{2}E^{2}/Q$ .

On the other hand, the upper limit for $E$ is the energy of an unstable spherical orbit, corresponding to a maximum of potential $U_{+}(r)$, eq. (2), or a minimum of function $R(r)$, eq. (6). To find the latter, we divide $R(r)$ by $Q$ and neglect terms $\propto 1/Q$. Then, the conditions for a nonstable spherical orbit read:

[TABLE]

where

[TABLE]

Solving the second equation for $\epsilon$ and substituting to the first equation, we find

[TABLE]

Note that the numerator is proportional to a derivative of the expression inside the square root. Also, recall that $y_{\pm}\equiv 1/r^{\pm}_{\rm ph}$ are solutions to equation $1-3y\pm 2ay^{3/2}=0$, i.e. they are roots of the denominator.

Other variables are:

[TABLE]

Notice that, if $0<a<1$, $\epsilon$ and $\lambda_{z}$ diverge as $y\to y_{\pm}$ while their ratio tends to a finite limit. The same property holds exactly in the case of a sufficiently small $q>0$, with the finite limit being independent of $q$ (see Appendix C for a proof). This allows us to evaluate exactly the maximal deviation from a circle for non-equatorial unstable spherical orbits with $Q>>1$. Namely,

[TABLE]

In this derivation we have used relations (92) and (93) from Appendix C to make sure that

[TABLE]

Finally, we are in a position to find the maximal deviation. Using (31) and (32), we obtain for corotating orbits:

[TABLE]

For counterrotating orbits, the largest deviation is not necessarily given by the respective limit taken at $y\to y_{-}$. Instead, orbits with $\lambda_{z}<0$ and $\lambda_{z}$ close to zero may give rise to larger deviations, at least for small $q$. For the sake of completeness, let us consider both cases.

At $y\to y_{-}$, provided that $\epsilon\to+\infty$ and $\lambda_{z}\to-\infty$, the answer reads:

[TABLE]

For $\lambda_{z}<0,|\lambda_{z}|<<1$, it is convenient to write the answer in a parametric form. That is, a spherical orbit of radius $1/y$ with $\lambda_{z}=0$ has

[TABLE]

where $y=y_{*}$ corresponds to $q=0$ and satisfies $1-3y+a^{2}y^{2}+a^{2}y^{3}=0$, cf. (32). Notice also that the denominator is positive and the numerator has exactly one root for $y\in(0,1)$, $a\neq 1$. The respective deviation is

[TABLE]

This deviation is maximal at $q=0$ and $a\to 1$, which implies $y_{*}=\sqrt{2}-1$. Hence,

[TABLE]

Concluding the case of unbounded orbits, their hodographs deviate from a circle by no more than $\approx 7\%$. Figure 3 shows such close-to-extreme hodographs.

5 Deviations from a circle

As we saw, the hodographs of the total angular momentum are nearly circular for both bounded and unbounded orbits that do not end up in the black hole. In any case the deviation from the circular shape does not exceed $\approx 10\%$ and $\approx 7\%$ for bounded and unbounded orbits, respectively. The circular shape implies precession while the deviations can be approximated by nutation (somewhat reminiscent of Ptolemaic epicycles). Interestingly, the more the deviation, the better nutation approximates it.

In the nutation approximation, the radius $\mathcal{L}$ in the $(L_{x},L_{y})$-plane is represented as a sum of two constant-length vectors (see Figure 4). Precession vector $\mathbf{D}$ rotates around the origin, whereas nutation vector $\mathbf{d}$ rotates around the end of the precession vector.

[TABLE]

where $\Omega_{\rm p}$ and $\Omega_{\rm n}$ are the precession and nutation frequencies, respectively; $\psi_{\rm p}$ and $\psi_{\rm n}$ are the corresponding phase shifts, and $\gamma$ is a parameter along the trajectory of a particle such as $\,{\rm d}\gamma=E\,{\rm d}\tau/\rho^{2}$. Thus, the squared radius is given by

[TABLE]

where $\Delta\Omega\equiv\Omega_{\rm p}-\Omega_{\rm n}$ and $\Delta\psi\equiv\psi_{\rm p}-\psi_{\rm n}$.

Now, let us show that (44) does provide an approximation for the actual $\mathcal{L}$. Indeed, consider equation of motion (7) in the $\theta$-direction and rewrite it in terms of $\epsilon$, $\lambda_{z}$, $q$, $\mu\equiv\cos{\theta}$, and the new parameter $\gamma$:

[TABLE]

In essence, this is the total energy of an anharmonic oscillator. Its degree of anharmonicity is given by the ratio of the third and second terms on the left-hand side:

[TABLE]

For both bounded and unbounded orbits with maximal deviation, the denominator is large while the numerator is on the order of one. Thus, the oscillator is approximately harmonic with the frequency

[TABLE]

Since the degree of anharmonicity is small, it can be taken into account by correcting the frequency. The answer is easy to find [20]:

[TABLE]

Assuming that a particle starts from the equatorial plane, we obtain $\mu(\gamma)=\mu_{\rm max}\sin{\omega\gamma}$, which results in

[TABLE]

That is, the hodograph’s motion is a sum of precession and nutation with $\Delta\Omega=2\omega$, cf. (44). The relative phase shift $\Delta\psi$ is fixed to [math] by choosing the equatorial plane as a starting point.

Therefore, either the deviation is small and the hodograph is approximated by a circle or the deviation is large and nutation almost perfectly accounts for a noncircular shape. Figures 5 and 6 show actual hodographs and their approximations by a sum of precession and nutation for bounded and unbounded orbits, respectively. The r.m.s. of the approximation’s relative error is $\lesssim 10^{-4}$.

6 The order of magnitude for the nutation

The Lense–Thirring precession frequency $\Omega_{\rm p}\sim GJ/(c^{2}r^{3})$, where $J$ is the angular momentum of the gravitating center and $r$ is the characteristic size of the orbit around that center, e.g. [3, § 40.7]. For Earth, the magnitude of the precession $\sim 100$ milliarcseconds per year, which was confirmed through measurements of the effect on Earth-orbiting satellites [22, 23].

To estimate the nutation frequency $\Omega_{\rm n}$, consider the angular velocity of vector $\mathcal{L}$ in the $(L_{x},L_{y})$-plane. Let $f\equiv\arctan{\frac{L_{y}}{L_{x}}}$. Then, by virtue of (42) and (43),

[TABLE]

Now, when $\mathcal{L}^{2}$ is maximal or minimal (the upper or lower sign in what follows),

[TABLE]

whence

[TABLE]

On the other hand, a direct differentiation of (1) and (2) yields

[TABLE]

which, after appropriate substitutions from (5), (7), and (8), results in

[TABLE]

At the starting point $\theta=\pi/2$, derivative $df/d\gamma$ is maximal (unbounded orbits) or minimal (bounded ones). In either case,

[TABLE]

Note that these relations give an order-of-magnitude estimate of the frequencies, where $r$ is a characteristic size of the orbit. For the unbounded orbits, the frequencies vary significantly as the particle approaches the black hole from infinity. This is because the orbital plane is approximately the same when the particle is far from the black hole. Figure 6 illustrates the fact: the start and the end of the trajectory manifest themselves by portions of the graphs with many integration points (in black).

In order to estimate the frequencies for Earth-orbiting satellites, we retain the first nonvanishing terms with respect to $a$ and express the frequencies with respect to proper time $\tau$:

[TABLE]

where we have used the fact that the size of the orbit $r>>1$, total energy $E\approx 1$, and $L_{z}^{2}+Q$ equals the squared total angular momentum $L^{2}$ up to terms $\propto a^{2}$. The precession frequency is obviously the Lense–Thirring frequency in a dimensionless form.

As for the nutation frequency, it is negative and its (dimensional) magnitude is

[TABLE]

where $T$ is the orbital period of a satellite on a circular orbit. That is, this frequency is twice the orbital frequency, and the nutation occurs in the direction opposite to that of the precession. For orbits with moderate inclinations $i$, the nutation effect would lead to slight oscillations of the orbital plane by an angle

[TABLE]

where we have substituted $|E^{2}-1|=1/L^{2}$ for circular orbits (in the dimensionless units). Or, in the dimensional form:

[TABLE]

where $l$ is the specific orbital angular momentum (per mass of the orbiting body) and $R_{0},T_{0}$ are the radius and rotation period of the gravitating center, respectively.

The linear displacement of the orbiting body due to the nutation is in turn

[TABLE]

where NS and WD stand for “neutron star” and “white dwarf”, respectively. For the latter, we took white dwarf G29-38 as a representative with $R_{0}=0.01R_{\astrosun}$ and $T_{0}=0.014$ days [24, 25] and we adopted that, first, a solar-mass neutron star with radius of $10$ km rotates at the maximum possible speed ($T_{0}\sim 1$ ms) and, second, the dimensionless moments of inertia of Earth, Jupiter, the neutron star, and the white dwarf are $0.35$.

7 Conclusion

We have studied the evolution of the orbital angular momentum of a massive particle following a geodesic path in the Kerr space–time. In doing so, we restricted ourselves to orbits, both bounded and unbounded, that do not end up in the black hole.

We have found that the end of the angular momentum approximately describes a circle, such that it deviates from the circle by no more than $\approx 10\%$ and $\approx 7\%$ for bounded and unbounded orbits, respectively. Curiously, for corotating orbits the relative maximal deviation takes a universal value independent of the type of orbit, $|\Delta Q|_{\rm max^{+}}/Q=1/3$. Second, we have demonstrated that nutation (precession around a precessing axis) accounts for this deviation within $0.01\%$ for the orbits with maximal deviation. These results imply that the Lense–Thirring precession of the orbital angular momentum, originally found in the weak-field limit, continues to be a valid description in the general case of (almost) arbitrary exact orbits. Also, introducing nutation makes this description highly accurate.

We have also estimated the nutation frequency and magnitude in the limit $a/r<<1$, where $r$ is a characteristic size of the orbit, equations (60)–(68). It turns out that the frequency is twice the orbital frequency in this limit. As for the magnitudes, the corresponding linear displacement for a satellite orbiting the Earth is $\sim 10^{-6}$ m, which eliminates any prospect of detecting the effect in the near future. For Jupiter and a white dwarf, this displacement is $\sim 10^{-2}$ m, which is, in principle, on the order of the available precision of a few millimeters [26]. However, the practical prospects are vague, given that there is no Jupiter’s gravity model similar to GRACE. For a neutron star, however, the linear displacement is considerable ($\sim 100$ m), which could be detected through the Doppler shift $\sim\Delta r_{\perp}^{\rm{NS}}/T_{0}$. The latter is comparable to the expected precision of third-generation spectohraphs ($\lesssim 10$ cm/s for ESPRESSO [27]), provided that the test body orbits the neutron star at a distance $\sim 500$ star radii. We leave a detailed study of radial velocity profiles for this case to future research.

Although this paper deals with the orbital Lense–Thirring precession, the effect of nutation as found in this work is somewhat reminiscent of the relativistic nutation of a gyroscope. Recall that the angular momentum of the Earth–Moon system plays the role of a gyroscope moving in the gravitational field of the Sun, which undergoes the geodetic [28, 29, 30] as well as a gravitomagnetic precession (see, for example, [31]). In the motion of the gyroscope, there are also the effects of relativistic nutation which are long- and short-period. The former combines with the precession at short timescales (see [32] and references therein), whereas the latter has the frequency twice that of the orbital motion and is effectively averaged out [33, p. 550].

The amplitude of these kinds of nutation is linear in $a$. Another kind of nutation with an amplitude proportional to the quadrupole moment of the central gravitating body (and, hence, to $a^{2}$) may result from a contribution to the tidal matrix which is periodic with a frequency twice that of the orbital motion [33, Appendix C]. Thus, at least the frequency of certain kinds of the relativistic nutation of a gyroscope is the same as the nutation frequency in the orbital Lense–Thirring effect (see Sec. 6). However, it remains unclear whether there is a deeper connection between the nutation of a gyroscope and the orbital nutation 444The precession of a gyroscope on essentially relativistic bounded orbits around a Kerr black hole was studied in [34]. However, the explicit expressions are quite cumbersome (even when restricted to the case of equatorial plane), and the properties of the nutation in this fully relativistic case are not evident..

Finally, as a by-product, we have also derived the parameters of unstable spherical timelike orbits for arbitrary rotation parameter $a$ and Carter’s constant $Q=1/q$ in an alternative form, i.e. as a function of $a$, $q$, and the radius of the unstable circular orbit, see (86) and (91)–(93). These relations provide an analytical way to separate unbounded orbits that will end in the black from those that will not. Previously, such relations were either restricted to extremely rotating black holes [9, 35] or given as a function of the radius and energy $E$ of the particle [15]. The alternative form of those relations has allowed us to notice a remarkable fact. Namely, if $q=1/Q\leq ay^{5/2}_{-}$, the orbit of a particle abruptly becomes confined to a plane as $r\to r_{\rm ph}^{\pm}$, even though $Q$ is large. This is not so if $a=1$ or $q>ay^{5/2}_{-}$ at either $y_{-}$ or $y_{+}$. This latter case requires further investigation.

Acknowledgements

The analytic part of this research was carried out by VNS and supported by a grant of the Russian Science Foundation (project no. 17-71-10260). The numeric simulations were carried out by SKh.

We are grateful to V.S. Beskin and Yu.Yu. Kovalev for providing a work environment which made completion of this Paper possible. We thank S.V. Repin for useful comments regarding the numerical simulation of geodesics and a presentation of this work before a committee. Last but not least, we thank an anonymous reviewer for pointing out (among other constructive comments) the analogy between the orbital Lense–Thirring precession and the precession of a gyroscope.

Symbolic computations were partially performed with SymPy [36]. We also made use of IPython [37], SciPy [38], Matplotlib [39], NumPy [40], and this preprint was typeset in arxiv-style [41].

Appendix A Numeric solution of geodesic equations

The geodesic equations can be recast into a form that is more convenient for numeric simulations [42, 43, 44]. Although the number of equations increases, now there is no need to track signs of the square roots present in the initial form. The equations to be solved are

[TABLE]

where $d\tau/d\gamma=\rho^{2}/E$, $\xi=L_{z}/E$, and $\eta=Q/E^{2}$. Also, the following relations hold between the functions:

[TABLE]

This set of differential equations was numerically solved with either the classical or variable step [21] fourth-order Runge–Kutta method. That numeric solution was then used to draw the hodographs in Figures 2, 3, 5, and 6. For convenience and reproducibility, we summarize the parameters for those hodographs in Table 2. Together with $r_{1}$ and $\theta_{1}$ evaluated through relations (74) and (75), they yield initial conditions for equations of motion (69)–(73). Also, we always choose $r_{1}<0$ and $\theta_{1}>0$ at the start.

Appendix B Behavior of functions $E(y,a)$ and $L_{z}(y,a)$

(the parameters of stable circular equatorial orbits)

Consider functions $E(y,a)$ and $L_{z}(y,a)$, which are given by (16) and (17), in the range $y\in(0,1/r_{\rm ph}^{\pm})$. Their derivatives are

[TABLE]

They vanish simultaneously at $y=1/r_{\rm isco}^{\pm}$, which satisfies [18]

[TABLE]

Besides, $E\to 1$ and $L_{z}\to\pm\infty$ as $y\to 0$, and $E\to+\infty$ and $L_{z}\to\pm\infty$ as $y\to y_{\pm}=1/r_{\rm ph}^{\pm}$. Therefore, at $y=1/r_{\rm isco}^{\pm}$, $E(y,a)$ must have a minimum, and $L_{z}(y,a)$ must have a minimum (maximum) at $y=1/r_{\rm isco}^{+}$ ($y=1/r_{\rm isco}^{-}$). In other words, both $E(y,a)$ and $L_{z}^{2}(y,a)$ have a minimum at $y=1/r_{\rm isco}^{\pm}$.

Appendix C Parameters of unstable spherical orbits with $E>1$ and arbitrary $Q$

Consider conditions (26) and (27) without omitting terms $\propto q\equiv 1/Q$:

[TABLE]

where $\Delta_{y}\equiv a^{2}y^{2}-2y+1>0$.

Solving the second of the equations for $\epsilon$ and subtituting the solution to the first equation leads to a quadratic equation with respect to $x^{2}$:

[TABLE]

where

[TABLE]

Let us now study how solutions to the quadratic equation behave at $y_{\pm}=1/r_{\rm ph}^{\pm}$. First, note that

[TABLE]

Since $C>0$ and $A\to 0^{-}$ as $y\to y_{\pm}$, the roots of (81) have opposite signs. One of them tends to infinity and the other tends to a finite number. Their specific signs depend on the sign of $K$. If $K>0$, the infinite root is positive and the finite root is negative. Otherwise, the infinite root is negative and the finite root is positive. If $K=0$, both roots tend to infinities of opposite signs.

If $q=0$ or, by continuity, $q<<1$, then $\lim\limits_{y\to y_{\pm}}{K}>0$, and at both $y_{-}$ and $y_{+}$ there is a positive infinite root. This root is the one that leads to solution (30)–(32) obtained under assumption $q<<1$.

For arbitrary $q>0$, it proves helpful to actually solve equation (81). A positive root is found from relation

[TABLE]

At $y_{+}$, which satisfies $3y-1=2ay^{3/2}>0$, the coefficient that multiplies $q$ is reduced as follows:

[TABLE]

Hence, as soon as $0<a<1$,

[TABLE]

At $y_{-}$, which satisfies $3y-1=-2ay^{3/2}<0$, the result depends on the magnitude of $q$. If $q\geq-y(3y-1)/2=ay^{5/2}$, a similar calculation leads to a right-hand side of (86) that vanishes as $A\to 0^{-}$. Thus, in this case the root tends to a finite positive value. Otherwise, the root tends to an infinite positive value.

The function $a[y_{-}(a)]^{5/2}$ grows in the range $0\leq a\leq 1$ reaching its maximal value $1/32$ at $a=1$. Consequently, for a given $a\in(0,1)$, there is a critical value of $q$ that is below $1/32$, which separates the cases of bounded and unbounded growth of $x^{2}$ in a neighborhood of $y_{-}$. Also, if $q$ is sufficiently small, the growth is always unbounded, which corresponds to the limiting case represented by (30)–(32).

Another feature of $x^{2}$ is that, for sufficiently small $q$, it vanishes at a point $y_{0}\in(y_{-},y_{+})$ where it also reaches a minimum. Indeed, for $q=0$ it is trivial to verify this. If $q<<1,q\neq 0$, the minimum value must simultaneously satisfy equations

[TABLE]

Since $C$ and $\partial C/\partial y$ have a common root, $x^{2}=0$ satisfies this last set of equations at that common root and, by continuity, is a minimum of the non-negative function $x^{2}$. The case of $q$ that is large enough for the common root to become complex is beyond the scope of this article.

To ensure that $\epsilon>0$, the following signs must be chosen for $x$:

[TABLE]

Finally, the $x$ found is used to evaluate $\epsilon$ and $\lambda_{z}$:

[TABLE]

Note that $\epsilon(y_{0})$ is a well-defined number, because term $1/x$ multiplies an expression that vanishes at that point. From (79), $\epsilon(y_{0})=\sqrt{(y_{0}^{2}+q)\Delta_{y}(y_{0})}$.

As it was pointed out, for $q$ that is small enough and $0<a<1$, $\lim\limits_{y\to y_{\pm}}{x^{2}}=+\infty$. Therefore,

[TABLE]

[TABLE]

[TABLE]

In the last equation we have used the fact that the ratio does not depend on $q$. Indeed, if we divide (79) by $x^{2}$, terms that contain $q$ vanish as $x^{2}\to+\infty$ and have no impact on the value of the limit.