Pruefung4.mw

Allgemeine Relativitätstheorie mit dem Computer 

General Theory of Relativity on the Computer 

Vorlesung gehalten an der J.W.Goethe-Universität in Frankfurt am Main (Sommersemester 2016)  

von Dr.phil.nat. Dr.rer.pol. Matthias Hanauske 

Frankfurt am Main 11.04.2016 

 

Erster Vorlesungsteil: Allgemeine Relativitätstheorie mit Maple  

Probe Prüfung 

Probe Prüfung 

Aufgabe 4. Maple Worksheet Vorlesung3.mw und Definition ds^2 aus Vorlesung2.mw verwendet. 

Bewegung eines Probekörpers um ein schwarzes Loch in der Ebene 

Im folgenden wird die Geodätengleichung in vorgegebener Schwarzschild Raumzeit betrachtet. Die Geodätengleichung beschreibt wie sich ein Probekörper (Masse = 0) im Raum bewegt und sagt voraus, dass diese Bewegung sich stehts entlang der kürzesten Kurve, in der durch die Metrik beschriebenen gekrümmten Raumzeit, vollzieht. 

Die Geodätengleichung und das effektive Potential V(r) 

> restart:
with( tensor ):
with(plots):
with(plottools):
 

Definition der kovarianten Raumzeit-Metrik eines schwarzen Lochs der Masse M in Schwarzschildkoordinaten: 

> coord := [t, r, theta, phi]:
g_compts := array(symmetric,sparse, 1..4, 1..4):
g_compts[1,1] := 1-2*M/r: g_compts[2,2] := -1/g_compts[1,1]:
g_compts[3,3] := -r^2:    g_compts[4,4] := -r^2*sin(theta)^2:
g := create( [-1,-1], eval(g_compts));
 

`:=`(g, table([compts = Matrix(%id = 25547560), index_char = [-1, -1]]))
`:=`(g, table([compts = Matrix(%id = 25547560), index_char = [-1, -1]]))
(2.1.1)
 

Berechnung der kontravarianten Metrik und der Christoffel Symbole: 

> ginv := invert( g, 'detg' ):
D1g := d1metric( g, coord ):  
D2g := d2metric( D1g, coord ):
Cf1 := Christoffel1( D1g ):
Cf2:= Christoffel2( ginv, Cf1 ):
 

Berechnung der Geodätengleichung als Funktion des affinen Parameters lambda. Die Geodätengleichung ist ein System gekoppelter Differentialgleichungen:   

> eqns:=geodesic_eqns( coord, lambda, Cf2 );
 

{`+`(diff(diff(t(lambda), lambda), lambda), `-`(`/`(`*`(2, `*`(M, `*`(diff(t(lambda), lambda), `*`(diff(r(lambda), lambda))))), `*`(r, `*`(`+`(`-`(r), `*`(2, `*`(M)))))))) = 0, `+`(diff(diff(phi(lambd...
{`+`(diff(diff(t(lambda), lambda), lambda), `-`(`/`(`*`(2, `*`(M, `*`(diff(t(lambda), lambda), `*`(diff(r(lambda), lambda))))), `*`(r, `*`(`+`(`-`(r), `*`(2, `*`(M)))))))) = 0, `+`(diff(diff(phi(lambd...
{`+`(diff(diff(t(lambda), lambda), lambda), `-`(`/`(`*`(2, `*`(M, `*`(diff(t(lambda), lambda), `*`(diff(r(lambda), lambda))))), `*`(r, `*`(`+`(`-`(r), `*`(2, `*`(M)))))))) = 0, `+`(diff(diff(phi(lambd...
{`+`(diff(diff(t(lambda), lambda), lambda), `-`(`/`(`*`(2, `*`(M, `*`(diff(t(lambda), lambda), `*`(diff(r(lambda), lambda))))), `*`(r, `*`(`+`(`-`(r), `*`(2, `*`(M)))))))) = 0, `+`(diff(diff(phi(lambd...
{`+`(diff(diff(t(lambda), lambda), lambda), `-`(`/`(`*`(2, `*`(M, `*`(diff(t(lambda), lambda), `*`(diff(r(lambda), lambda))))), `*`(r, `*`(`+`(`-`(r), `*`(2, `*`(M)))))))) = 0, `+`(diff(diff(phi(lambd...
(2.1.2)
 

Wir lassen nur ebene Bewegungen zu (theta=Pi/2, dtheta=0) und setzen M=2.5:Hier wurde das Maple-File verändert 

> setM:=2.5:
eq1:=subs({r=r(lambda),sin(theta)=1,cos(theta)=0,diff(theta(lambda),lambda)=0,theta=Pi/2,M=setM},eqns[1]):
eq2:=subs({r=r(lambda),sin(theta)=1,cos(theta)=0,diff(theta(lambda),lambda)=0,theta=Pi/2,M=setM},eqns[2]):
eq3:=subs({r=r(lambda),sin(theta)=1,cos(theta)=0,diff(theta(lambda),lambda)=0,theta=Pi/2,M=setM},eqns[3]):
eq4:=subs({r=r(lambda),sin(theta)=1,cos(theta)=0,diff(theta(lambda),lambda)=0,theta=Pi/2,M=setM},eqns[4]):
eq1:=subs({r(lambda)(lambda)=r(lambda)},eq1);
eq2:=subs({r(lambda)(lambda)=r(lambda)},eq2);
eq3:=subs({r(lambda)(lambda)=r(lambda)},eq3);
eq4:=subs({r(lambda)(lambda)=r(lambda)},eq4);
 

 

 

 

`+`(diff(diff(t(lambda), lambda), lambda), `-`(`/`(`*`(5.0, `*`(diff(t(lambda), lambda), `*`(diff(r(lambda), lambda)))), `*`(r(lambda), `*`(`+`(`-`(r(lambda)), 5.0)))))) = 0
`+`(diff(diff(phi(lambda), lambda), lambda), `/`(`*`(2, `*`(diff(r(lambda), lambda), `*`(diff(phi(lambda), lambda)))), `*`(r(lambda)))) = 0
0 = 0
`+`(diff(diff(r(lambda), lambda), lambda), `-`(`/`(`*`(2.5, `*`(`+`(`-`(r(lambda)), 5.0), `*`(`^`(diff(t(lambda), lambda), 2)))), `*`(`^`(r(lambda), 3)))), `/`(`*`(2.5, `*`(`^`(diff(r(lambda), lambda)...
`+`(diff(diff(r(lambda), lambda), lambda), `-`(`/`(`*`(2.5, `*`(`+`(`-`(r(lambda)), 5.0), `*`(`^`(diff(t(lambda), lambda), 2)))), `*`(`^`(r(lambda), 3)))), `/`(`*`(2.5, `*`(`^`(diff(r(lambda), lambda)...
(2.1.3)
 

Hier wurde das Maple-File verändert 

> VeffFb:=(r,M,l)->-M/r+l^2/(2*r^2)-M*l^2/r^3;
 

proc (r, M, l) options operator, arrow; `+`(`-`(`/`(`*`(M), `*`(r))), `/`(`*`(`/`(1, 2), `*`(`^`(l, 2))), `*`(`^`(r, 2))), `-`(`/`(`*`(M, `*`(`^`(l, 2))), `*`(`^`(r, 3))))) end proc (2.1.4)
 

 

Definition der Anfangswerte der innersten stabilen Kreisbahn und zweier weiterer stabilen kreisförmigen Bahnbewegungen.Hier wurde das Maple-File verändert 

Bemerkung: Der Anfangswert dt0 ergibt sich hierbei aus dem infinitesimalen Weglängenelements ds²=1 eines massiven Probekörpers: 

> dx:=create([1], array([dt,dr,dtheta,dphi])):
ds2:=get_compts(prod(dx,lower(g,dx,1),[1,1])):
ds2:=collect( simplify(ds2), [dt,dr,dtheta,dphi]);
Gl1:=subs(sin(theta)=1,cos(theta)=0,dtheta=0,dr=0,ds2=1);
dt=solve(Gl1,dt);
 

 

 

`+`(`-`(`/`(`*`(`+`(`*`(`^`(r, 2)), `-`(`*`(4, `*`(r, `*`(M)))), `*`(4, `*`(`^`(M, 2)))), `*`(`^`(dt, 2))), `*`(r, `*`(`+`(`-`(r), `*`(2, `*`(M))))))), `/`(`*`(`^`(dr, 2), `*`(r)), `*`(`+`(`-`(r), `*`...
`+`(`-`(`/`(`*`(`+`(`*`(`^`(r, 2)), `-`(`*`(4, `*`(r, `*`(M)))), `*`(4, `*`(`^`(M, 2)))), `*`(`^`(dt, 2))), `*`(r, `*`(`+`(`-`(r), `*`(2, `*`(M))))))), `/`(`*`(`^`(dr, 2), `*`(r)), `*`(`+`(`-`(r), `*`...
`+`(`-`(`/`(`*`(`+`(`*`(`^`(r, 2)), `-`(`*`(4, `*`(r, `*`(M)))), `*`(4, `*`(`^`(M, 2)))), `*`(`^`(dt, 2))), `*`(r, `*`(`+`(`-`(r), `*`(2, `*`(M))))))), `-`(`/`(`*`(`+`(`-`(`*`(`^`(r, 4))), `*`(2, `*`(...
dt = (`/`(`*`(`^`(`+`(`-`(`*`(`+`(`-`(r), `*`(2, `*`(M))), `*`(r, `*`(`+`(1, `*`(`^`(dphi, 2), `*`(`^`(r, 2))))))))), `/`(1, 2))), `*`(`+`(`-`(r), `*`(2, `*`(M))))), `+`(`-`(`/`(`*`(`^`(`+`(`-`(`*`(`+... (2.2.1)
 

 

> setM:=2.5:
t0:=0:
phi0:=0:
theta0:=Pi/2:
dtheta0:=0:
dr0:=0:
setlA:=2*sqrt(3)*setM;
r0A:=subs({M=setM,l=setlA},solve(simplify(diff(VeffFb(r,M,l),r)=0),r)[1]);
dphi0A:=setlA/r0A^2:
dt0A:=evalf(-sqrt(-(-r0A+2*setM)*r0A*(1+dphi0A^2*r0A^2))/(-r0A+2*setM)):
setEA:=dt0A*evalf((1-2*setM/r0A));
print("------------------------------------------");
setlB:=1.02*setlA;
r0B:=evalf(subs({M=setM,l=setlB},solve(simplify(diff(VeffFb(r,M,l),r)=0),r)[1]));
dphi0B:=setlB/r0B^2:
dt0B:=evalf(-sqrt(-(-r0B+2*setM)*r0B*(1+dphi0B^2*r0B^2))/(-r0B+2*setM)):
setEB:=dt0B*evalf((1-2*setM/r0B));
print("------------------------------------------");
setlC:=1.10*setlA;
r0C:=evalf(subs({M=setM,l=setlC},solve(simplify(diff(VeffFb(r,M,l),r)=0),r)[1]));
dphi0C:=setlC/r0C^2:
dt0C:=evalf(-sqrt(-(-r0C+2*setM)*r0C*(1+dphi0C^2*r0C^2))/(-r0C+2*setM)):
setEC:=dt0C*evalf((1-2/r0C));
 

 

 

 

 

 

 

 

 

 

 

`+`(`*`(5.0, `*`(`^`(3, `/`(1, 2)))))
15.00000000
.9428090412
------------------------------------------
`+`(`*`(5.100, `*`(`^`(3, `/`(1, 2)))))
18.68126195
.9466242861
------------------------------------------
`+`(`*`(5.500, `*`(`^`(3, `/`(1, 2)))))
25.71124990
1.095779007 (2.2.2)
 

Struktur des effektiven Potentials der innersten stabilen Kreisbahn (blauer Kurve) und zweier weiterer stabilen kreisförmigen Bahnbewegungen. 

> PotA:=plot(VeffFb(r,setM,setlA),r=2.85..52,view=[2.85..50,-0.08..-0.02],color=blue):
PotB:=plot(VeffFb(r,setM,setlB),r=2.85..52,view=[2.85..50,-0.08..-0.02],color=red):
PotC:=plot(VeffFb(r,setM,setlC),r=2.85..52,view=[2.85..50,-0.08..-0.02],color=green):
PlotEA:=pointplot({[r0A, VeffFb(r0A,setM,setlA)]}, symbol=solidcircle,symbolsize=20,color=blue):
PlotEB:=pointplot({[r0B, VeffFb(r0B,setM,setlB)]}, symbol=solidcircle,symbolsize=20,color=red):
PlotEC:=pointplot({[r0C, VeffFb(r0C,setM,setlC)]}, symbol=solidcircle,symbolsize=20,color=green):

display(PotA,PotB,PotC,PlotEA,PlotEB,PlotEC,title="Effektives Potential V(r) (Def. nach Ref. 1-3) für drei verschiedene Drehimpulse");
 

Plot_2d
 

Numerische Lösung der drei Bahnbewegungen: 

> LoesA:=dsolve({eq1,eq2,eq4,t(0)=t0,r(0)=r0A,D(r)(0)=dr0,D(t)(0)=dt0A,phi(0)=phi0,D(phi)(0)=dphi0A},{r(lambda),t(lambda),phi(lambda)},type=numeric,output=listprocedure):
LoesB:=dsolve({eq1,eq2,eq4,t(0)=t0,r(0)=r0B,D(r)(0)=dr0,D(t)(0)=dt0B,phi(0)=phi0,D(phi)(0)=dphi0B},{r(lambda),t(lambda),phi(lambda)},type=numeric,output=listprocedure):
LoesC:=dsolve({eq1,eq2,eq4,t(0)=t0,r(0)=r0C,D(r)(0)=dr0,D(t)(0)=dt0C,phi(0)=phi0,D(phi)(0)=dphi0C},{r(lambda),t(lambda),phi(lambda)},type=numeric,output=listprocedure):
 

Animation der Bewegungen: Der ISCO ist die blaue Kurve 

> frames:=80:
lend:=200:
BH:=display(disk([0,0],2,color=black,transparency=0.1)):
for i from 0 by 1 to frames do
KoerperA[i]:=display(disk([rhs(LoesA(i*lend/frames)[4])*cos(rhs(LoesA(i*lend/frames)[2])),rhs(LoesA(i*lend/frames)[4])*sin(rhs(LoesA(i*lend/frames)[2]))],0.6,color=blue)):
KoerperB[i]:=display(disk([rhs(LoesB(i*lend/frames)[4])*cos(rhs(LoesB(i*lend/frames)[2])),rhs(LoesB(i*lend/frames)[4])*sin(rhs(LoesB(i*lend/frames)[2]))],0.6,color=red)):
KoerperC[i]:=display(disk([rhs(LoesC(i*lend/frames)[4])*cos(rhs(LoesC(i*lend/frames)[2])),rhs(LoesC(i*lend/frames)[4])*sin(rhs(LoesC(i*lend/frames)[2]))],0.6,color=green)):
TrajA[i]:=listplot([seq([rhs(LoesA(j*lend/frames)[4])*cos(rhs(LoesA(j*lend/frames)[2])),rhs(LoesA(j*lend/frames)[4])*sin(rhs(LoesA(j*lend/frames)[2]))], j = 0 .. i)],color=blue,linestyle=dot):
TrajB[i]:=listplot([seq([rhs(LoesB(j*lend/frames)[4])*cos(rhs(LoesB(j*lend/frames)[2])),rhs(LoesB(j*lend/frames)[4])*sin(rhs(LoesB(j*lend/frames)[2]))], j = 0 .. i)],color=red,linestyle=dot):
TrajC[i]:=listplot([seq([rhs(LoesC(j*lend/frames)[4])*cos(rhs(LoesC(j*lend/frames)[2])),rhs(LoesC(j*lend/frames)[4])*sin(rhs(LoesC(j*lend/frames)[2]))], j = 0 .. i)],color=green,linestyle=dot):
Ani[i]:=display({KoerperA[i],KoerperB[i],KoerperC[i],TrajA[i],TrajB[i],TrajC[i],BH});
od:
 

> display([seq(Ani[i],i=0..frames)],insequence=true,scaling=constrained,title="Kreisförmige Bahnbewegungen und ISCO");
 

Plot_2d
 

 

>