Locally nonrotating frames and NEC.ipynb

Kernel: SageMath 8.1

Locally nonrotating frames and NEC

This Jupyter/SageMath notebook is related to the article Lamy et al, arXiv:1802.01635.

version()

'SageMath version 8.1, Release Date: 2017-12-07'

%display latex

Metric

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

4-dimensional differentiable manifold M

M1 = M.open_subset('M_1')
XBL1.<t,r,th,ph> = M1.chart(r't r:(0,+oo) th:(0,pi):\theta ph:(0,2*pi):\phi')
XBL1

M2 = M.open_subset('M_2')
forget(r>0)
XBL2.<t,r,th,ph> = M2.chart(r't r:(-oo,0) th:(0,pi):\theta ph:(0,2*pi):\phi')
M._top_charts.append(XBL2)
XBL2

forget(r<0)
assumptions()

g = M.lorentzian_metric('g')

a, b= var('a b')
Sigma = r^2 + a^2*cos(th)^2

m=r^3/(r^3+2*b^2)
m1 = m
Delta1 = r^2 - 2*m*r + a^2
g1 = g.restrict(M1)
g1[0,0] = -(1 - 2*r*m/Sigma)
g1[0,3] = -2*a*r*sin(th)^2*m/Sigma
g1[1,1] = Sigma/Delta1
g1[2,2] = Sigma
g1[3,3] = (r^2 + a^2 + 2*r*(a*sin(th))^2*m/Sigma)*sin(th)^2
g.display(XBL1.frame())

g.display_comp(XBL1.frame())

m = -r^3/(-r^3 + 2*b^2)
m2 = m
Delta2 = r^2 - 2*m*r + a^2
g2 = g.restrict(M2)
g2[0,0] = -(1 - 2*r*m/Sigma)
g2[0,3] = -2*a*r*sin(th)^2*m/Sigma
g2[1,1] = Sigma/Delta2
g2[2,2] = Sigma
g2[3,3] = (r^2 + a^2 + 2*r*(a*sin(th))^2*m/Sigma)*sin(th)^2
g.display(XBL2.frame())

g.display_comp(XBL2.frame())

Locally nonrotating frames

A1=(r^2+a^2)^2-a^2*Delta1*sin(th)^2
A2=(r^2+a^2)^2-a^2*Delta2*sin(th)^2

b1=XBL1.frame()
b1

change_frame = M1.automorphism_field()
change_frame[:] = [[sqrt(A1/(Sigma*Delta1)),0,0,0],
                   [0,sqrt(Delta1/Sigma),0,0],[0,0,1/sqrt(Sigma),0],
                   [2*m1*a*r/sqrt(A1*Sigma*Delta1),0,0,sqrt(Sigma/A1)/sin(th)]]; #old covariant basis in terms of the new covariant basis
e1 = b1.new_frame(change_frame, 'e1') ; e1

b2=XBL2.frame();
b2

change_frame = M2.automorphism_field()
change_frame[:] = [[sqrt(A2/(Sigma*Delta2)),0,0,0],
                   [0,sqrt(Delta2/Sigma),0,0],[0,0,1/sqrt(Sigma),0],
                   [2*m2*a*r/sqrt(A2*Sigma*Delta2),0,0,sqrt(Sigma/A2)/sin(th)]]; #old covariant basis in terms of the new covariant basis
e2 = b2.new_frame(change_frame, 'e2') ; e2

g.display_comp(e1)

g.display_comp(e2)

Energy momentum tensor

Ricci = g.ricci()
Ricci_scal = g.ricci_scalar()
G = Ricci-1/2*g*Ricci_scal; print G

Field of symmetric bilinear forms Ric(g)-unnamed metric on the 4-dimensional differentiable manifold M

G1=G.display(e1)

G2=G.display(e2)

G.display()

Null Energy Condition

k1 = M1.vector_field('k1')
k1[e1,:] = [1,1,0,0]
k1.display(e1)

k1d=k1.down(g)
k1k1=k1d['_a']*k1['^a']
k1k1.display()

k2 = M2.vector_field('k2')
k2[e2,:] = [1,1,0,0]
k2.display(e2)

k2d=k2.down(g)
k2k2=k2d['_a']*k2['^a']
k2k2.display()

NN=G['_{ab}']*k1['^a']
NEC1=8*pi*NN['_a']*k1['^a']
NEC1.display()

NNN=G['_{ab}']*k2['^a']
NEC2=8*pi*NNN['_a']*k2['^a']
NEC2.display()

nec1(a,b,r,th)=NEC1.expr()
nec2(a,b,r,th)=NEC2.expr()

a=0.9
b=1
plot1=plot((nec1(a,b,r,pi/6), nec1(a,b,r,pi/3), nec1(a,b,r,pi/2)), (r,0,4), 
           ymin=-70,ymax=70,
           legend_label=[r'$\theta=\frac{\pi}{6}$',r'$\theta=\frac{\pi}{3}$',
                         r'$\theta=\frac{\pi}{2}$'],
           axes_labels=[r'$r/m$',''], plot_points=400)

a=0.9
b=1
plot2=plot((nec2(a,b,r,pi/6), nec2(a,b,r,pi/3), nec2(a,b,r,pi/2)),(r,-4,0),
           axes_labels=[r'$r/m$',r'$T_{\hat{\mu}\hat{\nu}}k^{\hat{\mu}}k^{\hat{\nu}}$'],
           ymin=-70, ymax=70, plot_points=400)

graph = plot1 + plot2
show(graph)

graph.save('NEC_a=09_b=1.pdf')

ww1(a,b,r,th)=-12*((a^6*b^2*r^8 + a^4*b^2*r^10 - 2*a^6*b^4*r^5 - 2*a^4*b^4*r^7 - 2*a^4*b^2*r^9 - 8*a^6*b^6*r^2 - 8*a^4*b^6*r^4 + 8*a^4*b^4*r^6)*sin(th)^4 - (a^6*b^2*r^8 + 4*a^4*b^2*r^10 + 3*a^2*b^2*r^12 - 2*a^6*b^4*r^5 + 4*a^4*b^4*r^7 - 6*a^2*b^2*r^11 - 8*a^6*b^6*r^2 - 8*a^4*b^6*r^4 + 8*a^4*b^4*r^6 - 2*(a^4*b^2 - 3*a^2*b^4)*r^9)*sin(th)^2)/(a^2*r^18 + r^20 + 8*a^2*b^2*r^15 + 24*a^2*b^4*r^12 + 2*(a^2 + 4*b^2)*r^17 + 32*a^2*b^6*r^9 + 16*a^2*b^8*r^6 + 12*(a^2*b^2 + 2*b^4)*r^14 + 8*(3*a^2*b^4 + 4*b^6)*r^11 + 16*(a^2*b^6 + b^8)*r^8 + (a^8*r^12 + a^6*r^14 + 8*a^8*b^2*r^9 + 8*a^6*b^2*r^11 - 2*a^6*r^13 + 24*a^8*b^4*r^6 + 24*a^6*b^4*r^8 - 12*a^6*b^2*r^10 + 32*a^8*b^6*r^3 + 32*a^6*b^6*r^5 - 24*a^6*b^4*r^7 + 16*a^8*b^8 + 16*a^6*b^8*r^2 - 16*a^6*b^6*r^4)*cos(th)^6 + (3*a^6*r^14 + 3*a^4*r^16 + 24*a^6*b^2*r^11 - 4*a^4*r^15 + 72*a^6*b^4*r^8 - 24*a^4*b^2*r^12 + 96*a^6*b^6*r^5 - 48*a^4*b^4*r^9 + 48*a^6*b^8*r^2 - 32*a^4*b^6*r^6 + 2*(a^6 + 12*a^4*b^2)*r^13 + 12*(a^6*b^2 + 6*a^4*b^4)*r^10 + 24*(a^6*b^4 + 4*a^4*b^6)*r^7 + 16*(a^6*b^6 + 3*a^4*b^8)*r^4)*cos(th)^4 + (3*a^4*r^16 + 3*a^2*r^18 + 24*a^4*b^2*r^13 - 2*a^2*r^17 + 72*a^4*b^4*r^10 - 12*a^2*b^2*r^14 + 96*a^4*b^6*r^7 - 24*a^2*b^4*r^11 + 48*a^4*b^8*r^4 - 16*a^2*b^6*r^8 + 4*(a^4 + 6*a^2*b^2)*r^15 + 24*(a^4*b^2 + 3*a^2*b^4)*r^12 + 48*(a^4*b^4 + 2*a^2*b^6)*r^9 + 16*(2*a^4*b^6 + 3*a^2*b^8)*r^6)*cos(th)^2)

ww2(a,b,r,th)=12*((a^6*b^2*r^8 + a^4*b^2*r^10 + 2*a^6*b^4*r^5 + 2*a^4*b^4*r^7 - 2*a^4*b^2*r^9 - 8*a^6*b^6*r^2 - 8*a^4*b^6*r^4 - 8*a^4*b^4*r^6)*sin(th)^4 - (a^6*b^2*r^8 + 4*a^4*b^2*r^10 + 3*a^2*b^2*r^12 + 2*a^6*b^4*r^5 - 4*a^4*b^4*r^7 - 6*a^2*b^2*r^11 - 8*a^6*b^6*r^2 - 8*a^4*b^6*r^4 - 8*a^4*b^4*r^6 - 2*(a^4*b^2 + 3*a^2*b^4)*r^9)*sin(th)^2)/(a^2*r^18 + r^20 - 8*a^2*b^2*r^15 + 24*a^2*b^4*r^12 + 2*(a^2 - 4*b^2)*r^17 - 32*a^2*b^6*r^9 + 16*a^2*b^8*r^6 - 12*(a^2*b^2 - 2*b^4)*r^14 + 8*(3*a^2*b^4 - 4*b^6)*r^11 - 16*(a^2*b^6 - b^8)*r^8 + (a^8*r^12 + a^6*r^14 - 8*a^8*b^2*r^9 - 8*a^6*b^2*r^11 - 2*a^6*r^13 + 24*a^8*b^4*r^6 + 24*a^6*b^4*r^8 + 12*a^6*b^2*r^10 - 32*a^8*b^6*r^3 - 32*a^6*b^6*r^5 - 24*a^6*b^4*r^7 + 16*a^8*b^8 + 16*a^6*b^8*r^2 + 16*a^6*b^6*r^4)*cos(th)^6 + (3*a^6*r^14 + 3*a^4*r^16 - 24*a^6*b^2*r^11 - 4*a^4*r^15 + 72*a^6*b^4*r^8 + 24*a^4*b^2*r^12 - 96*a^6*b^6*r^5 - 48*a^4*b^4*r^9 + 48*a^6*b^8*r^2 + 32*a^4*b^6*r^6 + 2*(a^6 - 12*a^4*b^2)*r^13 - 12*(a^6*b^2 - 6*a^4*b^4)*r^10 + 24*(a^6*b^4 - 4*a^4*b^6)*r^7 - 16*(a^6*b^6 - 3*a^4*b^8)*r^4)*cos(th)^4 + (3*a^4*r^16 + 3*a^2*r^18 - 24*a^4*b^2*r^13 - 2*a^2*r^17 + 72*a^4*b^4*r^10 + 12*a^2*b^2*r^14 - 96*a^4*b^6*r^7 - 24*a^2*b^4*r^11 + 48*a^4*b^8*r^4 + 16*a^2*b^6*r^8 + 4*(a^4 - 6*a^2*b^2)*r^15 - 24*(a^4*b^2 - 3*a^2*b^4)*r^12 + 48*(a^4*b^4 - 2*a^2*b^6)*r^9 - 16*(2*a^4*b^6 - 3*a^2*b^8)*r^6)*cos(th)^2)