Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- (*||| Schwarzschild Simulator |||||||| yukterez.net ||||||||||| Simplified Version |||||||| Syntax: Mathematica |||*)
- (*|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||*)
- ClearAll["Global`*"] (* Variablen frei machen *)
- G = 1; M = 1; c = 1; rs = 2 G M/c^2; (* Einheiten *)
- wp = MachinePrecision; (* numerische Genauigkeit *)
- set= {"GlobalAdaptive", "MaxErrorIncreases" -> 100, Method -> "GaussKronrodRule"};
- τMax = Min[τM, plunge-0.000001]; (* Integrationsende *)
- para = 20; pstep = 1/3; (* Paraboloid Grid *)
- (*|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||*)
- j[v_] := Sqrt[1 - v^2/c^2]; J = j[v0]; (* Gammafaktor *)
- k[r_] := Sqrt[1 - rs/r]; κ = k[r0]; (* Schwarzschildfaktor *)
- (*|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||*)
- r0 = 3/2 rs; (* Startradius *)
- v0 = 7/10 c; (* lokale Startgeschwindigkeit *)
- φ = Pi/4; (* Abschusswinkel *)
- θ0 = 0; (* Startwinkel *)
- τM = 18 G M/c^3; (* Simulationsdauer *)
- (*|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||*)
- vr0 = v0 Sin[φ] κ/J; vθ0 = v0/r0 Cos[φ]/J; (* Längenkontraktion und Tiefenexpansion *)
- d1 = τM; d2 = d1; f = 1; (* Schweifdauer und Frameanzahl *)
- (*|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||*)
- sol = NDSolve[{ (* Differentialgleichung *)
- 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, τM},
- MaxSteps -> Infinity,
- Method -> Automatic,
- WorkingPrecision -> wp,
- InterpolationOrder -> All,
- StepMonitor :> (laststep=plunge; plunge=t;
- stepsize=plunge-laststep;), Method->{"EventLocator",
- "Event" :> (If[stepsize<1*^-5, 0, 1])}];
- (*|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||*)
- grr[r_]:=(1-2/Abs[r])^-1;
- R[t_] := Evaluate[r[t] /. sol][[1]];
- Ȓ[я_]:=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]]
- (*|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||*)
- s[text_] := Style[text, FontSize -> font]; font = 11; (* Stil der Anzeigetafel *)
- PR = 2 r0; (* Plot Range *)
- (*|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||*)
- Quiet[Do[Print[ (* Animation nach Eigenzeit *)
- Rasterize[Grid[{{Show[
- Graphics[{Table[{LightGray, Circle[{0, 0}, grid[n]]}, {n, pstep, para, pstep}]},
- Frame -> True, ImageSize -> 360,
- PlotRange -> PR, ImagePadding -> Automatic],
- Graphics[{Magenta, Opacity[0.1], Annulus[{0, 0}, {2, para}]}],
- Graphics[{White, Opacity[0.7], Disk[{0, 0}, 2]}],
- Graphics[{Cyan, Opacity[0.1], Disk[{0, 0}, 2]}],
- Graphics[{Black, Dashed, Circle[{0, 0}, r0]}],
- Graphics[{PointSize[0.01], Red, Point[{x[т], y[т]}]}],
- ParametricPlot[{x[η], y[η]}, {η, Max[1*^-16, т - d2], т},
- ColorFunction -> Function[{x, y, η},
- Hue[0, 1, 0.5, Max[Min[(-т + (η + d2))/d2, 1], 0]]],
- ColorFunctionScaling -> False]]},
- {Grid[{
- {s[" Eigenzeit"], " = ", s[N[т, 8]], s[" GM/c³"]},
- {s[" Radialkoordinate"], " = ", s[N[R[т] , 8]], s[" GM/c²"]},
- {s[" Winkelkoordinate"], " = ", s[N[Evaluate[(θ[т] /. sol) 180/Pi][[1]], 8]], s[" Grad°"]},
- {s[" "], " ", s[" "], s[" "]}
- }, Alignment -> Left, Spacings -> {0, 1/2}]}}, Alignment -> Left]]
- ], {т, τMax/f, τMax, τMax/f}]]
- (*|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||*)
- (*||||| Syntax: Mathematica ||||| schwarzschild.yukterez.net |||| kerr.yukterez.net ||||| Simon Tyran, Vienna |||||*)
RAW Paste Data