Two-body problem in general relativity

The Kepler problem in general relativity involves solving for the motion of two spherical bodies interacting with one another by gravitation, as described by the theory of general relativity.

Typically, and in this article, one body is assumed to have a mass m that is negligible compared to the mass M of the other body; this is a good approximation for the case of a planet revolving around the Sun, or a photon passing by a star. In such cases, one may assume that only the heavier body contributes to the curvature of space-time and that it is fixed in space. This curved space-time is described by the Schwarzschild solution to the vacuum Einstein equations of general relativity. The motion of the lighter body (called the "particle" below) is described by the space-time geodesics of the Schwarzschild solution. It is here assumed that the lighter body is point-like so that tidal forces can be ignored.

These geodesic solutions account for the anomalous precession of the planet Mercury, which is a key piece of evidence supporting the theory of general relativity. They also describe the deflection of light in a gravitational field, another prediction famously used as evidence for general relativity.

Orbit decay due to emission of gravitational radiation is not described by the Schwarzschild solution.

Historical context and intuitive understanding

See also: Introduction to general relativity and Tests of general relativity

In the absence of any other forces, a particle orbiting another under the influence of Newtonian gravity follows the same perfect ellipse eternally. The presence of other forces (such as the gravitation of other planets), causes this ellipse to rotate gradually. The rate of this rotation (called orbital precession) can be measured very accurately. The rate can also be predicted knowing the magnitudes and directions of the other forces. However, the predictions of Newtonian gravity do not match the observations, as discovered in 1859 from observations of Mercury.

In 1859, Urbain Le Verrier discovered that the orbital precession of the planet Mercury was not quite what it should be; the ellipse of its orbit was rotating (precessing) slightly faster than predicted by the traditional theory of Newtonian gravity, even after all the effects of the other planets had been accounted for. The effect is small (roughly 43 arcseconds of rotation per century), but well above the measurement error (roughly 0.1 arcseconds per century). Le Verrier realized the importance of his discovery immediately, and challenged astronomers and physicists alike to account for it. Several classical explanations were proposed, such as interplanetary dust, unobserved oblateness of the Sun, an undetected moon of Mercury, or a new planet named Vulcan. After these explanations were discounted, some physicists were driven to the more radical hypothesis that Newton's inverse-square law of gravitation was incorrect. For example, some physicists proposed a power law with an exponent that was slightly different from 2.

Others argued that Newton's law should be supplemented with a velocity-dependent potential. However, this implied a conflict with newtonian celestial dynamics. In his treatise on celestial mechanics, Laplace had shown that if the gravitational influence does not act instantaneously, then the motions of the planets themselves will not exactly conserve momentum (some of the momentum must then be ascribed to the mediator of the gravitational interaction, analogous to ascribing momentum to the mediator of the electromagnetic interaction.) As seen from a newtonian point of view, if gravitational influence does propagate at a finite speed, then at all points in time a planet is attracted to a point where the Sun was some time before, and not towards the instanteneous position of the Sun. On the assumption of the classical fundamentals, Laplace had shown that if gravity would propagate at a velocity in the order of the speed of light then the solar system would be unstable, and would not exist for a long time. The observation that the solar system is old allows one to put a lower limit on the speed of gravity that is many orders of magnitude faster than the speed of light.

To avoid those problems, between 1870-1900 many scientists used the electrodynamic laws of Wilhelm Eduard Weber, Carl Friedrich Gauß, Bernhard Riemann to produce stable orbits and to explain the Perihelion shift of Mercury's orbit. In 1890 Lévy succeeded in doing so by combining the laws of Weber and Riemann, whereby the speed of gravity is equal to the speed of light in his theory. And in another attempt Paul Gerber (1898) even succeeded in deriving the correct formula for the Perihelion shift (which was identical to that formula later used by Einstein). However, because the basic laws of Weber and others were wrong (for example, Weber's law was superseded by Maxwell's theory), those hypothesis were rejected. Another attempt by Hendrik Lorentz (1900), who already used Maxwell's theory, produced a Perihelion shift which was too low.

Around 1904-1905, the works of Hendrik Lorentz, Henri Poincaré and finally Albert Einstein's special theory of relativity, exclude the possibility of propagation of any effects faster than the speed of light. It followed that Newton's law of gravitation would have to be replaced with another law, compatible with the principle of relativity, while still obtaining the newtonian limit for circumstances where relativistic effects are negligible. Such attempts were made by Henri Poincaré (1905), Hermann Minkowski (1907) and Arnold Sommerfeld (1910). In 1907 Einstein came to the conclusion that to achieve this a successor to special relativity was needed. From 1907 to 1915, Einstein gradually groped towards a new theory, using his Equivalence Principle as a key concept to guide his way. According to this principle, a uniform gravitational field acts equally on everything within it and, therefore, cannot be detected by a free-falling observer. Conversely, all local gravitational effects should be reproducible in a linearly accelerating reference frame, and vice versa. Thus, gravity acts like a fictitious force such as the centrifugal force or the Coriolis force, which result from being in an accelerated reference frame; all fictitious forces are proportional to the inertial mass, just as gravity is. To effect the reconciliation of gravity and special relativity and to incorporate the equivalence principle, something had to be sacrificed; that something was the long-held classical assumption that our space obeys the laws of Euclidean geometry, e.g., that the Pythagorean theorem is true experimentally. Einstein used a more general geometry, Riemannian geometry, to allow for the curvature of space and time that was necessary for the reconciliation; after eight years of work (1907-1915), he succeeded in discovering the precise way in which space-time should be curved in order to reproduce the physical laws observed in Nature, particularly gravitation. Gravity is distinct from the fictitious forces centrifugal force and coriolis force in the sense that the curvature of spacetime is regarded as physically real, whereas the fictitious forces are not regarded as forces. The very first solutions of his field equations explained the anomalous precession of Mercury and predicted an unusual bending of light, which was confirmed after his theory was published. These solutions are explained below.

Geometrical background

In the normal Euclidean geometry, triangles obey the Pythagorean theorem, which states that the square distance ds² between two points in space is the sum of the squares of its perpendicular components

$ds^{2}=dx^{2}+dy^{2}+dz^{2}\,\!$

where dx, dy and dz represent the infinitesimal differences between the two points along the x, y and z axes of a Cartesian coordinate system (add Figure here). Now imagine a world in which this is not quite true; a world where the distance is instead given by

$ds^{2}=F(x,y,z)dx^{2}+G(x,y,z)dy^{2}+H(x,y,z)dz^{2}\,\!$

where F, G and H are arbitrary functions of position. It is not hard to imagine such a world; we live on one. The surface of the world is curved, which is why it's impossible to make a perfectly accurate flat map of the world. Non-Cartesian coordinate systems illustrate this well; for example, in the spherical coordinates (r, θ, φ), the Euclidean distance can be written

$ds^{2}=dr^{2}+r^{2}d\theta ^{2}+r^{2}\sin ^{2}\theta d\phi ^{2}\,\!$

Another illustration would be a world in which the rulers used to measure length were untrustworthy, rulers that changed their length with their position and even their orientation. In the most general case, one must allow for cross-terms when calculating the distance ds

$ds^{2}=g_{xx}dx^{2}+g_{xy}dxdy+g_{xz}dxdz+\cdots +g_{zy}dzdy+g_{zz}dz^{2}\,\!$

where the nine functions g_xx, g_xy constitute the metric tensor, which defines the geometry of the space in Riemannian geometry. In the spherical-coordinates example above, there are no cross-terms; the only nonzero metric tensor components are g_rr = 1, g_θθ = r² and g_φφ = r² sin² θ.

In his special theory of relativity, Albert Einstein showed that the distance ds between two spatial points is not constant, but depends on the motion of the observer. However, there is a measure of separation between two points in space-time — called "proper time" and denoted with the symbol dτ — that is invariant; in other words, it doesn't depend on the motion of the observer.

$c^{2}d\tau ^{2}=c^{2}dt^{2}-dx^{2}-dy^{2}-dz^{2}\,\!$

which may be written in spherical coordinates as

$c^{2}d\tau ^{2}=c^{2}dt^{2}-dr^{2}-r^{2}d\theta ^{2}-r^{2}\sin ^{2}\theta d\varphi ^{2}\,\!$

This formula is the natural extension of the Pythagorean theorem and similarly holds only when there is no curvature in space-time. In general relativity, however, space and time may have curvature, so this distance formula must be modified to a more general form

$c^{2}d\tau ^{2}=g_{\mu \nu }dx^{\mu }dx^{\nu }\,\!$

just as we generalized the formula to measure distance on the surface of the Earth. The exact form of the metric g_μν depends on the gravitating mass, momentum and energy, as described by the Einstein field equations. Einstein developed those field equations to match the then known laws of Nature; however, they predicted never-before-seen phenomena (such as the bending of light by gravity) that were confirmed later.

Schwarzschild metric

One solution to the Einstein field equations is the Schwarzschild metric, which corresponds to the gravitational field of an uncharged, non-rotating, spherically symmetric body of mass M. The Schwarzschild solution can be written as

$c^{2}{d\tau }^{2}=\left(1-{\frac {r_{s}}{r}}\right)c^{2}dt^{2}-{\frac {dr^{2}}{1-{\frac {r_{s}}{r}}}}-r^{2}d\theta ^{2}-r^{2}\sin ^{2}\theta \,d\varphi ^{2}$

where

τ is the proper time (time measured by a clock moving with the particle) in seconds,

c is the speed of light in meters per second,

t is the time coordinate (measured by a stationary clock at infinity) in seconds,

r is the radial coordinate (circumference of a circle centered on the star divided by 2π) in meters,

θ is the colatitude (angle from North) in radians,

φ is the longitude in radians, and

r_s is the Schwarzschild radius (in meters) of the massive body, which is related to its mass M by

$r_{s}={\frac {2GM}{c^{2}}}$

where G is the gravitational constant.

The classical Newtonian theory of gravity is recovered in the limit as the ratio r_s/r goes to zero. In that limit, the metric returns to the form given above for special relativity. In practice, this ratio is almost always extremely small. For example, the Schwarzschild radius r_s of the Earth is roughly 9 mm (3⁄8 inch), whereas a satellite in a geosynchronous orbit has a radius r that is roughly four billion times larger, at 42,164 km (26,200 miles). Even at the surface of the Earth, the corrections to Newtonian gravity are only one part in a billion. The ratio only becomes large close to black holes and other ultra-dense objects such as neutron stars.

Geodesic equation

According to Einstein's theory of general relativity, particles of negligible mass travel along geodesics in the space-time. In uncurved space-time, far from a source of gravity, these geodesics correspond to straight lines; however, they may deviate from straight lines when the space-time is curved. The equation for the geodesic lines is

${\frac {d^{2}x^{\mu }}{dq^{2}}}+\Gamma _{\nu \lambda }^{\mu }{\frac {dx^{\nu }}{dq}}{\frac {dx^{\lambda }}{dq}}=0$

where Γ represents the Christoffel symbol and the variable q parametrizes the particle's path through space-time, its so-called world line. The Christoffel symbol depends only on the metric tensor g_μν, or rather on how it changes with position. The variable q is a constant multiple of the proper time τ for timelike orbits (which are traveled by massive particles), and is usually taken to be equal to it. For lightlike orbits (which are traveled by massless particles such as the photon), the proper time is zero and, strictly speaking, cannot be used as the variable q. Nevertheless, lightlike orbits can be derived as the ultrarelativistic limit of timelike orbits, that is, the limit as the particle mass m goes to zero while holding its total energy fixed.

We may simplify the problem by using symmetry to eliminate one variable from consideration. Since the Schwarzschild metric is symmetrical about θ = π/2, any geodesic that begins moving in that plane will remain in that plane indefinitely (the plane is totally geodesic). Therefore, we orient the coordinate system so that the orbit of the particle lies in that plane, and fix the θ coordinate to be π/2 so that the metric (of this plane) simplifies to

$c^{2}d\tau ^{2}=\left(1-{\frac {r_{s}}{r}}\right)c^{2}dt^{2}-{\frac {dr^{2}}{1-{\frac {r_{s}}{r}}}}-r^{2}d\varphi ^{2}.$

Two constants of motion can then be identified (cf. the derivation given below)

$r^{2}{\frac {d\varphi }{d\tau }}={\frac {L}{m}},$

$\left(1-{\frac {r_{s}}{r}}\right){\frac {dt}{d\tau }}={\frac {E}{mc^{2}}}.$

Substituting these constants into the definition of the Schwarzschild metric

$c^{2}=\left(1-{\frac {r_{s}}{r}}\right)c^{2}\left({\frac {dt}{d\tau }}\right)^{2}-{\frac {1}{1-{\frac {r_{s}}{r}}}}\left({\frac {dr}{d\tau }}\right)^{2}-r^{2}\left({\frac {d\varphi }{d\tau }}\right)^{2},$

yields the equation of motion for the particle

$\left({\frac {dr}{d\tau }}\right)^{2}={\frac {E^{2}}{m^{2}c^{2}}}-\left(1-{\frac {r_{s}}{r}}\right)\left(c^{2}+{\frac {L^{2}}{m^{2}r^{2}}}\right).$

The dependence on proper time can be eliminated using the definition of L

$\left({\frac {dr}{d\varphi }}\right)^{2}=\left({\frac {dr}{d\tau }}\right)^{2}\left({\frac {d\tau }{d\varphi }}\right)^{2}=\left({\frac {dr}{d\tau }}\right)^{2}\left({\frac {mr^{2}}{L}}\right)^{2},$

which yields the equation for the orbit

$\left({\frac {dr}{d\varphi }}\right)^{2}={\frac {r^{4}}{b^{2}}}-\left(1-{\frac {r_{s}}{r}}\right)\left({\frac {r^{4}}{a^{2}}}+r^{2}\right)$

where, for brevity, two length-scales, a and b, have been defined by

$a={\frac {L}{mc}},$

$b={\frac {cL}{E}}.$

The same equation can also be derived using a Lagrangian approach or the Hamilton–Jacobi equation (see below). The solution of the orbit equation is

$\varphi =\int {\frac {dr}{r^{2}{\sqrt {{\frac {1}{b^{2}}}-\left(1-{\frac {r_{s}}{r}}\right)\left({\frac {1}{a^{2}}}+{\frac {1}{r^{2}}}\right)}}}}.$

Approximate formula for the bending of light

Eddington's 1919 measurements of the bending of star-light by the Sun's gravity led to the acceptance of general relativity worldwide.

In the limit as the particle mass m goes to zero (or, equivalently, as the length-scale a goes to infinity), the equation for the orbit becomes

$\varphi =\int {\frac {dr}{r^{2}{\sqrt {{\frac {1}{b^{2}}}-\left(1-{\frac {r_{s}}{r}}\right){\frac {1}{r^{2}}}}}}}$

Expanding in powers of r_s/r, the leading order term in this formula gives the approximate angular deflection δφ for a massless particle coming in from infinity and going back out to infinity:

$\delta \varphi \approx {\frac {2r_{s}}{b}}={\frac {4GM}{c^{2}b}}.$

Here, b can be interpreted as the distance of closest approach. Although this formula is approximate, it is accurate for most measurements of gravitational lensing, due to the smallness of the ratio r_s/r. For light grazing the surface of the sun, the approximate angular deflection is roughly 1.75 arcseconds.

Relation to classical mechanics and precession of elliptical orbits

The equation of motion for the particle derived above

$\left({\frac {dr}{d\tau }}\right)^{2}={\frac {E^{2}}{m^{2}c^{2}}}-c^{2}+{\frac {r_{s}c^{2}}{r}}-{\frac {L^{2}}{m^{2}r^{2}}}+{\frac {r_{s}L^{2}}{m^{2}r^{3}}}$

can be re-written using the definition of the Schwarzschild radius r_s as

${\frac {1}{2}}m\left({\frac {dr}{d\tau }}\right)^{2}=\left[{\frac {E^{2}}{2mc^{2}}}-{\frac {1}{2}}mc^{2}\right]+{\frac {GMm}{r}}-{\frac {L^{2}}{2mr^{2}}}+{\frac {GML^{2}}{c^{2}mr^{3}}}$

which is equivalent to a particle moving in a one-dimensional effective potential

$V(r)=-{\frac {GMm}{r}}+{\frac {L^{2}}{2mr^{2}}}-{\frac {GML^{2}}{c^{2}mr^{3}}}$

The first two terms are well-known classical energies, the first being the attractive Newtonian gravitational potential energy and the second corresponding to the repulsive "centrifugal" potential energy; however, the third term is an attractive energy unique to general relativity. As shown below and elsewhere, this inverse-cubic energy causes elliptical orbits to precess gradually by an angle δφ per revolution

$\delta \varphi \approx {\frac {6\pi GM}{c^{2}A\left(1-e^{2}\right)}}$

where A is the semi-major axis and e is the eccentricity.

The third term is attractive and dominates at small r values, giving a critical inner radius r_inner at which a particle is drawn inexorably inwards to r=0; this inner radius is a function of the particle's angular momentum per unit mass or, equivalently, the a length-scale defined above.

Circular orbits and their stability

Effective radial potential for various angular momenta. At small radii, the energy drops precipitously, causing the particle to be pulled inexorably inwards to r=0. However, when the normalized angular momentum a/r_s = L/mcr_s equals the square root of three, a metastable circular orbit is possible at the radius highlighted with a green circle. At higher angular momenta, there is a significant centrifugal barrier (orange curve) and an unstable inner radius, highlighted in red.

The effective potential V can be re-written in terms of the lengths a and b

$V(r)={\frac {mc^{2}}{2}}\left[-{\frac {r_{s}}{r}}+{\frac {a^{2}}{r^{2}}}-{\frac {r_{s}a^{2}}{r^{3}}}\right]$

Circular orbits are possible when the effective force is zero

$F=-{\frac {dV}{dr}}=-{\frac {mc^{2}}{2r^{4}}}\left[r_{s}r^{2}-2a^{2}r+3r_{s}a^{2}\right]=0$

i.e., when the two attractive forces — Newtonian gravity (first term) and the attraction unique to general relativity (third term) — are exactly balanced by the repulsive centrifugal force (second term). There are two radii at which this balancing can occur, denoted here as r_inner and r_outer

$r_{\mathrm {outer} }={\frac {a^{2}}{r_{s}}}\left(1+{\sqrt {1-{\frac {3r_{s}^{2}}{a^{2}}}}}\right)$

$r_{\mathrm {inner} }={\frac {a^{2}}{r_{s}}}\left(1-{\sqrt {1-{\frac {3r_{s}^{2}}{a^{2}}}}}\right)={\frac {3a^{2}}{r_{\mathrm {outer} }}}$

which are obtained using the quadratic formula. The inner radius r_inner is unstable, because the attractive third force strengthens much faster than the other two forces when r becomes small; if the particle slips slightly inwards from r_inner (where all three forces are in balance), the third force dominates the other two and draws the particle inexorably inwards to r=0. At the outer radius, however, the circular orbits are stable; the third term is less important and the system behaves more like the non-relativistic Kepler problem.

When a is much greater than r_s (the classical case), these formulae become approximately

$r_{\mathrm {outer} }\approx {\frac {2a^{2}}{r_{s}}}$

$r_{\mathrm {inner} }\approx {\frac {3}{2}}r_{s}$

The stable and unstable radii are plotted versus the normalized angular momentum a/r_s = L/mcr_s in blue and red, respectively. These curves meet at a unique circular orbit (green circle) when the normalized angular momentum equals the square root of three. For comparison, the classical radius predicted from the centripetal acceleration and Newton's law of gravity is plotted in black.

Substituting the definitions of a and r_s into r_outer yields the classical formula for a particle orbiting a mass M in a circle

$r_{\mathrm {outer} }^{3}\approx {\frac {GM}{\omega _{\varphi }^{2}}}$

where ω_φ is the orbital angular speed of the particle. This formula is obtained in non-relativistic mechanics by setting the centrifugal force equal to the Newtonian gravitational force

$m\omega _{\varphi }^{2}r={\frac {GMm}{r^{2}}}$

In our notation, the classical orbital angular speed equals

$\omega _{\varphi }^{2}\approx {\frac {GM}{r_{\mathrm {outer} }^{3}}}=\left({\frac {r_{s}c^{2}}{2r_{\mathrm {outer} }^{3}}}\right)=\left({\frac {r_{s}c^{2}}{2}}\right)\left({\frac {r_{s}^{3}}{8a^{6}}}\right)={\frac {c^{2}r_{s}^{4}}{16a^{6}}}$

At the other extreme, when a² approaches 3r_s² from above, the two radii converge to a single value

$r_{\mathrm {outer} }\approx r_{\mathrm {inner} }\approx 3r_{s}$

The quadratic solutions above ensure that r_outer is always greater than 3r_s, whereas r_inner lies between 3⁄2 r_s and 3r_s. Circular orbits smaller than 3⁄2 r_s are not possible. For massless particles, a goes to infinity, implying that there is a circular orbit for photons at r_inner = 3⁄2 r_s. The sphere of this radius is sometimes known as the photon sphere.

Precession of elliptical orbits

In the non-relativistic Kepler problem, a particle follows the same perfect ellipse (red orbit) eternally. General relativity introduces a third force that attracts the particle slightly more strongly than Newtonian gravity, especially at small radii. This third force causes the particle's elliptical orbit to precess (cyan orbit) in the direction of its rotation; this effect has been measured in Mercury, Venus and Earth. The yellow dot within the orbits represents the center of attraction, such as the Sun.

The orbital precession rate may be derived using this radial effective potential V. A small radial deviation from a circular orbit of radius r_outer will oscillate stably with an angular frequency

$\omega _{r}^{2}={\frac {1}{m}}\left[{\frac {d^{2}V}{dr^{2}}}\right]_{r=r_{\mathrm {outer} }}$

which equals

$\omega _{r}^{2}=\left({\frac {c^{2}r_{s}}{2r_{\mathrm {outer} }^{4}}}\right)\left(r_{\mathrm {outer} }-r_{\mathrm {inner} }\right)=\omega _{\varphi }^{2}{\sqrt {1-{\frac {3r_{s}^{2}}{a^{2}}}}}$

Taking the square root of both sides and expanding using the binomial theorem yields the formula

$\omega _{r}=\omega _{\varphi }\left(1-{\frac {3r_{s}^{2}}{4a^{2}}}+\cdots \right)$

Multiplying by the period T of one revolution gives the precession of the orbit per revolution

$\delta \varphi =T\left(\omega _{\varphi }-\omega _{r}\right)\approx 2\pi \left({\frac {3r_{s}^{2}}{4a^{2}}}\right)={\frac {3\pi m^{2}c^{2}}{2L^{2}}}r_{s}^{2}$

where we have used ω_φT = 2п and the definition of the length-scale a. Substituting the definition of the Schwarzschild radius r_s gives

$\delta \varphi \approx {\frac {3\pi m^{2}c^{2}}{2L^{2}}}\left({\frac {4G^{2}M^{2}}{c^{4}}}\right)={\frac {6\pi G^{2}M^{2}m^{2}}{c^{2}L^{2}}}$

This may be simplified using the elliptical orbit's semiaxis A and eccentricity e related by the formula

${\frac {L^{2}}{GMm^{2}}}=A\left(1-e^{2}\right)$

to give the most common form of the precession angle

$\delta \varphi \approx {\frac {6\pi GM}{c^{2}A\left(1-e^{2}\right)}}$

Orbital solution using elliptic functions

The equation for the orbit

$\left({\frac {dr}{d\varphi }}\right)^{2}={\frac {r^{4}}{b^{2}}}-\left(1-{\frac {r_{s}}{r}}\right)\left({\frac {r^{4}}{a^{2}}}+r^{2}\right)$

can be simplified by introducing a dimensionless variable

$\zeta ={\frac {r_{s}}{4r}}-{\frac {1}{12}}$

so that it reduces to

$\left({\frac {d\zeta }{d\varphi }}\right)^{2}=4\zeta ^{3}-g_{2}\zeta -g_{3},$

where the constant, dimensionless coefficients g₂ and g₃ are defined by

${\begin{aligned}g_{2}&={\frac {1}{12}}-{\frac {r_{s}^{2}}{4a^{2}}},\\g_{3}&={\frac {1}{216}}+{\frac {r_{s}^{2}}{24a^{2}}}-{\frac {r_{s}^{2}}{16b^{2}}}.\end{aligned}}$

The solution of this orbital equation is given by

$\varphi -\varphi _{0}=\int {\frac {d\zeta }{\sqrt {4\zeta ^{3}-g_{2}\zeta -g_{3}}}}.$

It follows that, up to a phase-shift, $\zeta =\wp (\varphi -\varphi _{0})$, where $\wp$ is Weierstrass's elliptic function with parameters g₂ and g₃, and φ₀ is a constant of integration (possibly complex).

Qualitative character of possible orbits

The orbits are described by the equation of motion

$\left({\frac {d\zeta }{d\varphi }}\right)^{2}=4\zeta ^{3}-g_{2}\zeta -g_{3},$

If the discriminant $\Delta =g_{2}^{3}-27g_{3}^{2}$ is greater than zero, the cubic equation

$G(\zeta )=4\zeta ^{3}-g_{2}\zeta -g_{3}=0\,$

has three distinct real roots, e₁, e₂, and e₃, which may be listed in decreasing order

$e_{1}>e_{2}>e_{3}.$

In such cases, the solution $\zeta =\wp (\varphi -\varphi _{0})$ is an elliptic function with two half-periods, one completely real

$\omega _{1}=\int _{e_{1}}^{\infty }{\frac {dz}{\sqrt {4z^{3}-g_{2}z-g_{3}}}}$

and the other completely imaginary

$\omega _{3}=i\int _{-e_{3}}^{\infty }{\frac {dz}{\sqrt {4z^{3}-g_{2}z-g_{3}}}}$

The remaining root defines a complex half-period ω₂ = -ω₁ - ω₃. These three half-periods are related to the corresponding roots by the equation $\wp (\omega _{i})=e_{i}$, where i can equal 1, 2, or 3. Therefore, if φ₀ is set equal to any of these half-periods, the derivative of ζ is zero, indicating a periapse or an apoapse, the points of closest or furthest approach, respectively.

${\frac {d\zeta }{d\phi }}=0\ \mathrm {when} \ \zeta =\wp (-\omega _{i})=e_{i}$

since

$\left({\frac {d\zeta }{d\varphi }}\right)^{2}=G(\zeta )=4\zeta ^{3}-g_{2}\zeta -g_{3}=4\left(\zeta -e_{1}\right)\left(\zeta -e_{2}\right)\left(\zeta -e_{3}\right),$

The roots e_i are the critical values of ζ with respect to φ, i.e., the values at which the derivative is zero.

The qualitative character of the orbits depends on the choice of φ₀. Solutions with φ₀ equal to ω₂ correspond to oscillatory orbits that vary between ζ=e₂ and ζ=e₃, or else diverge to infinity (ζ=-1/12). By contrast, solutions with φ₀ equal to ω₁ (or any other real number) correspond to orbits that decay to zero radius, since (to be a real number) ζ cannot be less than e₁ and hence increases inexorably to infinity.

Quasi-elliptical orbits

Solutions $\zeta =\wp (\phi -\phi _{0})$ in which φ₀ equals ω₂ give a real value of ζ provided that the energy E satisfies the inequality E² < m² c⁴. For such solutions, the variable ζ is confined between e₃ ≤ ζ ≤ e₂. If both roots are greater than -1/12, then ζ never becomes -1/12, at which point the radius r goes to infinity. Hence, such solutions correspond to a gradually precessing elliptical orbit. As the particle (or planet) revolves around the origin, its radius oscillates between a minimum radius

$r_{min}={\frac {3r_{s}}{1+12e_{2}}}$

and a maximum radius

$r_{max}={\frac {3r_{s}}{1+12e_{1}}}$

corresponding to the value of ζ at the two extrema of radius. The real periodicity of Weierstrass' elliptic function is 2ω₁; thus the particle returns to the same radius after revolving an angle 2ω₁, which may not equal 2π. Hence, the orbit is in general precessing. In general, the amount of precession per orbit (2ω₁ - 2π) is quite small.

Stable circular orbits

A special case occurs when 2e₂ = 2e₃ = −e₃, i.e., two of the roots of G(ζ) are equal and negative, while the third is positive. In this case there is a solution with ζ constant, equal to the repeated root, e = e₂ = e₃. This corresponds to circular orbits, specifically, the classical r_outer solution derived above; as shown there, the radii of these orbits must be greater than 3r_s. Such circular orbits are stable, because a small perturbation of the parameters will separate the repeated roots, resulting in a quasi-elliptical orbit. For example, giving a small radial "kick" to a particle in a classical circular orbit pushes it into an elliptical orbit that gradually precesses, as derived above.

Unbounded orbits

An unbounded orbit occurs when r goes to infinity, corresponding to ζ = −¹⁄₁₂. Unbounded orbits correspond to an oscillatory orbit in which -1/12 falls between the two limiting roots, i.e., when e₃ ≤ −¹⁄₁₂ ≤ ζ ≤ e₂.

Asymptotically circular orbits

Another special case occurs when −e₃ = 2e₂ = 2e₁, i.e., two of the roots of G(ζ) are equal and positive, whereas the third root e₃ is negative. Denoting the repeated root by e = n²/3, the orbits are asymptotically circular at positive and negative infinite φ:

$\zeta ={\frac {r_{s}}{4r}}-{\frac {1}{12}}=e-{\frac {n^{2}}{\cosh ^{2}n\varphi }}.$

as may be verified by substitution. As φ goes to positive or negative infinity, the orbit approaches asymptotically to the circle

${\frac {r_{s}}{4r}}-{\frac {1}{12}}=e.$

In such cases, the radius of the orbit must remain between 2r_s and 3r_s.

The asymptotic formula may also be derived from the expression for Weierstrass' elliptic function in terms of Jacobi's elliptic functions

$\zeta =\wp (\phi -\phi _{0})=e_{1}+\left(e_{1}-e_{3}\right){\frac {\mathrm {cn} ^{2}w}{\mathrm {sn} ^{2}w}}$

where $w=(\phi -\phi _{0}){\sqrt {e_{1}-e_{3}}}$ and the modulus equals

$k={\sqrt {\frac {e_{2}-e_{3}}{e_{1}-e_{3}}}}$

In the limit as e₂ approaches e₁, the modulus goes to one and w goes to n(φ-φ₀). Finally choosing φ₀ to be the imaginary number $iK^{\prime }$ (a quarter-period) gives the asymptotic formula above.

Decaying orbits

Real solutions for $\zeta =\wp (\phi -\phi _{0})$ in which φ₀ equals ω₁ (or some other real number) have the property that ζ is never less than e₁. Since the equation of motion

$\left({\frac {d\zeta }{d\varphi }}\right)^{2}=4\zeta ^{3}-g_{2}\zeta -g_{3}=4\left(\zeta -e_{1}\right)\left(\zeta -e_{2}\right)\left(\zeta -e_{3}\right)$

is positive for all values of ζ>e₁, ζ increases without bound, corresponding to the particle falling inexorably to the origin r = 0.

Experimentally observed decreases of the orbital period of the binary pulsar PSR B1913+16 (blue dots) match the predictions of general relativity (black curve) almost exactly.

Corrections to the geodesic solutions

See also: PSR B1913+16 and Gravitational radiation

According to general relativity, two bodies rotating about one another will emit gravitational radiation, causing the orbits to differ slightly from the geodesics calculated above. This has been observed indirectly in a binary star system known as PSR B1913+16, for which Russell Alan Hulse and Joseph Hooton Taylor, Jr. were awarded the 1993 Nobel Prize in Physics. The two neutron stars of this system are extremely close and rotate about one another very quickly, completing a revolution in roughly 465 minutes. Their orbit is highly elliptical, with an eccentricity of 0.62 (62%). According to general relativity, the short orbital period and high eccentricity should make the system an excellent emitter of gravitational radiation, thereby losing energy and decreasing the orbital period still further. The observed decrease in the orbital period over thirty years matches the predictions of general relativity within even the most precise measurements. General relativity predicts that, in another 300 million years, these two stars will spiral into one another.

Two neutron stars rotating rapidly around one another gradually lose energy by emitting gravitational radiation. As they lose energy, they revolve about each other more quickly and more closely to one another.

The formulae describing the loss of energy and angular momentum due to gravitational radiation from the two bodies of the Kepler problem have been calculated. The rate of losing energy (averaged over a complete orbit) is given by

$-\langle {\frac {dE}{dt}}\rangle ={\frac {32G^{4}m_{1}^{2}m_{2}^{2}\left(m_{1}+m_{2}\right)}{5c^{5}a^{5}\left(1-e^{2}\right)^{7/2}}}\left(1+{\frac {73}{24}}e^{2}+{\frac {37}{96}}e^{4}\right)$

where e is the orbital eccentricity and a is the semimajor axis of the elliptical orbit. The angular brackets on the left-hand side of the equation represent the averaging over a single orbit. Similarly, the average rate of losing angular momentum equals

$-\langle {\frac {dL_{z}}{dt}}\rangle ={\frac {32G^{7/2}m_{1}^{2}m_{2}^{2}{\sqrt {m_{1}+m_{2}}}}{5c^{5}a^{7/2}\left(1-e^{2}\right)^{2}}}\left(1+{\frac {7}{8}}e^{2}\right)$

The losses in energy and angular momentum increase significantly as the eccentricity approaches one, i.e., as the ellipse of the orbit becomes ever more elongated. The radiation losses also increase significantly with a decreasing size a of the orbit.

Mathematical derivations of the orbit equation

Hamilton–Jacobi approach

The orbital equation can be derived from the Hamilton–Jacobi equation. The advantage of this approach is that it equates the motion of the particle with the propagation of a wave, and leads neatly into the derivation of the deflection of light by gravity in general relativity, through Fermat's principle. The basic idea is that, due to gravitational slowing of time, parts of a wave-front closer to a gravitating mass move more slowly than those further away, thus bending the direction of the wave-front's propagation (add Figure).

Using general covariance, the Hamilton–Jacobi equation for a single particle in arbitrary coordinates can be expressed as

$g^{\mu \nu }{\frac {\partial S}{\partial x^{\mu }}}{\frac {\partial S}{\partial x^{\nu }}}=m^{2}c^{2}.$

Using the Schwarzschild metric g^μν, this equation becomes

${\frac {1}{c^{2}\left(1-{\frac {r_{s}}{r}}\right)}}\left({\frac {\partial S}{\partial t}}\right)^{2}-\left(1-{\frac {r_{s}}{r}}\right)\left({\frac {\partial S}{\partial r}}\right)^{2}-{\frac {1}{r^{2}}}\left({\frac {\partial S}{\partial \varphi }}\right)^{2}=m^{2}c^{2}$

where we again orient the spherical coordinate system with the plane of the orbit. The time t and longitude φ are cyclic coordinates, so that the solution for Hamilton's principal function S can be written

$S=-Et+L\varphi +S_{r}(r)\,$

where E and L again represent the particle's energy and angular momentum, respectively. The Hamilton–Jacobi equation gives an integral solution for the radial part S_r(r)

$S_{r}(r)=\int {\frac {Ldr}{1-{\frac {r_{s}}{r}}}}{\sqrt {{\frac {1}{b^{2}}}-\left(1-{\frac {r_{s}}{r}}\right)\left({\frac {1}{a^{2}}}+{\frac {1}{r^{2}}}\right)}}.$

Taking the derivative of Hamilton's principal function S in the usual way

${\frac {\partial S}{\partial L}}=\varphi +{\frac {\partial S_{r}}{\partial L}}=\mathrm {constant}$

yields the orbital equation derived earlier

$\left({\frac {dr}{d\varphi }}\right)^{2}={\frac {r^{4}}{b^{2}}}-\left(1-{\frac {r_{s}}{r}}\right)\left({\frac {r^{4}}{a^{2}}}+r^{2}\right).$

This approach can also be used to derive the rate of orbital precession in an elegant way.

In the limit of zero mass m (or, equivalently, infinite a), the radial part of Hamilton's principal function S becomes

$S_{r}(r)={\frac {E}{c}}\int dr{\sqrt {{\frac {r^{2}}{\left(r-r_{s}\right)^{2}}}-{\frac {b^{2}}{r\left(r-r_{s}\right)}}}}$

from which the equation for the deflection of light can be derived.

Lagrangian approach

In general relativity, free particles of negligible mass m follow geodesics in space-time, owing to the equivalence principle. Geodesics in space-time are defined as curves for which small local variations in their coordinates (while holding their endpoints events fixed) make no significant change in their overall length s. This may be expressed mathematically using the calculus of variations

$0=\delta s=\delta \int ds=\delta \int {\sqrt {g_{\mu \nu }{\frac {dx^{\mu }}{d\tau }}{\frac {dx^{\nu }}{d\tau }}}}d\tau =\delta \int {\sqrt {2T}}d\tau$

where τ is the proper time, s=cτ is the arc-length in space-time and T is defined as

$2T=c^{2}=\left({\frac {ds}{d\tau }}\right)^{2}=g_{\mu \nu }{\frac {dx^{\mu }}{d\tau }}{\frac {dx^{\nu }}{d\tau }}=\left(1-{\frac {r_{s}}{r}}\right)c^{2}\left({\frac {dt}{d\tau }}\right)^{2}-{\frac {1}{1-{\frac {r_{s}}{r}}}}\left({\frac {dr}{d\tau }}\right)^{2}-r^{2}\left({\frac {d\varphi }{d\tau }}\right)^{2}$

in analogy with kinetic energy. If the derivative with respect to proper time is represented by a dot for brevity

${\dot {x}}^{\mu }={\frac {dx^{\mu }}{d\tau }}$

T may be written as

$2T=c^{2}=\left(1-{\frac {r_{s}}{r}}\right)c^{2}\left({\dot {t}}\right)^{2}-{\frac {1}{1-{\frac {r_{s}}{r}}}}\left({\dot {r}}\right)^{2}-r^{2}\left({\dot {\varphi }}\right)^{2}$

Constant factors (such as c or the square root of two) don't affect the answer to the variational problem; therefore, taking the variation inside the integral yields Hamilton's principle

$0=\delta \int {\sqrt {2T}}d\tau =\int {\frac {\delta T}{\sqrt {2T}}}d\tau ={\frac {1}{c}}\delta \int Td\tau .$

The solution of the variational problem is given by Lagrange's equations

${\frac {d}{d\tau }}\left({\frac {\partial T}{\partial {\dot {x}}^{\sigma }}}\right)={\frac {\partial T}{\partial x^{\sigma }}}.$

When applied to t and φ, these equations reveal two constants of motion

${\frac {d}{d\tau }}\left[r^{2}{\frac {d\varphi }{d\tau }}\right]=0,$

${\frac {d}{d\tau }}\left[\left(1-{\frac {r_{s}}{r}}\right){\frac {dt}{d\tau }}\right]=0,$

which may be written as the equations for L and E

$r^{2}{\frac {d\varphi }{d\tau }}={\frac {L}{m}},$

$\left(1-{\frac {r_{s}}{r}}\right){\frac {dt}{d\tau }}={\frac {E}{mc^{2}}}.$

As shown above, substitution of these equations into the definition of the Schwarzschild metric yields the equation for the orbit.

Hamilton's principle

The action integral for a particle affected only by gravity is

$S=\int {-mc^{2}d\tau }=-mc\int {c{\frac {d\tau }{dq}}dq}=-mc\int {{\sqrt {g_{\mu \nu }{\frac {dx^{\mu }}{dq}}{\frac {dx^{\nu }}{dq}}}}dq}$

where τ is the proper time and q is any smooth parameterization of the particle's world line. If one applies the calculus of variations to this, one again gets the equations for a geodesic. The calculations are simplified, if we first take the variation of the square of the integrand. For the metric and coordinates of this case, that square is

$\left(c{\frac {d\tau }{dq}}\right)^{2}=g_{\mu \nu }{\frac {dx^{\mu }}{dq}}{\frac {dx^{\nu }}{dq}}=\left(1-{\frac {r_{s}}{r}}\right)c^{2}\left({\frac {dt}{dq}}\right)^{2}-{\frac {1}{1-{\frac {r_{s}}{r}}}}\left({\frac {dr}{dq}}\right)^{2}-r^{2}\left({\frac {d\varphi }{dq}}\right)^{2}$

Taking variation of this, we get

$\delta \left(c{\frac {d\tau }{dq}}\right)^{2}=2c^{2}{\frac {d\tau }{dq}}\delta {\frac {d\tau }{dq}}=\delta \left[\left(1-{\frac {r_{s}}{r}}\right)c^{2}\left({\frac {dt}{dq}}\right)^{2}-{\frac {1}{1-{\frac {r_{s}}{r}}}}\left({\frac {dr}{dq}}\right)^{2}-r^{2}\left({\frac {d\varphi }{dq}}\right)^{2}\right]$

If we vary with respect to longitude φ only, we get

$2c^{2}{\frac {d\tau }{dq}}\delta {\frac {d\tau }{dq}}=-2r^{2}{\frac {d\varphi }{dq}}\delta {\frac {d\varphi }{dq}}$

we divide by $2c{\frac {d\tau }{dq}}$ to get the variation of the integrand itself

$c\,\delta {\frac {d\tau }{dq}}=-{\frac {r^{2}}{c}}{\frac {d\varphi }{d\tau }}\delta {\frac {d\varphi }{dq}}=-{\frac {r^{2}}{c}}{\frac {d\varphi }{d\tau }}{\frac {d\delta \varphi }{dq}}$

Thus we have

$0=\delta \int {c{\frac {d\tau }{dq}}dq}=\int {c\delta {\frac {d\tau }{dq}}dq}=\int {-{\frac {r^{2}}{c}}{\frac {d\varphi }{d\tau }}{\frac {d\delta \varphi }{dq}}dq}$

integrating by parts gives

$0=-{\frac {r^{2}}{c}}{\frac {d\varphi }{d\tau }}\delta \varphi -\int {{\frac {d}{dq}}\left[-{\frac {r^{2}}{c}}{\frac {d\varphi }{d\tau }}\right]\delta \varphi dq}$

The variation of the longitude is assumed to be zero at the end points, so the first term disappears. The integral can be made nonzero by a perverse choice of δφ unless the other factor inside is zero everywhere. So we get the equation of motion

${\frac {d}{dq}}\left[-{\frac {r^{2}}{c}}{\frac {d\varphi }{d\tau }}\right]=0$

If we vary with respect to time t only, we get

$2c^{2}{\frac {d\tau }{dq}}\delta {\frac {d\tau }{dq}}=2\left(1-{\frac {r_{s}}{r}}\right)c^{2}{\frac {dt}{dq}}\delta {\frac {dt}{dq}}$

we divide by $2c{\frac {d\tau }{dq}}$ to get the variation of the integrand itself

$c\delta {\frac {d\tau }{dq}}=c\left(1-{\frac {r_{s}}{r}}\right){\frac {dt}{d\tau }}\delta {\frac {dt}{dq}}=c\left(1-{\frac {r_{s}}{r}}\right){\frac {dt}{d\tau }}{\frac {d\delta t}{dq}}$

Thus we have

$0=\delta \int {c{\frac {d\tau }{dq}}dq}=\int {c\left(1-{\frac {r_{s}}{r}}\right){\frac {dt}{d\tau }}{\frac {d\delta t}{dq}}dq}$

integrating by parts gives

$0=c\left(1-{\frac {r_{s}}{r}}\right){\frac {dt}{d\tau }}\delta t-\int {{\frac {d}{dq}}\left[c\left(1-{\frac {r_{s}}{r}}\right){\frac {dt}{d\tau }}\right]\delta tdq}$

giving the equation of motion

${\frac {d}{dq}}\left[c\left(1-{\frac {r_{s}}{r}}\right){\frac {dt}{d\tau }}\right]=0$

If we integrate these equations of motion and determine the constants of integration, we again get the equations

$r^{2}{\frac {d\varphi }{d\tau }}={\frac {L}{m}},$

$\left(1-{\frac {r_{s}}{r}}\right){\frac {dt}{d\tau }}={\frac {E}{mc^{2}}}.$

These two equations for the constants of motion L and E can be combined to form one equation that is true even for photons and other massless particles for which the proper time along a geodesic is zero.

${\frac {r^{2}}{bc}}{\frac {d\varphi }{dt}}=1-{\frac {r_{s}}{r}}$

See also

• Kepler problem
• Schwarzschild metric
• Tests of general relativity
• Laplace–Runge–Lenz vector
• Deriving the Schwarzschild solution

References

• Adler, R (1965). Introduction to General Relativity. New York: McGraw-Hill Book Company. pp. pp. 177–193. ISBN 978-0-07-000420-7. {{cite book}}: |pages= has extra text (help); Unknown parameter |coauthors= ignored (|author= suggested) (help)

• Einstein, A (1956). The Meaning of Relativity (5th. ed. ed.). Princeton, NJ: Princeton University Press. pp. pp. 92–97. ISBN 978-0-691-02352-6. {{cite book}}: |edition= has extra text (help); |pages= has extra text (help)

• Hagihara, Y (1931). "Theory of the relativistic trajectories in a gravitational field of Schwarzschild". Japanese Journal of Astronomy and Geophysics. 8: 67–176. ISSN 0368-346X.

• Lanczos, C. The Variational Principles of Mechanics (4th ed. ed.). New York: Dover Publications. pp. pp. 330–338. ISBN 978-0-486-65067-8. {{cite book}}: |edition= has extra text (help); |pages= has extra text (help)

• Landau, LD (1975). The Classical Theory of Fields (Course of Theoretical Physics, Vol. 2) (revised 4th English ed. ed.). New York: Pergamon Press. pp. pp. 299–309. ISBN 978-0-08-018176-9. {{cite book}}: |edition= has extra text (help); |pages= has extra text (help); Unknown parameter |coauthors= ignored (|author= suggested) (help)

• Misner, CW (1973). Gravitation. San Francisco: W. H. Freeman. pp. Chapter 25 (pp. 636–687), §33.5 (pp. 897–901), and §40.5 (pp. 1110–1116). ISBN 978-0-7167-0344-0. {{cite book}}: Unknown parameter |coauthors= ignored (|author= suggested) (help) (See Gravitation (book).)

• Pais, A. (1982). Subtle is the Lord: The Science and the Life of Albert Einstein. Oxford University Press. pp. pp. 253-256. ISBN 0-19-520438-7. {{cite book}}: |pages= has extra text (help)

• Pauli, W (1958). Theory of Relativity (translated by G. Field ed.). New York: Dover Publications. pp. pp. 40–41, 166–169. ISBN 978-0-486-64152-2. {{cite book}}: |pages= has extra text (help)

• Rindler, W (1977). Essential Relativity: Special, General, and Cosmological (revised 2nd ed. ed.). New York: Springer Verlag. pp. pp. 143–149. ISBN 978-0-387-10090-6. {{cite book}}: |edition= has extra text (help); |pages= has extra text (help)

• Roseveare, N. T (1982). Mercury's perihelion, from Leverrier to Einstein. Oxford: University Press. ISBN 0198581742.

• Synge, JL (1960). Relativity: The General Theory. Amsterdam: North-Holland Publishing. pp. pp. 289–298. ISBN 978-0-7204-0066-3. {{cite book}}: |pages= has extra text (help)

• Wald, RM (1984). General Relativity. Chicago: The University of Chicago Press. pp. pp. 136–146. ISBN 978-0-226-87032-8. {{cite book}}: |pages= has extra text (help); Cite has empty unknown parameter: |1= (help)

• Walter, S. (2007). "Breaking in the 4-vectors: the four-dimensional movement in gravitation, 1905–1910". In Renn, J. (ed.). The Genesis of General Relativity. Vol. 3. Berlin: Springer. pp. 193–252. {{cite book}}: External link in |chapterurl= (help); Unknown parameter |chapterurl= ignored (|chapter-url= suggested) (help)

• Weinberg, S (1972). Gravitation and Cosmology. New York: John Wiley and Sons. pp. pp. 185–201. ISBN 978-0-471-92567-5. {{cite book}}: |pages= has extra text (help)

• Whittaker, ET (1937). A Treatise on the Analytical Dynamics of Particles and Rigid Bodies, with an Introduction to the Problem of Three Bodies (4th ed. ed.). New York: Dover Publications. pp. pp. 389–393. ISBN 978-1-114-28944-4. {{cite book}}: |edition= has extra text (help); |pages= has extra text (help)