Kerr_princ_null_geod.ipynb

Kernel: SageMath 9.3.rc2

Plot of principal null geodesics in Kerr spacetime

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

The computations make use of tools developed through the SageManifolds project.

version()

'SageMath version 9.3.rc2, Release Date: 2021-04-06'

First we set up the notebook to display mathematical objects using LaTeX rendering:

%display latex

Spacetime manifold

We declare the Kerr spacetime as a 4-dimensional Lorentzian manifold $M$:

M = Manifold(4, 'M', structure='Lorentzian')
print(M)

4-dimensional Lorentzian manifold M

and introduce the (3+1 version of) Kerr coordinates $(\tilde{t},r,\theta,\tilde{\varphi})$ as a chart KC on $M$, via the method chart(). The argument of the latter is a string (delimited by r"..." because of the backslash symbols) expressing the coordinates names, their ranges (the default is $(-\infty,+\infty)$) and their LaTeX symbols:

KC.<tt,r,th,tph> = M.chart(r"tt:\tilde{t} r th:(0,pi):\theta tph:(0,2*pi):periodic:\tilde{\varphi}") 
print(KC); KC

Chart (M, (tt, r, th, tph))

$$\renewcommand{\Bold}[1]{\mathbf{#1}}\left(M,({\tilde{t}}, r, {\theta}, {\tilde{\varphi}})\right)$$

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

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

Plot of the principal null geodesics in the $(\tilde{t}, r)$ plane

tt_in(r, v) = v - r
tt_in

$$\renewcommand{\Bold}[1]{\mathbf{#1}}\left( r, v \right) \ {\mapsto} \ -r + v$$

tt_out(r, u) = u + r + 2*m/sqrt(m^2 - a^2)*(rp*ln(abs((r - rp)/(2*m))) - rm*ln(abs((r - rm)/(2*m))))
tt_out

$$\renewcommand{\Bold}[1]{\mathbf{#1}}\left( r, u \right) \ {\mapsto} \ r + u - 2.58831467741124 \, \log\left({\left| \frac{1}{2} \, r - 0.282055052822966 \right|}\right) + 6.58831467741124 \, \log\left({\left| \frac{1}{2} \, r - 0.717944947177034 \right|}\right)$$

rmin = -8; rmax = 8
graph = Graphics()
for u0 in range(-20, 20, 2):
    graph += plot(tt_out(r, u0), (r, rmin, rmax), color='green', ticks=2,
                  axes_labels=[r"$r/m$", r"$\tilde{t}/m$"])

for v0 in range(-20, 20, 2):
    graph += plot(tt_in(r, v0), (r, rmin, rmax), color='green', linestyle='--')

H = line(((rp, -8), (rp, 8)), color='black', thickness=3) + \
    text(r'$\mathscr{H}$', (rp+0.5, 4.7), color='black', fontsize=20)
Hin = line(((rm, -8), (rm, 8)), color='peru', thickness=3) + \
      text(r'$\mathscr{H}_{\rm in}$', (rm-0.6, 4.7), color='peru', fontsize=20)
graph += H + Hin
show(graph, aspect_ratio=1, ymin=-5, ymax=5, figsize=8)

Adding the vectors $k$ and $\ell$ to the plot

k = M.vector_field(1, -1, 0, 0, name='k')
k.display()

$$\renewcommand{\Bold}[1]{\mathbf{#1}}k = \frac{\partial}{\partial {\tilde{t}} }-\frac{\partial}{\partial r }$$

el = M.vector_field(1/2 + m*r/(r^2 + a^2),
                    1/2 - m*r/(r^2 + a^2),
                    0,
                    a/(r^2 + a^2),
                    name='el', latex_name=r'\ell')
el.display()

$$\renewcommand{\Bold}[1]{\mathbf{#1}}\ell = \left( \frac{r}{r^{2} + 0.810000000000000} + \frac{1}{2} \right) \frac{\partial}{\partial {\tilde{t}} } + \left( -\frac{r}{r^{2} + 0.810000000000000} + \frac{1}{2} \right) \frac{\partial}{\partial r } + \left( \frac{0.900000000000000}{r^{2} + 0.810000000000000} \right) \frac{\partial}{\partial {\tilde{\varphi}} }$$

We add the vectors $k$ and $\ell$ at the intersection of the $v=6m$ ingoing geodesic with the $u=-6m$ outgoing one:

u0, v0 = -6, 6
r0 = RDF(find_root(tt_in(r, v0) == tt_out(r, u0), 2, 6))
tt0 = tt_in(r0, v0)
tt0, r0

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

p0 = M((tt0, r0, pi/2, 0), name='p_0')
k0 = k.at(p0)
print(k0)

Tangent vector k at Point p_0 on the 4-dimensional Lorentzian manifold M

el0 = el.at(p0)
print(el0)

Tangent vector el at Point p_0 on the 4-dimensional Lorentzian manifold M

graph +=  k0.plot(chart=KC, ambient_coords=(r, tt), color='green', 
                  scale=1.5, fontsize=18, label_offset=0.3)  \
          + el0.plot(chart=KC, ambient_coords=(r, tt), color='green', 
                     parameters={m: 1}, scale=1.5, fontsize=18,
                     label_offset=0.25)
graph.save("ker_princ_null_geod.pdf", aspect_ratio=1, ymin=-5, ymax=5,
           figsize=8)
show(graph, aspect_ratio=1, ymin=-5, ymax=5, figsize=8)

Plot of $u - \tilde{t}$ and $\tilde{\tilde{\varphi}} - \tilde{\varphi}$ as functions of $r$

U(r) = u - tt_out(r, u)
U(r)

$$\renewcommand{\Bold}[1]{\mathbf{#1}}-r + 2.58831467741124 \, \log\left({\left| \frac{1}{2} \, r - 0.282055052822966 \right|}\right) - 6.58831467741124 \, \log\left({\left| \frac{1}{2} \, r - 0.717944947177034 \right|}\right)$$

graph = plot(U(r), (r, -8, 8), color='red', thickness=2, 
             axes_labels=[r"$r/m$", r"$(u - \tilde{t})/m$"], frame=True, 
             gridlines=True)
graph += line(((rp, -10), (rp, 10)), color='black', thickness=1.5, linestyle='--') + \
         line(((rm, -10), (rm, 10)), color='peru', thickness=1.5, linestyle='--')
graph.save('ker_u_r.pdf', aspect_ratio=1, ymin=-8, ymax=8)
show(graph, aspect_ratio=1, ymin=-8, ymax=8)

ttphi(r) = - a/sqrt(m^2 - a^2)*ln(abs((r - rp)/(r - rm)))
ttphi(r)

$$\renewcommand{\Bold}[1]{\mathbf{#1}}-2.06474160483506 \, \log\left({\left| \frac{r - 1.43588989435407}{r - 0.564110105645933} \right|}\right)$$

graph = plot(ttphi(r), (r, -8, 8), color='red', thickness=2, 
             axes_labels=[r"$r/m$", r"$\tilde{\tilde{\varphi}} - \tilde{\varphi}$"], 
             frame=True, gridlines=True)
graph += line(((rp, -10), (rp, 10)), color='black', thickness=1.5, linestyle='--') + \
         line(((rm, -10), (rm, 10)), color='peru', thickness=1.5, linestyle='--')
graph.save('ker_ttphi_r.pdf', aspect_ratio=1, ymin=-8, ymax=8)
show(graph, aspect_ratio=1, ymin=-8, ymax=8)