Yukterez

Gravity of a Planet with a Disk or Ring, Vector Field

Apr 28th, 2020
34
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. (* ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  2. (* ||||| Mathematica Syntax | yukterez.net | Planet, Ring, Disk Fieldlines | Version 2 ||||||| *)
  3. (* ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
  4.  
  5. G = 1; (* Gravitationskonstante *)
  6. M = 1; (* Masse der Scheibe *)
  7. Ḿ = 1; (* Masse der zentralen Kugel *)
  8.  
  9. я1 = 15; (* Scheibeninnenradius *)
  10. я2 = 20; (* Scheibenaußenradius *)
  11. я3 = 10; (* Kugelradius *)
  12.  
  13. PR = 25; (* Plot Range *)
  14. IS = 400; (* Bildgröße *)
  15.  
  16. m[x_, y_, z_] := If[я3==0, Ḿ, Ḿ Sqrt[x^2+y^2+z^2]^3/я3^3] (* innere Kugelrestmasse *)
  17. ρ->M/(π я2^2-π я1^2)
  18. ε = 1/1000;
  19.  
  20. g[{x_, y_, z_}]:=-{
  21.  
  22. If[Abs[x]<ε, 0, ((2 Sqrt[я2] G M x (-EllipticE[(4 я2 Sqrt[x^2+y^2])/(я2^2+x^2+y^2+
  23. 2 я2 Sqrt[x^2+y^2]+z^2)]+(1-(2 я2 Sqrt[x^2+y^2])/(я2^2+x^2+y^2+2 я2 Sqrt[x^2+y^2]+
  24. z^2)) EllipticK[(4 я2 Sqrt[x^2+y^2])/(я2^2+x^2+y^2+2 я2 Sqrt[x^2+y^2]+z^2)]))/((-я1^2 π+
  25. я2^2 π) (x^2+y^2)^(3/4) Sqrt[(я2 Sqrt[x^2+y^2])/(я2^2+x^2+y^2+2 я2 Sqrt[x^2+y^2]+z^2)]))-
  26. ((2 Sqrt[я1] G M x (-EllipticE[(4 я1 Sqrt[x^2+y^2])/(я1^2+x^2+y^2+
  27. 2 я1 Sqrt[x^2+y^2]+z^2)]+(1-(2 я1 Sqrt[x^2+y^2])/(я1^2+x^2+y^2+2 я1 Sqrt[x^2+y^2]+
  28. z^2)) EllipticK[(4 я1 Sqrt[x^2+y^2])/(я1^2+x^2+y^2+2 я1 Sqrt[x^2+y^2]+z^2)]))/((-я1^2 π+
  29. я2^2 π) (x^2+y^2)^(3/4) Sqrt[(я1 Sqrt[x^2+y^2])/(я1^2+x^2+y^2+2 я1 Sqrt[x^2+y^2]+z^2)]))+
  30. If[Ḿ==0, 0, G Min[m[x, y, z], Ḿ] x/Sqrt[(x^2+y^2+z^2)^3]]], (* x Beschleunigung *)
  31.  
  32. If[Abs[y]<ε, 0, ((2 Sqrt[я2] G M y (-EllipticE[(4 я2 Sqrt[x^2+y^2])/(я2^2+x^2+y^2+
  33. 2 я2 Sqrt[x^2+y^2]+z^2)]+(1-(2 я2 Sqrt[x^2+y^2])/(я2^2+x^2+y^2+2 я2 Sqrt[x^2+y^2]+
  34. z^2)) EllipticK[(4 я2 Sqrt[x^2+y^2])/(я2^2+x^2+y^2+2 я2 Sqrt[x^2+y^2]+z^2)]))/((-я1^2 π+
  35. я2^2 π) (x^2+y^2)^(3/4) Sqrt[(я2 Sqrt[x^2+y^2])/(я2^2+x^2+y^2+2 я2 Sqrt[x^2+y^2]+z^2)]))-
  36. ((2 Sqrt[я1] G M y (-EllipticE[(4 я1 Sqrt[x^2+y^2])/(я1^2+x^2+y^2+2 я1 Sqrt[x^2+
  37. y^2]+z^2)]+(1-(2 я1 Sqrt[x^2+y^2])/(я1^2+x^2+y^2+2 я1 Sqrt[x^2+y^2]+
  38. z^2)) EllipticK[(4 я1 Sqrt[x^2+y^2])/(я1^2+x^2+y^2+2 я1 Sqrt[x^2+y^2]+z^2)]))/((-я1^2 π+
  39. я2^2 π) (x^2+y^2)^(3/4) Sqrt[(я1 Sqrt[x^2+y^2])/(я1^2+x^2+y^2+2 я1 Sqrt[x^2+y^2]+z^2)]))+
  40. If[Ḿ==0, 0, G Min[m[x, y, z], Ḿ] y/Sqrt[(x^2+y^2+z^2)^3]]], (* y Beschleunigung *)
  41.  
  42. If[Abs[z]<ε, 0, -1/((-я1^2 π+я2^2 π) Abs[z]) 2 G M z ((-π+((я2+Sqrt[x^2+y^2+z^2]) ((Sqrt[x^2+
  43. y^2]+Sqrt[x^2+y^2+z^2])/(-Sqrt[x^2+y^2]+Sqrt[x^2+y^2+z^2]))^(1/2) EllipticPi[-((2 Sqrt[x^2+
  44. y^2])/(-Sqrt[x^2+y^2]+Sqrt[x^2+y^2+z^2])), (4 я2 Sqrt[x^2+y^2])/(я2^2+x^2+y^2+2 я2 Sqrt[x^2+
  45. y^2]+z^2)]+(-я2+Sqrt[x^2+y^2+z^2]) Sqrt[(-Sqrt[x^2+y^2]+Sqrt[x^2+y^2+z^2])/(Sqrt[x^2+y^2]+
  46. Sqrt[x^2+y^2+z^2])] EllipticPi[(2 Sqrt[x^2+y^2])/(Sqrt[x^2+y^2]+Sqrt[x^2+y^2+z^2]),
  47. (4 я2 Sqrt[x^2+y^2])/(я2^2+x^2+y^2+2 я2 Sqrt[x^2+y^2]+z^2)])/(Sqrt[я2^2+x^2+y^2+2 я2 Sqrt[x^2+
  48. y^2]+z^2]))-(-π+((я1+Sqrt[x^2+y^2+z^2]) Sqrt[(Sqrt[x^2+y^2]+Sqrt[x^2+y^2+
  49. z^2])/(-Sqrt[x^2+y^2]+Sqrt[x^2+y^2+z^2])] EllipticPi[-((2 Sqrt[x^2+y^2])/(-Sqrt[x^2+y^2]+
  50. Sqrt[x^2+y^2+z^2])), (4 я1 Sqrt[x^2+y^2])/(я1^2+x^2+y^2+2 я1 Sqrt[x^2+y^2]+z^2)]+(-я1+Sqrt[x^2+
  51. y^2+z^2]) Sqrt[(-Sqrt[x^2+y^2]+Sqrt[x^2+y^2+z^2])/(Sqrt[x^2+y^2]+Sqrt[x^2+y^2+
  52. z^2])] EllipticPi[(2 Sqrt[x^2+y^2])/(Sqrt[x^2+y^2]+Sqrt[x^2+y^2+z^2]), (4 я1 Sqrt[x^2+
  53. y^2])/(я1^2+x^2+y^2+2 я1 Sqrt[x^2+y^2]+z^2)])/(Sqrt[я1^2+x^2+y^2+2 я1 Sqrt[x^2+y^2]+z^2])))+
  54. If[Ḿ==0, 0, G Min[m[x, y, z], Ḿ] z/Sqrt[(x^2+y^2+z^2)^3]]]}; (* z Beschleunigung *)
  55.  
  56. annulus3dF[color_: GrayLevel[2/5], o : OptionsPattern[]]:= (* Scheibe *)
  57. Graphics3D[{Glow[color], Opacity[0.3], EdgeForm[None], Black,
  58. ChartElementData["CylindricalSector3D"][{##}, 1]}, o,
  59. Boxed->False] &;
  60.  
  61. Show[ (* 3D Vektorplot *)
  62. VectorPlot3D[{g[{x, y, z}][[1]], g[{x, y, z}][[2]], g[{x, y, z}][[3]]},
  63. {x, -PR, PR}, {y, -PR, PR}, {z, -PR, PR},
  64. ImageSize->IS, VectorPoints->Fine],
  65. Graphics3D[{Glow[GrayLevel[2/5]], Black, Opacity[0.3], Sphere[{0, 0, 0}, я3]},
  66. PlotRange->PR,
  67. SphericalRegion->False,
  68. ImagePadding->1],
  69. annulus3dF[][{0, 2 π}, {я1, я2}, {-PR/150, PR/150}]]
  70.  
  71. vcp1 = Show[ (* x,z-Ebene *)
  72. VectorPlot[{g[{x, ε, z}][[1]], g[{x, ε, z}][[3]]}, {x, -PR, PR}, {z, -PR, PR},
  73. ImageSize->IS, VectorPoints->31, VectorScale->0.05, PlotRange->PR],
  74. Graphics[{Glow[GrayLevel[2/5]], Opacity[0.3], Disk[{0, 0}, я3]}],
  75. Graphics[{Glow[GrayLevel[2/5]], Opacity[0.3], Thickness[0.01], Line[{{-я2, 0}, {я2, 0}}]}]];
  76.  
  77. vcp2 = Show[ (* x,y-Ebene *)
  78. VectorPlot[{g[{x, y, ε}][[1]], g[{x, y, ε}][[2]]}, {x, -PR, PR}, {y, -PR, PR},
  79. ImageSize->IS, VectorPoints->31, VectorScale->0.05, PlotRange->PR],
  80. Graphics[{Glow[GrayLevel[2/5]], Opacity[0.3], Disk[{0, 0}, я3]}],
  81. Graphics[{Glow[GrayLevel[2/5]], Opacity[0.3],
  82. If[я1==0, Disk[{0, 0}, я2], Annulus[{0, 0}, {я1, я2}]]}]];
  83.  
  84. Grid[{{vcp1, vcp2}}] (* 2D Vektorplot *)
  85.  
  86. ctp1 = Show[Table[ (* x,z-Ebene *)
  87. ContourPlot[{Norm@g[{x, ε, z}]==n}, {x, -PR, PR}, {z, -PR, PR},
  88. ImageSize->IS, MaxRecursion->8, PlotPoints->50, ClippingStyle->Automatic, PlotRange->PR],
  89. {n, 0.001, 0.01, 0.001}], (* Plot Range und Intervall für die Linien konstanter Gravitation *)
  90. Graphics[{Glow[GrayLevel[2/5]], Opacity[0.3], Disk[{0, 0}, я3]}],
  91. Graphics[{Glow[GrayLevel[2/5]], Opacity[0.3], Thickness[0.01], Line[{{-я2, 0}, {я2, 0}}]}]];
  92.  
  93. ctp2 = Show[Table[ (* x,z-Ebene *)
  94. ContourPlot[{Norm@g[{x, y, ε}]==n}, {x, -PR, PR}, {y, -PR, PR},
  95. ImageSize->IS, MaxRecursion->8, PlotPoints->50, ClippingStyle->Automatic, PlotRange->PR],
  96. {n, 0.001, 0.01, 0.001}], (* Plot Range und Intervall für die Linien konstanter Gravitation *)
  97. Graphics[{Glow[GrayLevel[2/5]], Opacity[0.3], Disk[{0, 0}, я3]}],
  98. Graphics[{Glow[GrayLevel[2/5]], Opacity[0.3],
  99. If[я1==0, Disk[{0, 0}, я2], Annulus[{0, 0}, {я1, я2}]]}]];
  100.  
  101. Grid[{{ctp1, ctp2}}] (* Konturplot *)
RAW Paste Data