Finding the Path of Maximum Proper Time
In general relativity, a freely falling particle follows a geodesic — the path through spacetime that extremizes (maximizes) proper time. This is the relativistic generalization of “straight line” motion.
The proper time along a worldline is:
$$\tau = \int \sqrt{-g_{\mu\nu} \frac{dx^\mu}{d\lambda} \frac{dx^\nu}{d\lambda}} , d\lambda$$
The geodesic equation is:
$$\frac{d^2 x^\mu}{d\tau^2} + \Gamma^\mu_{\nu\sigma} \frac{dx^\nu}{d\tau} \frac{dx^\sigma}{d\tau} = 0$$
where the Christoffel symbols are:
$$\Gamma^\mu_{\nu\sigma} = \frac{1}{2} g^{\mu\lambda} \left( \partial_\nu g_{\lambda\sigma} + \partial_\sigma g_{\lambda\nu} - \partial_\lambda g_{\nu\sigma} \right)$$
Example Problem: Geodesics in Schwarzschild Spacetime
We consider a particle orbiting a non-rotating black hole. The Schwarzschild metric is:
$$ds^2 = -\left(1 - \frac{r_s}{r}\right)c^2 dt^2 + \left(1 - \frac{r_s}{r}\right)^{-1} dr^2 + r^2 d\Omega^2$$
where $r_s = 2GM/c^2$ is the Schwarzschild radius. Restricting to the equatorial plane ($\theta = \pi/2$), the conserved quantities are:
$$E = \left(1 - \frac{r_s}{r}\right) \frac{dt}{d\tau}, \quad L = r^2 \frac{d\phi}{d\tau}$$
The radial equation becomes:
$$\left(\frac{dr}{d\tau}\right)^2 = E^2 - \left(1 - \frac{r_s}{r}\right)\left(1 + \frac{L^2}{r^2}\right)$$
The effective potential is:
$$V_{\text{eff}}(r) = \left(1 - \frac{r_s}{r}\right)\left(1 + \frac{L^2}{r^2}\right)$$
The Innermost Stable Circular Orbit (ISCO) is at $r = 3r_s = 6GM/c^2$.
Python Code
1 | import numpy as np |
Code Walkthrough
Units & Setup. We use geometrized units $G = c = 1$, so the Schwarzschild radius is simply $r_s = 2M$.
geodesic_rhs. This function encodes the four geodesic equations for the Schwarzschild metric on the equatorial plane. The non-zero Christoffel symbols give three second-order ODEs:
$$\frac{d^2t}{d\tau^2} = -\frac{r_s}{r^2 f} \dot{t}\dot{r}, \quad \frac{d^2r}{d\tau^2} = -\frac{r_s f}{2r^2}\dot{t}^2 + \frac{r_s}{2r^2 f}\dot{r}^2 + rf\dot{\phi}^2, \quad \frac{d^2\phi}{d\tau^2} = -\frac{2}{r}\dot{r}\dot{\phi}$$
where $f = 1 - r_s/r$ and dots denote $d/d\tau$.
circular_IC. For a circular orbit, $\dot{r}=0$. The exact formula $\dot{t} = \sqrt{r/(r-\frac{3}{2}r_s)}$ follows from the geodesic condition and ensures the normalization $g_{\mu\nu}\dot{x}^\mu\dot{x}^\nu = -1$.
Three test geodesics. The circular ISCO orbit at $r=6M$ is the marginally stable case. A slight radial kick creates an eccentric orbit. A plunging orbit starts at $r=8M$ with inward velocity and reduced angular momentum — it spirals into the horizon.
Integration. solve_ivp with the DOP853 (8th-order Dormand–Prince) solver is used with tight tolerances (rtol=1e-10) to preserve the geodesic constraint. An event function terminates integration when $r \to r_s$.
Effective potential. $V_{\rm eff}(r) = (1-r_s/r)(1+L^2/r^2)$ encodes the competition between gravity and centrifugal repulsion. The ISCO sits at the inflection point where both the first and second derivatives vanish.
Figures & Results
Figure 1 — Orbital Trajectories
The circular ISCO orbit is a perfect ring at $r=6M$. The eccentric orbit oscillates radially around $6M$. The plunging orbit spirals inward and crosses the horizon.
Figure 2 — Effective Potential
The potential barrier shrinks for smaller $L$ (plunging case). At the ISCO, the barrier and well merge into a flat inflection point — the slightest perturbation sends the particle in.
Figure 3 — Proper Time vs Coordinate Time
The slope $d\tau/dt < 1$ everywhere, reflecting time dilation. The plunging orbit shows a dramatic slowdown in $\tau/t$ near the horizon: a distant observer sees the infalling particle “freeze” at $r_s$, while the particle’s own clock keeps ticking.
Figure 4 — 3D Spacetime Worldlines
Each worldline is a helix in $(x, y, t)$. The circular geodesic is a perfect helix with constant pitch. The plunging worldline curves sharply upward in $t$ (infinite coordinate time) as it asymptotically approaches the horizon.
Figure 5 — Time Dilation Factor $d\tau/dt$
For the circular ISCO orbit, the ratio is constant at $\sqrt{1 - 3M/r_{\rm ISCO}} = \sqrt{1/2} \approx 0.707$, meaning the orbiting clock ticks at 70.7% the rate of a clock at infinity — a combination of gravitational redshift and special-relativistic time dilation.
Execution Log — Geodesic Summary
1 | ======================================================= |
Reading the Numbers
Circular (ISCO). The ratio $\tau/t = 0.7071 = 1/\sqrt{2}$ is exact. It matches the analytical result for a circular orbit at $r = 6M$:
$$\frac{d\tau}{dt}\bigg|_{\rm ISCO} = \sqrt{1 - \frac{3M}{r}} = \sqrt{1 - \frac{1}{2}} = \frac{1}{\sqrt{2}} \approx 0.7071$$
This is a perfect sanity check — the numerical integrator reproduces the exact analytic value.
Eccentric. The ratio $\tau/t = 0.7174$ is slightly higher than the ISCO value. This seems counterintuitive at first, but makes sense: the eccentric orbit spends part of its time at larger $r$ where time dilation is weaker, pulling the average ratio up. However, $r_{\rm min} = 2.106 M$ tells us it dipped dangerously close to the horizon before being flung back out.
Plunging. The ratio drops sharply to $\tau/t = 0.5257$, and the orbit terminates in only $\tau = 14.3 M$ of proper time. From a distant observer’s perspective, coordinate time $t = 27.3 M$ passes — and near the horizon the particle appears to freeze entirely. The particle itself crosses the horizon in finite proper time with no drama.
Key Takeaways
The geodesic principle — maximizing proper time — is the spacetime analogue of the straight-line principle in flat space. In Schwarzschild spacetime it predicts precessing orbits, an innermost stable orbit at $r=6M$, and the dramatic time dilation experienced by observers near a black hole. All three effects emerge naturally from integrating a single set of second-order ODEs derived from the metric.