CoCalc Public FilesBHLectures / sage / Schwarzschild_horizon.ipynbOpen with one click!
Author: Eric Gourgoulhon
Views : 223
Compute Environment: Ubuntu 18.04 (Deprecated)

Schwarzschild horizon in 3+1 Eddington-Finkelstein coordinates

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

The involved computations are based on tools developed through the SageManifolds project.

NB: a version of SageMath at least equal to 8.2 is required to run this worksheet:

In [1]:
version()
'SageMath version 9.1.beta8, Release Date: 2020-03-18'

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

In [2]:
%display latex

Spacetime

We declare the spacetime manifold MM:

In [3]:
M = Manifold(4, 'M', structure='Lorentzian') print(M)
4-dimensional Lorentzian manifold M

and the 3+1 Eddington-Finkelstein coordinates (t,r,θ,ϕ)(t,r,\theta,\phi) as a chart on MM:

In [4]:
X.<t,r,th,ph> = M.chart(r't r:(0,+oo) th:(0,pi):\theta ph:(0,2*pi):\phi') X
(M,(t,r,θ,ϕ))\left(M,(t, r, {\theta}, {\phi})\right)

The mass parameter and the metric:

In [5]:
var('m', domain='real') assume(m>=0)
In [6]:
g = M.metric() g[0,0] = -(1-2*m/r) g[0,1] = 2*m/r g[1,1] = 1+2*m/r g[2,2] = r^2 g[3,3] = (r*sin(th))^2 g.display()
g=(2mr1)dtdt+2mrdtdr+2mrdrdt+(2mr+1)drdr+r2dθdθ+r2sin(θ)2dϕdϕg = \left( \frac{2 \, m}{r} - 1 \right) \mathrm{d} t\otimes \mathrm{d} t + \frac{2 \, m}{r} \mathrm{d} t\otimes \mathrm{d} r + \frac{2 \, m}{r} \mathrm{d} r\otimes \mathrm{d} t + \left( \frac{2 \, m}{r} + 1 \right) \mathrm{d} r\otimes \mathrm{d} r + r^{2} \mathrm{d} {\theta}\otimes \mathrm{d} {\theta} + r^{2} \sin\left({\theta}\right)^{2} \mathrm{d} {\phi}\otimes \mathrm{d} {\phi}
In [7]:
g[:]
(2mr12mr002mr2mr+10000r20000r2sin(θ)2)\left(\begin{array}{rrrr} \frac{2 \, m}{r} - 1 & \frac{2 \, m}{r} & 0 & 0 \\ \frac{2 \, m}{r} & \frac{2 \, m}{r} + 1 & 0 & 0 \\ 0 & 0 & r^{2} & 0 \\ 0 & 0 & 0 & r^{2} \sin\left({\theta}\right)^{2} \end{array}\right)
In [8]:
g.inverse()[:]
(2m+rr2mr002mr2mrr00001r200001r2sin(θ)2)\left(\begin{array}{rrrr} -\frac{2 \, m + r}{r} & \frac{2 \, m}{r} & 0 & 0 \\ \frac{2 \, m}{r} & -\frac{2 \, m - r}{r} & 0 & 0 \\ 0 & 0 & \frac{1}{r^{2}} & 0 \\ 0 & 0 & 0 & \frac{1}{r^{2} \sin\left({\theta}\right)^{2}} \end{array}\right)
In [9]:
g.christoffel_symbols_display()
Γtttttt=2m2r3Γttrttr=2m2+mrr3Γtrrtrr=2(m2+mr)r3Γtθθtθθ=2mΓtϕϕtϕϕ=2msin(θ)2Γrttrtt=2m2mrr3Γrtrrtr=2m2r3Γrrrrrr=2m2+mrr3Γrθθrθθ=2mrΓrϕϕrϕϕ=(2mr)sin(θ)2Γθrθθrθ=1rΓθϕϕθϕϕ=cos(θ)sin(θ)Γϕrϕϕrϕ=1rΓϕθϕϕθϕ=cos(θ)sin(θ)\begin{array}{lcl} \Gamma_{ \phantom{\, t} \, t \, t }^{ \, t \phantom{\, t} \phantom{\, t} } & = & \frac{2 \, m^{2}}{r^{3}} \\ \Gamma_{ \phantom{\, t} \, t \, r }^{ \, t \phantom{\, t} \phantom{\, r} } & = & \frac{2 \, m^{2} + m r}{r^{3}} \\ \Gamma_{ \phantom{\, t} \, r \, r }^{ \, t \phantom{\, r} \phantom{\, r} } & = & \frac{2 \, {\left(m^{2} + m r\right)}}{r^{3}} \\ \Gamma_{ \phantom{\, t} \, {\theta} \, {\theta} }^{ \, t \phantom{\, {\theta}} \phantom{\, {\theta}} } & = & -2 \, m \\ \Gamma_{ \phantom{\, t} \, {\phi} \, {\phi} }^{ \, t \phantom{\, {\phi}} \phantom{\, {\phi}} } & = & -2 \, m \sin\left({\theta}\right)^{2} \\ \Gamma_{ \phantom{\, r} \, t \, t }^{ \, r \phantom{\, t} \phantom{\, t} } & = & -\frac{2 \, m^{2} - m r}{r^{3}} \\ \Gamma_{ \phantom{\, r} \, t \, r }^{ \, r \phantom{\, t} \phantom{\, r} } & = & -\frac{2 \, m^{2}}{r^{3}} \\ \Gamma_{ \phantom{\, r} \, r \, r }^{ \, r \phantom{\, r} \phantom{\, r} } & = & -\frac{2 \, m^{2} + m r}{r^{3}} \\ \Gamma_{ \phantom{\, r} \, {\theta} \, {\theta} }^{ \, r \phantom{\, {\theta}} \phantom{\, {\theta}} } & = & 2 \, m - r \\ \Gamma_{ \phantom{\, r} \, {\phi} \, {\phi} }^{ \, r \phantom{\, {\phi}} \phantom{\, {\phi}} } & = & {\left(2 \, m - r\right)} \sin\left({\theta}\right)^{2} \\ \Gamma_{ \phantom{\, {\theta}} \, r \, {\theta} }^{ \, {\theta} \phantom{\, r} \phantom{\, {\theta}} } & = & \frac{1}{r} \\ \Gamma_{ \phantom{\, {\theta}} \, {\phi} \, {\phi} }^{ \, {\theta} \phantom{\, {\phi}} \phantom{\, {\phi}} } & = & -\cos\left({\theta}\right) \sin\left({\theta}\right) \\ \Gamma_{ \phantom{\, {\phi}} \, r \, {\phi} }^{ \, {\phi} \phantom{\, r} \phantom{\, {\phi}} } & = & \frac{1}{r} \\ \Gamma_{ \phantom{\, {\phi}} \, {\theta} \, {\phi} }^{ \, {\phi} \phantom{\, {\theta}} \phantom{\, {\phi}} } & = & \frac{\cos\left({\theta}\right)}{\sin\left({\theta}\right)} \end{array}

Let us check that we are dealing with a solution of Einstein equation in vacuum:

In [10]:
g.ricci().display()
Ric(g)=0\mathrm{Ric}\left(g\right) = 0

The scalar field uu defining the horizon

In [11]:
u = M.scalar_field(coord_expression={X: (1-r/(2*m))*exp((r-t)/(4*m))}, name='u') u.display()
u:MR(t,r,θ,ϕ)12(rm2)e(rt4m)\begin{array}{llcl} u:& M & \longrightarrow & \mathbb{R} \\ & \left(t, r, {\theta}, {\phi}\right) & \longmapsto & -\frac{1}{2} \, {\left(\frac{r}{m} - 2\right)} e^{\left(\frac{r - t}{4 \, m}\right)} \end{array}
In [12]:
plm2 = implicit_plot(u.expr().subs(m=1)+2, (r, 0, 15), (t, -10, 15), color='green') + \ text('$u=-2$', (11,12), color='green') plm1 = implicit_plot(u.expr().subs(m=1)+1, (r, 0, 15), (t, -10, 15), color='green') + \ text('$u=-1$', (6.3,12), color='green') pl0 = implicit_plot(u.expr().subs(m=1), (r, 0, 15), (t, -10, 15), color='green') + \ text('$u=0$', (3.1, 12), color='green') pl1 = implicit_plot(u.expr().subs(m=1)-1, (r, 0, 15), (t, -10, 15), color='green') + \ text('$u=1$', (0.8, 0.5), color='green') pl2 = implicit_plot(u.expr().subs(m=1)-2, (r, 0, 15), (t, -10, 15), color='green') + \ text('$u=2$', (-0.2, -2.5), color='green', background_color='white') graph = plm2+plm1+pl0+pl1+pl2 show(graph, aspect_ratio=True, axes_labels=[r'$r/m$', r'$t/m$'],)
In [13]:
graph.save('def_Schwarz_Hu.pdf', aspect_ratio=True, axes_labels=[r'$r/m$', r'$t/m$'])
In [14]:
du = u.differential() print(du) du.display()
1-form du on the 4-dimensional Lorentzian manifold M
du=(2mr)e(r4mt4m)8m2dt(2m+r)e(r4mt4m)8m2dr\mathrm{d}u = -\frac{{\left(2 \, m - r\right)} e^{\left(\frac{r}{4 \, m} - \frac{t}{4 \, m}\right)}}{8 \, m^{2}} \mathrm{d} t -\frac{{\left(2 \, m + r\right)} e^{\left(\frac{r}{4 \, m} - \frac{t}{4 \, m}\right)}}{8 \, m^{2}} \mathrm{d} r
In [15]:
grad_u = du.up(g) print(grad_u) grad_u.display()
Vector field on the 4-dimensional Lorentzian manifold M
(2m+r)e(r4mt4m)8m2t+(2mr)e(r4mt4m)8m2r-\frac{{\left(2 \, m + r\right)} e^{\left(\frac{r}{4 \, m} - \frac{t}{4 \, m}\right)}}{8 \, m^{2}} \frac{\partial}{\partial t } + \frac{{\left(2 \, m - r\right)} e^{\left(\frac{r}{4 \, m} - \frac{t}{4 \, m}\right)}}{8 \, m^{2}} \frac{\partial}{\partial r }

Let us check that each hypersurface u=constu={\rm const} is a null hypersurface:

In [16]:
g(grad_u, grad_u).expr()
00

The null normal \ell

In [17]:
rho = - log(-grad_u[[0]]) print(rho) rho.display()
Scalar field on the 4-dimensional Lorentzian manifold M
MR(t,r,θ,ϕ)12mlog(2)4mlog(2m+r)+8mlog(m)r+t4m\begin{array}{llcl} & M & \longrightarrow & \mathbb{R} \\ & \left(t, r, {\theta}, {\phi}\right) & \longmapsto & \frac{12 \, m \log\left(2\right) - 4 \, m \log\left(2 \, m + r\right) + 8 \, m \log\left(m\right) - r + t}{4 \, m} \end{array}
In [18]:
l = - exp(rho) * grad_u l.set_name('l', latex_name=r'\ell') print(l) l.display()
Vector field l on the 4-dimensional Lorentzian manifold M
=t+(2mr2m+r)r\ell = \frac{\partial}{\partial t } + \left( -\frac{2 \, m - r}{2 \, m + r} \right) \frac{\partial}{\partial r }
In [19]:
graph_l = l.plot(ambient_coords=(r,t), ranges={r:(0.01,8), t:(0,8)}, fixed_coords={th:pi/2, ph:pi}, parameters={m:1}, color='green', scale=0.8, aspect_ratio=1) show(graph_l)

Let us check that \ell is a null vector everywhere:

In [20]:
g(l,l).expr()
00
In [21]:
l_form = l.down(g) l_form.set_name('lf', latex_name=r'\underline{\ell}') print(l_form) l_form.display()
1-form lf on the 4-dimensional Lorentzian manifold M
=(2mr2m+r)dt+dr\underline{\ell} = \left( \frac{2 \, m - r}{2 \, m + r} \right) \mathrm{d} t +\mathrm{d} r
In [22]:
nab = g.connection() print(nab) nab
Levi-Civita connection nabla_g associated with the Lorentzian metric g on the 4-dimensional Lorentzian manifold M
g\nabla_{g}
In [23]:
nab_l_form = nab(l_form) nab_l_form.display()
g=(2m2mr2mr2+r3)dtdt+(4m3+4m2r3mr24m2r2+4mr3+r4)dtdr+mr2drdt+(2m2+3mr2mr2+r3)drdr+(2mrr22m+r)dθdθ+((2mrr2)sin(θ)22m+r)dϕdϕ\nabla_{g} \underline{\ell} = \left( \frac{2 \, m^{2} - m r}{2 \, m r^{2} + r^{3}} \right) \mathrm{d} t\otimes \mathrm{d} t + \left( \frac{4 \, m^{3} + 4 \, m^{2} r - 3 \, m r^{2}}{4 \, m^{2} r^{2} + 4 \, m r^{3} + r^{4}} \right) \mathrm{d} t\otimes \mathrm{d} r + \frac{m}{r^{2}} \mathrm{d} r\otimes \mathrm{d} t + \left( \frac{2 \, m^{2} + 3 \, m r}{2 \, m r^{2} + r^{3}} \right) \mathrm{d} r\otimes \mathrm{d} r + \left( -\frac{2 \, m r - r^{2}}{2 \, m + r} \right) \mathrm{d} {\theta}\otimes \mathrm{d} {\theta} + \left( -\frac{{\left(2 \, m r - r^{2}\right)} \sin\left({\theta}\right)^{2}}{2 \, m + r} \right) \mathrm{d} {\phi}\otimes \mathrm{d} {\phi}
In [24]:
nab_l_form.symmetrize().display()
(2m2mr2mr2+r3)dtdt+(4m3+4m2rmr24m2r2+4mr3+r4)dtdr+(4m3+4m2rmr24m2r2+4mr3+r4)drdt+(2m2+3mr2mr2+r3)drdr+(2mrr22m+r)dθdθ+((2mrr2)sin(θ)22m+r)dϕdϕ\left( \frac{2 \, m^{2} - m r}{2 \, m r^{2} + r^{3}} \right) \mathrm{d} t\otimes \mathrm{d} t + \left( \frac{4 \, m^{3} + 4 \, m^{2} r - m r^{2}}{4 \, m^{2} r^{2} + 4 \, m r^{3} + r^{4}} \right) \mathrm{d} t\otimes \mathrm{d} r + \left( \frac{4 \, m^{3} + 4 \, m^{2} r - m r^{2}}{4 \, m^{2} r^{2} + 4 \, m r^{3} + r^{4}} \right) \mathrm{d} r\otimes \mathrm{d} t + \left( \frac{2 \, m^{2} + 3 \, m r}{2 \, m r^{2} + r^{3}} \right) \mathrm{d} r\otimes \mathrm{d} r + \left( -\frac{2 \, m r - r^{2}}{2 \, m + r} \right) \mathrm{d} {\theta}\otimes \mathrm{d} {\theta} + \left( -\frac{{\left(2 \, m r - r^{2}\right)} \sin\left({\theta}\right)^{2}}{2 \, m + r} \right) \mathrm{d} {\phi}\otimes \mathrm{d} {\phi}

Check of the identity μαμ=0\ell^\mu \nabla_\alpha \ell_\mu=0:

In [25]:
v = l.contract(nab_l_form, 0) v.display()
00

The null normal as a pregeodesic vector field

In [26]:
nab_l = nab(l) nab_l[:]
(mr22m2+3mr2mr2+r3002m2mr2mr2+r34m3+4m2r3mr24m2r2+4mr3+r400002mr2mr+r200002mr2mr+r2)\left(\begin{array}{rrrr} \frac{m}{r^{2}} & \frac{2 \, m^{2} + 3 \, m r}{2 \, m r^{2} + r^{3}} & 0 & 0 \\ -\frac{2 \, m^{2} - m r}{2 \, m r^{2} + r^{3}} & -\frac{4 \, m^{3} + 4 \, m^{2} r - 3 \, m r^{2}}{4 \, m^{2} r^{2} + 4 \, m r^{3} + r^{4}} & 0 & 0 \\ 0 & 0 & -\frac{2 \, m - r}{2 \, m r + r^{2}} & 0 \\ 0 & 0 & 0 & -\frac{2 \, m - r}{2 \, m r + r^{2}} \end{array}\right)
In [27]:
div_l = nab_l.trace() print(div_l) div_l.display()
Scalar field on the 4-dimensional Lorentzian manifold M
MR(t,r,θ,ϕ)2(4m22mrr2)4m2r+4mr2+r3\begin{array}{llcl} & M & \longrightarrow & \mathbb{R} \\ & \left(t, r, {\theta}, {\phi}\right) & \longmapsto & -\frac{2 \, {\left(4 \, m^{2} - 2 \, m r - r^{2}\right)}}{4 \, m^{2} r + 4 \, m r^{2} + r^{3}} \end{array}
In [28]:
div_l.expr().factor()
2(4m22mrr2)(2m+r)2r-\frac{2 \, {\left(4 \, m^{2} - 2 \, m r - r^{2}\right)}}{{\left(2 \, m + r\right)}^{2} r}
In [29]:
div_l.expr().subs(r=2*m)
14m\frac{1}{4 \, m}
In [30]:
acc_l = l.contract(0,nab_l,1) print(acc_l) acc_l.display()
Vector field on the 4-dimensional Lorentzian manifold M
(4m4m2+4mr+r2)t+(4(2m2mr)8m3+12m2r+6mr2+r3)r\left( \frac{4 \, m}{4 \, m^{2} + 4 \, m r + r^{2}} \right) \frac{\partial}{\partial t } + \left( -\frac{4 \, {\left(2 \, m^{2} - m r\right)}}{8 \, m^{3} + 12 \, m^{2} r + 6 \, m r^{2} + r^{3}} \right) \frac{\partial}{\partial r }

The non-affinity parameter κ\kappa:

In [31]:
kappa = l(rho) kappa.display()
MR(t,r,θ,ϕ)4m4m2+4mr+r2\begin{array}{llcl} & M & \longrightarrow & \mathbb{R} \\ & \left(t, r, {\theta}, {\phi}\right) & \longmapsto & \frac{4 \, m}{4 \, m^{2} + 4 \, m r + r^{2}} \end{array}

Check of the pregeodesic equation =κ\nabla_{\ell} \ell = \kappa \ell:

In [32]:
acc_l == kappa * l
True\mathrm{True}
In [33]:
kappa.expr().factor()
4m(2m+r)2\frac{4 \, m}{{\left(2 \, m + r\right)}^{2}}

Value of κ\kappa on the horizon:

In [34]:
kappaH = kappa.expr().subs(r=2*m) kappaH
14m\frac{1}{4 \, m}

The complementary null vector field kk

In [35]:
k = M.vector_field(name='k') k[0] = 1/2 + m/r k[1] = -1/2 - m/r k.display()
k=(mr+12)t+(mr12)rk = \left( \frac{m}{r} + \frac{1}{2} \right) \frac{\partial}{\partial t } + \left( -\frac{m}{r} - \frac{1}{2} \right) \frac{\partial}{\partial r }
In [36]:
g(k,k).expr()
00
In [37]:
g(k,l).expr()
1-1
In [38]:
graph_k = k.plot(ambient_coords=(r,t), ranges={r:(1,8), t:(0,8)}, number_values={r:8, t:9}, fixed_coords={th:pi/2, ph:pi}, parameters={m:1}, color='red', scale=0.8, aspect_ratio=1) graph_lk = graph_l+graph_k show(graph_lk) graph_lk.save('def_plot_lk.pdf')
In [39]:
k_form = k.down(g) k_form.set_name('kf', latex_name=r'\underline{k}') k_form.display()
k=(2m+r2r)dt+(2m+r2r)dr\underline{k} = \left( -\frac{2 \, m + r}{2 \, r} \right) \mathrm{d} t + \left( -\frac{2 \, m + r}{2 \, r} \right) \mathrm{d} r

The 2-metric qq

We define q=g+k+kq = g + \underline{\ell}\otimes \underline{k} + \underline{k}\otimes \underline{\ell}:

In [40]:
q = g + l_form*k_form + k_form*l_form q.set_name('q') q.display()
q=r2dθdθ+r2sin(θ)2dϕdϕq = r^{2} \mathrm{d} {\theta}\otimes \mathrm{d} {\theta} + r^{2} \sin\left({\theta}\right)^{2} \mathrm{d} {\phi}\otimes \mathrm{d} {\phi}
In [41]:
q_up = q.up(g) print(q_up) q_up.display()
Tensor field of type (2,0) on the 4-dimensional Lorentzian manifold M
1r2θθ+1r2sin(θ)2ϕϕ\frac{1}{r^{2}} \frac{\partial}{\partial {\theta} }\otimes \frac{\partial}{\partial {\theta} } + \frac{1}{r^{2} \sin\left({\theta}\right)^{2}} \frac{\partial}{\partial {\phi} }\otimes \frac{\partial}{\partial {\phi} }

Expansion along the null normal

We compute θ()\theta_{(\ell)} as θ()=qμνμν\theta_{(\ell)} = q^{\mu\nu}\nabla_\mu \ell_\nu:

In [42]:
theta = q_up.contract(0,1,nab(l_form),0,1) theta.expr()
2(2mr)2mr+r2-\frac{2 \, {\left(2 \, m - r\right)}}{2 \, m r + r^{2}}

Check of the formula θ()=κ\theta_{(\ell)} = \nabla\cdot\ell - \kappa:

In [43]:
theta == div_l - kappa
True\mathrm{True}

Check of the forumla θ()=12Llndetq\theta_{(\ell)} = \frac{1}{2} \mathcal{L}_{\ell} \ln \det q:

In [44]:
detq = M.scalar_field({X: r^4*sin(th)^2}) theta == 1/2*ln(detq).lie_der(l)
True\mathrm{True}
In [45]:
theta == 1/2*l(ln(detq))
True\mathrm{True}

Deformation rate tensor of the cross-sections

We compute Θ\Theta as Θ=12Lq\Theta = \frac{1}{2} \mathcal{L}_{\ell} q:

In [46]:
Theta = 1/2 * q.lie_der(l) Theta.set_name('Theta', latex_name=r'\Theta') print(Theta) Theta.display()
Tensor field Theta of type (0,2) on the 4-dimensional Lorentzian manifold M
Θ=(2mrr22m+r)dθdθ+((2mrr2)sin(θ)22m+r)dϕdϕ\Theta = \left( -\frac{2 \, m r - r^{2}}{2 \, m + r} \right) \mathrm{d} {\theta}\otimes \mathrm{d} {\theta} + \left( -\frac{{\left(2 \, m r - r^{2}\right)} \sin\left({\theta}\right)^{2}}{2 \, m + r} \right) \mathrm{d} {\phi}\otimes \mathrm{d} {\phi}

Expansion of the cross-sections along the null normal kk:

We compute θ(k)\theta_{(k)} as θ(k)=qμνμkν\theta_{(k)} = q^{\mu\nu}\nabla_\mu k_\nu:

In [47]:
theta_k = q_up.contract(0,1,nab(k_form),0,1) theta_k.expr()
2m+rr2-\frac{2 \, m + r}{r^{2}}

Value of θ(k)\theta_{(k)} at the horizon:

In [48]:
theta_k.expr().subs(r=2*m)
1m-\frac{1}{m}