This Jupyter/SageMath notebook is relative to the lectures Geometry and physics of black holes.
version()
'SageMath version 9.4.beta3, Release Date: 2021-06-21'
%display latex
r = var('r')
mu = var('mu', latex_name=r'\mu')
rho(r, mu) = sqrt(1 - mu^2)*sqrt(r^2 + 1 + 2*r*(1 - mu^2)/(r^2 + mu^2))
rho
diff(rho(r, mu), r).simplify_full()
zp(r, mu) = sqrt( (r^2 + mu^2)/(r - 1)^2 - (diff(rho(r), r))^2 ).simplify_full()
zp(r, mu)
zp2 = (2*r^7 + 4*r^5 - 4*r^4 + 2*r^3 - r^2 + 2*r - 1)/(r^3*(r^3 + r + 2)*(r - 1)^2)
zp2
bool(zp(r, 0)^2 == zp2)
s = (zp(r, mu)^2).numerator()
s.collect(r)
s = (zp(r, mu)^2).denominator()
s.factor()
#mu0 = sqrt(3)/2 # theta = pi/6
mu0 = 0 # theta = pi/2
zpf = fast_callable(zp(r, mu0), vars=[r])
def zz(r0, r1, verbose=False):
numint = numerical_integral(zpf, r0, r1, algorithm='qags')
error = numint[1]
if verbose:
print("Error = {}".format(error))
if error > 1e-3:
print("Warning: error = {}".format(error))
return numint[0]
zz(2, 3, verbose=True)
Error = 1.6403995233359683e-14
zz(1.5, 3, verbose=True)
Error = 5.314212959789082e-11
zz(1.001, 3, verbose=True)
Error = 2.628708897735641e-09
zz(1.00001, 3, verbose=True)
Error = 8.104079286898108e-07
r_0 = 2
rmin = 1.0001
plot(lambda r: zz(r_0, r), (r, rmin, 4),
axes_labels = [r'$r/m$', r'$Z(r)/m$'])
rhof = fast_callable(rho(r, mu0), vars=[r])
plot(rhof, (1, 4),
axes_labels = [r'$r/m$', r'$P(r)/m$'])
rmin = 1.00001
r_0 = 2
def xe(r, ph):
return rhof(r)*cos(ph)
def ye(r, ph):
return rhof(r)*sin(ph)
def ze(r, ph):
return zz(r_0, r)
rmax = 10
graph = parametric_plot3d([xe, ye, ze], (rmin, rmax), (0, 2*pi),
color='lightsteelblue')
graph
Constant $r$ curve:
def const_r(r0, color='blue', thickness=5):
return parametric_plot3d([lambda ph: xe(r0, ph), lambda ph: ye(r0, ph),
lambda ph: ze(r0, ph)], (0, 2*pi),
color=color, thickness=thickness)
Adding the ergosurface, i.e. the curve $r=m(1 + \sin\theta)$:
graph += const_r(1 + sqrt(1 - mu0^2), thickness=7)
Adding the curves $r/m = 1 + 10^{-5},\ldots, 1 + 10^{-1}$:
for i in range(1, 6):
graph += const_r(1 + 10^(-i), color='black')
show(graph, axes_labels_style=dict(fontsize=12, fontfamily='Liberation Sans'))