This Jupyter/SageMath worksheet is relative to the lectures Geometry and physics of black holes.
Click here to download the worksheet file (ipynb format). To run it, you must start SageMath with the Jupyter notebook, with the command sage -n jupyter
version()
'SageMath version 8.1.rc3, Release Date: 2017-11-23'
%display latex
var('eta r0')
assume(eta>0, eta<pi, r0>2)
The function $\frac{\mathrm{d}t}{\mathrm{d}\eta}$:
dtdeta = r0/2*sqrt(r0/2-1)*(1+cos(eta))^2/(1+cos(eta)-4/r0)
dtdeta
The primitive:
dtdeta.integrate(eta).simplify_full()
Since this is a complicated expression, let us try to find something simpler by taking $u=\eta/2$ as a variable; we have then $\frac{\mathrm{d}t}{\mathrm{d}u}= 2 \frac{\mathrm{d}t}{\mathrm{d}\eta}$:
var('u')
assume(u>0, u<pi/2)
dtdu = 2*dtdeta.subs({eta: 2*u}).simplify_full()
dtdu
This time, integrate
yields a simpler expression for $t$:
t = dtdu.integrate(u).simplify_full()
t
We force some obvious simplification by
t = (t*sqrt(2*r0-4)/2/sqrt(r0/2-1)).simplify_full()
t
We note that, given the range of $u$, the argument of the logarithm is negative; let us enforce $\ln|x|$ as the primitve of $1/x$ instead of $\ln x$ or $\ln(-x)$:
x = t.operands()[1].operands()[0].operands()[0]
x
t = t.subs({x: abs(x)})
t
We move back to $\eta$, with some extra simplification enforced:
t = t.subs({sqrt(2*r0-4): 2*sqrt(r0/2-1), u: eta/2})
t
Finally we invoke trig_reduce()
to let appear $\sin\eta$:
t = t.trig_reduce().simplify_full()
t
tau = sqrt(r0^3/8)*(eta + sin(eta))
tau
r = r0/2*(1 + cos(eta))
r
r0_values = [2.1, 3, 4, 5, 6]
graph = Graphics()
for r0_val in r0_values:
tau1 = tau.subs({r0: r0_val})
r1 = r.subs({r0: r0_val})
graph += parametric_plot((tau1, r1), (eta, 0, pi),
color=hue(r0_val/6), thickness=2,
legend_label="$r_0 = {:.2f}\\; m$".format(float(r0_val)),
axes_labels=[r'$\tau/m$', r'$r/m$'])
show(graph, aspect_ratio=1, gridlines=['minor', 'major'])
graph.save("ges_radial_infall_tau.pdf", aspect_ratio=1,
gridlines=['minor', 'major'])
graph = Graphics()
for r0_val in r0_values:
t1 = t.subs({r0: r0_val})
r1 = r.subs({r0: r0_val})
graph += parametric_plot((r1, t1), (eta, 0, pi), plot_points=300,
color=hue(r0_val/6), thickness=2,
legend_label="$r_0 = {:.2f}\\; m$".format(float(r0_val)),
axes_labels=[r'$r/m$', r'$t/m$'])
graph += polygon([(0,0), (2,0), (2,26), (0,26)], color='lightgrey', alpha=0.7)
graph += line([(2,0), (2, 26)], thickness=1.5, color='black')
show(graph, xmin=0, xmax=8, ymin=0, ymax=25,
aspect_ratio=1, legend_loc=(0.8,0.8))
graph.save("ges_radial_infall_t.pdf", xmin=0, xmax=8, ymin=0, ymax=25,
aspect_ratio=1, legend_loc=(0.8,0.8))
ti = 2*(sqrt(r0/2-1)*(eta + r0/4*(eta + sin(eta))) +
2*ln(cos(eta/2) + 1/sqrt(r0/2-1)*sin(eta/2)))
ti
graph = Graphics()
for r0_val in r0_values:
ti1 = ti.subs({r0: r0_val})
r1 = r.subs({r0: r0_val})
graph += parametric_plot((r1, ti1), (eta, 0, pi),
color=hue(r0_val/6), thickness=2,
legend_label="$r_0 = {:.2f}\\; m$".format(float(r0_val)),
axes_labels=[r'$r/m$', r'$\tilde{t}/m$'])
graph += polygon([(0,0), (2,0), (2,26), (0,26)], color='lightgrey', alpha=0.7)
graph += line([(2,0), (2, 26)], thickness=1.5, color='black')
show(graph, xmin=0, xmax=8, ymin=0, ymax=25,
aspect_ratio=1, legend_loc=(0.8,0.8))
graph.save("ges_radial_infall_IEF.pdf", xmin=0, xmax=8, ymin=0, ymax=25,
aspect_ratio=1, legend_loc=(0.8,0.8))
tif = 2*(pi*sqrt(r0/2-1)*(r0/4+1) - ln(r0/2-1))
tif
tih = 2*(2*(1+r0/4)*sqrt(r0/2-1)*atan(sqrt(r0/2-1)) + r0/2-1 -ln(r0/2))
tih
graph = plot(tif, (r0, 2.001, 6), thickness=2,
legend_label=r'$\tilde{t}_{\rm f}$',
axes_labels=[r'$r_0/m$', r'$\tilde{t}/m$'], gridlines=True)
graph += plot(tih, (r0, 2.001, 6), linestyle='--', thickness=2,
legend_label=r'$\tilde{t}_{\rm h}$', color='purple')
graph
tauf = pi*sqrt(r0^3/8)
tauf
tauh = sqrt(r0^3/2)*(atan(sqrt(r0/2-1)) + sqrt(2/r0*(1-2/r0)))
tauh
graph = plot(tauf, (r0, 2.001, 6), thickness=2,
legend_label=r'$\tau_{\rm f}$',
axes_labels=[r'$r_0/m$', r'$\tau/m$'], gridlines='minor')
graph += plot(tauh, (r0, 2.001, 6), linestyle='--', thickness=2,
legend_label=r'$\tau_{\rm h}$', color='purple')
graph
graph.save("ges_infall_time.pdf")
taui = (tauf - tauh).simplify_full()
taui
lim(taui, r0 = 2)
lim(taui, r0 = +oo)
We may check the above limit by asking for a Taylor expansion in terms of $u=1/r_0$:
var('u')
s = taui.subs({r0: 1/u}).simplify_full()
s
s.taylor(u, 0, 2)
graph = plot(tauf - tauh, (r0, 2.001, 20), thickness=2,
axes_labels=[r'$r_0/m$', r'$\Delta\tau_{\rm in}/m$'],
ymin=1.3, ymax = 3.2, gridlines='minor')
graph
graph.save("ges_time_inside.pdf")