Kerr_null_geod_zero_ener.ipynb

Kernel: SageMath 9.2.rc3

Zero-energy null geodesics in Kerr spacetime

This Jupyter/SageMath notebook is relative to the lectures Geometry and physics of black holes.

version()

'SageMath version 9.2.rc3, Release Date: 2020-10-18'

%display latex

a = 0.9

r_min(L) = 1 - sqrt(1 - a^2/(1 + L^2))
r_min

$$\renewcommand{\Bold}[1]{\mathbf{#1}}L \ {\mapsto}\ -\sqrt{-\frac{0.810000000000000}{L^{2} + 1} + 1} + 1$$

r_max(L) = 1 + sqrt(1 - a^2/(1 + L^2))
r_max

$$\renewcommand{\Bold}[1]{\mathbf{#1}}L \ {\mapsto}\ \sqrt{-\frac{0.810000000000000}{L^{2} + 1} + 1} + 1$$

th_min(L) = atan(-L)
th_max(L) = pi - th_min(L)

$r$ as a function of $\theta$

def rf(th, K, L, s):
    b = sqrt(1 - a^2/(1 + L^2))
    return 1 + b*sin(K - s*arcsin(sqrt(1 + L^2)*cos(th)))

def Kf(rm, L, s):
    return n(arcsin((rm - 1)/sqrt(1 - a^2/(1 + L^2))) + s*pi/2)

L1 = -1
r_min(L1), r_max(L1)

$$\renewcommand{\Bold}[1]{\mathbf{#1}}\left(0.228637568972924, 1.77136243102708\right)$$

n(th_min(L1)), n(th_max(L1))

$$\renewcommand{\Bold}[1]{\mathbf{#1}}\left(0.785398163397448, 2.35619449019234\right)$$

Plot in O'Neill coordinates

$$x = e^{r/m}\sin\theta$$$$y = e^{r/m}\cos\theta$$

def plot_geod(L1, rm):
    K1 = Kf(rm, L1, 1)
    K2 = Kf(rm, L1, -1)
    xx = lambda th: exp(rf(th, K1, L1, 1))*sin(th)
    yy = lambda th: exp(rf(th, K1, L1, 1))*cos(th)
    g = parametric_plot((xx, yy), (th_min(L1), th_max(L1)), 
                        thickness=2., color='green')
    xx = lambda th: exp(rf(th, K2, L1, -1))*sin(th)
    yy = lambda th: exp(rf(th, K2, L1, -1))*cos(th)
    g += parametric_plot((xx, yy), (th_min(L1), th_max(L1)), 
                         thickness=2., color='green')
    return g

def plot_arrow(L1, rm, th, dth, s):
    K1 = Kf(rm, L1, s)
    x1 = n(exp(rf(th, K1, L1, s))*sin(th))
    y1 = n(exp(rf(th, K1, L1, s))*cos(th))
    th2 = th + dth
    x2 = n(exp(rf(th2, K1, L1, s))*sin(th2))
    y2 = n(exp(rf(th2, K1, L1, s))*cos(th2))
    v = vector([x2 - x1, y2 - y1])
    g = v.plot(start=(x1, y1), color='green')
    return g

g = plot_geod(-1, 1.5)

rp = 1 + sqrt(1 - a^2)
rm = 1 - sqrt(1 - a^2)
rp, rm

$$\renewcommand{\Bold}[1]{\mathbf{#1}}\left(1.43588989435407, 0.564110105645933\right)$$

theta = var('theta')
Rp = exp(rp)
Rm = exp(rm)
Hp = parametric_plot((Rp*sin(theta), Rp*cos(theta)), (theta, 0, pi), color='black', thickness=3) + \
     text(r'$r\!=\!r_+$', (2.5, 4.), color='black', fontsize=16)
Hm = parametric_plot((Rm*sin(theta), Rm*cos(theta)), (theta, 0, pi), color='peru', thickness=3) + \
     text(r'$r\!=\!r_-$', (1.2, 2.), color='peru', fontsize=16)
R0 = parametric_plot((sin(theta), cos(theta)), (theta, 0, pi), color='darkorange', 
                     linestyle='dotted', thickness=3) + \
      text(r'$r\!=\!0$', (0,-1.3), color='darkorange', fontsize=16)
sing = circle((1,0), 0.1, color='red', fill=True) 
rminf = circle((0,0), 0.1, edgecolor='black', facecolor='white', fill=True)
g += Hp + Hm + R0 + sing + rminf 
g

Arrows:

g += plot_arrow(-1, 1.5, 1.4, -0.01, 1)
g += plot_arrow(-1, 1.5, 1.1, 0.01, -1)

Point $A$:

L1 = -1
rm = 1.5
K1 = Kf(rm, L1, 1)
xA = exp(rf(pi/2, K1, L1, 1))
g += circle((xA,0), 0.1, color='green', fill=True)
g += text(r'$A$', (4.7, 0.4), color='green', fontsize=16)
g

Plot of the ergoregion

x,y = var('x y')
r = 1/2*ln(x^2 + y^2)
costh2 = y^2/(x^2+y^2)
f_ergo = r^2 - 2*r + a^2*costh2
ergo = region_plot(f_ergo < 0, (x, 0, 8), (y, -5, 5), incol='whitesmoke', 
                   bordercol='grey',
                   axes_labels=[r'$\mathrm{e}^{r/m}\sin\theta$', 
                                r'$\mathrm{e}^{r/m}\cos\theta$'])
g += ergo
g

show(g, figsize=7.5)

g.save("gik_zero_ener_merid.pdf", figsize=7.5)