(* schwarzschild.yukterez.net *)

ClearAll["Global`*"] 
G=1;
M=1;
c=1;
rs=2 G M/(c^2);
wp=MachinePrecision;
set={"GlobalAdaptive","MaxErrorIncreases"->100,Method->"GaussKronrodRule"};
para=20;
pstep=1/3;
j[v_]:=Sqrt[1-(v^2)/(c^2)];
J=j[v0];
k[r_]:=Sqrt[1-rs/r];
κ=k[r0];
r0=3;
v0=1/3;
φ=Pi/4;
θ0=0;
υ=20;
vr0=v0 Sin[φ] κ/J;vθ0=v0/r0 Cos[φ]/J;
d1=υ/10;d2=d1;f=3;
Needs["DifferentialEquations`NDSolveProblems`"];Needs["DifferentialEquations`NDSolveUtilities`"];
i[x_]:=If[r[t]==2,1,If[r[t]<2,Im[x],x]];
sol=NDSolve[{r''[t]==-((G M)/((r[t])^2))+r[t] ((θ'[t])^2)-(3 G M)/(c^2) ((θ'[t])^2),r'[0]==vr0,r[0]==r0,θ''[t]==-((2 r'[t] θ'[t])/r[t]),θ'[0]==vθ0,θ[0]==θ0},{r,θ},{t,0,υ},MaxSteps->Infinity,Method->Automatic,WorkingPrecision->wp,InterpolationOrder->All];
ṫ[t_]:=Evaluate[(Sqrt[(c^2) r[t]+r[t] ((r'[t])^2)-(c^2) rs+((r[t])^3) ((θ'[t])^2)-((r[t])^2) rs ((θ'[t])^2)]/(c Sqrt[r[t]-rs] Sqrt[1-rs/r[t]]))/.sol][[1]];
тi[t_]:=NIntegrate[ṫ[ttt],{ttt,0,t},Method->set,MaxRecursion->10000];
grr[r_]:=If[Abs[r]>2,(1-2/Abs[r])^-1,(1-(2((Abs[r])^2))/(2^3))^-1];
gtt[r_]:=If[Abs[r]>2,(1-2/Abs[r])^-1,4(3Sqrt[1-2/2]-Sqrt[1-(2((Abs[r])^2))/(2^3)])^-2];
Ȓ[я_]:=Quiet[NIntegrate[Sqrt[Abs[grr[r]]],{r,0,Abs[я]},Method->set,MaxRecursion->100]];
Ř[я_]:=Quiet[NIntegrate[Sqrt[Abs[grr[r]-1]],{r,Abs[я],para},Method->set,MaxRecursion->100]];
ř[я_]:=Quiet[Ř[0]-Ř[я]];
grid[я_]:=Quiet[x/.FindRoot[Ȓ[x]-я,{x,1}]];
x[t_]:=(Sin[Evaluate[θ[t]/.sol]] Evaluate[r[t]/.sol])[[1]]
y[t_]:=(Cos[Evaluate[θ[t]/.sol]] Evaluate[r[t]/.sol])[[1]]
R[t_]:=Evaluate[r[t]/.sol][[1]];
ς[t_]:=1/Sqrt[1-2/R[t]];
g[t_]:=ṫ[t]/ς[t];
v[t_]:=Sqrt[-1+((g[t])^2)]/g[t];
grph=Graphics[{Table[{LightGray,Circle[{0,0},grid[n]]},{n,pstep,para,pstep}]},Frame->True,ImageSize->440,PlotRange->PR,ImagePadding->Automatic];
s[text_]:=Style[text,FontSize->font];font=11;
PR=4;

t3=t/.FindRoot[R[t]-3,{t,5}]

t2=t/.FindRoot[R[t]-2,{t,5}]

t0=7.107284;

Quiet[Do[Print[(*Animation nach Eigenzeit*)Rasterize[Grid[{{Show[grph,Graphics[{Magenta,Opacity[0.05],Annulus[{0,0},{2,para}]}],Graphics[{White,Opacity[0.7],Disk[{0,0},2]}],Graphics[{Cyan,Opacity[0.1],Disk[{0,0},2]}],Graphics[{PointSize[0.01],Red,Point[{x[т],y[т]}]}],
grph,ParametricPlot[{x[η],y[η]},{η,Max[1*^-16,т-d2],т},ColorFunction->Function[{x,y,η},Hue[0,1,0.5,Max[Min[(-т+(η+d2))/d2,1],0]]],ColorFunctionScaling->False],ParametricPlot[{x[η],y[η]},{η,0,Max[1*^-16,т-d1]},PlotStyle->{Thickness[0.004],LightGray}]]},{Grid[{{s["  τ"]," = ",s[N[т,8]],s[" GM/c³"]},{s["  t"]," = ",s[N[тi[т]]],s[" GM/c³"]},{s["  ṫ"]," = ",s[N[ṫ[т],8]],s[" dτ/dt"]},
{s["  ς"]," = ",s[N[ς[т],8]],s[" dτ/dt"]},
{s["  γ"]," = ",s[N[g[т],8]],s[" dτ/dt"]},{s["  φ"]," = ",s[N[Evaluate[(θ[т]/.sol) 180/Pi][[1]],8]],s[" grad"]},{s["  r"]," = ",s[N[R[т],8]],s[" GM/c²"]},{s["  x"]," = ",s[N[x[т],8]],s[" GM/c²"]},{s["  y"]," = ",s[N[y[т],8]],s[" GM/c²"]},
{s["  v"]," = ",s[N[v[т],8]],s[" c"]},{s[" "],"   ",s["                 "],s[" "]}},Alignment->Left,Spacings->{0,1/2}]}},Alignment->Left]]],{т,{0,t3,t2,t0}}]]

Plot[R[t],{t,0,t0},ImageSize->260,Frame->True,ImagePadding->32,GridLines->{{t3,t2},{0,2,3}}]

Plot[тi[t],{t,0,t0},ImageSize->260,Frame->True,ImagePadding->32,GridLines->{{t3,t2},{0,10.4355}}]

Export["tau_vs_t_interior.html",EvaluationNotebook[],"GraphicsOutput"->"PNG"]