Yukterez

Kerr Newman 3D Simulator, Boyer Lindquist

Aug 20th, 2017
265
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. (* Simulator-Code für Photonen, geladene und neutrale Teilchen *)
  2. (* in Boyer Lindquist Koordinaten *)
  3.  
  4. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  5. (* |||| Mathematica | kerr.newman.yukterez.net | 06.08.2017 - 13.06.2020, Version 25 |||| *)
  6. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  7.  
  8. wp=MachinePrecision;
  9. mt1=Automatic;
  10. mt2={"EquationSimplification"-> "Residual"};
  11. mt3={"ImplicitRungeKutta", "DifferenceOrder"-> 20};
  12. mt4={"StiffnessSwitching", Method-> {"ExplicitRungeKutta", Automatic}};
  13. mt5={"EventLocator", "Event"-> (r[τ]-1001/1000 rA)};
  14. mta=mt1; (* mt1: Speed, mt3: Accuracy *)
  15. dgl=1; (* 1: Full d²x/dτ², 2: Mixed dx/dτ, 3: Testmodul, 4: Weak Field, 5: Newton *)
  16.  
  17. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  18. (* |||||||| 1) STARTBEDINGUNGEN EINGEBEN |||||||||||||||||||||||||||||||||||||||||||||||| *)
  19. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  20.  
  21. A=a; (* pseudosphärisch [BL]: A=0, kartesisch [KS]: A=a *)
  22.  
  23. tmax=300; (* Eigenzeit *)
  24. Tmax=300; (* Koordinatenzeit *)
  25. TMax=Min[Tmax, т[plunge-1/100]]; tMax=Min[tmax, plunge-1/100]; (* Integrationsende *)
  26.  
  27. r0 = Sqrt[7^2-a^2]; (* Startradius *)
  28. r1 = r0+2; (* Endradius wenn v0=vr0=vr1 *)
  29. θ0 = π/2; (* Breitengrad *)
  30. φ0 = 0; (* Längengrad *)
  31. a = 9/10; (* Spinparameter *)
  32. ℧ = 2/5; (* spezifische Ladung des schwarzen Lochs *)
  33. q = 0; (* spezifische Ladung des Testkörpers *)
  34.  
  35. v0 = 2/5; (* Anfangsgeschwindigkeit *)
  36. α0 = 0; (* vertikaler Abschusswinkel *)
  37. i0 = ArcTan[5/6]; (* Bahninklinationswinkel *)
  38.  
  39. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  40. (* |||||||| 2) GESCHWINDIGKEITS-, ENERGIE UND DREHIMPULSKOMPONENTEN ||||||||||||||||||||| *)
  41. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  42.  
  43. vr0=v0 Sin[α0]; (* radiale Geschwindigkeitskomponente *)
  44. vθ0=v0 Cos[α0] Sin[i0]; (* longitudinale Geschwindigkeitskomponente *)
  45. vφ0=v0 Cos[α0] Cos[i0]; (* latitudinale Geschwindigkeitskomponente *)
  46.  
  47. vrj[τ_]:=R'[τ]/Sqrt[Δi[τ]] Sqrt[Σi[τ] (1+μ v[τ]^2)];
  48. vθj[τ_]:=Θ'[τ] Sqrt[Σi[τ] (1+μ v[τ]^2)];
  49. vφj[τ_]:=Evaluate[(-(((a^2 Cos[θ[τ]]^2+r[τ]^2) (a^2+℧^2-2 r[τ]+r[τ]^2) Sin[θ[τ]] Sqrt[1-
  50. μ^2 v[τ]^2] (-φ'[τ]-(a q ℧ r[τ])/((a^2 Cos[θ[τ]]^2+r[τ]^2) (a^2+℧^2-2 r[τ]+r[τ]^2))+
  51. (ε Csc[θ[τ]]^2 (a (-a^2-℧^2+2 r[τ]-r[τ]^2) Sin[θ[τ]]^2+a (a^2+
  52. r[τ]^2) Sin[θ[τ]]^2))/((a^2 Cos[θ[τ]]^2+r[τ]^2) (a^2+℧^2-2 r[τ]+r[τ]^2))+(a q ℧ r[τ] (a^2+
  53. ℧^2-2 r[τ]+r[τ]^2-a^2 Sin[θ[τ]]^2))/((a^2 Cos[θ[τ]]^2+r[τ]^2)^2 (a^2+℧^2-2 r[τ]+
  54. r[τ]^2) (1-μ^2 v[τ]^2))))/((a^2+℧^2-2 r[τ]+r[τ]^2-a^2 Sin[θ[τ]]^2) Sqrt[((a^2+r[τ]^2)^2-
  55. a^2 (a^2+℧^2-2 r[τ]+r[τ]^2) Sin[θ[τ]]^2)/(a^2 Cos[θ[τ]]^2+r[τ]^2)]))) /. sol][[1]]
  56. vtj[τ_]:=Sqrt[vrj[τ]^2+vθj[τ]^2+vφj[τ]^2];
  57. vr[τ_]:=vrj[τ]/vtj[τ]*v[τ];
  58. vθ[τ_]:=vθj[τ]/vtj[τ]*v[τ];
  59. vφ[τ_]:=vφj[τ]/vtj[τ]*v[τ];
  60. VΦ[τ_]:=Sqrt[v[τ]^2-vθ[τ]^2-vr[τ]^2];
  61. Vφ[τ_]:=If[q==0, Vφ[τ], VΦ[τ]];
  62.  
  63. x0[A_]:=Sqrt[r0^2+A^2] Sin[θ0] Cos[φ0]; (* Anfangskoordinaten *)
  64. y0[A_]:=Sqrt[r0^2+A^2] Sin[θ0] Sin[φ0];
  65. z0[A_]:=r0 Cos[θ0];
  66.  
  67. ε0=Sqrt[δ Ξ/χ]/j[v0]+Lz ω0;
  68. ε=ε0+((q r0 ℧)/(r0^2+a^2 Cos[θ0]^2));
  69. εζ:=Sqrt[Δ Σ/Χ]/j[ν]+Lz ωζ+((q r[τ] ℧)/(r[τ]^2+a^2 Cos[θ[τ]]^2));
  70. LZ=vφ0 Ы/j[v0];
  71. Lz=LZ+((q a r0 ℧ Sin[θ0]^2)/(r0^2+a^2 Cos[θ0]^2)) j[v0]^2;
  72. Lζ:=vφ0 я/j[ν]+0((q a r[τ] ℧ Sin[θ[τ]]^2)/(r[τ]^2+a^2 Cos[θ[τ]]^2));
  73. pθ0=vθ0 Sqrt[Ξ]/j[v0]; pθζ:=θ'[τ] Σ;
  74. pr0=vr0 Sqrt[(Ξ/δ)/j[v0]^2];
  75. Qk=Limit[pθ0^2+(Lz^2 Csc[θ1]^2-a^2 (ε^2+μ)) Cos[θ1]^2, θ1->θ0]; (* Carter Konstante *)
  76. Q=Limit[pθ0^2+(Lz^2 Csc[θ1]^2-a^2 (ε^2+μ)) Cos[θ1]^2, θ1->θ0];
  77. Qζ:=pθζ^2+(Lz^2 Csc[θ[τ]]^2-a^2 (εζ^2+μ)) Cos[θ[τ]]^2;
  78. k=Q+Lz^2+a^2 (ε^2+μ); kζ:=Qζ+Lz^2+a^2 (εζ^2+μ);
  79. (* ISCO *)
  80. isco = rISCO/.Solve[0 == rISCO (6 rISCO-rISCO^2-9 ℧^2+3 a^2)+4 ℧^2 (℧^2-a^2)-
  81. 8 a (rISCO-℧^2)^(3/2) && rISCO>=rA, rISCO][[1]];
  82. μ=If[Abs[v0]==1, 0, If[Abs[v0]<1, -1, 1]]; (* Baryon: μ=-1, Photon: μ=0 *)
  83.  
  84. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  85. (* |||||||| 3) FLUCHTGESCHWINDIGKEIT UND BENÖTIGTER ABSCHUSSWINKEL |||||||||||||||||||||| *)
  86. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  87.  
  88. vEsc=If[q==0, ж0, Abs[(\[Sqrt](r0^2 (r0^2 (δ Ξ-χ)+2 q r0 χ ℧-q^2 χ ℧^2)+
  89. 2 a^2 r0 (r0 δ Ξ-r0 χ+q χ ℧) Cos[θ0]^2+a^4 (δ Ξ-
  90. χ) Cos[θ0]^4))/(Sqrt[χ] (r0 (r0-q ℧)+a^2 Cos[θ0]^2))]];
  91. (* horizontaler Photonenkreiswinkel, i0 *)
  92. iP[r0_, a_]:=Normal[iPh/.NSolve[1/(8 (r0^2+a^2 Cos[θ0]^2)^3) (a^2+(-2+r0) r0+
  93. ℧^2) (8 r0 (r0^2+a^2 Cos[θ0]^2) Sin[iPh]^2+1/((a^2-2 r0+r0^2+℧^2) (r0^2+
  94. a^2 Cos[θ0]^2)) 8 a (Cos[iPh] Sin[θ0] (a^2-2 r0+r0^2+℧^2-a^2 Sin[θ0]^2) Sqrt[((a^2+
  95. r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)/(r0^2+a^2 Cos[θ0]^2)]+(a (a^2+r0^2) Sin[θ0]^2+
  96. a (-a^2+2 r0-r0^2-℧^2) Sin[θ0]^2) (Sqrt[((a^2-2 r0+r0^2+℧^2) (r0^2+a^2 Cos[θ0]^2))/((a^2+
  97. r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)]+(a (2 r0-℧^2) Cos[iPh] Sin[θ0] Sqrt[((a^2+
  98. r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)/(r0^2+a^2 Cos[θ0]^2)])/((a^2+r0^2)^2-a^2 (a^2-
  99. 2 r0+r0^2+℧^2) Sin[θ0]^2))) (-(1/((a^2-2 r0+r0^2+℧^2) (r0^2+
  100. a^2 Cos[θ0]^2)))2 a^2 Cot[θ0]^2 (Cos[iPh] Sin[θ0] (-a (a^2+r0^2) Sin[θ0]^2+a (a^2-2 r0+
  101. r0^2+℧^2) Sin[θ0]^2) Sqrt[((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)/(r0^2+
  102. a^2 Cos[θ0]^2)]+((a^2+r0^2)^2 Sin[θ0]^2+a^2 (-a^2+2 r0-r0^2-
  103. ℧^2) Sin[θ0]^4) (Sqrt[((a^2-2 r0+r0^2+℧^2) (r0^2+a^2 Cos[θ0]^2))/((a^2+r0^2)^2-a^2 (a^2-
  104. 2 r0+r0^2+℧^2) Sin[θ0]^2)]+(a (2 r0-℧^2) Cos[iPh] Sin[θ0] Sqrt[((a^2+r0^2)^2-a^2 (a^2-
  105. 2 r0+r0^2+℧^2) Sin[θ0]^2)/(r0^2+a^2 Cos[θ0]^2)])/((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+
  106. ℧^2) Sin[θ0]^2)))+1/((a^2-2 r0+r0^2+℧^2) (r0^2+a^2 Cos[θ0]^2)) 2 r0 (r0-
  107. ℧^2) Csc[θ0]^2 (Cos[iPh] Sin[θ0] (-a (a^2+r0^2) Sin[θ0]^2+a (a^2-2 r0+r0^2+
  108. ℧^2) Sin[θ0]^2) Sqrt[((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)/(r0^2+
  109. a^2 Cos[θ0]^2)]+((a^2+r0^2)^2 Sin[θ0]^2+a^2 (-a^2+2 r0-r0^2-℧^2) Sin[θ0]^4) (Sqrt[((a^2-
  110. 2 r0+r0^2+℧^2) (r0^2+a^2 Cos[θ0]^2))/((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)]+
  111. (a (2 r0-℧^2) Cos[iPh] Sin[θ0] Sqrt[((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)/(r0^2+
  112. a^2 Cos[θ0]^2)])/((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2))))+1/((a^2-2 r0+r0^2+
  113. ℧^2) (r0^2+a^2 Cos[θ0]^2)) 8 Csc[θ0]^2 (Cos[iPh] Sin[θ0] (-a (a^2+r0^2) Sin[θ0]^2+a (a^2-
  114. 2 r0+r0^2+℧^2) Sin[θ0]^2) Sqrt[((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)/(r0^2+
  115. a^2 Cos[θ0]^2)]+((a^2+r0^2)^2 Sin[θ0]^2+a^2 (-a^2+2 r0-r0^2-℧^2) Sin[θ0]^4) (Sqrt[((a^2-
  116. 2 r0+r0^2+℧^2) (r0^2+a^2 Cos[θ0]^2))/((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)]+
  117. (a (2 r0-℧^2) Cos[iPh] Sin[θ0] Sqrt[((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)/(r0^2+
  118. a^2 Cos[θ0]^2)])/((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2))) (1/((a^2-2 r0+r0^2+
  119. ℧^2) (r0^2+a^2 Cos[θ0]^2)) a^2 Cot[θ0]^2 (Cos[iPh] Sin[θ0] (-a (a^2+r0^2) Sin[θ0]^2+
  120. a (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2) Sqrt[((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+
  121. ℧^2) Sin[θ0]^2)/(r0^2+a^2 Cos[θ0]^2)]+((a^2+r0^2)^2 Sin[θ0]^2+a^2 (-a^2+2 r0-r0^2-
  122. ℧^2) Sin[θ0]^4) (Sqrt[((a^2-2 r0+r0^2+℧^2) (r0^2+a^2 Cos[θ0]^2))/((a^2+r0^2)^2-
  123. a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)]+(a (2 r0-℧^2) Cos[iPh] Sin[θ0] Sqrt[((a^2+
  124. r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)/(r0^2+a^2 Cos[θ0]^2)])/((a^2+r0^2)^2-a^2 (a^2-
  125. 2 r0+r0^2+℧^2) Sin[θ0]^2)))+1/((a^2-2 r0+r0^2+℧^2) (r0^2+a^2 Cos[θ0]^2)) r0 (-r0+
  126. ℧^2) Csc[θ0]^2 (Cos[iPh] Sin[θ0] (-a (a^2+r0^2) Sin[θ0]^2+a (a^2-2 r0+r0^2+
  127. ℧^2) Sin[θ0]^2) Sqrt[((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)/(r0^2+a^2 Cos[θ0]^2)]+
  128. ((a^2+r0^2)^2 Sin[θ0]^2+a^2 (-a^2+2 r0-r0^2-℧^2) Sin[θ0]^4) (Sqrt[((a^2-2 r0+r0^2+
  129. ℧^2) (r0^2+a^2 Cos[θ0]^2))/((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)]+(a (2 r0-
  130. ℧^2) Cos[iPh] Sin[θ0] Sqrt[((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)/(r0^2+
  131. a^2 Cos[θ0]^2)])/((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2))))+1/((a^2-2 r0+r0^2+
  132. ℧^2)^2 (r0^2+a^2 Cos[θ0]^2)^2) Csc[θ0]^2 (Cos[iPh] Sin[θ0] (a^2-2 r0+r0^2+℧^2-
  133. a^2 Sin[θ0]^2) Sqrt[((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)/(r0^2+a^2 Cos[θ0]^2)]+
  134. (a (a^2+r0^2) Sin[θ0]^2+a (-a^2+2 r0-r0^2-℧^2) Sin[θ0]^2) (Sqrt[((a^2-2 r0+r0^2+℧^2) (r0^2+
  135. a^2 Cos[θ0]^2))/((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)]+(a (2 r0-
  136. ℧^2) Cos[iPh] Sin[θ0] Sqrt[((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)/(r0^2+
  137. a^2 Cos[θ0]^2)])/((a^2+r0^2)^2-a^2 (a^2-2 r0+r0^2+℧^2) Sin[θ0]^2)))^2 (r0 (a^2 (3 a^2+
  138. 4 ℧^2+4 (a-℧) (a+℧) Cos[2 θ0]+a^2 Cos[4 θ0])+8 r0 (r0^3+2 a^2 r0 Cos[θ0]^2-
  139. a^2 Sin[θ0]^2))+2 a^4 Sin[2 θ0]^2))==0,iPh,Reals]][[1]]/.C[1]->0
  140.  
  141. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  142. (* |||||||| 4) HORIZONTE UND ERGOSPHÄREN RADIEN ||||||||||||||||||||||||||||||||||||||||| *)
  143. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  144.  
  145. rE=1+Sqrt[1-a^2 Cos[θ]^2-℧^2]; (* äußere Ergosphäre *)
  146. RE[A_, w1_, w2_]:=Xyz[xyZ[
  147. {Sqrt[rE^2+A^2] Sin[θ] Cos[φ], Sqrt[rE^2+A^2] Sin[θ] Sin[φ], rE Cos[θ]}, w1], w2];
  148. rG=1-Sqrt[1-a^2 Cos[θ]^2-℧^2]; (* innere Ergosphäre *)
  149. RG[A_, w1_, w2_]:=Xyz[xyZ[
  150. {Sqrt[rG^2+A^2] Sin[θ] Cos[φ], Sqrt[rG^2+A^2] Sin[θ] Sin[φ], rG Cos[θ]}, w1], w2];
  151. rA=1+Sqrt[1-a^2-℧^2]; (* äußerer Horizont *)
  152. RA[A_, w1_, w2_]:=Xyz[xyZ[
  153. {Sqrt[rA^2+A^2] Sin[θ] Cos[φ], Sqrt[rA^2+A^2] Sin[θ] Sin[φ], rA Cos[θ]}, w1], w2];
  154. rI=1-Sqrt[1-a^2-℧^2]; (* innerer Horizont *)
  155. RI[A_, w1_, w2_]:=Xyz[xyZ[
  156. {Sqrt[rI^2+A^2] Sin[θ] Cos[φ], Sqrt[rI^2+A^2] Sin[θ] Sin[φ], rI Cos[θ]}, w1], w2];
  157.  
  158. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  159. (* |||||||| 5) HORIZONTE UND ERGOSPHÄREN PLOT ||||||||||||||||||||||||||||||||||||||||||| *)
  160. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  161.  
  162. horizons[A_, mesh_, w1_, w2_]:=Show[
  163. ParametricPlot3D[RE[A, w1, w2], {φ, 0, 2 π}, {θ, 0, π},
  164. Mesh -> mesh, PlotPoints -> plp, PlotStyle -> Directive[Blue, Opacity[0.10]]],
  165. ParametricPlot3D[RA[A, w1, w2], {φ, 0, 2 π}, {θ, 0, π},
  166. Mesh -> None, PlotPoints -> plp, PlotStyle -> Directive[Cyan, Opacity[0.15]]],
  167. ParametricPlot3D[RI[A, w1, w2], {φ, 0, 2 π}, {θ, 0, π},
  168. Mesh -> None, PlotPoints -> plp, PlotStyle -> Directive[Red, Opacity[0.25]]],
  169. ParametricPlot3D[RG[A, w1, w2], {φ, 0, 2 π}, {θ, 0, π},
  170. Mesh -> None, PlotPoints -> plp, PlotStyle -> Directive[Red, Opacity[0.35]]]];
  171. BLKS:=Grid[{{horizons[a, 35, 0, 0], horizons[0, 35, 0, 0]}}];
  172.  
  173. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  174. (* |||||||| 6) FUNKTIONEN ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  175. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  176.  
  177. j[v_]:=Sqrt[1-μ^2 v^2]; (* Lorentzfaktor *)
  178. mirr=Sqrt[2-℧^2+2 Sqrt[1-a^2-℧^2]]/2; (* irreduzible Masse *)
  179. я=Sqrt[Χ/Σ]Sin[θ[τ]]; (* axialer Umfangsradius *)
  180. яi[τ_]:=Sqrt[Χi[τ]/Σi[τ]]Sin[Θ[τ]];
  181. Ы=Sqrt[χ/Ξ]Sin[θ0];
  182.  
  183. Σ=r[τ]^2+a^2 Cos[θ[τ]]^2; (* poloidialer Umfangsradius *)
  184. Σi[τ_]:=R[τ]^2+a^2 Cos[Θ[τ]]^2;
  185. Ξ=r0^2+a^2 Cos[θ0]^2;
  186.  
  187. Δ=r[τ]^2-2r[τ]+a^2+℧^2;
  188. Δr[r_]:=r^2-2r+a^2+℧^2;
  189. Δi[τ_]:=R[τ]^2-2R[τ]+a^2+℧^2;
  190. δ=r0^2-2r0+a^2+℧^2;
  191.  
  192. Χ=(r[τ]^2+a^2)^2-a^2 Sin[θ[τ]]^2 Δ;
  193. Χi[τ_]:=(R[τ]^2+a^2)^2-a^2 Sin[Θ[τ]]^2 Δi[τ];
  194. χ=(r0^2+a^2)^2-a^2 Sin[θ0]^2 δ;
  195.  
  196. Xj=a Sin[θ0]^2;
  197. xJ[τ_]:=a Sin[Θ[τ]]^2;
  198. XJ=a Sin[θ[τ]]^2;
  199.  
  200. Pr[r_]:=ε(r^2+a^2)+℧ q r-a Lz;
  201. Pt[τ_]:=ε(R[τ]^2+a^2)+℧ q R[τ]-a Lz;
  202. Pτ=ε(r[τ]^2+a^2)+℧ q r[τ]-a Lz;
  203. pτ=ε(r0^2+a^2)+℧ q r0-a Lz;
  204.  
  205. Vr[r_]:=Pr[r]^2-Δr[r](μ^2 r^2+(Lz-a ε)^2+Q); (* effektives radiales Potential *)
  206. Vτ=Pτ^2-Δ(μ^2 r[τ]^2+(Lz-a ε)^2+Q);
  207. Vθ[θ_]:=Q-Cos[θ]^2(a^2 (μ^2-ε^2)+Lz^2 Sin[θ]^(-2)); (* effektives latitudinales Potential *)
  208.  
  209. т[τ_]:=Evaluate[t[τ]/.sol][[1]]; (* Koordinatenzeit nach Eigenzeit *)
  210. д[ξ_]:=Quiet[zt /.FindRoot[т[zt]-ξ, {zt, 0}]]; (* Eigenzeit nach Koordinatenzeit *)
  211. T :=Quiet[д[tk]];
  212.  
  213. ю[τ_]:=Evaluate[t'[τ]/.sol][[1]];
  214. γ[τ_]:=If[μ==0, "Infinity", ю[τ]]; (* totale ZD *)
  215. R[τ_]:=Evaluate[r[τ]/.sol][[1]]; (* Boyer-Lindquist Radius *)
  216. Φ[τ_]:=Evaluate[φ[τ]/.sol][[1]];
  217. Θ[τ_]:=Evaluate[θ[τ]/.sol][[1]];
  218. ß[τ_]:=Sqrt[X'[τ]^2+Y'[τ]^2+Z'[τ]^2 ]/ю[τ];
  219.  
  220. ς[τ_]:=Sqrt[Χi[τ]/Δi[τ]/Σi[τ]]; ς0=Sqrt[χ/δ/Ξ]; (* gravitative ZD *)
  221. ω[τ_]:=(a(2R[τ]-℧^2))/Χi[τ]; ω0=(a(2r0-℧^2))/χ; ωζ=(a(2r[τ]-℧^2))/Χ; (* F-Drag Winkelg *)
  222. Ω[τ_]:=ω[τ] Sqrt[X[τ]^2+Y[τ]^2]; (* Frame Dragging beobachtete Geschwindigkeit *)
  223. й[τ_]:=ω[τ] яi[τ] ς[τ]; й0=ω0 Ы ς0; (* Frame Dragging lokale Geschwindigkeit *)
  224.  
  225. ж[τ_]:=Sqrt[ς[τ]^2-1]/ς[τ]; ж0=Sqrt[ς0^2-1]/ς0; (* Fluchtgeschwindigkeit *)
  226. V[τ_]:=If[μ==0, 1, Re[Sqrt[-ς[τ]^2+ю[τ]^2]/ю[τ]]];
  227. (* Fluchtgeschwindigkeit von r0 nach r1 *)
  228. vd1:=v1/.NSolve[Sqrt[δ Ξ/χ]/Sqrt[1-v1^2]+((q r0 ℧)/(r0^2+a^2 Cos[θ0]^2))==Sqrt[((a^2+
  229. (-2+r1) r1+℧^2) (r1^2+a^2 Cos[θ0]^2))/((a^2+r1^2)^2-a^2 (a^2+(-2+r1) r1+℧^2) Sin[θ0]^2)]+
  230. (a^2 q r1 ℧ (2 r1-℧^2) Sin[θ0]^2)/((r1^2+a^2 Cos[θ0]^2) ((a^2+r1^2)^2-a^2 (a^2+(-2+r1) r1+
  231. ℧^2) Sin[θ0]^2))+((q r1 ℧)/(r1^2+a^2 Cos[θ0]^2))&&v1>0,v1][[1]];
  232. (* lokale Dreiergeschwindigkeit *)
  233. vd[τ_]:=Abs[(Sqrt[Δ Σ^3 Χ-ε^2 Σ^2 Χ^2-2 a Lz ε Σ^2 Χ ℧^2-a^2 Lz^2 Σ^2 ℧^4+
  234. 4 a Lz ε Σ^2 Χ r[τ]+2 q ε Σ Χ^2 ℧ r[τ]+4 a^2 Lz^2 Σ^2 ℧^2 r[τ]+2 a Lz q Σ Χ ℧^3 r[τ]-
  235. 4 a^2 Lz^2 Σ^2 r[τ]^2-4 a Lz q Σ Χ ℧ r[τ]^2-q^2 Χ^2 ℧^2 r[τ]^2])/(ε Σ Χ+
  236. a Lz Σ ℧^2-2 a Lz Σ r[τ]-q Χ ℧ r[τ])];
  237.  
  238. v[τ_]:=If[μ==0, 1, Evaluate[vlt'[τ]/.sol][[1]]];
  239.  
  240. vnt[τ_]:=Evaluate[Sqrt[(φ'[τ] r[τ]/Csc[θ[τ]])^2+(θ'[τ] r[τ])^2+r'[τ]^2]/.sol][[1]]
  241. ν:=If[μ==0, 1, Re[Sqrt[(Δ Σ-Χ(εζ-Lζ ωζ)^2)/(μ Χ (εζ-Lζ ωζ)^2)]]];
  242. vesc[τ_]:=Abs[(\[Sqrt](R[τ]^2 (R[τ]^2 (Δi[τ] Σi[τ]-Χi[τ])+2 q R[τ] Χi[τ] ℧-q^2 Χi[τ] ℧^2)+
  243. 2 a^2 R[τ] (R[τ] Δi[τ] Σi[τ]-R[τ] Χi[τ]+q Χi[τ] ℧) Cos[Θ[τ]]^2+a^4 (Δi[τ] Σi[τ]-
  244. Χi[τ]) Cos[Θ[τ]]^4))/(Sqrt[Χi[τ]] (R[τ] (R[τ]-q ℧)+a^2 Cos[Θ[τ]]^2))];
  245.  
  246. dst[τ_]:=Evaluate[str[τ]/.sol][[1]]; (* Strecke *)
  247.  
  248. pΘ[τ_]:=Evaluate[Ξ θ'[τ] /. sol][[1]];
  249. pR[τ_]:=Evaluate[r'[τ] Ξ/δ /. sol][[1]];
  250.  
  251. epot[τ_]:=ε+μ-ekin[τ]; (* potentielle Energie *)
  252. ekin[τ_]:=If[μ==0, ς[τ], 1/Sqrt[1-v[τ]^2]-1]; (* kinetische Energie *)
  253.  
  254. drwf=vr0/(Sqrt[(2+r0)/r0] Sqrt[1-v0^2 μ^2]);
  255. duwf=vθ0/(Sqrt[r0^2] Sqrt[1-v0^2 μ^2]);
  256. dfwf=(vφ0 Csc[θ0]^4 (r0^2 Sin[θ0]^2)^(3/2))/(r0^4 Sqrt[1-v0^2 μ^2]);
  257. dtwf=1/(-2+r0) (-2 a Sin[θ0]^2 dfwf+r0 Sqrt[-μ+(2 μ)/r0+(1-4/r0^2) drwf^2+(-2+r0) r0 duwf^2-
  258. 2 r0 Sin[θ0]^2 dfwf^2+r0^2 Sin[θ0]^2 dfwf^2+(4 a^2 Sin[θ0]^4 dfwf^2)/r0^2]);
  259.  
  260. (* beobachtete Inklination *)
  261. ink0:=б/. Solve[Z'[0]/ю[0] Cos[б]==-Y'[0]/ю[0] Sin[б]&&б>0&&б<2π&&б<δp[r0, a], б][[1]];
  262.  
  263. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  264. (* |||||||| 7) DIFFERENTIALGLEICHUNG |||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  265. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  266.  
  267. dp= Style[\!\(\*SuperscriptBox[\(Y\),\(Y\)]\), White]; n0[z_] := Chop[Re[N[Simplify[z]]]];
  268.  
  269. initcon = NSolve[
  270. dr0 == pr0 δ/Ξ
  271. &&
  272. dθ0 == pθ0/Ξ
  273. &&
  274. dφ0 == 1/(δ Ξ Sin[θ0]^2) (ε (-δ Xj+a Sin[θ0]^2 (r0^2+a^2))+Lz (δ-a^2 Sin[θ0]^2)-
  275. q ℧ r0 a Sin[θ0]^2)
  276. &&
  277. -μ == -(((r0^2+a^2 Cos[θ0]^2) dr0^2)/(a^2-2 r0+r0^2+℧^2))+((a^2-2 r0+
  278. r0^2+℧^2-a^2 Sin[θ0]^2) (dt0)^2)/(r0^2+a^2 Cos[θ0]^2)+(-r0^2-
  279. a^2 Cos[θ0]^2) dθ0^2+(2 a (2 r0-℧^2) Sin[θ0]^2 dt0 dφ0)/(r0^2+
  280. a^2 Cos[θ0]^2)+((-(a^2+r0^2)^2 Sin[θ0]^2+a^2 (a^2-2 r0+r0^2+
  281. ℧^2) Sin[θ0]^4) dφ0^2)/(r0^2+a^2 Cos[θ0]^2)
  282. &&
  283. dt0 > 0,
  284. {dθ0, dr0, dt0, dφ0}, Reals];
  285.  
  286. initkon = NSolve[
  287. dr0 == pr0 δ/Ξ
  288. &&
  289. dθ0 == pθ0/Ξ
  290. &&
  291. dφ0 == 1/(δ Ξ Sin[θ0]^2) (ε (-δ Xj+a Sin[θ0]^2 (r0^2+a^2))+Lz (δ-a^2 Sin[θ0]^2)-
  292. q ℧ r0 a Sin[θ0]^2)
  293. &&
  294. dt0 == 1/(δ Ξ Sin[θ0]^2) (Lz (δ Xj-a Sin[θ0]^2 (r0^2+a^2))+ε (-δ Xj^2+
  295. Sin[θ0]^2 (r0^2+a^2)^2)-q ℧ r0 Sin[θ0]^2 (r0^2+a^2))
  296. &&
  297. dt0 > 0,
  298. {dθ0, dr0, dt0, dφ0}, Reals];
  299.  
  300. DG1={
  301.  
  302. t''[τ]==-(((r'[τ] ((a^2+r[τ]^2) (a^2 Cos[θ[τ]]^2 (q ℧-2 t'[τ])+r[τ] (-q ℧ r[τ]+
  303. 2 (-℧^2+r[τ]) t'[τ]))+a (2 a^4 Cos[θ[τ]]^2+a^2 ℧^2 (3+Cos[2 θ[τ]]) r[τ]-
  304. a^2 (3+Cos[2 θ[τ]]) r[τ]^2+4 ℧^2 r[τ]^3-6 r[τ]^4) Sin[θ[τ]]^2 φ'[τ]))/(a^2+℧^2+(-2+
  305. r[τ]) r[τ])+a^2 θ'[τ] (Sin[2 θ[τ]] (q ℧ r[τ]+(℧^2-2 r[τ]) t'[τ])-2 a Cos[θ[τ]] (℧^2-
  306. 2 r[τ]) Sin[θ[τ]]^3 φ'[τ]))/(a^2 Cos[θ[τ]]^2+r[τ]^2)^2),
  307.  
  308. t'[0]==dt0/.initcon[[1]],
  309. t[0]==0,
  310.  
  311. r''[τ]==((-1+r[τ])/(a^2+℧^2+(-2+r[τ]) r[τ])-r[τ]/(a^2 Cos[θ[τ]]^2+r[τ]^2)) r'[τ]^2+
  312. (a^2 Sin[2 θ[τ]] r'[τ] θ'[τ])/(a^2 Cos[θ[τ]]^2+r[τ]^2)+(1/(8 (a^2 Cos[θ[τ]]^2+
  313. r[τ]^2)^3))(a^2+℧^2+(-2+r[τ]) r[τ]) (8 t'[τ] (a^2 Cos[θ[τ]]^2 (-q ℧+t'[τ])+
  314. r[τ] (q ℧ r[τ]+(℧^2-r[τ]) t'[τ]))+8 r[τ] (a^2 Cos[θ[τ]]^2+r[τ]^2)^2 θ'[τ]^2+
  315. 8 a Sin[θ[τ]]^2 (a^2 Cos[θ[τ]]^2 (q ℧-2 t'[τ])+r[τ] (-q ℧ r[τ]+2 (-℧^2+r[τ]) t'[τ])) φ'[τ]+
  316. Sin[θ[τ]]^2 (r[τ] (a^2 (3 a^2+4 ℧^2+4 (a-℧) (a+℧) Cos[2 θ[τ]]+a^2 Cos[4 θ[τ]])+
  317. 8 r[τ] (2 a^2 Cos[θ[τ]]^2 r[τ]+r[τ]^3-a^2 Sin[θ[τ]]^2))+2 a^4 Sin[2 θ[τ]]^2) φ'[τ]^2),
  318.  
  319. r'[0]==dr0/.initcon[[1]],
  320. r[0]==r0,
  321.  
  322. θ''[τ]==-((a^2 Cos[θ[τ]] Sin[θ[τ]] r'[τ]^2)/((a^2+℧^2+(-2+r[τ]) r[τ]) (a^2 Cos[θ[τ]]^2+
  323. r[τ]^2)))-(2 r[τ] r'[τ] θ'[τ])/(a^2 Cos[θ[τ]]^2+r[τ]^2)+(1/(16 (a^2 Cos[θ[τ]]^2+
  324. r[τ]^2)^3))Sin[2 θ[τ]] (a^2 (-8 t'[τ] (2 q ℧ r[τ]+(℧^2-2 r[τ]) t'[τ])+8 (a^2 Cos[θ[τ]]^2+
  325. r[τ]^2)^2 θ'[τ]^2)+16 a (a^2+r[τ]^2) (q ℧ r[τ]+(℧^2-2 r[τ]) t'[τ]) φ'[τ]+(3 a^6-5 a^4 ℧^2+
  326. 10 a^4 r[τ]+11 a^4 r[τ]^2-8 a^2 ℧^2 r[τ]^2+16 a^2 r[τ]^3+16 a^2 r[τ]^4+8 r[τ]^6+
  327. a^4 Cos[4 θ[τ]] (a^2+℧^2+(-2+r[τ]) r[τ])+4 a^2 Cos[2 θ[τ]] (a^2+℧^2+(-2+
  328. r[τ]) r[τ]) (a^2+2 r[τ]^2)) φ'[τ]^2),
  329.  
  330. θ'[0]==dθ0/.initcon[[1]],
  331. θ[0]==θ0,
  332.  
  333. φ''[τ]==-(1/(4 (a^2 Cos[θ[τ]]^2+r[τ]^2)^2))((r'[τ] (4 a q ℧ (a^2 Cos[θ[τ]]^2-r[τ]^2)-
  334. 8 a (a^2 Cos[θ[τ]]^2+(℧^2-r[τ]) r[τ]) t'[τ]+(a^2 (3 a^2+8 ℧^2+a^2 (4 Cos[2 θ[τ]]+
  335. Cos[4 θ[τ]])) r[τ]-4 a^2 (3+Cos[2 θ[τ]]) r[τ]^2+8 (a^2+℧^2+a^2 Cos[2 θ[τ]]) r[τ]^3-
  336. 16 r[τ]^4+8 r[τ]^5+2 a^4 Sin[2 θ[τ]]^2) φ'[τ]))/(a^2+℧^2+(-2+r[τ]) r[τ])+
  337. θ'[τ] (8 a Cot[θ[τ]] (q ℧ r[τ]+(℧^2-2 r[τ]) t'[τ])+(8 Cot[θ[τ]] (a^2+r[τ]^2)^2-
  338. 2 a^2 (3 a^2+2 ℧^2+4 (-1+r[τ]) r[τ]) Sin[2 θ[τ]]-a^4 Sin[4 θ[τ]]) φ'[τ])),
  339.  
  340. φ'[0]==dφ0/.initcon[[1]],
  341. φ[0]==φ0,
  342.  
  343. str'[τ]==If[μ==0, 1, vd[τ]/Abs[Sqrt[1-vd[τ]^2]]],
  344. str[0]==0,
  345. vlt'[τ]==If[μ==0, 1, vd[τ]],
  346. vlt[0]==0
  347.  
  348. };
  349.  
  350. DG2={
  351.  
  352. t'[τ]==1/(Δ Σ Sin[θ[τ]]^2) (Lz (Δ XJ-a Sin[θ[τ]]^2 (r[τ]^2+a^2))+ε (-Δ XJ^2+
  353. Sin[θ[τ]]^2 (r[τ]^2+a^2)^2)-q ℧ r[τ] Sin[θ[τ]]^2 (r[τ]^2+a^2)),
  354. t[0]==0,
  355.  
  356. r''[τ]==((-1+r[τ])/(a^2+℧^2+(-2+r[τ]) r[τ])-r[τ]/(a^2 Cos[θ[τ]]^2+r[τ]^2)) r'[τ]^2+
  357. (a^2 Sin[2 θ[τ]] r'[τ] θ'[τ])/(a^2 Cos[θ[τ]]^2+r[τ]^2)+(1/(8 (a^2 Cos[θ[τ]]^2+
  358. r[τ]^2)^3))(a^2+℧^2+(-2+r[τ]) r[τ]) (8 t'[τ] (a^2 Cos[θ[τ]]^2 (-q ℧+t'[τ])+
  359. r[τ] (q ℧ r[τ]+(℧^2-r[τ]) t'[τ]))+8 r[τ] (a^2 Cos[θ[τ]]^2+r[τ]^2)^2 θ'[τ]^2+
  360. 8 a Sin[θ[τ]]^2 (a^2 Cos[θ[τ]]^2 (q ℧-2 t'[τ])+r[τ] (-q ℧ r[τ]+2 (-℧^2+r[τ]) t'[τ])) φ'[τ]+
  361. Sin[θ[τ]]^2 (r[τ] (a^2 (3 a^2+4 ℧^2+4 (a-℧) (a+℧) Cos[2 θ[τ]]+a^2 Cos[4 θ[τ]])+
  362. 8 r[τ] (2 a^2 Cos[θ[τ]]^2 r[τ]+r[τ]^3-a^2 Sin[θ[τ]]^2))+2 a^4 Sin[2 θ[τ]]^2) φ'[τ]^2),
  363.  
  364. r'[0]==(pr0 δ)/Ξ,
  365. r[0]==r0,
  366.  
  367. θ''[τ]==-((a^2 Cos[θ[τ]] Sin[θ[τ]] r'[τ]^2)/((a^2+℧^2+(-2+r[τ]) r[τ]) (a^2 Cos[θ[τ]]^2+
  368. r[τ]^2)))-(2 r[τ] r'[τ] θ'[τ])/(a^2 Cos[θ[τ]]^2+r[τ]^2)+(1/(16 (a^2 Cos[θ[τ]]^2+
  369. r[τ]^2)^3))Sin[2 θ[τ]] (a^2 (-8 t'[τ] (2 q ℧ r[τ]+(℧^2-2 r[τ]) t'[τ])+8 (a^2 Cos[θ[τ]]^2+
  370. r[τ]^2)^2 θ'[τ]^2)+16 a (a^2+r[τ]^2) (q ℧ r[τ]+(℧^2-2 r[τ]) t'[τ]) φ'[τ]+(3 a^6-5 a^4 ℧^2+
  371. 10 a^4 r[τ]+11 a^4 r[τ]^2-8 a^2 ℧^2 r[τ]^2+16 a^2 r[τ]^3+16 a^2 r[τ]^4+8 r[τ]^6+
  372. a^4 Cos[4 θ[τ]] (a^2+℧^2+(-2+r[τ]) r[τ])+4 a^2 Cos[2 θ[τ]] (a^2+℧^2+(-2+
  373. r[τ]) r[τ]) (a^2+2 r[τ]^2)) φ'[τ]^2),
  374.  
  375. θ'[0]==pθ0/Ξ,
  376. θ[0]==θ0,
  377.  
  378. φ'[τ]==1/(Δ Σ Sin[θ[τ]]^2) (ε (-Δ XJ+a Sin[θ[τ]]^2 (r[τ]^2+a^2))+Lz (Δ-a^2 Sin[θ[τ]]^2)-
  379. q ℧ r[τ] a Sin[θ[τ]]^2),
  380. φ[0]==φ0,
  381.  
  382. str'[τ]==If[μ==0, 1, vd[τ]/Abs[Sqrt[1-vd[τ]^2]]],
  383. str[0]==0,
  384. vlt'[τ]==If[μ==0, 1, vd[τ]],
  385. vlt[0]==0
  386.  
  387. };
  388.  
  389. DG3={
  390.  
  391. t'[τ]==1/(Δ Σ Sin[θ[τ]]^2) (Lz (Δ XJ-a Sin[θ[τ]]^2 (r[τ]^2+a^2))+ε (-Δ XJ^2+
  392. Sin[θ[τ]]^2 (r[τ]^2+a^2)^2)-q ℧ r[τ] Sin[θ[τ]]^2 (r[τ]^2+a^2)),
  393. t[0]==0,
  394.  
  395. r'[τ]==Sign[vr0] Sqrt[Vτ]/Σ,
  396. r[0]==r0,
  397.  
  398. θ'[τ]==Sign[vθ0] Sqrt[Q-Cos[θ[τ]]^2 (a^2 (μ^2-ε^2)+Lz^2/Sin[θ[τ]]^2)]/Σ,
  399. θ[0]==θ0,
  400.  
  401. φ'[τ]==1/(Δ Σ Sin[θ[τ]]^2) (ε (-Δ XJ+a Sin[θ[τ]]^2 (r[τ]^2+a^2))+Lz (Δ-a^2 Sin[θ[τ]]^2)-
  402. q ℧ r[τ] a Sin[θ[τ]]^2),
  403. φ[0]==φ0,
  404.  
  405. str'[τ]==If[μ==0, 1, vd[τ]/Abs[Sqrt[1-vd[τ]^2]]],
  406. str[0]==0,
  407. vlt'[τ]==If[μ==0, 1, vd[τ]],
  408. vlt[0]==0
  409.  
  410. };
  411.  
  412. DG4={
  413.  
  414. t''[τ]==((-4 a^2 r[τ] Sin[2 θ[τ]] t'[τ] θ'[τ]+r'[τ] (4 a^2 Sin[θ[τ]]^2 t'[τ]+
  415. r[τ]^3 (q ℧-2 t'[τ]+6 a Sin[θ[τ]]^2 φ'[τ])))/((-2+r[τ]) r[τ]^4+4 a^2 r[τ] Sin[θ[τ]]^2)),
  416. t'[0]==dtwf,
  417. t[0]==0,
  418.  
  419. r''[τ]==((r'[τ]^2-t'[τ]^2+t'[τ] (q ℧+2 a Sin[θ[τ]]^2 φ'[τ])+r[τ]^3 (θ'[τ]^2+
  420. Sin[θ[τ]]^2 φ'[τ]^2))/(r[τ] (2+r[τ]))),
  421. r'[0]==drwf,
  422. r[0]==r0,
  423.  
  424. θ''[τ]==((-2 r[τ]^2 r'[τ] θ'[τ]+Cos[θ[τ]] Sin[θ[τ]] φ'[τ] (-4 a t'[τ]+
  425. r[τ]^3 φ'[τ]))/(r[τ]^3)),
  426. θ'[0]==duwf,
  427. θ[0]==θ0,
  428.  
  429. φ''[τ]==-((2 (r'[τ] (-a q ℧+a r[τ] t'[τ]+((-2+r[τ]) r[τ]^3-2 a^2 Sin[θ[τ]]^2) φ'[τ])+
  430. r[τ] θ'[τ] (-2 a Cot[θ[τ]] (-2+r[τ]) t'[τ]+(Cot[θ[τ]] (-2+r[τ]) r[τ]^3+
  431. 2 a^2 Sin[2 θ[τ]]) φ'[τ])))/((-2+r[τ]) r[τ]^4+4 a^2 r[τ] Sin[θ[τ]]^2)),
  432. φ'[0]==dfwf,
  433. φ[0]==φ0,
  434.  
  435. str'[τ]==vd[τ],
  436. str[0]==0,
  437. vlt'[τ]==vd[τ],
  438. vlt[0]==0
  439.  
  440. };
  441.  
  442. DG5={
  443.  
  444. t''[τ]==0,
  445.  
  446. t'[0]==1,
  447. t[0]==0,
  448.  
  449. r''[τ]==-(1-℧ q)/r[τ]^2+(Sqrt[vφ0^2+vθ0^2] r0)^2/r[τ]^3,
  450.  
  451. r'[0]==vr0,
  452. r[0]==r0,
  453.  
  454. θ''[τ]==-2 θ'[τ] r'[τ]/r[τ]+Sin[θ[τ]] Cos[θ[τ]] φ'[τ]^2,
  455.  
  456. θ'[0]==vθ0/r0,
  457. θ[0]==θ0,
  458.  
  459. φ''[τ]==-2 φ'[τ] (r'[τ]+r[τ] θ'[τ] Cot[θ[τ]])/r[τ],
  460.  
  461. φ'[0]==vφ0/r0 Csc[θ0],
  462. φ[0]==φ0,
  463.  
  464. vlt'[τ]==Sqrt[r'[τ]^2+φ'[τ]^2 r[τ]^2],
  465. vlt[0]==0,
  466. str'[τ]==Sqrt[r'[τ]^2+φ'[τ]^2 r[τ]^2],
  467. str[0]==0
  468.  
  469. };
  470.  
  471. DGL = If[dgl==1, DG1, If[dgl==2, DG2, If[dgl==3, DG3, If[dgl==4, DG4, DG5]]]];
  472.  
  473. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  474. (* |||||||| 8) INTEGRATION |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  475. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  476.  
  477. sol=NDSolve[DGL, {t, r, θ, φ, vlt, str}, {τ, 0, tmax+1/1000},
  478. WorkingPrecision-> wp,
  479. MaxSteps-> Infinity,
  480. Method-> mta,
  481. InterpolationOrder-> All,
  482. StepMonitor :> (laststep=plunge; plunge=τ;
  483. stepsize=plunge-laststep;), Method->{"EventLocator",
  484. "Event" :> (If[stepsize<1*^-4, 0, 1])}];
  485.  
  486. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  487. (* |||||||| 9) KOORDINATEN |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  488. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  489.  
  490. X[τ_]:=Evaluate[Sqrt[r[τ]^2+a^2] Sin[θ[τ]] Cos[φ[τ]]/.sol][[1]]; (* kartesisch *)
  491. Y[τ_]:=Evaluate[Sqrt[r[τ]^2+a^2] Sin[θ[τ]] Sin[φ[τ]]/.sol][[1]];
  492. Z[τ_]:=Evaluate[r[τ] Cos[θ[τ]]/.sol][[1]];
  493.  
  494. x[τ_]:=Evaluate[Sqrt[r[τ]^2+A^2] Sin[θ[τ]] Cos[φ[τ]]/.sol][[1]]; (* Plotkoordinaten *)
  495. y[τ_]:=Evaluate[Sqrt[r[τ]^2+A^2] Sin[θ[τ]] Sin[φ[τ]]/.sol][[1]];
  496. z[τ_]:=Z[τ];
  497.  
  498. XYZ[τ_]:=Sqrt[X[τ]^2+Y[τ]^2+Z[τ]^2]; XY[τ_]:=Sqrt[X[τ]^2+Y[τ]^2]; (* kartesischer Radius *)
  499.  
  500. Xyz[{x_, y_, z_}, α_]:={x Cos[α]-y Sin[α], x Sin[α]+y Cos[α], z}; (* Rotationsmatrix *)
  501. xYz[{x_, y_, z_}, β_]:={x Cos[β]+z Sin[β], y, z Cos[β]-x Sin[β]};
  502. xyZ[{x_, y_, z_}, ψ_]:={x, y Cos[ψ]-z Sin[ψ], y Sin[ψ]+z Cos[ψ]};
  503.  
  504. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  505. (* |||||||| 10) PLOT EINSTELLUNGEN |||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  506. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  507.  
  508. PR=r1; (* Plot Range *)
  509. VP={r0, r0, r0}; (* Perspektive x,y,z *)
  510. d1=10; (* Schweiflänge *)
  511. plp=50; (* Flächenplot Details *)
  512. Plp=Automatic; (* Kurven Details *)
  513.  
  514. w1l=0; w2l=0; w1r=0; w2r=0; (* Startperspektiven *)
  515. Mrec=100; mrec=10; (* Parametric Plot Subdivisionen *)
  516. imgsize=380; (* Bildgröße *)
  517.  
  518. s[text_]:=Style[text, FontFamily->"Consolas", FontSize->11]; (* Anzeigestil *)
  519.  
  520. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  521. (* |||||||| 11) PLOT NACH EIGENZEIT ||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  522. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  523.  
  524. Plot[R[tt], {tt, 0, plunge},
  525. Frame->True, PlotStyle->Red, AspectRatio->1/5, ImagePadding->{{40, 10}, {20, 10}},
  526. ImageSize->600, PlotRange->{{0, plunge}, All}, GridLines->{{}, {rA, rI}},
  527. PlotLabel -> "r(τ)"]
  528.  
  529. Plot[Mod[180/Pi Θ[tt], 360], {tt, 0, plunge},
  530. Frame->True, PlotStyle->Cyan, AspectRatio->1/5, ImagePadding->{{40, 10}, {20, 10}},
  531. ImageSize->600, PlotRange->{{0, plunge}, {0, 360}}, GridLines->{{}, {90, 180, 270}},
  532. PlotLabel -> "θ(τ)"]
  533.  
  534. Plot[Mod[180/Pi Φ[tt], 360], {tt, 0, plunge},
  535. Frame->True, PlotStyle->Magenta, AspectRatio->1/5, ImagePadding->{{40, 10}, {20, 10}},
  536. ImageSize->600, PlotRange->{{0, plunge}, {0, 360}}, GridLines->{{}, {90, 180, 270}},
  537. PlotLabel -> "φ(τ)"]
  538.  
  539. Plot[v[tt], {tt, 0, plunge},
  540. Frame->True, PlotStyle->Orange, AspectRatio->1/5, ImagePadding->{{40, 10}, {20, 10}},
  541. ImageSize->600, PlotRange->{{0, plunge}, All}, GridLines->{{}, {0, 1}},
  542. PlotLabel -> "v(τ)"]
  543.  
  544. displayP[T_]:=Grid[{
  545. {If[μ==0, s[" affineP"], s[" τ propr"]], " = ", s[n0[tp]], s["GM/c³"], s[dp]},
  546. {s[" t coord"], " = ", s[n0[т[tp]]], s["GM/c³"], s[dp]},
  547. {s[" ṫ total"], " = ", s[n0[ю[tp]]], s["dt/dτ"], s[dp]},
  548. {s[" ς gravt"], " = ", s[n0[If[dgl==5, 1, ς[tp]]]], s["dt/dτ"], s[dp]},
  549. {s[" γ kinet"], " = ", If[μ==0, s[n0[ς[tp]]], s[n0[If[dgl==5,
  550. vnt[tp]^2/2, 1/Sqrt[1-v[tp]^2]]]]], s["dt/dτ"], s[dp]},
  551. {s[" R carts"], " = ", s[n0[XYZ[tp]]], s["GM/c²"], s[dp]},
  552. {s[" x carts"], " = ", s[n0[X[tp]]], s["GM/c²"], s[dp]},
  553. {s[" y carts"], " = ", s[n0[Y[tp]]], s["GM/c²"], s[dp]},
  554. {s[" z carts"], " = ", s[n0[Z[tp]]], s["GM/c²"], s[dp]},
  555. {s[" s dstnc"], " = ", s[n0[dst[tp]]], s["GM/c²"], s[dp]},
  556.  
  557. {s[" r coord"], " = ", s[n0[R[tp]]], s["GM/c²"], s[dp]},
  558. {s[" φ longd"], " = ", s[n0[Φ[tp] 180/π]], s["deg"], s[dp]},
  559. {s[" θ lattd"], " = ", s[n0[Θ[tp] 180/π]], s["deg"], s[dp]},
  560. {s[" d¹r/dτ¹"], " = ", s[n0[R'[tp]]], s["c"], s[dp]},
  561. {s[" d¹φ/dτ¹"], " = ", s[n0[Φ'[tp]]], s["c\.b3/G/M"], s[dp]},
  562. {s[" d¹θ/dτ¹"], " = ", s[n0[Θ'[tp]]], s["c\.b3/G/M"], s[dp]},
  563. {s[" d\.b2r/dτ\.b2"], " = ", s[n0[R''[tp]]], s["c⁴/G/M"], s[dp]},
  564. {s[" d\.b2φ/dτ\.b2"], " = ", s[n0[Φ''[tp]]], s["c⁶/G\.b2/M\.b2"], s[dp]},
  565. {s[" d\.b2θ/dτ\.b2"], " = ", s[n0[Θ''[tp]]], s["c⁶/G\.b2/M\.b2"], s[dp]},
  566. {s[" a SpinP"], " = ", s[n0[a]], s["GM²/c"], s[dp]},
  567.  
  568. {s[" ℧ cntrl"], " = ", s[n0[℧]], s["Q/M"], s[dp]},
  569. {s[" q prtcl"], " = ", s[n0[q]], s["q/m"], s[dp]},
  570. {s[" M irred"], " = ", s[N[mirr]], s["M"], s[dp]},
  571. {s[" E kinet"], " = ", s[n0[If[dgl==5, vnt[tp]^2/2, ekin[tp]]]], s["mc²"], s[dp]},
  572. {s[" E poten"], " = ", s[n0[If[dgl==5, -(1-℧ q)/R[tp], epot[tp]]]], s["mc²"], s[dp]},
  573. {s[" E total"], " = ", s[n0[If[dgl==5, v0^2/2-(1-℧ q)/r0+1, ε]]], s["mc²"], s[dp]},
  574. {s[" CarterQ"], " = ", s[n0[Qk]], s["GMm/c"], s[dp]},
  575. {s[" L axial"], " = ", s[n0[If[dgl==5, vφ0 r0, Lz]]], s["GMm/c"], s[dp]},
  576. {s[" L polar"], " = ", s[n0[If[dgl==5, vθ0 r0, pΘ[tp]]]], s["GMm/c"], s[dp]},
  577. {s[" g acclr"], " = ", s[n0[v''[tp]]], s["c⁴/G/M"], s[dp]},
  578.  
  579. {s[" ω fdrag"], " = ", s[n0[Abs[ω[tp]]]], s["c³/G/M"], s[dp]},
  580. {s[" v fdrag"], " = ", s[n0[Abs[й[tp]]]], s["c"], s[dp]},
  581. {s[" Ω fdrag"], " = ", s[n0[Abs[Ω[tp]]]], s["c"], s[dp]},
  582. {s[" v propr"], " = ", s[n0[If[dgl==5, vnt[tp], v[tp]/Sqrt[1-v[tp]^2]]]], s["c"], s[dp]},
  583. {s[" v escpe"], " = ", s[n0[vesc[tp]]], s["c"], s[dp]},
  584. {s[" v obsvd"], " = ", s[n0[ß[tp]]], s["c"], s[dp]},
  585. {s[" v r,loc"], " = ", s[n0[If[dgl==5, R'[tp], vr[tp]]]], s["c"], s[dp]},
  586. {s[" v θ,loc"], " = ", s[n0[If[dgl==5, Θ'[tp] R[tp], vθ[tp]]]], s["c"], s[dp]},
  587. {s[" v φ,loc"], " = ", s[n0[If[dgl==5, Φ'[tp] R[tp]/Csc[Θ[tp]], vφ[tp]]]], s["c"], s[dp]},
  588. {s[" v local"], " = ", s[n0[If[dgl==5, vnt[tp], v[tp]]]], s["c"], s[dp]},
  589. {s[" "], s[" "], s[" "], s[" "]}},
  590. Alignment-> Left, Spacings-> {0, 0}];
  591.  
  592. plot1b[{xx_, yy_, zz_, tk_, w1_, w2_}]:= (* Animation *)
  593. Show[
  594.  
  595. Graphics3D[{
  596. {PointSize[0.011], Red, Point[
  597. Xyz[xyZ[{x[tp], y[tp], z[tp]}, w1], w2]]}},
  598. ImageSize-> imgsize,
  599. PlotRange-> {
  600. {-(2 Sign[Abs[xx]]+1) PR, +(2 Sign[Abs[xx]]+1) PR},
  601. {-(2 Sign[Abs[yy]]+1) PR, +(2 Sign[Abs[yy]]+1) PR},
  602. {-(2 Sign[Abs[zz]]+1) PR, +(2 Sign[Abs[zz]]+1) PR}
  603. },
  604. SphericalRegion->False,
  605. ImagePadding-> 1],
  606.  
  607. horizons[A, None, w1, w2],
  608.  
  609. If[a==0, {}, ParametricPlot3D[
  610. Xyz[xyZ[{
  611. Sin[prm] a,
  612. Cos[prm] a,
  613. 0}, w1], w2],
  614. {prm, 0, 2π},
  615. PlotStyle -> {Thickness[0.005], Orange}]],
  616.  
  617. If[a==0, {},
  618. Graphics3D[{{PointSize[0.009], Purple, Point[
  619. Xyz[xyZ[{
  620. Sin[-φ0-ω0 т[tp]+π/2] Sqrt[x0[A]^2+y0[A]^2],
  621. Cos[-φ0-ω0 т[tp]+π/2] Sqrt[x0[A]^2+y0[A]^2],
  622. z0[A]}, w1], w2]]}}]],
  623.  
  624. If[tk==0, {}, If[a==0, {},
  625. ParametricPlot3D[
  626. Xyz[xyZ[{
  627. Sin[-φ0-ω0 т[tt]+π/2] Sqrt[x0[A]^2+y0[A]^2],
  628. Cos[-φ0-ω0 т[tt]+π/2] Sqrt[x0[A]^2+y0[A]^2],
  629. z0[A]}, w1], w2],
  630. {tt, Max[0, д[т[tp]-1/2 π/ω0]], tp},
  631. PlotStyle -> {Thickness[0.001], Dashed, Purple},
  632. PlotPoints-> Automatic,
  633. MaxRecursion-> 12]]],
  634.  
  635. If[tk==0, {},
  636. Block[{$RecursionLimit = Mrec},
  637. ParametricPlot3D[
  638. Xyz[xyZ[{x[tt], y[tt], z[tt]}, w1], w2], {tt, If[tp<0, Min[0, tp+d1], Max[0, tp-d1]], tp},
  639. PlotStyle-> {Thickness[0.004]},
  640. ColorFunction-> Function[{x, y, z, t},
  641. Hue[0, 1, 0.5, If[tp<0,
  642. Max[Min[(+tp+(-t+d1))/d1, 1], 0], Max[Min[(-tp+(t+d1))/d1, 1], 0]]]],
  643. ColorFunctionScaling-> False,
  644. PlotPoints-> Automatic,
  645. MaxRecursion-> mrec]]],
  646.  
  647. If[tk==0, {},
  648. Block[{$RecursionLimit = Mrec},
  649. ParametricPlot3D[
  650. Xyz[xyZ[{x[tt], y[tt], z[tt]}, w1], w2],
  651. {tt, 0, If[tp<0, Min[-1*^-16, tp+d1/3], Max[1*^-16, tp-d1/3]]},
  652. PlotStyle-> {Thickness[0.004], Opacity[0.6], Darker[Gray]},
  653. PlotPoints-> Plp,
  654. MaxRecursion-> mrec]]],
  655.  
  656. ViewPoint-> {xx, yy, zz}];
  657.  
  658. Do[
  659. Print[Rasterize[Grid[{{
  660. plot1b[{0, -Infinity, 0, tp, w1l, w2l}],
  661. plot1b[{0, 0, +Infinity, tp, w1r, w2r}],
  662. displayP[tp]
  663. }, {" ", " ", " "}
  664. }, Alignment->Left]]],
  665. {tp, 0, tMax, tMax/1}]
  666.  
  667. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  668. (* |||||||| 12) PLOT NACH KOORDINATENZEIT ||||||||||||||||||||||||||||||||||||||||||||||| *)
  669. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  670.  
  671. Plot[R[д[tt]], {tt, 0, TMax},
  672. Frame->True, PlotStyle->Red, AspectRatio->1/5, ImagePadding->{{40, 10}, {20, 10}},
  673. ImageSize->600, PlotRange->{{0, TMax}, All}, GridLines->{{}, {rA, rI}},
  674. PlotLabel -> "r(t)"]
  675.  
  676. Plot[Mod[180/Pi Θ[д[tt]], 360], {tt, 0, TMax},
  677. Frame->True, PlotStyle->Cyan, AspectRatio->1/5, ImagePadding->{{40, 10}, {20, 10}},
  678. ImageSize->600, PlotRange->{{0, TMax}, {0, 360}}, GridLines->{{}, {90, 180, 270}},
  679. PlotLabel -> "θ(t)"]
  680.  
  681. Plot[Mod[180/Pi Φ[д[tt]], 360], {tt, 0, TMax},
  682. Frame->True, PlotStyle->Magenta, AspectRatio->1/5, ImagePadding->{{40, 10}, {20, 10}},
  683. ImageSize->600, PlotRange->{{0, TMax}, {0, 360}}, GridLines->{{}, {90, 180, 270}},
  684. PlotLabel -> "φ(t)"]
  685.  
  686. Plot[v[д[tt]], {tt, 0, TMax},
  687. Frame->True, PlotStyle->Orange, AspectRatio->1/5, ImagePadding->{{40, 10}, {20, 10}},
  688. ImageSize->600, PlotRange->{{0, TMax}, All}, GridLines->{{}, {0, 1}},
  689. PlotLabel -> "v(t)"]
  690.  
  691. displayC[T_]:=Grid[{
  692. {s[" t coord"], " = ", s[n0[tk]], s["GM/c³"], s[dp]},
  693. {If[μ==0, s[" affineP"], s[" τ propr"]], " = ", s[n0[T]], s["GM/c³"], s[dp]},
  694. {s[" ṫ total"], " = ", s[n0[ю[T]]], s["dt/dτ"], s[dp]},
  695. {s[" ς gravt"], " = ", s[n0[If[dgl==5, 1, ς[T]]]], s["dt/dτ"], s[dp]},
  696. {s[" γ kinet"], " = ", If[μ==0, s[n0[ς[T]]], s[n0[If[dgl==5,
  697. vnt[T]^2/2, 1/Sqrt[1-v[T]^2]]]]], s["dt/dτ"], s[dp]},
  698. {s[" R carts"], " = ", s[n0[XYZ[T]]], s["GM/c²"], s[dp]},
  699. {s[" x carts"], " = ", s[n0[X[T]]], s["GM/c²"], s[dp]},
  700. {s[" y carts"], " = ", s[n0[Y[T]]], s["GM/c²"], s[dp]},
  701. {s[" z carts"], " = ", s[n0[Z[T]]], s["GM/c²"], s[dp]},
  702. {s[" s dstnc"], " = ", s[n0[dst[T]]], s["GM/c²"], s[dp]},
  703.  
  704. {s[" r coord"], " = ", s[n0[R[T]]], s["GM/c²"], s[dp]},
  705. {s[" φ longd"], " = ", s[n0[Φ[T] 180/π]], s["deg"], s[dp]},
  706. {s[" θ lattd"], " = ", s[n0[Θ[T] 180/π]], s["deg"], s[dp]},
  707. {s[" d¹r/dτ¹"], " = ", s[n0[R'[T]]], s["c"], s[dp]},
  708. {s[" d¹φ/dτ¹"], " = ", s[n0[Φ'[T]]], s["c\.b3/G/M"], s[dp]},
  709. {s[" d¹θ/dτ¹"], " = ", s[n0[Θ'[T]]], s["c\.b3/G/M"], s[dp]},
  710. {s[" d\.b2r/dτ\.b2"], " = ", s[n0[R''[T]]], s["c⁴/G/M"], s[dp]},
  711. {s[" d\.b2φ/dτ\.b2"], " = ", s[n0[Φ''[T]]], s["c⁶/G\.b2/M\.b2"], s[dp]},
  712. {s[" d\.b2θ/dτ\.b2"], " = ", s[n0[Θ''[T]]], s["c⁶/G\.b2/M\.b2"], s[dp]},
  713. {s[" a SpinP"], " = ", s[n0[a]], s["GM²/c"], s[dp]},
  714.  
  715. {s[" ℧ cntrl"], " = ", s[n0[℧]], s["Q/M"], s[dp]},
  716. {s[" q prtcl"], " = ", s[n0[q]], s["q/m"], s[dp]},
  717. {s[" M irred"], " = ", s[N[mirr]], s["M"], s[dp]},
  718. {s[" E kinet"], " = ", s[n0[If[dgl==5, vnt[T]^2/2, ekin[T]]]], s["mc²"], s[dp]},
  719. {s[" E poten"], " = ", s[n0[If[dgl==5, -(1-℧ q)/R[T], epot[T]]]], s["mc²"], s[dp]},
  720. {s[" E total"], " = ", s[n0[If[dgl==5, v0^2/2-(1-℧ q)/r0+1, ε]]], s["mc²"], s[dp]},
  721. {s[" CarterQ"], " = ", s[n0[Qk]], s["GMm/c"], s[dp]},
  722. {s[" L axial"], " = ", s[n0[If[dgl==5, vφ0 r0, Lz]]], s["GMm/c"], s[dp]},
  723. {s[" L polar"], " = ", s[n0[If[dgl==5, vθ0 r0, pΘ[T]]]], s["GMm/c"], s[dp]},
  724. {s[" g acclr"], " = ", s[n0[v''[T]]], s["c⁴/G/M"], s[dp]},
  725.  
  726. {s[" ω fdrag"], " = ", s[n0[Abs[ω[T]]]], s["c³/G/M"], s[dp]},
  727. {s[" v fdrag"], " = ", s[n0[Abs[й[T]]]], s["c"], s[dp]},
  728. {s[" Ω fdrag"], " = ", s[n0[Abs[Ω[T]]]], s["c"], s[dp]},
  729. {s[" v propr"], " = ", s[n0[If[dgl==5, vnt[T], v[T]/Sqrt[1-v[T]^2]]]], s["c"], s[dp]},
  730. {s[" v escpe"], " = ", s[n0[vesc[T]]], s["c"], s[dp]},
  731. {s[" v obsvd"], " = ", s[n0[ß[T]]], s["c"], s[dp]},
  732. {s[" v r,loc"], " = ", s[n0[If[dgl==5, R'[T], vr[T]]]], s["c"], s[dp]},
  733. {s[" v θ,loc"], " = ", s[n0[If[dgl==5, Θ'[T] R[T], vθ[T]]]], s["c"], s[dp]},
  734. {s[" v φ,loc"], " = ", s[n0[If[dgl==5, Φ'[T] R[T]/Csc[Θ[T]], vφ[T]]]], s["c"], s[dp]},
  735. {s[" v local"], " = ", s[n0[If[dgl==5, vnt[T], v[T]]]], s["c"], s[dp]},
  736. {s[" "], s[" "], s[" "], s[" "]}},
  737. Alignment-> Left, Spacings-> {0, 0}];
  738.  
  739. plot1a[{xx_, yy_, zz_, tk_, w1_, w2_}]:= (* Animation *)
  740. Show[
  741.  
  742. Graphics3D[{
  743. {PointSize[0.011], Red, Point[
  744. Xyz[xyZ[{x[T], y[T], z[T]}, w1], w2]]}},
  745. ImageSize-> imgsize,
  746. PlotRange-> {
  747. {-(2 Sign[Abs[xx]]+1) PR, +(2 Sign[Abs[xx]]+1) PR},
  748. {-(2 Sign[Abs[yy]]+1) PR, +(2 Sign[Abs[yy]]+1) PR},
  749. {-(2 Sign[Abs[zz]]+1) PR, +(2 Sign[Abs[zz]]+1) PR}
  750. },
  751. SphericalRegion->False,
  752. ImagePadding-> 1],
  753.  
  754. horizons[A, None, w1, w2],
  755.  
  756. If[a==0, {}, ParametricPlot3D[
  757. Xyz[xyZ[{
  758. Sin[prm] a,
  759. Cos[prm] a,
  760. 0}, w1], w2],
  761. {prm, 0, 2π},
  762. PlotStyle -> {Thickness[0.005], Orange}]],
  763.  
  764. If[a==0, {},
  765. Graphics3D[{{PointSize[0.009], Purple, Point[
  766. Xyz[xyZ[{
  767. Sin[-φ0-ω0 tk+π/2] Sqrt[x0[A]^2+y0[A]^2],
  768. Cos[-φ0-ω0 tk+π/2] Sqrt[x0[A]^2+y0[A]^2],
  769. z0[A]}, w1], w2]]}}]],
  770.  
  771. If[tk==0, {}, If[a==0, {},
  772. ParametricPlot3D[
  773. Xyz[xyZ[{
  774. Sin[-φ0-ω0 tt+π/2] Sqrt[x0[A]^2+y0[A]^2],
  775. Cos[-φ0-ω0 tt+π/2] Sqrt[x0[A]^2+y0[A]^2],
  776. z0[A]}, w1], w2],
  777. {tt, Max[0, tk-1/2 π/ω0], tk},
  778. PlotStyle -> {Thickness[0.001], Dashed, Purple},
  779. PlotPoints-> Automatic,
  780. MaxRecursion-> mrec]]],
  781.  
  782. Block[{$RecursionLimit = Mrec},
  783. If[tk==0, {},
  784. ParametricPlot3D[
  785. Xyz[xyZ[{x[tt], y[tt], z[tt]}, w1], w2], {tt, If[TMax<0, Min[0, T+d1], Max[0, T-d1]], T},
  786. PlotStyle-> {Thickness[0.004]},
  787. ColorFunction-> Function[{x, y, z, t},
  788. Hue[0, 1, 0.5, If[TMax<0, Max[Min[(+T+(-t+d1))/d1, 1], 0], Max[Min[(-T+(t+d1))/d1, 1], 0]]]],
  789. ColorFunctionScaling-> False,
  790. PlotPoints-> Automatic,
  791. MaxRecursion-> mrec]]],
  792.  
  793. If[tk==0, {},
  794. Block[{$RecursionLimit = Mrec},
  795. ParametricPlot3D[
  796. Xyz[xyZ[{x[tt], y[tt], z[tt]}, w1], w2],
  797. {tt, 0, If[Tmax<0, Min[-1*^-16, T+d1/3], Max[1*^-16, T-d1/3]]},
  798. PlotStyle-> {Thickness[0.004], Opacity[0.6], Darker[Gray]},
  799. PlotPoints-> Plp,
  800. MaxRecursion-> mrec]]],
  801.  
  802. ViewPoint-> {xx, yy, zz}];
  803.  
  804. Quiet[Do[
  805. Print[Rasterize[Grid[{{
  806. plot1a[{0, -Infinity, 0, tk, w1l, w2l}],
  807. plot1a[{0, 0, Infinity, tk, w1r, w2r}],
  808. displayC[Quiet[д[tk]]]
  809. }, {" ", " ", " "}
  810. }, Alignment->Left]]],
  811. {tk, 0, TMax, TMax/1}]]
  812.  
  813. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  814. (* |||||||| 13) EXPORTOPTIONEN |||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  815. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  816.  
  817. (* Export als HTML Dokument *)
  818. (* Export["Y:\\export\\dateiname.html", EvaluationNotebook[], "GraphicsOutput" -> "PNG"] *)
  819. (* Export direkt als Bildsequenz *)
  820. (* Do[Export["Y:\\export\\dateiname" <> ToString[tk] <> ".png", Rasterize[...] *)
  821.  
  822. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  823. (* |||||||||||| http://kerr.newman.yukerez.net ||||| Simon Tyran, Vienna |||||||||||||||| *)
  824. (* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
RAW Paste Data