Contact
CoCalc Logo Icon
StoreFeaturesDocsShareSupport News Sign UpSign In
| Download

SageMath notebooks associated to the Black Hole Lectures (https://luth.obspm.fr/~luthier/gourgoulhon/bh16)

Project: BHLectures
Views: 20094
Kernel: SageMath 8.1.rc3

Radial free fall in Schwarzschild spacetime

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

Schwarzschild-Droste coordinate tt as a function of η\eta

var('eta r0') assume(eta>0, eta<pi, r0>2)

The function dtdη\frac{\mathrm{d}t}{\mathrm{d}\eta}:

dtdeta = r0/2*sqrt(r0/2-1)*(1+cos(eta))^2/(1+cos(eta)-4/r0) dtdeta
r012r01(cos(η)+1)22(4r0cos(η)1)\renewcommand{\Bold}[1]{\mathbf{#1}}-\frac{r_{0} \sqrt{\frac{1}{2} \, r_{0} - 1} {\left(\cos\left(\eta\right) + 1\right)}^{2}}{2 \, {\left(\frac{4}{r_{0}} - \cos\left(\eta\right) - 1\right)}}

The primitive:

dtdeta.integrate(eta).simplify_full()
(r012r01(2arctan(sin(η)cos(η)+1)+sin(η))+812r01arctan(sin(η)cos(η)+1))2r04812r01log(2r04(cos(η)+1)2sin(η)2r04(cos(η)+1)+2sin(η))22r04\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{{\left(r_{0} \sqrt{\frac{1}{2} \, r_{0} - 1} {\left(2 \, \arctan\left(\frac{\sin\left(\eta\right)}{\cos\left(\eta\right) + 1}\right) + \sin\left(\eta\right)\right)} + 8 \, \sqrt{\frac{1}{2} \, r_{0} - 1} \arctan\left(\frac{\sin\left(\eta\right)}{\cos\left(\eta\right) + 1}\right)\right)} \sqrt{2 \, r_{0} - 4} - 8 \, \sqrt{\frac{1}{2} \, r_{0} - 1} \log\left(-\frac{\sqrt{2 \, r_{0} - 4} {\left(\cos\left(\eta\right) + 1\right)} - 2 \, \sin\left(\eta\right)}{\sqrt{2 \, r_{0} - 4} {\left(\cos\left(\eta\right) + 1\right)} + 2 \, \sin\left(\eta\right)}\right)}{2 \, \sqrt{2 \, r_{0} - 4}}

Since this is a complicated expression, let us try to find something simpler by taking u=η/2u=\eta/2 as a variable; we have then dtdu=2dtdη\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
2r0212r01cos(u)4r0cos(u)22\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{2 \, r_{0}^{2} \sqrt{\frac{1}{2} \, r_{0} - 1} \cos\left(u\right)^{4}}{r_{0} \cos\left(u\right)^{2} - 2}

This time, integrate yields a simpler expression for tt:

t = dtdu.integrate(u).simplify_full() t
2r04r012r01cos(u)sin(u)+(r012r01+412r01)2r04u412r01log(2r04cos(u)2sin(u)2r04cos(u)+2sin(u))2r04\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{\sqrt{2 \, r_{0} - 4} r_{0} \sqrt{\frac{1}{2} \, r_{0} - 1} \cos\left(u\right) \sin\left(u\right) + {\left(r_{0} \sqrt{\frac{1}{2} \, r_{0} - 1} + 4 \, \sqrt{\frac{1}{2} \, r_{0} - 1}\right)} \sqrt{2 \, r_{0} - 4} u - 4 \, \sqrt{\frac{1}{2} \, r_{0} - 1} \log\left(-\frac{\sqrt{2 \, r_{0} - 4} \cos\left(u\right) - 2 \, \sin\left(u\right)}{\sqrt{2 \, r_{0} - 4} \cos\left(u\right) + 2 \, \sin\left(u\right)}\right)}{\sqrt{2 \, r_{0} - 4}}

We force some obvious simplification by

t = (t*sqrt(2*r0-4)/2/sqrt(r0/2-1)).simplify_full() t
12(r0cos(u)sin(u)+(r0+4)u)2r042log(2r04cos(u)2sin(u)2r04cos(u)+2sin(u))\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{2} \, {\left(r_{0} \cos\left(u\right) \sin\left(u\right) + {\left(r_{0} + 4\right)} u\right)} \sqrt{2 \, r_{0} - 4} - 2 \, \log\left(-\frac{\sqrt{2 \, r_{0} - 4} \cos\left(u\right) - 2 \, \sin\left(u\right)}{\sqrt{2 \, r_{0} - 4} \cos\left(u\right) + 2 \, \sin\left(u\right)}\right)

We note that, given the range of uu, the argument of the logarithm is negative; let us enforce lnx\ln|x| as the primitve of 1/x1/x instead of lnx\ln x or ln(x)\ln(-x):

x = t.operands()[1].operands()[0].operands()[0] x
2r04cos(u)2sin(u)2r04cos(u)+2sin(u)\renewcommand{\Bold}[1]{\mathbf{#1}}-\frac{\sqrt{2 \, r_{0} - 4} \cos\left(u\right) - 2 \, \sin\left(u\right)}{\sqrt{2 \, r_{0} - 4} \cos\left(u\right) + 2 \, \sin\left(u\right)}
t = t.subs({x: abs(x)}) t
12(r0cos(u)sin(u)+(r0+4)u)2r042log(2r04cos(u)2sin(u)2r04cos(u)+2sin(u))\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{2} \, {\left(r_{0} \cos\left(u\right) \sin\left(u\right) + {\left(r_{0} + 4\right)} u\right)} \sqrt{2 \, r_{0} - 4} - 2 \, \log\left({\left| \frac{\sqrt{2 \, r_{0} - 4} \cos\left(u\right) - 2 \, \sin\left(u\right)}{\sqrt{2 \, r_{0} - 4} \cos\left(u\right) + 2 \, \sin\left(u\right)} \right|}\right)

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
12(2r0cos(12η)sin(12η)+η(r0+4))12r012log(12r01cos(12η)sin(12η)12r01cos(12η)+sin(12η))\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{2} \, {\left(2 \, r_{0} \cos\left(\frac{1}{2} \, \eta\right) \sin\left(\frac{1}{2} \, \eta\right) + \eta {\left(r_{0} + 4\right)}\right)} \sqrt{\frac{1}{2} \, r_{0} - 1} - 2 \, \log\left({\left| \frac{\sqrt{\frac{1}{2} \, r_{0} - 1} \cos\left(\frac{1}{2} \, \eta\right) - \sin\left(\frac{1}{2} \, \eta\right)}{\sqrt{\frac{1}{2} \, r_{0} - 1} \cos\left(\frac{1}{2} \, \eta\right) + \sin\left(\frac{1}{2} \, \eta\right)} \right|}\right)

Finally we invoke trig_reduce() to let appear sinη\sin\eta:

t = t.trig_reduce().simplify_full() t
12(η+sin(η))r012r01+2η12r012log(12r01cos(12η)sin(12η)12r01cos(12η)+sin(12η))\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{2} \, {\left(\eta + \sin\left(\eta\right)\right)} r_{0} \sqrt{\frac{1}{2} \, r_{0} - 1} + 2 \, \eta \sqrt{\frac{1}{2} \, r_{0} - 1} - 2 \, \log\left(\frac{{\left| \sqrt{\frac{1}{2} \, r_{0} - 1} \cos\left(\frac{1}{2} \, \eta\right) - \sin\left(\frac{1}{2} \, \eta\right) \right|}}{{\left| \sqrt{\frac{1}{2} \, r_{0} - 1} \cos\left(\frac{1}{2} \, \eta\right) + \sin\left(\frac{1}{2} \, \eta\right) \right|}}\right)

Plot of rr as a function of τ\tau

tau = sqrt(r0^3/8)*(eta + sin(eta)) tau
1212r03(η+sin(η))\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{2} \, \sqrt{\frac{1}{2}} \sqrt{r_{0}^{3}} {\left(\eta + \sin\left(\eta\right)\right)}
r = r0/2*(1 + cos(eta)) r
12r0(cos(η)+1)\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{2} \, r_{0} {\left(\cos\left(\eta\right) + 1\right)}
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'])
Image in a Jupyter notebook
graph.save("ges_radial_infall_tau.pdf", aspect_ratio=1, gridlines=['minor', 'major'])

Plot of rr as a function of tt

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))
Image in a Jupyter notebook
graph.save("ges_radial_infall_t.pdf", xmin=0, xmax=8, ymin=0, ymax=25, aspect_ratio=1, legend_loc=(0.8,0.8))

Spacetime diagram based on IEF coordinates

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
12((η+sin(η))r0+4η)12r01+4log(sin(12η)12r01+cos(12η))\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{2} \, {\left({\left(\eta + \sin\left(\eta\right)\right)} r_{0} + 4 \, \eta\right)} \sqrt{\frac{1}{2} \, r_{0} - 1} + 4 \, \log\left(\frac{\sin\left(\frac{1}{2} \, \eta\right)}{\sqrt{\frac{1}{2} \, r_{0} - 1}} + \cos\left(\frac{1}{2} \, \eta\right)\right)
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))
Image in a Jupyter notebook
graph.save("ges_radial_infall_IEF.pdf", xmin=0, xmax=8, ymin=0, ymax=25, aspect_ratio=1, legend_loc=(0.8,0.8))

Infall time

tif = 2*(pi*sqrt(r0/2-1)*(r0/4+1) - ln(r0/2-1)) tif
12π(r0+4)12r012log(12r01)\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{2} \, \pi {\left(r_{0} + 4\right)} \sqrt{\frac{1}{2} \, r_{0} - 1} - 2 \, \log\left(\frac{1}{2} \, r_{0} - 1\right)
tih = 2*(2*(1+r0/4)*sqrt(r0/2-1)*atan(sqrt(r0/2-1)) + r0/2-1 -ln(r0/2)) tih
(r0+4)12r01arctan(12r01)+r02log(12r0)2\renewcommand{\Bold}[1]{\mathbf{#1}}{\left(r_{0} + 4\right)} \sqrt{\frac{1}{2} \, r_{0} - 1} \arctan\left(\sqrt{\frac{1}{2} \, r_{0} - 1}\right) + r_{0} - 2 \, \log\left(\frac{1}{2} \, r_{0}\right) - 2
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
Image in a Jupyter notebook
tauf = pi*sqrt(r0^3/8) tauf
1212πr03\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{2} \, \sqrt{\frac{1}{2}} \pi \sqrt{r_{0}^{3}}
tauh = sqrt(r0^3/2)*(atan(sqrt(r0/2-1)) + sqrt(2/r0*(1-2/r0))) tauh
12r03(22r01r0+arctan(12r01))\renewcommand{\Bold}[1]{\mathbf{#1}}\sqrt{\frac{1}{2}} \sqrt{r_{0}^{3}} {\left(\sqrt{2} \sqrt{-\frac{\frac{2}{r_{0}} - 1}{r_{0}}} + \arctan\left(\sqrt{\frac{1}{2} \, r_{0} - 1}\right)\right)}
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
Image in a Jupyter notebook
graph.save("ges_infall_time.pdf")

Proper time spent inside the black hole

taui = (tauf - tauh).simplify_full() taui
14(2π22arctan(12r01)4r02r02)r032\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{4} \, {\left(\sqrt{2} \pi - 2 \, \sqrt{2} \arctan\left(\sqrt{\frac{1}{2} \, r_{0} - 1}\right) - 4 \, \sqrt{\frac{r_{0} - 2}{r_{0}^{2}}}\right)} r_{0}^{\frac{3}{2}}
lim(taui, r0 = 2)
π\renewcommand{\Bold}[1]{\mathbf{#1}}\pi
lim(taui, r0 = +oo)
43\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{4}{3}

We may check the above limit by asking for a Taylor expansion in terms of u=1/r0u=1/r_0:

var('u') s = taui.subs({r0: 1/u}).simplify_full() s
2π4u2u1u22arctan(1222u1u)4u32\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{\sqrt{2} \pi - 4 \, u \sqrt{-\frac{2 \, u - 1}{u}} - 2 \, \sqrt{2} \arctan\left(\frac{1}{2} \, \sqrt{2} \sqrt{-\frac{2 \, u - 1}{u}}\right)}{4 \, u^{\frac{3}{2}}}
s.taylor(u, 0, 2)
67u2+45u+43\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{6}{7} \, u^{2} + \frac{4}{5} \, u + \frac{4}{3}
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
Image in a Jupyter notebook
graph.save("ges_time_inside.pdf")