Important Black Hole Surfaces¶

This notebook will plot important surfaces of black holes, i.e. the location of the event horizon and innermost stable circular orbit as a function of black hole spin.

Units are $G=c=M=1$, with radii in $r_g=GM/c^2$. All surfaces are in the midplane ($\theta=\pi/2$). Retrograde orbits denoted by $a<0$.

Set-up¶

In [1]:

import numpy as np
import matplotlib.pyplot as plt
plt.style.use('lia.mplstyle')
from bh_accretion_toolkit.gr import kerr_bl as bl
from bh_accretion_toolkit import plotting
from bh_accretion_toolkit.gr import horizons as hz
from bh_accretion_toolkit.gr import orbits
from bh_accretion_toolkit.constants import cgs
from scipy.optimize import brentq
from importlib import reload

spin_range = np.linspace(-1., 1., 200)

Event Horizon¶

The most notorious! The event horizon is where no light can escape. Mathematically, it occurs where $g_{rr} = \rho^2/\Delta$ blows up, i.e. where $$\Delta = r^2 - 2r + a^2 = 0 \quad\Rightarrow\quad r_\pm = 1 \pm \sqrt{1 - a^2}.$$ The outer root $r_+$ is the event horizon. The inner (Cauchy) horizon $r_-$ is hidden inside it and never matters for accretion. For $|a| > 1$ there is no real root and so no horizon (a "naked singularity"), which is why spin is bounded by $|a| \le 1$. In practice, accretion spins a black hole up only to $a \approx 0.998$, because the disk photons it swallows preferentially carry angular momentum that opposes the spin (Thorne 1974).

The horizon depends only on $a^2$, so prograde and retrograde disks see the same one: it shrinks from $2\,r_g$ at $a = 0$ to $1\,r_g$ at $|a| = 1$.

Approaching the horizon, the lapse function $\alpha$ (the rate at which a local, zero-angular-momentum observer's clock ticks relative to one at infinity) drops to zero. Light emitted near the horizon is redshifted without limit on its way out, so infalling gas fades from view rather than being seen to cross.

For more info:

Andrew Hamilton's webpage on falling into a black hole
kerr_metric_properties notebook for more on the lapse function.

In [2]:

a_values = [0., 0.5, 0.9, 0.998]

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(10, 4), layout='constrained')
ax1.plot(spin_range, hz.get_EH_location(spin_range), color='k', label=r"$r_+$ (event horizon)")
ax1.plot(spin_range, hz.get_inner_horizon_location(spin_range), color='k', ls='--', label=r"$r_-$ (inner horizon)")
ax1.set_xlim(spin_range[0], spin_range[-1])
ax1.set_ylim(0, 2.1)
ax1.set_xlabel(r"Black Hole Spin $a$")
ax1.set_ylabel(r"Radial Coordinate $r/r_g$")
ax1.legend()

for a in a_values:
    r = np.linspace(hz.get_EH_location(a)*(1. + 1e-6), 10., 1000)
    ax2.plot(r, bl.lapse_function(r, np.pi/2, a), label=fr"$a = {a}$")
ax2.set_xlim(0, 10)
ax2.set_ylim(0, 1)
ax2.set_xlabel(r"Radial Coordinate $r/r_g$")
ax2.set_ylabel(r"Midplane lapse Function $\alpha$")
ax2.legend()
plotting.save_fig("event_horizon", "01_black_hole_surfaces")

Out[2]:

[PosixPath('figures/01_black_hole_surfaces/event_horizon.pdf'),
 PosixPath('figures/01_black_hole_surfaces/event_horizon.png')]

No description has been provided for this image

Ergosphere¶

The ergosphere is important for extracting energy from the black hole, through the Penrose process or through the corresponding Blandford-Znajek mechanism, which uses magnetic fields threading the horizon to tap the hole's spin. Blandford-Znajek is the leading model for powering relativistic jets (Blandford & Znajek 1977).

Both mechanisms rely on frame dragging. Inside the ergosphere, spacetime rotates so fast that not even light can stay at rest relative to the distant stars. Unlike the horizon, though, the ergosphere can still be escaped. In the Penrose process (Penrose 1969), a particle splits in two inside the ergosphere. One piece falls in on an orbit with negative energy (as measured from infinity), which is only possible inside the ergosphere. That means the other piece escapes with more energy than the original. Energy is conserved by taking some out of the black hole's rotational energy, spinning it down ever so slightly.

Mathematically, the ergosphere occurs where a static observer's worldline would have to become null, i.e. where $$g_{tt} = -\left(1 - \frac{2r}{\rho^2}\right) = 0 \quad\Rightarrow\quad r_{\rm erg} = 1 + \sqrt{1 - a^2\cos^2\theta}.$$ It touches the horizon at the poles and reaches $r = 2$ in the midplane for every spin, so it is oblate. For $a = 0$ it coincides with the horizon: a non-spinning hole has no ergosphere. There's no rotational energy to extract.

The below plot shows the cross-section of the ergosphere for two different spins (left panel), and the ergosphere location as a function of spin at three theta locations (right panel). The cross-section plots $x = r\sin\theta$, $z = r\cos\theta$ as a schematic, since the Boyer-Lindquist $r$ is not a Cartesian radius.

Notice that the midplane ergosphere radius is independent of spin. What does that mean?

For more:

Lowell et al. 2024 on magnetic fields noticeably changing black hole spin

In [ ]:

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(10, 4.5), layout='constrained')

# Left: meridional cross-section, a = 0.7 on the left half and a = 0.998 on the right half
theta = np.linspace(0., np.pi, 300)
for side, a, color in [(-1, 0.7, 'C0'), (1, 0.998, 'C1')]:
    r_erg, r_plus = hz.get_ergosphere_location(a, theta), hz.get_EH_location(a)
    ax1.fill(side*r_erg*np.sin(theta), r_erg*np.cos(theta), color=color, alpha=0.3, lw=0)
    ax1.plot(side*r_erg*np.sin(theta), r_erg*np.cos(theta), color=color)
    ax1.fill(side*r_plus*np.sin(theta), r_plus*np.cos(theta), color='k')
    ax1.text(side*1.25, 2.1, fr"$a = {a}$", ha='center', color=color)
ax1.set_aspect('equal')
ax1.set_xlim(-2.5, 2.5)
ax1.set_ylim(-2.5, 2.5)
ax1.set_xlabel(r"$r\sin\theta$")
ax1.set_ylabel(r"$r\cos\theta$")
ax1.set_title("Shaded: ergosphere. Black: horizon.")

# Right: ergosphere radius vs spin at three polar angles
for th, label in [(np.pi/2, r"$\theta = \pi/2$ (midplane)"), (np.pi/4, r"$\theta = \pi/4$"), (np.pi/8, r"$\theta = \pi/8$")]:
    ax2.plot(spin_range, hz.get_ergosphere_location(spin_range, th), label=label)
ax2.fill_between(spin_range, 1., hz.get_EH_location(spin_range), color='0.8', label="Inside event horizon")
ax2.plot(spin_range, hz.get_EH_location(spin_range), color='k')
ax2.set_xlim(spin_range[0], spin_range[-1])
ax2.set_ylim(1, 2.2)
ax2.set_xlabel(r"Black Hole Spin $a$")
ax2.set_ylabel(r"Ergosphere Radius $r_{\rm erg}/r_g$")
ax2.legend(loc='lower center')
plotting.save_fig("ergosphere", "01_black_hole_surfaces")

Photon orbit¶

Astrophysical relevance: the photon orbit is what the Event Horizon Telescope (EHT) sees. Light passing near it can loop around the hole before escaping, piling up into a thin, bright photon ring around a dark "shadow". For a non-spinning hole the shadow radius is $\sqrt{27}\,r_g \approx 5.2\,r_g$ (larger than $r_{\rm ph} = 3\,r_g$ because of lensing), and it changes only weakly with spin. So we measure a ring size of $GM/c^2D$, the black hole's mass over distance (see the M87* estimate below).

Mathematical origin: the energy per unit rest mass of a circular orbit (Bardeen, Press & Teukolsky 1972) is $$E_{\rm circ} = \frac{1 - 2/r + a\,r^{-3/2}}{\sqrt{1 - 3/r + 2a\,r^{-3/2}}}.$$ It diverges where the denominator vanishes. Because a massive particle would need infinite energy to orbit there, only light can orbit at that radius. Solving $1 - 3/r + 2a\,r^{-3/2} = 0$ gives $r_{\rm ph} = 2\left[1 + \cos\left(\tfrac{2}{3}\arccos(-a)\right)\right]$, from $3\,r_g$ at $a = 0$ down to $1\,r_g$ at $a = 1$ (prograde) and up to $4\,r_g$ at $a = -1$ (retrograde). The orbit is unstable: a slightly perturbed photon either falls in or escapes, which is why the ring is thin.

For more:

Bardeen, Press, Teukolsky 1972 Eq. 2.18 for derivation of photon orbit
Notebooks on 02_effective_potential for circular orbits

In [ ]:

fig, ax = plt.subplots()
for a in [-0.9, 0., 0.9]:
    r_ph = hz.get_photon_orbit_location(a)
    r = np.linspace(r_ph*(1. + 1e-4), 12., 1000)
    line, = ax.plot(r, orbits.circular_orbit_energy(r, a), label=fr"$a = {a}$")
    ax.axvline(r_ph, color=line.get_color(), ls=':')
ax.set_yscale('log')
ax.set_xlim(0, 12)
ax.set_ylim(0.8, 30)
ax.set_xlabel(r"Radial Coordinate $r/r_g$")
ax.set_ylabel(r"Circular Orbit Energy $E_{\rm circ}$")
ax.set_title(r"$E_{\rm circ} \to \infty$ at the photon orbit (dotted)")
ax.legend()
plotting.save_fig("photon_orbit_energy", "01_black_hole_surfaces")

# M87*: EHT mass and distance (EHT Collaboration 2019, Paper I)
M, D = 6.5e9*cgs['MSOLAR'], 16.8e6*cgs['PC']
theta_g = cgs['GNEWT']*M/cgs['CL']**2/D * 206264.806e6  # radians -> microarcseconds
print(f"GM/c^2D = {theta_g:.1f} uas, so the predicted shadow diameter is 2 sqrt(27) GM/c^2D = {2.*np.sqrt(27.)*theta_g:.0f} uas")
print("EHT measured a ring diameter of 42 +/- 3 uas")

GM/c^2D = 3.8 uas, so the predicted shadow diameter is 2 sqrt(27) GM/c^2D = 40 uas
EHT measured a ring diameter of 42 +/- 3 uas

Marginally bound orbit¶

Physical meaning: the marginally bound orbit is the circular orbit with $E_{\rm circ} = 1$, exactly the energy of a particle falling in from rest at infinity. Circular orbits between $r_{\rm ph}$ and $r_{\rm mb}$ have $E_{\rm circ} > 1$: they exist, but they are unbound, and an outward nudge sends the particle to infinity. Equivalently, $r_{\rm mb}$ is the closest a parabolic orbit can come to the hole without plunging in. From $r_{\rm mb} = 2 - a + 2\sqrt{1 - a}$, it is $4\,r_g$ at $a = 0$, $1\,r_g$ at $a = 1$, and $3 + 2\sqrt{2} \approx 5.8\,r_g$ at $a = -1$.

Astrophysical relevance:

Tidal disruption events: stars reach the black hole on nearly parabolic orbits. If a star's tidal radius lies inside $r_{\rm mb}$, it is swallowed whole instead of being torn apart, and no flare is seen. This sets the maximum black hole mass that can produce a flare (the Hills mass, $\sim 10^8\,M_\odot$ for a Sun-like star), which grows with prograde spin because $r_{\rm mb}$ shrinks.

Innermost stable circular orbit¶

Physical meaning: the ISCO is the smallest radius at which a circular orbit is stable. Outside it, gas in a thin disk spirals slowly inward through nearly circular orbits as turbulence removes its angular momentum, radiating away its binding energy as it goes. At the ISCO the orbits become unstable, and the gas plunges into the hole on a dynamical time. That's why the region inside the ISCO is often referred to as the "plunging region".

Mathematically, the ISCO is where $E_{\rm circ}(r)$ has its minimum ($dE_{\rm circ}/dr = 0$). Outside the minimum, a slightly perturbed orbit oscillates about its circular radius. Inside, the effective potential has no minimum left to trap it. This gives $r_{\rm ISCO} = 3 + Z_2 \mp \sqrt{(3 - Z_1)(3 + Z_1 + 2Z_2)}$ (Bardeen, Press & Teukolsky 1972), from $6\,r_g$ at $a = 0$ down to $1\,r_g$ at $a = 1$ and up to $9\,r_g$ at $a = -1$. For a prograde orbit, the hole's rotation acts like extra centrifugal support, so stable orbits survive closer in.

Curiosity: for $a > 2\sqrt{2}/3 \approx 0.943$, the prograde ISCO lies inside the ergosphere ($r_{\rm ISCO} < 2$). Why that value? Please contact me if you have an intuitive reason. The gas co-rotates anyway, so its orbit barely notices, but the inner disk, and the magnetic fields anchored in it, now sit in the frame-dragging region where energy can be extracted from the hole's spin. A retrograde ISCO never enters the ergosphere.

Astrophysical relevance:

Often considered inner edge of thin accretion disk. Can be used to measure black hole spin because of strong dependence onf $a$. Open questions about radiation from within the plunging region and actual stress at the ISCO

TLDR; Overview of Surfaces¶

See e.g. Reynolds 2021 review Fig. 1.

Event horizon: $r_+ = 1 + \sqrt{1 - a^2}$: the outer root of $\Delta = r^2 - 2r + a^2 = 0$. Nothing inside can escape. See Boyer & Lindquist 1967.
Ergosphere (static limit): $r_{\rm erg} = 1 + \sqrt{1 - a^2\cos^2\theta}$: where $g_{tt} = 0$. Inside it, frame dragging forces everything to co-rotate with the hole. In the midplane, $r_{\rm erg} = 2$ for every spin. See Boyer & Lindquist 1967.
Photon orbit: $r_{\rm ph} = 2\left[1 + \cos\left(\tfrac{2}{3}\arccos(-a)\right)\right]$: the innermost circular orbit, where $E_{\rm circ} \to \infty$, so only light can follow it. See Bardeen, Press & Teukolsky 1972, Eq. 2.18.
Marginally bound orbit: $r_{\rm mb} = 2 - a + 2\sqrt{1 - a}$: the circular orbit with $E_{\rm circ} = 1$. Circular orbits inside it are unbound. See Bardeen, Press & Teukolsky 1972, Eq. 2.19.
Innermost stable circular orbit (ISCO): $r_{\rm ISCO} = 3 + Z_2 \mp \sqrt{(3 - Z_1)(3 + Z_1 + 2Z_2)}$ (upper sign prograde): circular orbits inside it are unstable, so gas plunges in if perturbed. It is the usual inner edge of a thin disk. See Bardeen, Press & Teukolsky 1972, Eq. 2.21.

In [ ]:

reload(hz)
fig, ax = plt.subplots()

def label_line(ax, x, f, text, color, offset=3):
    # Write text just above the curve r = f(a) at a = x, rotated to follow it on screen
    y, h = f(x), 1e-3
    data_angle = np.degrees(np.arctan2(f(x + h) - f(x - h), 2*h))
    screen_angle = ax.transData.transform_angles(np.array([data_angle]), np.array([[x, y]]))[0]
    ax.annotate(text, (x, y), xytext=(0, offset), textcoords='offset points', color=color,
                rotation=screen_angle, rotation_mode='anchor', ha='center', va='bottom')

surfaces = [(hz.get_photon_orbit_location, "Photon Orbit", -0.55),
            (hz.get_marginally_bound_location, "Marginally Bound Orbit", -0.52),
            (hz.get_ISCO_location, "ISCO", -0.6),
            (hz.get_ergosphere_location, "Ergosphere", -0.6)]
for f, name, x_label in surfaces:
    line, = ax.plot(spin_range, f(spin_range))
    ax.set_xlim(spin_range[0], spin_range[-1])
    ax.set_ylim(1, 9.0)
    label_line(ax, x_label, f, name, line.get_color())

ax.set_xlabel(r"Black Hole Spin $a$")
ax.set_ylabel(r"Radial Coordinate $r/r_g$")

ax.set_axisbelow(False)
grid_kw = dict(color=plt.rcParams['grid.color'], ls=plt.rcParams['grid.linestyle'],
               lw=plt.rcParams['grid.linewidth'], zorder=0.5)
for x in ax.get_xticks():
    ax.axvline(x, **grid_kw)
for y in ax.get_yticks():
    ax.axhline(y, **grid_kw)
ax.axvline(0., color='k', ls='--', lw=1, zorder=1)
ax.fill_between(spin_range, 1., hz.get_EH_location(spin_range), color='0.8', zorder=1.5)
ax.plot(spin_range, hz.get_EH_location(spin_range), color='k', zorder=2)
ax.text(0., 1.45, "Inside Event Horizon", ha='center', va='center', zorder=2)
plotting.save_fig("black_hole_surfaces", "01_black_hole_surfaces", dpi=300)

Out[ ]:

[PosixPath('figures/01_black_hole_surfaces/black_hole_surfaces.pdf'),
 PosixPath('figures/01_black_hole_surfaces/black_hole_surfaces.png')]

In [ ]: