Kepler orbits and orbital mechanics

Kepler found his three laws by fitting curves to Tycho Brahe's observations of Mars, without knowing why they held. Newton showed that all three follow from an inverse-square attraction. This lesson completes that derivation, using the Binet equation from lesson 7, and then turns it into the working tools of orbital mechanics: the relation between energy and orbit size, the vis-viva equation for speed at any point, Kepler's equation for position as a function of time, and the Hohmann transfer, the most fuel-efficient way to move between circular orbits. The same mathematics plans the paths of satellites, interplanetary probes and the James Webb Space Telescope.

The orbit is a conic section

For an attractive inverse-square force \(F = -k/r^2\), the Binet equation of lesson 7 reads

This is the equation of a harmonic oscillator in the angle \(\phi\) with a constant term. The general solution is the constant particular solution plus an oscillation: \(u = (\mu k/\ell^2)\big(1 + e\cos(\phi - \phi_0)\big)\). Choosing the angle so that \(\phi_0 = 0\) and writing \(p = \ell^2/\mu k\),

This is the polar equation of a conic section with one focus at the center of force. The constant \(e\) is the eccentricity and \(p\), the semi-latus rectum, is the distance from the focus to the orbit measured perpendicular to the major axis. For \(e = 0\) the orbit is a circle; for \(0 < e < 1\) an ellipse; for \(e = 1\) a parabola; for \(e > 1\) one branch of a hyperbola. The orbit closes for \(e < 1\) because \(r\) depends on \(\phi\) through \(\cos\phi\), with period exactly \(2\pi\), as lesson 7's ratio \(\omega_r/\omega_\phi = 1\) required. This is Kepler's first law: planets move on ellipses with the Sun at one focus.

Energy, size and eccentricity

The semi-major axis \(a\) is half the sum of the perihelion and aphelion distances. Setting \(\phi = 0\) and \(\pi\) gives \(r_{\min} = p/(1 + e)\) and \(r_{\max} = p/(1 - e)\), so \(a = p/(1 - e^2)\), and the orbit is often written \(r = a(1 - e^2)/(1 + e\cos\phi)\). The turning points also satisfy \(Er^2 + kr - \ell^2/2\mu = 0\) from lesson 7, so their sum is \(-k/E\). Therefore

The energy depends only on the semi-major axis, not on the eccentricity. A circular orbit and a long thin ellipse with the same \(a\) have the same energy and, as we will see, the same period; angular momentum decides the shape. The second formula confirms the classification: \(e < 1\) exactly when \(E < 0\) (bound), \(e = 1\) when \(E = 0\) (parabolic escape), and \(e > 1\) when \(E > 0\) (hyperbolic).

Combining the energy result with \(E = \tfrac12\mu v^2 - k/r\) gives the speed at any distance, the vis-viva equation (from the Latin for living force):

written here for a small body orbiting a large mass \(M\), so that \(k/\mu = GM\). For a circular orbit \(r = a\) and \(v = \sqrt{GM/r}\). Escape from a circular orbit requires \(a \to \infty\), giving \(v = \sqrt{2GM/r}\), exactly \(\sqrt2\) times the circular speed. Vis-viva is the single most used formula in astrodynamics: given where a spacecraft is and how big an orbit you want, it tells you how fast it must go.

Kepler's third law

By Kepler's second law the area swept per unit time is \(\ell/2\mu\), and over one period \(T\) the whole ellipse, of area \(\pi ab\), is swept. So \(T = 2\pi ab\mu/\ell\). With \(b = a\sqrt{1 - e^2}\) and \(\ell^2 = \mu kp = \mu ka(1 - e^2)\), the eccentricity cancels:

The square of the period is proportional to the cube of the semi-major axis. Kepler stated this in 1619 as the same ratio for every planet; Newton's version adds the total mass, which differs slightly from planet to planet because each planet adds its own mass to the Sun's. The correction is small for the solar system (Jupiter adds 0.1%) but essential for binary stars, where the third law measured from the orbit gives the total mass of the pair. It is how almost all stellar masses are known, and how the mass of the black hole at the center of our galaxy, about four million solar masses, was found from the orbits of nearby stars (Nobel Prize 2020).

Position as a function of time

The orbit equation gives \(r\) as a function of angle, but navigation needs position as a function of time. Kepler solved this with a geometric trick. Introduce the eccentric anomaly \(E\) (not to be confused with energy), defined so that \(r = a(1 - e\cos E)\). The mean anomaly \(M = 2\pi t/T\), measured from perihelion, increases uniformly in time. They are related by Kepler's equation,

and the true angle \(\theta\) follows from \(\tan(\theta/2) = \sqrt{(1 + e)/(1 - e)}\tan(E/2)\). Kepler's equation cannot be solved for \(E\) in closed form, which made it one of the most studied equations in mathematics. Newton's method solves it in a few iterations: each step replaces \(E\) by \(E - (E - e\sin E - M)/(1 - e\cos E)\).

import math

GM_earth, R_earth = 3.986004e14, 6.371e6
GM_sun, AU = 1.32712e20, 1.495979e11

# 1. Kepler's equation M = E - e sin E solved by Newton's method, then the position on the orbit
def kepler_E(M, e):
    E = M if e < 0.8 else math.pi
    for _ in range(50):
        dE = (E - e * math.sin(E) - M) / (1 - e * math.cos(E))
        E -= dE
        if abs(dE) < 1e-14: break
    return E

a, e = 1.0, 0.6                                   # an orbit with period 1 in these units
for frac in (0.1, 0.25, 0.5):
    M = 2 * math.pi * frac
    E = kepler_E(M, e)
    nu = 2 * math.atan2(math.sqrt(1 + e) * math.sin(E / 2), math.sqrt(1 - e) * math.cos(E / 2))
    r = a * (1 - e * math.cos(E))
    print(f"t = {frac:.2f} T: eccentric anomaly {math.degrees(E):7.2f} deg, true anomaly {math.degrees(nu) % 360:7.2f} deg, r = {r:.4f} a")

# 2. Hohmann transfers: two burns, tangent to both circles
def hohmann(GM, r1, r2):
    at = (r1 + r2) / 2
    v1, v2 = math.sqrt(GM / r1), math.sqrt(GM / r2)
    vp, va = math.sqrt(GM * (2 / r1 - 1 / at)), math.sqrt(GM * (2 / r2 - 1 / at))   # vis-viva
    return vp - v1, v2 - va, math.pi * math.sqrt(at ** 3 / GM)

r_leo, r_geo = R_earth + 300e3, 42164e3
d1, d2, t = hohmann(GM_earth, r_leo, r_geo)
print(f"LEO (300 km) to GEO: burns {d1 / 1000:.3f} + {d2 / 1000:.3f} = {(d1 + d2) / 1000:.3f} km/s, coast {t / 3600:.2f} h")
d1, d2, t = hohmann(GM_sun, 1.000 * AU, 1.524 * AU)
print(f"Earth to Mars (heliocentric): burns {d1 / 1000:.3f} + {d2 / 1000:.3f} = {(d1 + d2) / 1000:.3f} km/s, coast {t / 86400:.0f} days")
print(f"escape speed from Earth's surface: {math.sqrt(2 * GM_earth / R_earth) / 1000:.2f} km/s")
print(f"geostationary radius from Kepler's third law: {(GM_earth * (86164.1 / (2 * math.pi)) ** 2) ** (1 / 3) / 1000:.0f} km")

t = 0.10 T: eccentric anomaly   67.84 deg, true anomaly  106.73 deg, r = 0.7737 a
t = 0.25 T: eccentric anomaly  119.82 deg, true anomaly  147.69 deg, r = 1.2984 a
t = 0.50 T: eccentric anomaly  180.00 deg, true anomaly  180.00 deg, r = 1.6000 a
LEO (300 km) to GEO: burns 2.428 + 1.468 = 3.895 km/s, coast 5.27 h
Earth to Mars (heliocentric): burns 2.946 + 2.650 = 5.596 km/s, coast 259 days
escape speed from Earth's surface: 11.19 km/s
geostationary radius from Kepler's third law: 42164 km

For an orbit with \(e = 0.6\), a tenth of the period after perihelion the body has already swept through \(106.7^\circ\) of true anomaly, while at a quarter period it has gone \(147.7^\circ\). Most of the angle is covered near perihelion, where the body moves fastest, exactly as the equal-areas law demands. At half the period it is at aphelion, \(r = 1.6a = a(1 + e)\).

Hohmann transfers

To move a satellite from a circular orbit of radius \(r_1\) to one of radius \(r_2\), the cheapest two-burn route in fuel is the transfer ellipse that just touches both circles, proposed by the German engineer Walter Hohmann in 1925. Its perihelion is \(r_1\) and its aphelion \(r_2\), so its semi-major axis is \(a_t = (r_1 + r_2)/2\). The first burn, at \(r_1\), raises the speed from the circular value \(\sqrt{GM/r_1}\) to the transfer-orbit speed given by vis-viva. After coasting half the transfer ellipse, the satellite arrives at \(r_2\) moving too slowly for a circular orbit there, and the second burn raises it to \(\sqrt{GM/r_2}\). The coast takes half the transfer orbit's period, \(\pi\sqrt{a_t^3/GM}\).

The program gives 3.90 km/s in total for a transfer from a 300 km orbit to geostationary orbit, split 2.43 then 1.47 km/s, with a coast of 5.3 hours. In practice a satellite launched from Cape Canaveral, at latitude \(28.5^\circ\) north, must also tilt its orbit into the equatorial plane, which is usually combined with the second burn and adds a few hundred meters per second. For a transfer from Earth's orbit to Mars's, treating both as circles around the Sun, the total is 5.60 km/s and the coast takes 259 days. A Hohmann transfer to Mars works only if Mars is in the right place when the spacecraft arrives, which happens once every synodic period of about 26 months. This is why missions to Mars launch in narrow windows roughly two years apart, as NASA's Perseverance rover did in July 2020.

Gravity assists and the limits of the two-body picture

A spacecraft passing close to a planet follows a hyperbola relative to the planet and leaves with the same speed relative to the planet that it arrived with. But the planet is moving around the Sun, and in the Sun's frame the spacecraft's velocity has been rotated relative to a moving body, which can add or remove up to twice the planet's orbital speed. Voyager 2 used assists at Jupiter, Saturn, Uranus and Neptune to tour all four giant planets between 1979 and 1989. The method shows the strength of the two-body solution: a complicated trajectory is built from pieces, each a conic section around one body, joined where one body's gravity takes over from another's.

The patched-conic method is an approximation. With three or more bodies interacting at comparable strengths the problem has no general closed-form solution, as Poincaré proved, and real mission design uses numerical integration of all the significant forces, including the Earth's equatorial bulge, the Moon and Sun, solar radiation pressure and atmospheric drag for low orbits. The bulge alone rotates the plane of a low orbit by several degrees per day, which engineers exploit to build Sun-synchronous orbits that keep the same lighting on the ground below.