(* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
(*  Mathematica Syntax | EINSTEIN-MAXWELL TENSOR+GEODESIC SOLVER | geodesics.yukterez.net *)
(* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)

ClearAll["Local`*"];
smp[y_]:=FullSimplify[y,
X∈Reals&&Y∈Reals&&Z∈Reals&&β>0&&r>0];

(* kovariante metrische Komponenten *)
g11=gtt=1-2/r;
g22=gXX=-1-2(X^2)/(r^2)/(r-2);
g33=gYY=-1-2(Y^2)/(r^2)/(r-2);
g44=gZZ=-1-2(Z^2)/(r^2)/(r-2);
g12=gtX=0;
g13=gtY=0;
g14=gtZ=0;
g23=gXY=-2 X Y/(r^2)/(r-2);
g24=gXZ=-2 X Z/(r^2)/(r-2);
g34=gYZ=-2 Y Z/(r^2)/(r-2);

(* Radialkoordinate *)
r=(Sqrt[(X^2)+(Y^2)+(Z^2)]);

(* Dimensionen *)
x={t, X, Y, Z}; n=4; 

"Metrischer Tensor"
mt={
{g11, g12, g13, g14},
{g12, g22, g23, g24},
{g13, g23, g33, g34},
{g14, g24, g34, g44}
};
Subscript["g", μσ] -> MatrixForm[mt]
it=smp[Inverse[mt]];
"g"^μσ -> MatrixForm[it]
mx=ParallelTable[smp[Sum[
it[[i, k]] mt[[k, j]], {k, 1, n}]],
{i, 1, n}, {j, 1, n}];
Subsuperscript["g","μ","σ"] -> MatrixForm[mx]

"Christoffelsymbole"
chr=ParallelTable[smp[(1/2)Sum[(it[[i, s]])
(D[mt[[s, j]], x[[k]]]+D[mt[[s, k]], x[[j]]] -D[mt[[j, k]], x[[s]]]), {s, 1, n}]],
{i, 1, n}, {j, 1, n}, {k, 1, n}];
crs=ParallelTable[If[UnsameQ[chr[[i, j, k]], 0],
{ToString[Γ[i, j, k]]"->", chr[[i, j, k]]}], {i, 1, n}, {j, 1, n}, {k, 1, j}];
TableForm[Partition[DeleteCases[Flatten[crs], Null], 2]]

"4-Geschwindigkeit"
γ=1/(Sqrt[1-(β^2)]);ut=γ(Sqrt[it[[1,1]]]);uX=β γ/(Sqrt[-mt[[2,2]]]);
u=smp[{ut,uX,0,0}];"(u^μ)"->u

"4-Beschleunigung"
Г=smp[Table[-Sum[chr[[i,j,k]]u[[j]]u[[k]], {j,1,n}, {k,1,n}], {i,1,n}]]; "\!\(\*SuperscriptBox[\(Г\), \(μ\)]\)"->Г
ẍ=smp[{D[ut,X] uX,D[uX,X]uX,0,0}];  "\!\(\*SuperscriptBox[\(ẍ\), \(μ\)]\)"->ẍ
α=smp[ẍ-Г];  "(α^μ)"->α

"Kraft"
F=Simplify[(Sqrt[Sum[-mt[[j,k]]α[[j]]α[[k]], {j,1,n}, {k,1,n}]])];"F"->F

(* Beispiele *)

c=299792458m/sek;
G=667384*^-16(m^3)(kg^(-1))(sek^(-2));  
M=6*^24kg;
R=6370000m;

(* Radiale Bewegung, gleicher Schwerpunktabstand im System der Erde *)

Clear[r];

X=r;
Y=0;
Z=0;

"(u^μ)"->smp[u]
"\!\(\*SuperscriptBox[\(ẍ\), \(μ\)]\)"->smp[ẍ]
"\!\(\*SuperscriptBox[\(Г\), \(μ\)]\)"->smp[Г]
"(α^μ)"->smp[α]
"F"->smp[F]

r=R (c^2)/G/M;
Grid[Table[{"β"->β,", ","F"->N[Evaluate[F]((c^4)/(G M))]},{β,0,1,0.1}],Alignment->Left]
Plot[F((c^4)/(G M))((sek^2)/m),{β,0,1},Frame->True, ImagePadding->{{40,1},{20,1}}]

(* Radiale Bewegung, gleicher Schwerpunktabstand im System der Waage *)

Clear[r];

X=r γ;
Y=0;
Z=0;

"(u^μ)"->smp[u]
"\!\(\*SuperscriptBox[\(ẍ\), \(μ\)]\)"->smp[ẍ]
"\!\(\*SuperscriptBox[\(Г\), \(μ\)]\)"->smp[Г]
"(α^μ)"->smp[α]
"F"->smp[F]

r=R (c^2)/G/M;
Grid[Table[{"β"->β,", ","F"->N[Evaluate[F]((c^4)/(G M))]},{β,0,1,0.1}],Alignment->Left]
Plot[F((c^4)/(G M))((sek^2)/m),{β,0,1},Frame->True, ImagePadding->{{40,1},{20,1}}]

(* Transversale Bewegung *)

Clear[r];

X=0;
Y=r;
Z=0;

"(u^μ)"->smp[u]
"\!\(\*SuperscriptBox[\(ẍ\), \(μ\)]\)"->smp[ẍ]
"\!\(\*SuperscriptBox[\(Г\), \(μ\)]\)"->smp[Г]
"(α^μ)"->smp[α]
"F"->smp[F]

r=R (c^2)/G/M;
Grid[Table[{"β"->β,", ","F"->N[Evaluate[F]((c^4)/(G M))]},{β,0,1,0.1}],Alignment->Left]
Plot[F((c^4)/(G M))((sek^2)/m),{β,0,1},Frame->True, ImagePadding->{{40,1},{20,1}}]

Clear[X];

Y=r;
Z=0;

kg=1;
m=1;
sek=1;

β=0.999;

dp= Style[(J^J), White];
Manipulate[Grid[{{Show[
Graphics[{Black,Line[{{-6R,R},{+2R,R}}]},
PlotRange->{{-6R,+2R},{-7/5R,+7/5R}},
Frame->True,ImageSize->720, 
AspectRatio->2 7/5/8],
Graphics[{Darker[Cyan],Disk[{c β t,0},{R/γ,R}]}],
Graphics[{Red,PointSize->0.01,Point[{0,R}]}]]},
{"      t"->N@t "sek" ,dp},
{"      v"->N@β "c" ,dp},
{"      v"->N@β c "m/sek" ,dp},
{"      x"->N@c β t "m" ,dp},
{"      F"->N@Block[{X=β t/(G M/(c^3)) γ},F] (c^4)/G/M "m/sek²" ,dp}},
Alignment->Left,Spacings->{0,0}],
{{t,0},-205/1000,85/1000,5/10000}]