Orbital Mechanics Math, Step by Step

From your first orbit equation to Lambert’s problem: derivations, worked examples and practice, freshman to PhD.

Constants used on this page

μ⊕=398 600.4418 km3/s2\mu_\oplus=398\,600.4418\ \text{km}^3/\text{s}^2
R⊕=6378.137 kmR_\oplus=6378.137\ \text{km}
J2=1.08262668×10−3J_2=1.08262668\times10^{-3}
g0=9.80665 m/s2g_0=9.80665\ \text{m/s}^2
Tsid=86 164.0905 sT_{sid}=86\,164.0905\ \text{s}

Every number on this page is computed by the site’s build script (scripts/learn_calc.py) and checked: LEO 300 km to GEO gives about 3.89 km/s, matching standard textbook results; the Lambert example is checked by integrating the orbit numerically.

Try it in the Orbit Lab

Circular orbits, period and escape speed

  • Teens
  • Community college
  • Undergraduate (freshman–senior)

Track: Orbital Mechanics · About 15 min · Numbers verified in build:

You will be able to

Derive the circular orbit speed from Newton’s laws; compute an orbital period; explain why escape speed is √2 times circular speed.

Before you start

Algebra, Newton’s second law, centripetal acceleration v²/r. Orbits Explained

Intuition first

In a circular orbit gravity supplies exactly the centripetal force needed to keep turning. The constant μ = GM (the gravitational parameter) is known far more precisely than G or M separately, so we always use μ.

🧒 Explain like I’m 10

A satellite is like a ball on a string you swing around your head. Gravity is the string. Swing too slowly and the ball drops; swing just right and it goes round and round.

Step by step

  1. GMmr2=mv2r\frac{GMm}{r^2}=\frac{mv^2}{r}
    Gravity = mass × centripetal acceleration.
  2. vc=μrv_c=\sqrt{\frac{\mu}{r}}
    Cancel m and solve for v.
  3. T=2πrvc=2πr3μT=\frac{2\pi r}{v_c}=2\pi\sqrt{\frac{r^3}{\mu}}
    Period = circumference / speed.
  4. 12vesc2−μr=0  ⇒  vesc=2μr=2 vc\tfrac12 v_{esc}^2-\frac{\mu}{r}=0\;\Rightarrow\; v_{esc}=\sqrt{\frac{2\mu}{r}}=\sqrt2\,v_c
    Escape: total energy just reaches zero.

Worked example

Given
Circular orbit at 400 km altitude.
Find
Speed v_c and period T.
Strategy
Use r = R⊕ + h, then v_c = √(μ/r) and T = 2πr/v_c.
Steps
  1. r = 6378.137 + 400 = 6778.137 km
  2. v_c = √(398600.4418 / 6778.137) = 7.6686 km/s
  3. T = 2π · 6778.137 / 7.6686 s = 92.56 min
Answer
v_c = 7.669 km/s, T = 92.6 min.
Sanity check
NASA lists the space station at roughly 7.7 km/s and about 90 minutes per orbit, so the numbers agree.
Significance
About 15–16 orbits a day: the basis for ground-station pass planning.
Show the Python
import math
mu, RE = 398600.4418, 6378.137   # km^3/s^2, km
r = RE + 400
v = math.sqrt(mu / r)
T = 2 * math.pi * r / v
print(v, T / 60)   # km/s, min

Auto-checked practice

  1. km/s

    Full solution

    r = 7178.137 km; v = √(398600.4418/7178.137) = 7.4518 km/s.

  2. min

    Full solution

    T = 2π√(7178.137³/398600.4418) s = 100.87 min.

Practice

  1. Find the radius of an orbit whose period is one sidereal day (86 164.09 s).

    Show solution

    r = (μ T²/4π²)^(1/3) = 42164 km (geostationary), where v_c = 3.075 km/s.

Common mistake

Using altitude instead of radius. r is measured from Earth’s centre: r = R⊕ + h.

References

  • NASA, Basics of Space Flight.
  • H. D. Curtis, Orbital Mechanics for Engineering Students, 4th ed., Butterworth-Heinemann (Elsevier), 2020.

Kepler’s three laws

  • Teens
  • Community college
  • Undergraduate (freshman–senior)

Track: Orbital Mechanics · About 15 min · Numbers verified in build:

You will be able to

State Kepler’s laws; show that equal areas in equal times is conservation of angular momentum; derive the third law.

Before you start

Ellipse geometry (a, e), cross product, derivative of an area. Orbits Explained

Intuition first

1) Orbits are conic sections with the central body at one focus. 2) The line from the focus sweeps equal areas in equal times. 3) T² is proportional to a³.

🧒 Explain like I’m 10

Planets speed up when they come close to the Sun and slow down when far away, like a skateboarder in a half-pipe: fastest at the bottom, slowest at the top.

Step by step

  1. dA=12 ∣r×dr∣  ⇒  dAdt=12∣r×v∣=h2dA=\tfrac12\,|\mathbf r\times d\mathbf r|\;\Rightarrow\;\frac{dA}{dt}=\tfrac12|\mathbf r\times\mathbf v|=\frac{h}{2}
    Area swept per unit time is half the specific angular momentum h.
  2. h˙=r˙×v+r×(−μr3r)=0\dot{\mathbf h}=\dot{\mathbf r}\times\mathbf v+\mathbf r\times\left(-\frac{\mu}{r^3}\mathbf r\right)=\mathbf 0
    For a central force h is constant, so dA/dt is constant (2nd law).
  3. T=πabh/2,h=μa(1−e2),b=a1−e2T=\frac{\pi ab}{h/2},\quad h=\sqrt{\mu a(1-e^2)},\quad b=a\sqrt{1-e^2}
    Total area of the ellipse divided by the areal rate.
  4. T=2πa3μT=2\pi\sqrt{\frac{a^3}{\mu}}
    Third law: depends only on a, not on e.

Worked example

Given
Geostationary radius r = 42 164.17 km (circular).
Find
The orbital period T.
Strategy
Kepler’s third law with a = r.
Steps
  1. T = 2π √(a³/μ) = 2π √(42164.17³ / 398600.4418)
  2. T = 23.9345 h
Answer
T = 23.9345 h (23 h 56 min).
Sanity check
It equals one sidereal day (86 164 s), exactly what "geostationary" requires.
Significance
Kepler’s law, found from planet data in 1619, still sets where every TV and weather satellite is parked.
Show the Python
import math
mu, RE = 398600.4418, 6378.137   # km^3/s^2, km
a = 42164.17
T = 2 * math.pi * math.sqrt(a**3 / mu)
print(T / 3600)   # h

Practice

  1. Two orbits have the same a but e = 0 and e = 0.7. Which has the longer period?

    Show solution

    Neither: T depends only on a.

Common mistake

Using 24 h for geostationary orbit. Earth turns once relative to the stars in 23 h 56 min 4 s.

References

  • H. D. Curtis, Orbital Mechanics for Engineering Students, 4th ed., Butterworth-Heinemann (Elsevier), 2020.
  • R. R. Bate, D. D. Mueller, J. E. White, W. W. Saylor, Fundamentals of Astrodynamics, 2nd ed., Dover, 2020.

The two-body problem and the orbit equation

  • Undergraduate (freshman–senior)

Track: Orbital Mechanics · About 30 min · Numbers verified in build:

You will be able to

Reduce two bodies to one relative equation; integrate it with vector calculus to get the conic r(ν); identify the eccentricity vector.

Before you start

Vector calculus (dot/cross products, BAC-CAB identity), first-order ODEs. Orbits Explained

Intuition first

Subtracting the two equations of motion leaves one ODE for the relative position. Two vector integrals, h and e, then give the whole shape of the orbit without solving for time.

🧒 Explain like I’m 10

Earth and a satellite both pull on each other. We can pretend Earth stands still and only watch the satellite, and the path always comes out as a circle, an oval, or a curve that flies away.

Step by step

  1. r¨=−μr3r,μ=G(m1+m2)\ddot{\mathbf r}=-\frac{\mu}{r^3}\mathbf r,\qquad \mu=G(m_1+m_2)
    Relative equation of motion.
  2. r¨×h=−μr3 r×(r×r˙)=μddt ⁣(rr)\ddot{\mathbf r}\times\mathbf h=-\frac{\mu}{r^3}\,\mathbf r\times(\mathbf r\times\dot{\mathbf r})=\mu\frac{d}{dt}\!\left(\frac{\mathbf r}{r}\right)
    Cross with constant h and use BAC-CAB.
  3. r˙×h−μrr=μ e\dot{\mathbf r}\times\mathbf h-\mu\frac{\mathbf r}{r}=\mu\,\mathbf e
    Integrate once: e is a constant vector pointing to periapsis.
  4. r⋅(r˙×h)=h2=μr+μrecos⁡ν  ⇒  r=h2/μ1+ecos⁡ν\mathbf r\cdot(\dot{\mathbf r}\times\mathbf h)=h^2=\mu r+\mu r e\cos\nu\;\Rightarrow\; r=\frac{h^2/\mu}{1+e\cos\nu}
    Dot with r: the conic equation.

Worked example

Given
r = [7000, 1000, 500] km, v = [−1.0, 7.4, 1.0] km/s.
Find
h, e and the orbit size a.
Strategy
Use the two vector integrals h = r × v and μe = v × h − μ r/r, then the energy.
Steps
  1. |h| = |r × v| = 53398.3 km²/s
  2. |e| = 0.019306
  3. ε = v²/2 − μ/r = -27.850 km²/s²
  4. a = −μ/(2ε) = 7156.15 km
Answer
Nearly circular (e = 0.0193) with a = 7156.1 km.
Sanity check
Check h²/μ = a(1 − e²): both sides agree to rounding.
Significance
Two conserved vectors give the full shape of the orbit from a single radar fix.
Show the Python
import numpy as np
mu = 398600.4418
r = np.array([7000., 1000., 500.]); v = np.array([-1.0, 7.4, 1.0])
h = np.cross(r, v)
e = np.cross(v, h) / mu - r / np.linalg.norm(r)
eps = v @ v / 2 - mu / np.linalg.norm(r)
print(np.linalg.norm(h), np.linalg.norm(e), -mu / (2 * eps))

Practice

  1. Show that e < 1 gives a closed orbit and e > 1 an open one.

    Show solution

    The denominator 1 + e cosν stays positive for every ν only if e < 1; for e > 1 it vanishes at cosν = -1/e, so r → ∞.

Common mistake

Forgetting that μ = G(m₁+m₂). For a satellite m₂ is negligible, but not for the Earth–Moon pair.

References

  • H. D. Curtis, Orbital Mechanics for Engineering Students, 4th ed., Butterworth-Heinemann (Elsevier), 2020.
  • D. A. Vallado, Fundamentals of Astrodynamics and Applications, 4th ed., Microcosm Press, 2013.
  • R. R. Bate, D. D. Mueller, J. E. White, W. W. Saylor, Fundamentals of Astrodynamics, 2nd ed., Dover, 2020.

Energy and the vis-viva equation

  • Community college
  • Undergraduate (freshman–senior)

Track: Orbital Mechanics · About 20 min · Numbers verified in build:

You will be able to

Show that specific orbital energy is conserved and equals -μ/2a; use vis-viva to find speed anywhere on an orbit.

Before you start

Kinetic and potential energy; dot product; the conic r(ν). Orbits Explained

Intuition first

Vis-viva ("living force") is the single most useful equation in mission design: it links speed, distance and orbit size.

🧒 Explain like I’m 10

Going up a hill you slow down; coming down you speed up. An orbit is the same: close to Earth you are fast, far away you are slow, and the total energy stays the same.

Step by step

  1. r˙⋅r¨+μr3 r⋅r˙=0  ⇒  ddt ⁣(v22−μr)=0\dot{\mathbf r}\cdot\ddot{\mathbf r}+\frac{\mu}{r^3}\,\mathbf r\cdot\dot{\mathbf r}=0\;\Rightarrow\;\frac{d}{dt}\!\left(\frac{v^2}{2}-\frac{\mu}{r}\right)=0
    Dot the equation of motion with velocity.
  2. ε=vp22−μrp,rp=a(1−e),vp=hrp\varepsilon=\frac{v_p^2}{2}-\frac{\mu}{r_p},\quad r_p=a(1-e),\quad v_p=\frac{h}{r_p}
    Evaluate at periapsis, where v is perpendicular to r.
  3. ε=−μ2a\varepsilon=-\frac{\mu}{2a}
    Simplify with h² = μa(1-e²).
  4. v2=μ(2r−1a)v^2=\mu\left(\frac{2}{r}-\frac{1}{a}\right)
    Vis-viva.

Worked example

Given
Transfer ellipse from r_p = 6678.137 km to r_a = 42 164.17 km.
Find
Speeds at perigee and apogee.
Strategy
a = (r_p + r_a)/2, then vis-viva at each end.
Steps
  1. a = 24421.2 km
  2. v_p = √(μ(2/r_p − 1/a)) = 10.1515 km/s
  3. v_a = √(μ(2/r_a − 1/a)) = 1.6078 km/s
Answer
v_p = 10.151 km/s, v_a = 1.608 km/s.
Sanity check
Angular momentum check: r_p v_p = r_a v_a (both about 67 793 km²/s).
Significance
These two speeds are exactly what the Hohmann burns are measured against.
Show the Python
import math
mu, RE = 398600.4418, 6378.137   # km^3/s^2, km
rp, ra = RE + 300, 42164.17
a = (rp + ra) / 2
vv = lambda r: math.sqrt(mu * (2/r - 1/a))
print(vv(rp), vv(ra), rp*vv(rp), ra*vv(ra))

Auto-checked practice

  1. km/s

    Full solution

    v = √(398600.4418 × (2/8000 − 1/10000)) = 7.7324 km/s.

Practice

  1. What is the sign of ε for an escape (parabolic) trajectory?

    Show solution

    Zero (a → ∞). Bound orbits have ε < 0, hyperbolas ε > 0.

Common mistake

Mixing units: with μ in km³/s², r and a must be in km and v comes out in km/s.

References

  • H. D. Curtis, Orbital Mechanics for Engineering Students, 4th ed., Butterworth-Heinemann (Elsevier), 2020.
  • R. R. Bate, D. D. Mueller, J. E. White, W. W. Saylor, Fundamentals of Astrodynamics, 2nd ed., Dover, 2020.

Orbital elements from a state vector

  • Undergraduate (freshman–senior)
  • Graduate

Track: Orbital Mechanics · About 30 min · Numbers verified in build:

You will be able to

Convert position and velocity (r, v) into the six classical elements a, e, i, Ω, ω, ν using linear algebra.

Before you start

Two-body problem, vectors, arccos quadrant checks. Orbits Explained

Intuition first

Six numbers fix an orbit: size (a), shape (e), tilt (i), where it crosses the equator going north (Ω), where periapsis lies (ω) and where the satellite is now (ν).

🧒 Explain like I’m 10

Six numbers are like an address for an orbit: how big, how stretched, how tilted, which way it faces, where its closest point is, and where the satellite is right now.

Step by step

  1. h=r×v,n=k^×h,e=v×hμ−rr\mathbf h=\mathbf r\times\mathbf v,\quad \mathbf n=\hat{\mathbf k}\times\mathbf h,\quad \mathbf e=\frac{\mathbf v\times\mathbf h}{\mu}-\frac{\mathbf r}{r}
    Three key vectors.
  2. a=−μ2ε,cos⁡i=hzh,cos⁡Ω=nxn  (ny<0⇒Ω=360∘−Ω)a=-\frac{\mu}{2\varepsilon},\quad \cos i=\frac{h_z}{h},\quad \cos\Omega=\frac{n_x}{n}\;(n_y<0\Rightarrow\Omega=360^\circ-\Omega)
    Size, tilt, node (with quadrant check).
  3. cos⁡ω=n⋅ene  (ez<0⇒360∘−ω),cos⁡ν=e⋅rer  (r⋅v<0⇒360∘−ν)\cos\omega=\frac{\mathbf n\cdot\mathbf e}{ne}\;(e_z<0\Rightarrow 360^\circ-\omega),\quad \cos\nu=\frac{\mathbf e\cdot\mathbf r}{er}\;(\mathbf r\cdot\mathbf v<0\Rightarrow 360^\circ-\nu)
    Periapsis and current position.

Worked example

Given
r = [7000, 1000, 500] km, v = [−1.0, 7.4, 1.0] km/s.
Find
a, e, i, Ω, ω, ν.
Strategy
Build h, n, e; read off each angle with its quadrant check (Vallado, Algorithm 9).
Steps
  1. i = arccos(h_z/h) = 8.585°
  2. Ω = 340.20° (n_y < 0)
  3. ω = 326.44° (e_z < 0)
  4. ν = 61.76° (r·v > 0)
  5. a = 7156.15 km, e = 0.019306
Answer
a = 7156.1 km, e = 0.0193, i = 8.59°, Ω = 340.2°, ω = 326.4°, ν = 61.8°.
Sanity check
Converting the elements back to r and v returns the inputs (round trip).
Significance
Elements are how orbits are published (e.g. two-line element sets) and compared.
Show the Python
import numpy as np
mu = 398600.4418
r = np.array([7000., 1000., 500.]); v = np.array([-1.0, 7.4, 1.0])
h = np.cross(r, v); n = np.cross([0, 0, 1], h)
e = np.cross(v, h) / mu - r / np.linalg.norm(r)
i = np.degrees(np.arccos(h[2] / np.linalg.norm(h)))
W = np.degrees(np.arccos(n[0] / np.linalg.norm(n))); W = 360 - W if n[1] < 0 else W
print(i, W)

Practice

  1. Why are Ω and ω undefined for an equatorial circular orbit?

    Show solution

    n = 0 (no node line) and e = 0 (no periapsis); use true longitude instead.

Common mistake

Skipping the quadrant checks: arccos only returns 0–180°.

References

  • D. A. Vallado, Fundamentals of Astrodynamics and Applications, 4th ed., Microcosm Press, 2013.
  • H. D. Curtis, Orbital Mechanics for Engineering Students, 4th ed., Butterworth-Heinemann (Elsevier), 2020.

Kepler’s equation and Newton–Raphson

  • Undergraduate (freshman–senior)
  • Graduate

Track: Orbital Mechanics · About 30 min · Numbers verified in build:

You will be able to

Relate time to position with mean, eccentric and true anomaly; solve the transcendental equation with Newton–Raphson.

Before you start

Derivatives, Taylor series to first order, radians. Orbits Explained

Intuition first

There is no closed-form x(t) for an ellipse. Kepler’s equation must be solved numerically; Newton–Raphson converges quadratically from E₀ = M (or M + e for large e).

🧒 Explain like I’m 10

Some puzzles can’t be solved in one step, so you guess, check how wrong you are, and fix the guess. Each try gets much closer, like a game of hot and cold.

Step by step

  1. M=n(t−tp),n=μ/a3M=n(t-t_p),\qquad n=\sqrt{\mu/a^3}
    Mean anomaly grows uniformly in time.
  2. M=E−esin⁡EM=E-e\sin E
    Kepler’s equation.
  3. Ek+1=Ek−Ek−esin⁡Ek−M1−ecos⁡EkE_{k+1}=E_k-\frac{E_k-e\sin E_k-M}{1-e\cos E_k}
    Newton–Raphson: f/f′.
  4. tan⁡ν2=1+e1−etan⁡E2\tan\frac{\nu}{2}=\sqrt{\frac{1+e}{1-e}}\tan\frac{E}{2}
    Convert to true anomaly.

Worked example

Given
e = 0.37255, M = 3.6029 rad (the same inputs as Curtis, Example 3.2; our own computation).
Find
The eccentric anomaly E to 10⁻¹⁰ rad.
Strategy
Newton–Raphson on f(E) = E − e sin E − M with f′(E) = 1 − e cos E, starting at E₀ = M (e ≤ 0.8); stop when |ΔE| < 10⁻¹⁰.
Steps
Iteration table
nE_n (rad)|ΔE|
13.47855373871.24e-01
23.47942201008.68e-04
33.47942204433.44e-08
43.47942204430.00e+00
Answer
E = 3.479422 rad after 4 iterations.
Sanity check
Substitute back: E − e sin E returns M = 3.6029 rad. Since M > π the body is past apoapsis, and E < M as expected when sin E < 0. Curtis reports E = 3.479 rad.
Significance
The error roughly squares each step (quadratic convergence), which is why flight software solves Kepler’s equation thousands of times per second.

Try it: Kepler solver, one step at a time

nE_n (rad)|ΔE|
Show the Python
import math
def kepler_E(M, e, tol=1e-10):
    E = math.pi if e > 0.8 else M
    while True:
        dE = (E - e*math.sin(E) - M) / (1 - e*math.cos(E))
        E -= dE
        print(E, abs(dE))
        if abs(dE) < tol:
            return E
kepler_E(3.6029, 0.37255)

Auto-checked practice

  1. rad

    Full solution

    Newton–Raphson converges in a few steps to E = 1.185324 rad.

Practice

  1. Why must M and E be in radians inside the iteration?

    Show solution

    Because the E term and the e sin E term must share units; e sin E is dimensionless (radians).

Common mistake

Mixing degrees and radians in M = E - e sin E.

References

  • D. A. Vallado, Fundamentals of Astrodynamics and Applications, 4th ed., Microcosm Press, 2013.
  • R. H. Battin, An Introduction to the Mathematics and Methods of Astrodynamics, rev. ed., AIAA, 1999.

Numerical propagation with RK4

  • Undergraduate (freshman–senior)
  • Graduate

Track: Orbital Mechanics · About 35 min · Numbers verified in build:

You will be able to

Write the orbit as a first-order ODE system; integrate it with the classical fourth-order Runge–Kutta method; check accuracy with energy conservation.

Before you start

ODEs, Taylor series, basic programming. Orbits Explained

Intuition first

Once perturbations are added (J2, drag, third bodies) there is no exact solution, so we integrate. RK4 has local error O(Δt⁵) and global error O(Δt⁴).

🧒 Explain like I’m 10

A computer moves the satellite forward in tiny steps, like a flip-book: each page is a fraction of a second later. Smarter steps (RK4) peek ahead before moving.

Step by step

  1. y=[rv],y˙=f(y)=[v−μr/r3]\mathbf y=\begin{bmatrix}\mathbf r\\ \mathbf v\end{bmatrix},\qquad \dot{\mathbf y}=\mathbf f(\mathbf y)=\begin{bmatrix}\mathbf v\\ -\mu\mathbf r/r^3\end{bmatrix}
    Second-order ODE becomes a 6-state first-order system.
  2. k1=f(yn), k2=f(yn+Δt2k1), k3=f(yn+Δt2k2), k4=f(yn+Δt k3)\mathbf k_1=\mathbf f(\mathbf y_n),\ \mathbf k_2=\mathbf f(\mathbf y_n+\tfrac{\Delta t}{2}\mathbf k_1),\ \mathbf k_3=\mathbf f(\mathbf y_n+\tfrac{\Delta t}{2}\mathbf k_2),\ \mathbf k_4=\mathbf f(\mathbf y_n+\Delta t\,\mathbf k_3)
    Four slope samples.
  3. yn+1=yn+Δt6(k1+2k2+2k3+k4)\mathbf y_{n+1}=\mathbf y_n+\frac{\Delta t}{6}\left(\mathbf k_1+2\mathbf k_2+2\mathbf k_3+\mathbf k_4\right)
    Weighted average (Simpson-like).

Worked example

Given
Initial state from the Lambert example: r₁ = [7000, 0, 0] km and its v₁; 3600 s; 20 000 RK4 steps.
Find
Final position and its distance from the target r₂ = [0, 8000, 0] km.
Strategy
Integrate ẏ = f(y) with RK4, Δt = 0.18 s.
Steps
  1. v₁ = (4.6069, 5.8532, 0) km/s
  2. Δt = 3600 / 20000 = 0.18 s
  3. |r(3600) − r₂| = 1.7×10⁻¹⁰ km
Answer
The orbit lands on the target to within 1.7×10⁻¹⁰ km.
Sanity check
Two independent methods (Lambert solver, RK4) agree, so both are working.
Significance
The same integrator, with J2 and drag added to f(y), is how real orbits are predicted.
Show the Python
import numpy as np
mu = 398600.4418
def f(y): r = y[:3]; return np.r_[y[3:], -mu * r / np.linalg.norm(r)**3]
def rk4(y, dt):
    k1 = f(y); k2 = f(y + dt/2*k1); k3 = f(y + dt/2*k2); k4 = f(y + dt*k3)
    return y + dt/6*(k1 + 2*k2 + 2*k3 + k4)

Practice

  1. Halving Δt reduces RK4 global error by roughly what factor?

    Show solution

    About 2⁴ = 16.

Common mistake

Trusting a long run without checking a conserved quantity (energy, h). RK4 is not symplectic and drifts slowly.

References

  • D. A. Vallado, Fundamentals of Astrodynamics and Applications, 4th ed., Microcosm Press, 2013.
  • MIT OpenCourseWare, 16.346 Astrodynamics (CC BY-NC-SA 4.0).

Hohmann transfer

  • Community college
  • Undergraduate (freshman–senior)

Track: Orbital Mechanics · About 25 min · Numbers verified in build:

You will be able to

Design a two-burn transfer between circular coplanar orbits; compute each Δv and the transfer time.

Before you start

Vis-viva, Kepler’s third law. Orbits Explained

Intuition first

Burn forward at the inner orbit to enter an ellipse that just touches the outer orbit, coast half an orbit, then burn forward again to circularise.

🧒 Explain like I’m 10

To move to a higher lane on a round track, you speed up once to swing outward, then speed up again when you reach the new lane so you stay there.

Step by step

  1. at=r1+r22a_t=\frac{r_1+r_2}{2}
    Transfer ellipse touches both circles.
  2. Δv1=μ(2r1−1at)−μr1\Delta v_1=\sqrt{\mu\left(\frac{2}{r_1}-\frac{1}{a_t}\right)}-\sqrt{\frac{\mu}{r_1}}
    Perigee burn (vis-viva minus circular).
  3. Δv2=μr2−μ(2r2−1at)\Delta v_2=\sqrt{\frac{\mu}{r_2}}-\sqrt{\mu\left(\frac{2}{r_2}-\frac{1}{a_t}\right)}
    Apogee burn.
  4. tt=πat3μt_{t}=\pi\sqrt{\frac{a_t^3}{\mu}}
    Half the transfer-orbit period.

Worked example

Given
μ = 398 600.4418 km³/s², start: circular 300 km (r₁ = 6678.1 km), target: GEO (r₂ = 42164.2 km), same plane.
Find
Δv₁, Δv₂, total Δv and transfer time.
Strategy
Circular speeds from v = √(μ/r); speeds on the transfer ellipse from vis-viva; subtract; half the ellipse period.
Steps
  1. a_t = (r₁ + r₂)/2 = 24421.2 km
  2. v_c1 = 7.7258 km/s, v_p = 10.1515 km/s → Δv₁ = 2.4257 km/s
  3. v_a = 1.6078 km/s, v_c2 = 3.0747 km/s → Δv₂ = 1.4668 km/s
  4. t = π√(a_t³/μ) = 18990 s
Answer
Δv = 3.893 km/s total; transfer time 5.28 h.
Sanity check
Matches the textbook ≈ 3.9 km/s. The transfer orbit’s period (5.28 h × 2) lies between the LEO and GEO periods, as it must.
Significance
Launching from 28.5° latitude adds a plane change; doing it at apogee with the second burn gives 4.26 km/s in total, the classic geostationary transfer budget.

Try it: Hohmann calculator

Show the Python
import math
mu, RE = 398600.4418, 6378.137          # km^3/s^2, km
r1, r2 = RE + 300, 42164.17
at = (r1 + r2) / 2
vv = lambda r, a: math.sqrt(mu * (2/r - 1/a))
dv1 = vv(r1, at) - math.sqrt(mu/r1)
dv2 = math.sqrt(mu/r2) - vv(r2, at)
tof = math.pi * math.sqrt(at**3 / mu)
print(dv1, dv2, dv1 + dv2, tof/3600)   # km/s, km/s, km/s, h

Auto-checked practice

  1. km/s

    Full solution

    Δv₁ = 0.1608, Δv₂ = 0.1575; total 0.3183 km/s.

Practice

  1. Why do real GEO missions need more than this total?

    Show solution

    Launch sites are not on the equator, so a plane change is also needed (see Plane changes).

Common mistake

Adding the circular speeds instead of differences: Δv is the change in speed at each burn.

References

  • H. D. Curtis, Orbital Mechanics for Engineering Students, 4th ed., Butterworth-Heinemann (Elsevier), 2020.
  • R. R. Bate, D. D. Mueller, J. E. White, W. W. Saylor, Fundamentals of Astrodynamics, 2nd ed., Dover, 2020.
  • NASA, Basics of Space Flight.

Bi-elliptic transfer

  • Undergraduate (freshman–senior)
  • Graduate

Track: Orbital Mechanics · About 20 min · Numbers verified in build:

You will be able to

Compute a three-burn bi-elliptic transfer and decide when it beats Hohmann.

Before you start

Hohmann transfer, vis-viva. Orbits Explained

Intuition first

Go far out to r_b, where speeds are tiny, then drop to the target. For r₂/r₁ above about 11.94 some choice of r_b beats Hohmann; above about 15.58 every r_b > r₂ does, at the cost of a much longer flight.

🧒 Explain like I’m 10

Sometimes the long way round is cheaper: climb very high where moving is easy, adjust there, then come back down to where you want to be.

Step by step

  1. a1=r1+rb2,a2=r2+rb2a_1=\frac{r_1+r_b}{2},\qquad a_2=\frac{r_2+r_b}{2}
    Two transfer ellipses.
  2. Δv=∣vp1−vc1∣+∣va2−va1∣+∣vc2−vp2∣\Delta v=\big|v_{p1}-v_{c1}\big|+\big|v_{a2}-v_{a1}\big|+\big|v_{c2}-v_{p2}\big|
    Burn at r₁, at r_b, and at r₂, each from vis-viva.

Worked example

Given
r₁ = 6678.137 km, r₂ = 20 r₁, intermediate r_b = 40 r₁.
Find
Total Δv of bi-elliptic vs Hohmann.
Strategy
Three vis-viva burns at r₁, r_b, r₂; compare with the two-burn Hohmann.
Steps
  1. a₁ = (r₁ + r_b)/2, a₂ = (r₂ + r_b)/2
  2. Δv_bi = 4.0609 km/s
  3. Δv_H = 4.1312 km/s
Answer
Bi-elliptic 4.061 km/s beats Hohmann 4.131 km/s.
Sanity check
r₂/r₁ = 20 is above the 15.58 threshold, so any r_b > r₂ must win: consistent.
Significance
The saving (about 70 m/s here) costs a far longer flight; mission design is a trade-off.
Show the Python
import math
mu, RE = 398600.4418, 6378.137   # km^3/s^2, km
vv = lambda r, a: math.sqrt(mu * (2/r - 1/a)); vc = lambda r: math.sqrt(mu / r)
r1 = RE + 300; r2, rb = 20*r1, 40*r1
a1, a2 = (r1+rb)/2, (r2+rb)/2
dv = (vv(r1,a1)-vc(r1)) + (vv(rb,a2)-vv(rb,a1)) + (vv(r2,a2)-vc(r2))
print(abs(dv))

Practice

  1. For LEO to GEO (r₂/r₁ ≈ 6.3), should you use a bi-elliptic transfer?

    Show solution

    No: below a ratio of about 11.94 Hohmann is always cheaper.

Common mistake

Ignoring flight time: the bi-elliptic saving can cost days or weeks.

References

  • D. A. Vallado, Fundamentals of Astrodynamics and Applications, 4th ed., Microcosm Press, 2013.
  • H. D. Curtis, Orbital Mechanics for Engineering Students, 4th ed., Butterworth-Heinemann (Elsevier), 2020.

Plane changes

  • Undergraduate (freshman–senior)

Track: Orbital Mechanics · About 20 min · Numbers verified in build:

You will be able to

Derive the plane-change Δv from the law of cosines; show why combining it with an apogee burn saves fuel.

Before you start

Vectors, law of cosines, Hohmann transfer. Orbits Explained

Intuition first

Turning the velocity vector costs speed proportional to how fast you are going, so plane changes are cheapest where the orbit is slowest.

🧒 Explain like I’m 10

Turning a fast bike is hard; turning a slow one is easy. Spacecraft tilt their orbits where they move slowest.

Step by step

  1. Δv=v12+v22−2v1v2cos⁡Δi\Delta v=\sqrt{v_1^2+v_2^2-2v_1v_2\cos\Delta i}
    Law of cosines on the two velocity vectors.
  2. v1=v2=v  ⇒  Δv=2vsin⁡Δi2v_1=v_2=v\;\Rightarrow\;\Delta v=2v\sin\frac{\Delta i}{2}
    Pure plane change at constant speed.

Worked example

Given
GEO insertion from a 28.5° transfer orbit: v_a = 1.6078 km/s, v_GEO = 3.0747 km/s.
Find
Δv for combined vs separate circularise-and-turn.
Strategy
Law of cosines on the two velocity vectors.
Steps
  1. Δv_comb = √(v_a² + v_GEO² − 2 v_a v_GEO cos 28.5°) = 1.8302 km/s
  2. Δv_sep = (v_GEO − v_a) + 2 v_GEO sin 14.25° = 2.9805 km/s
Answer
Combined: 1.830 km/s; separate: 2.980 km/s.
Sanity check
The combined burn is less than the sum, as the triangle inequality requires.
Significance
Saving about 1.15 km/s is why launches to GEO always combine the turn with the apogee burn.
Show the Python
import math
va, vg, di = 1.6078, 3.0747, math.radians(28.5)
comb = math.sqrt(va**2 + vg**2 - 2*va*vg*math.cos(di))
sep = (vg - va) + 2*vg*math.sin(di/2)
print(comb, sep)

Auto-checked practice

  1. km/s

    Full solution

    Δv = 2 × 7.5 × sin 5° = 1.3073 km/s.

Practice

  1. Where on an elliptical orbit is a plane change cheapest?

    Show solution

    At apoapsis (lowest speed), and it must be done at a node.

Common mistake

Forgetting that a plane change must be made where the two planes intersect (a node).

References

  • H. D. Curtis, Orbital Mechanics for Engineering Students, 4th ed., Butterworth-Heinemann (Elsevier), 2020.
  • D. A. Vallado, Fundamentals of Astrodynamics and Applications, 4th ed., Microcosm Press, 2013.

The rocket equation

  • Teens
  • Community college
  • Undergraduate (freshman–senior)

Track: Orbital Mechanics · About 15 min · Numbers verified in build:

You will be able to

Derive Tsiolkovsky’s equation from conservation of momentum by separating variables; turn a Δv budget into propellant mass.

Before you start

Momentum, integrals of 1/m (natural log). Orbits Explained

Intuition first

Δv is the currency of spaceflight. Because the rocket must carry its own propellant, mass grows exponentially with Δv.

🧒 Explain like I’m 10

A rocket pushes stuff out the back really fast, and that pushes the rocket forward, like letting go of a balloon. But it must carry all that stuff, so going faster needs a lot more fuel.

Step by step

  1. m dv=−ve dmm\,dv=-v_e\,dm
    Momentum balance over a short burn (no gravity or drag).
  2. ∫0Δvdv=−ve∫m0mfdmm\int_{0}^{\Delta v}dv=-v_e\int_{m_0}^{m_f}\frac{dm}{m}
    Separate variables and integrate.
  3. Δv=veln⁡m0mf,ve=Isp g0\Delta v=v_e\ln\frac{m_0}{m_f},\qquad v_e=I_{sp}\,g_0
    Tsiolkovsky; g₀ = 9.80665 m/s².

Worked example

Given
Δv = 3.8926 km/s (LEO→GEO Hohmann), I_sp = 450 s, g₀ = 9.80665 m/s².
Find
Mass ratio and propellant fraction.
Strategy
v_e = I_sp g₀, then m₀/m_f = exp(Δv/v_e).
Steps
  1. v_e = 450 · 9.80665 m/s = 4.413 km/s
  2. m₀/m_f = exp(3.8926 / 4.413) = 2.4159
  3. 1 − 1/2.4159 = 58.6 %
Answer
m₀/m_f = 2.416; 58.6% of the starting mass is propellant.
Sanity check
Δv is smaller than v_e, so the ratio must be below e ≈ 2.72: it is.
Significance
More than half the stage is fuel just to reach GEO from LEO, which is why upper stages use high-I_sp engines.
Show the Python
import math
dv, isp, g0 = 3.8926, 450, 9.80665e-3   # km/s, s, km/s^2
mr = math.exp(dv / (isp * g0))
print(mr, 1 - 1/mr)

Auto-checked practice

    Full solution

    m₀/m_f = exp(2/2.942) = 1.974.

Practice

  1. Doubling Δv at fixed I_sp does what to the mass ratio?

    Show solution

    It is squared, not doubled: e^(2x) = (e^x)².

Common mistake

Using g at altitude for g₀. g₀ in I_sp is a fixed conversion constant, 9.80665 m/s².

References

  • H. D. Curtis, Orbital Mechanics for Engineering Students, 4th ed., Butterworth-Heinemann (Elsevier), 2020.
  • NASA, Basics of Space Flight.

J2 perturbations and sun-synchronous orbits

  • Graduate

Track: Orbital Mechanics · About 40 min · Numbers verified in build:

You will be able to

Model Earth’s equatorial bulge with the J2 term of the gravity potential; derive secular node drift by averaging; design a sun-synchronous orbit.

Before you start

Gradients of a potential, Legendre polynomials, Lagrange planetary equations (statement). Orbits Explained

Intuition first

The extra mass at the equator pulls on the orbit and makes the orbital plane precess. Averaged over one orbit, the short-period effects cancel and a steady drift of Ω and ω remains.

🧒 Explain like I’m 10

Earth is a little fat around the middle. That extra bulge tugs on satellites and slowly turns their orbits, and engineers use this to keep some satellites in step with the Sun.

Step by step

  1. U=−μr[1−J2(R⊕r)23sin⁡2ϕ−12]U=-\frac{\mu}{r}\left[1-J_2\left(\frac{R_\oplus}{r}\right)^2\frac{3\sin^2\phi-1}{2}\right]
    Potential with the J2 zonal term (φ = latitude).
  2. Ω˙=−32 n J2(R⊕p)2cos⁡i,p=a(1−e2)\dot\Omega=-\frac{3}{2}\,n\,J_2\left(\frac{R_\oplus}{p}\right)^2\cos i,\qquad p=a(1-e^2)
    Secular node rate after orbit-averaging (Lagrange equations).
  3. Ω˙SSO=360∘365.2422 d≈0.9856∘/d  ⇒  cos⁡i=Ω˙SSO−32nJ2(R⊕/p)2\dot\Omega_{SSO}=\frac{360^\circ}{365.2422\ \text{d}}\approx 0.9856^\circ/\text{d}\;\Rightarrow\;\cos i=\frac{\dot\Omega_{SSO}}{-\tfrac32 nJ_2(R_\oplus/p)^2}
    Sun-synchronous: the node follows the Sun once a year.

Worked example

Given
Circular orbit, 700 km altitude; J2 = 1.08263×10⁻³; target node rate 360°/365.2422 d.
Find
The sun-synchronous inclination i.
Strategy
Set the secular node rate equal to the target and solve for cos i.
Steps
  1. a = 7078.137 km, n = √(μ/a³)
  2. cos i = −Ω̇_SSO / (1.5 n J2 (R⊕/a)²)
  3. i = 98.188°
Answer
i = 98.19°.
Sanity check
Published sun-synchronous missions near 700 km fly at about 98.2°, matching.
Significance
A small natural perturbation, used on purpose, gives Earth-imaging satellites constant lighting for free.
Show the Python
import math
mu, RE = 398600.4418, 6378.137   # km^3/s^2, km
J2 = 1.08262668e-3
a = RE + 700; n = math.sqrt(mu / a**3)
target = 2*math.pi / (365.2422*86400)
print(math.degrees(math.acos(-target / (1.5*n*J2*(RE/a)**2))))

Auto-checked practice

  1. °

    Full solution

    cos i is negative and i = 98.60°, slightly higher than at 700 km.

Practice

  1. Why are sun-synchronous orbits slightly retrograde (i > 90°)?

    Show solution

    The node must move east (positive rate); the formula has -cos i, so cos i < 0.

Common mistake

Using the solar day: the node must turn 360° per tropical year, about 0.9856°/day.

References

  • D. A. Vallado, Fundamentals of Astrodynamics and Applications, 4th ed., Microcosm Press, 2013.
  • N. K. Pavlis et al., "The development and evaluation of the Earth Gravitational Model 2008 (EGM2008)", J. Geophys. Res. 117, B04406 (2012). doi:10.1029/2011JB008916
  • R. H. Battin, An Introduction to the Mathematics and Methods of Astrodynamics, rev. ed., AIAA, 1999.

Lambert’s problem

  • Graduate
  • PhD / research

Track: Orbital Mechanics · About 60 min · Numbers verified in build:

You will be able to

Given two positions and a flight time, find the connecting orbit using universal variables and a root finder; verify by numerical propagation; know the research literature.

Before you start

Universal variables, Stumpff functions, root finding (bisection, Newton, Halley), RK4. Orbits Explained

Intuition first

Lambert’s problem is the core of rendezvous, targeting and interplanetary "porkchop" plots. It is a two-point boundary-value problem; Lambert’s theorem says the flight time depends only on a, r₁ + r₂ and the chord c.

🧒 Explain like I’m 10

You want to throw a ball from here to there and have it arrive in exactly one hour. Lambert’s problem finds the exact throw.

Step by step

  1. A=sin⁡Δνr1r21−cos⁡Δν,y(z)=r1+r2+A zS(z)−1C(z)A=\sin\Delta\nu\sqrt{\frac{r_1r_2}{1-\cos\Delta\nu}},\qquad y(z)=r_1+r_2+A\,\frac{zS(z)-1}{\sqrt{C(z)}}
    Geometry constant A and auxiliary y(z) (Stumpff C, S).
  2. μ Δt=(yC)3/2S+Ay\sqrt{\mu}\,\Delta t=\left(\frac{y}{C}\right)^{3/2}S+A\sqrt{y}
    Time-of-flight equation: solve for z by root finding.
  3. f=1−yr1,g=Ayμ,g˙=1−yr2,v1=r2−fr1g,v2=g˙ r2−r1gf=1-\frac{y}{r_1},\quad g=A\sqrt{\frac{y}{\mu}},\quad \dot g=1-\frac{y}{r_2},\quad \mathbf v_1=\frac{\mathbf r_2-f\mathbf r_1}{g},\quad \mathbf v_2=\frac{\dot g\,\mathbf r_2-\mathbf r_1}{g}
    Lagrange coefficients give the velocities.

Worked example

Given
r₁ = [7000, 0, 0] km, r₂ = [0, 8000, 0] km, Δt = 3600 s, prograde.
Find
v₁ and v₂.
Strategy
Universal-variable time-of-flight equation; find z by bisection; Lagrange f, g, ġ give the velocities (Curtis, Algorithm 5.2).
Steps
  1. Δν = 90°, A = √(r₁r₂)
  2. v₁ = (4.6069, 5.8532, 0) km/s
  3. v₂ = (-5.1216, -3.8753, 0) km/s
  4. RK4 check: 1.7×10⁻¹⁰ km
Answer
v₁ = (4.607, 5.853, 0) km/s; v₂ = (-5.122, -3.875, 0) km/s.
Sanity check
Energy and angular momentum computed at both ends are equal, and RK4 from (r₁, v₁) hits r₂.
Significance
Repeating this over a grid of departure and arrival dates gives the "porkchop" plots used to plan Mars missions.
Show the Python
# universal-variable Lambert (Curtis Alg. 5.2): see scripts/learn_calc.py lambert()
from learn_calc import lambert, propagate
v1, v2 = lambert([7000, 0, 0], [0, 8000, 0], 3600)
print(v1, v2)

Practice

  1. Research: what does Izzo (2015) change compared with the universal-variable method?

    Show solution

    A single non-dimensional variable x with a good initial guess and Householder (third-order) iterations, handling multi-revolution cases robustly in very few iterations.

Common mistake

Picking the wrong transfer direction: Δν and 360° - Δν give two different orbits.

References

  • D. Izzo, "Revisiting Lambert’s problem", Celestial Mechanics and Dynamical Astronomy 121, 1–15 (2015). doi:10.1007/s10569-014-9587-y
  • R. H. Battin, An Introduction to the Mathematics and Methods of Astrodynamics, rev. ed., AIAA, 1999.
  • D. A. Vallado, Fundamentals of Astrodynamics and Applications, 4th ed., Microcosm Press, 2013.
  • H. D. Curtis, Orbital Mechanics for Engineering Students, 4th ed., Butterworth-Heinemann (Elsevier), 2020.

Check your understanding

Doubling a satellite’s orbit radius makes its circular speed…

Where is a plane change cheapest on an elliptical orbit?

Which anomaly grows uniformly with time?

Go deeper

Next: the aeronautical engineering and physics tracks (coming), Vallado’s algorithms and Izzo’s Lambert solver for research-level work.

For teachers and instructors

Suggested use: one 50-minute session per two sections. Undergraduate path: circular → Kepler → vis-viva → Hohmann → rocket equation, with the Orbit Lab as a live demo of the Hohmann worked example. Graduate path: elements → Kepler’s equation → RK4 → J2 → Lambert (assign the Python snippets as labs). Every numeric answer is computed by the site’s build script, so the practice answers double as an answer key. Aligned with MIT OCW 16.346 Astrodynamics topics.

References

Help / Donate