n-Body, Gravitation & Ladung

Deutschsprachige Version
Benutzeravatar
Yukterez
Administrator
Beiträge: 194
Registriert: Mi 21. Okt 2015, 02:16

n-Body, Gravitation & Ladung

Beitragvon Yukterez » Di 9. Feb 2016, 07:33

Mehrkörpersimulator, Gravitation und Ladung (yukterez.net/nbody)
Bild

Startpositionen links: M1 @ {½, 0, √¾}, M2 @ {0, 0, 0}, M3 @ {1, 0, 0} (gleichseitiges ◬)
Startpositionen rechts: M1 @ {½, 0, √¾}, M2 @ {0.1, 0, 0}, M3 @ {0.8, 0, 0.2}
Anfangsgeschwindigkeiten: v=0

Bild

Diskussion: Uwudl + Stackexchange
Bild
Simon Tyran aka Симон Тыран @ vk || wikipedia || stackexchange || wolfram

Benutzeravatar
Yukterez
Administrator
Beiträge: 194
Registriert: Mi 21. Okt 2015, 02:16

n-Body, geschlossene Loops

Beitragvon Yukterez » Sa 16. Feb 2019, 06:14

Geschlossene, aber instabile Orbits mit 5 Körpern; Periodendauer: 2π sek, Schweiflänge: 1 sek:

Bild

Die Startbedingungen stammen aus dem Werk von Carles Simó.
Bild
Simon Tyran aka Симон Тыран @ vk || wikipedia || stackexchange || wolfram

Benutzeravatar
Yukterez
Administrator
Beiträge: 194
Registriert: Mi 21. Okt 2015, 02:16

n-Body, Code

Beitragvon Yukterez » Sa 16. Feb 2019, 06:15

Mehrkörpersimulatoren für Masse und Ladung
Bild

2 Körper:

Code: Alles auswählen

(* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
(* ||| Mathematica Syntax || yukterez.net || 2 Body Newtonian Mass & Charge Simulator ||| *)
(* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)

ClearAll["Global`*"]
Amp = 1; kg = 1; m = 1; sek = 1; km = 1000 m; (* SI Einheiten *)
 
mt1 = {"StiffnessSwitching", Method-> {"ExplicitRungeKutta", Automatic}};
mt2 = {"ImplicitRungeKutta", "DifferenceOrder"-> 20};
mt3 = {"EquationSimplification"-> "Residual"};
mt0 = Automatic;
mta = mt2;
wp  = MachinePrecision;
 
(* Plot Optionen *)
 
Tmax      = 24 sek;
trail     = 12 sek;
point     = 0.015;
thk       = 0.004;
plotrange = {{-1.2 m, +1.2 m}, {-1.2 m, +1.2 m}, {-1.2 m, +1.2 m}};
imagesize = 430;
 
(* Masse 1 *)
 
m1  = 1000000000 kg;
q1  = 0;

x1x = 1/2 m;
y1y = 0 m;
z1z = 0 m;

v1x = 0 m/sek;
v1y = 0 m/sek;
v1z = Sqrt[G m2] Sqrt[1/2];
 
(* Masse 2 *)
 
m2  = m1/2;
q2  = 0 Amp sek;
 
x2x = -1/2 m;
y2y = 0 m;
z2z = 0 m;
 
v2x = 0 m/sek;
v2y = 0 m/sek;
v2z = -Sqrt[G m1] Sqrt[1/2];
 
(* Konstanten *)
 
G  = 667384/10^16 m^3/kg/sek^2;
ε0 = 8.854187817*^-12 Amp^2 sek^4/kg/m^3;
Au = 149597870690 m;
m0 = 1988*10^27 kg;
yr = 36525*24*36 sek;
 
(* Differentialgleichung *)
 
nds=NDSolve[{
 
x1'[t] == vx1[t], y1'[t] == vy1[t], z1'[t] == vz1[t],
x2'[t] == vx2[t], y2'[t] == vy2[t], z2'[t] == vz2[t],
 
vx1'[t] ==
 
(G m2 (x2[t]-x1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(-q1*q2/(4Pi ε0 )/m1 (x2[t]-x1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3],
 
 vy1'[t] ==
 
(G m2 (y2[t]-y1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(-q1*q2/(4Pi ε0 )/m1 (y2[t]-y1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3],
 
vz1'[t] ==
 
(G m2 (z2[t]-z1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(-q1*q2/(4Pi ε0 )/m1 (z2[t]-z1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3],
 
vx2'[t] ==
 
(G m1 (x1[t]-x2[t]))/Sqrt[((x1[t]-x2[t])^2+(y1[t]-y2[t])^2+(z1[t]-z2[t])^2)^3]+
(-q2*q1/(4Pi ε0 )/m2 (x1[t]-x2[t]))/Sqrt[((x1[t]-x2[t])^2+(y1[t]-y2[t])^2+(z1[t]-z2[t])^2)^3],
 
 vy2'[t] ==
 
(G m1 (y1[t]-y2[t]))/Sqrt[((x1[t]-x2[t])^2+(y1[t]-y2[t])^2+(z1[t]-z2[t])^2)^3]+
(-q2*q1/(4Pi ε0 )/m2 (y1[t]-y2[t]))/Sqrt[((x1[t]-x2[t])^2+(y1[t]-y2[t])^2+(z1[t]-z2[t])^2)^3],
 
vz2'[t] ==
 
(G m1 (z1[t]-z2[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(-q2*q1/(4Pi ε0 )/m2 (z1[t]-z2[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3],
 
x1[0] == x1x, y1[0] == y1y, z1[0] == z1z,
x2[0] == x2x, y2[0] == y2y, z2[0] == z2z,
 
vx1[0] == v1x, vy1[0] == v1y, vz1[0] == v1z,
vx2[0] == v2x, vy2[0] == v2y, vz2[0] == v2z},
 
{x1, x2, y1, y2, z1, z2,
vx1, vx2, vy1, vy2, vz1, vz2},
 
{t, 0, Tmax},

WorkingPrecision-> wp,
MaxSteps-> Infinity,
Method-> mta,
InterpolationOrder-> All];
 
(* Position, Geschwindigkeit *)
 
f2p={{x1[t], y1[t], z1[t]}, {x2[t], y2[t], z2[t]}}/.nds[[1]];
f2v={{vx1[t], vy1[t], vz1[t]}, {vx2[t], vy2[t], vz2[t]}}/.nds[[1]];
 
(* Formatierung *)
 
s[text_]:=Style[text, FontSize->11];
sw[text_]:=Style[text, White, FontSize->11];
colorfunc[n_]:=Function[{x, y, z, t},
Hue[0, n, 0.5,
If[Tmax<0, Max[Min[(+T+(-t+trail))/trail, 1], 0],
Max[Min[(-T+(t+trail))/d1, 1], 0]]]];
 
(* Animation *)
 
Do[Print[Rasterize[
Grid[{{
Show[

If[T == 0, {},

ParametricPlot3D[Evaluate[f2p],
{t, Max[0, T-trail], T},

PlotStyle->{
{Thickness[thk], Red},
{Thickness[thk], Blue}},

PlotRange->plotrange, AspectRatio->1, MaxRecursion->15, Axes->True, ImageSize->imagesize]],
 
Graphics3D[{
 
{PointSize[2point/3], Lighter[Red],     Point[{x1x, y1y, z1z}]},
{PointSize[2point/3], Lighter[Blue],    Point[{x2x, y2y, z2z}]}},

PlotRange->plotrange, AspectRatio->1, Axes->True, ImageSize->imagesize],
 
Graphics3D[{PointSize[point], Red,      Point[With[{t=T}, Evaluate[f2p]][[1]]]}],
Graphics3D[{PointSize[point], Blue,     Point[With[{t=T}, Evaluate[f2p]][[2]]]}],
 
ViewPoint->{0, Infinity, 0}]},
 
{ },
{s["t"->N[T]]},
{ },
{s["p1{x,y,z}"->With[{t=T}, Evaluate[f2p[[1]]]]], sw[1/2]},
{s["v1{x,y,z}"->With[{t=T}, Evaluate[f2v[[1]]]]], sw[1/2]},
{s["v1{total}"->{With[{t=T}, Evaluate[Norm[f2v[[1]]]]]}], sw[1/2]},
{ },
{s["p2{x,y,z}"->With[{t=T}, Evaluate[f2p[[2]]]]], sw[1/2]},
{s["v2{x,y,z}"->With[{t=T}, Evaluate[f2v[[2]]]]], sw[1/2]},
{s["v2{total}"->{With[{t=T}, Evaluate[Norm[f2v[[2]]]]]}], sw[1/2]}
}, Alignment->Left]]],
 
(* Zeitregler *)
 
{T, 0, Tmax, Tmax/5}]

(* Export als HTML Dokument *)
(* Export["dateiname.html", EvaluationNotebook[], "GraphicsOutput" -> "PNG"] *)
(* Export direkt als Bildsequenz *)
(* Do[Export["dateiname" <> ToString[T] <> ".png", Rasterize[...] ], {T, 0, 10, 5}] *)










3 Körper:

Code: Alles auswählen

(* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
(* ||| Mathematica Syntax || yukterez.net || 3 Body Newtonian Mass & Charge Simulator ||| *)
(* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)

ClearAll["Global`*"]
Amp = 1; kg = 1; m = 1; sek = 1; km = 1000 m; (* SI Einheiten *)
 
mt1 = {"StiffnessSwitching", Method-> {"ExplicitRungeKutta", Automatic}};
mt2 = {"ImplicitRungeKutta", "DifferenceOrder"-> 20};
mt3 = {"EquationSimplification"-> "Residual"};
mt0 = Automatic;
mta = mt2;
wp  = MachinePrecision;
 
(* Plot Optionen *)
 
Tmax      = 10000 sek;
trail     = 3000 sek;
point     = 0.015;
thk       = 0.004;
plotrange = {{-0.2 m, +1.2 m}, {-0.6 m, +0.6 m}, {-0.2 m, +1.2 m}};
imagesize = 430;
 
(* Masse 1 *)
 
m1  = 1000 kg;
q1  = 0;

x1x = 1/2 m;
y1y = 0 m;
z1z = Sqrt[3]/2 m;

v1x = 0 m/sek;
v1y = 0 m/sek;
v1z = 0 m/sek;
 
(* Masse 2 *)
 
m2  = 2000/3 kg;
q2  = 0 Amp sek;
 
x2x = 1/10 m;
y2y = 0 m;
z2z = 0 m;
 
v2x = 0 m/sek;
v2y = 0 m/sek;
v2z = 0 m/sek;
 
(* Masse 3 *)
 
m3  = 500 kg;
q3  = 0 Amp sek;
 
x3x = 4/5 m;
y3y = 0 m;
z3z = 1/5 m;
 
v3x = 0 m/sek;
v3y = 0 m/sek;
v3z = 0 m/sek;
 
(* Konstanten *)
 
G  = 667384/10^16 m^3/kg/sek^2;
ε0 = 8.854187817*^-12 Amp^2 sek^4/kg/m^3;
Au = 149597870690 m;
m0 = 1988*10^27 kg;
yr = 36525*24*36 sek;
 
(* Differentialgleichung *)
 
nds=NDSolve[{
 
x1'[t] == vx1[t], y1'[t] == vy1[t], z1'[t] == vz1[t],
x2'[t] == vx2[t], y2'[t] == vy2[t], z2'[t] == vz2[t],
x3'[t] == vx3[t], y3'[t] == vy3[t], z3'[t] == vz3[t],
 
vx1'[t] ==
 
(G m2 (x2[t]-x1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(G m3 (x3[t]-x1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3]+
 
(-q1*q2/(4Pi ε0 )/m1 (x2[t]-x1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(-q1*q3/(4Pi ε0 )/m1 (x3[t]-x1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3],
 
 vy1'[t] ==
 
(G m2 (y2[t]-y1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(G m3 (y3[t]-y1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3]+
 
(-q1*q2/(4Pi ε0 )/m1 (y2[t]-y1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(-q1*q3/(4Pi ε0 )/m1 (y3[t]-y1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3],
 
vz1'[t] ==
 
(G m2 (z2[t]-z1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(G m3 (z3[t]-z1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3]+
 
(-q1*q2/(4Pi ε0 )/m1 (z2[t]-z1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(-q1*q3/(4Pi ε0 )/m1 (z3[t]-z1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3],
 
vx2'[t] ==
 
(G m1 (x1[t]-x2[t]))/Sqrt[((x1[t]-x2[t])^2+(y1[t]-y2[t])^2+(z1[t]-z2[t])^2)^3]+
(G m3 (x3[t]-x2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3]+
 
(-q2*q1/(4Pi ε0 )/m2 (x1[t]-x2[t]))/Sqrt[((x1[t]-x2[t])^2+(y1[t]-y2[t])^2+(z1[t]-z2[t])^2)^3]+
(-q2*q3/(4Pi ε0 )/m2 (x3[t]-x2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3],
 
 vy2'[t] ==
 
(G m1 (y1[t]-y2[t]))/Sqrt[((x1[t]-x2[t])^2+(y1[t]-y2[t])^2+(z1[t]-z2[t])^2)^3]+
(G m3 (y3[t]-y2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3]+
 
(-q2*q1/(4Pi ε0 )/m2 (y1[t]-y2[t]))/Sqrt[((x1[t]-x2[t])^2+(y1[t]-y2[t])^2+(z1[t]-z2[t])^2)^3]+
(-q2*q3/(4Pi ε0 )/m2 (y3[t]-y2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3],
 
vz2'[t] ==
 
(G m1 (z1[t]-z2[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(G m3 (z3[t]-z2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3]+
 
(-q2*q1/(4Pi ε0 )/m2 (z1[t]-z2[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(-q2*q3/(4Pi ε0 )/m2 (z3[t]-z2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3],
 
vx3'[t] ==
 
(G m1 (x1[t]-x3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(G m2 (x2[t]-x3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3]+
 
(-q3*q1/(4Pi ε0 )/m3 (x1[t]-x3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(-q3*q2/(4Pi ε0 )/m3 (x2[t]-x3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3],
 
vy3'[t] ==
 
(G m1 (y1[t]-y3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(G m2 (y2[t]-y3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3]+
 
(-q3*q1/(4Pi ε0 )/m3 (y1[t]-y3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(-q3*q2/(4Pi ε0 )/m3 (y2[t]-y3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3],
 
vz3'[t] ==
 
(G m1 (z1[t]-z3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(G m2 (z2[t]-z3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3]+
 
(-q3*q1/(4Pi ε0 )/m3 (z1[t]-z3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(-q3*q2/(4Pi ε0 )/m3 (z2[t]-z3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3],
 
x1[0] == x1x, y1[0] == y1y, z1[0] == z1z,
x2[0] == x2x, y2[0] == y2y, z2[0] == z2z,
x3[0] == x3x, y3[0] == y3y, z3[0] == z3z,
 
vx1[0] == v1x, vy1[0] == v1y, vz1[0] == v1z,
vx2[0] == v2x, vy2[0] == v2y, vz2[0] == v2z,
vx3[0] == v3x, vy3[0] == v3y, vz3[0] == v3z},
 
{x1, x2, x3, y1, y2, y3, z1, z2, z3,
vx1, vx2, vx3, vy1, vy2, vy3, vz1, vz2, vz3},
 
{t, 0, Tmax},

WorkingPrecision-> wp,
MaxSteps-> Infinity,
Method-> mta,
InterpolationOrder-> All];
 
(* Position, Geschwindigkeit *)
 
f2p={{x1[t], y1[t], z1[t]}, {x2[t], y2[t], z2[t]}, {x3[t], y3[t], z3[t]}}/.nds[[1]];
f2v={{vx1[t], vy1[t], vz1[t]}, {vx2[t], vy2[t], vz2[t]}, {vx3[t], vy3[t], vz3[t]}}/.nds[[1]];
 
(* Formatierung *)
 
s[text_]:=Style[text, FontSize->11];
sw[text_]:=Style[text, White, FontSize->11];
colorfunc[n_]:=Function[{x, y, z, t},
Hue[0, n, 0.5,
If[Tmax<0, Max[Min[(+T+(-t+trail))/trail, 1], 0],
Max[Min[(-T+(t+trail))/d1, 1], 0]]]];
 
(* Animation *)
 
Do[Print[Rasterize[
Grid[{{
Show[

If[T == 0, {},

ParametricPlot3D[Evaluate[f2p],
{t, Max[0, T-trail], T},

PlotStyle->{
{Thickness[thk], Red},
{Thickness[thk], Blue},
{Thickness[thk], Green}},

PlotRange->plotrange, AspectRatio->1, MaxRecursion->15, Axes->True, ImageSize->imagesize]],
 
Graphics3D[{
 
{PointSize[2point/3], Lighter[Red],     Point[{x1x, y1y, z1z}]},
{PointSize[2point/3], Lighter[Blue],    Point[{x2x, y2y, z2z}]},
{PointSize[2point/3], Lighter[Green],   Point[{x3x, y3y, z3z}]}},

PlotRange->plotrange, AspectRatio->1, Axes->True, ImageSize->imagesize],
 
Graphics3D[{PointSize[point], Red,      Point[With[{t=T}, Evaluate[f2p]][[1]]]}],
Graphics3D[{PointSize[point], Blue,     Point[With[{t=T}, Evaluate[f2p]][[2]]]}],
Graphics3D[{PointSize[point], Green,    Point[With[{t=T}, Evaluate[f2p]][[3]]]}],
 
ViewPoint->{0, Infinity, 0}]},
 
{ },
{s["t"->N[T]]},
{ },
{s["p1{x,y,z}"->With[{t=T}, Evaluate[f2p[[1]]]]], sw[1/2]},
{s["v1{x,y,z}"->With[{t=T}, Evaluate[f2v[[1]]]]], sw[1/2]},
{s["v1{total}"->{With[{t=T}, Evaluate[Norm[f2v[[1]]]]]}], sw[1/2]},
{ },
{s["p2{x,y,z}"->With[{t=T}, Evaluate[f2p[[2]]]]], sw[1/2]},
{s["v2{x,y,z}"->With[{t=T}, Evaluate[f2v[[2]]]]], sw[1/2]},
{s["v2{total}"->{With[{t=T}, Evaluate[Norm[f2v[[2]]]]]}], sw[1/2]},
{ },
{s["p3{x,y,z}"->With[{t=T}, Evaluate[f2p[[3]]]]], sw[1/2]},
{s["v3{x,y,z}"->With[{t=T}, Evaluate[f2v[[3]]]]], sw[1/2]},
{s["v3{total}"->{With[{t=T}, Evaluate[Norm[f2v[[3]]]]]}], sw[1/2]}
}, Alignment->Left]]],
 
(* Zeitregler *)
 
{T, 0, Tmax, Tmax/5}]

(* Export als HTML Dokument *)
(* Export["dateiname.html", EvaluationNotebook[], "GraphicsOutput" -> "PNG"] *)
(* Export direkt als Bildsequenz *)
(* Do[Export["dateiname" <> ToString[T] <> ".png", Rasterize[...] ], {T, 0, 10, 5}] *)










4 Körper:

Code: Alles auswählen

(* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
(* ||| Mathematica Syntax || yukterez.net || 4 Body Newtonian Mass & Charge Simulator ||| *)
(* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)

ClearAll["Global`*"]
Amp = 1; kg = 1; m = 1; sek = 1; km = 1000 m; (* SI Einheiten *)
 
mt1 = {"StiffnessSwitching", Method-> {"ExplicitRungeKutta", Automatic}};
mt2 = {"ImplicitRungeKutta", "DifferenceOrder"-> 20};
mt3 = {"EquationSimplification"-> "Residual"};
mt0 = Automatic;
mta = mt2;
wp  = MachinePrecision;
 
(* Plot Optionen *)
 
Tmax      = 2π sek;
trail     = π/3 sek;
point     = 0.015;
thk       = 0.004;
plotrange = {{-1.2 m, +1.2 m}, {-1.2 m, +1.2 m}, {-1.2 m, +1.2 m}};
imagesize = 430;
 
(* Masse 1 *)
 
m1  = 1 m^3/sek^2/G;
q1  = 0 Amp sek;
 
x1x = 1.00231488346205 m;
y1y = 0 m;
z1z = 0 m;
 
v1x = 0 m/sek;
v1y = 0 m/sek;
v1z = -0.293790277732029 m/sek;
 
(* Masse 2 *)
 
m2  = 1 m^3/sek^2/G;
q2  = 0 Amp sek;
 
x2x = -0.52869409402363 m;
y2y = 0 m;
z2z = 0.567125954067238 m;
 
v2x = -0.175826619093916 m/sek;
v2y = 0 m/sek;
v2z = 1.02361310165052 m/sek;
 
(* Masse 3 *)
 
m3  = 1 m^3/sek^2/G;
q3  = 0 Amp sek;
 
x3x = 0.0550733045852099 m;
y3y = 0 m;
z3z = 0 m;
 
v3x = 0 m/sek;
v3y = 0 m/sek;
v3z = -1.75343592556901 m/sek;
 
(* Masse 4 *)
 
m4  = 1 m^3/sek^2/G;
q4  = 0 Amp sek;
 
x4x = -0.528694094023634 m;
y4y = 0 m;
z4z = -0.567125954067235 m;
 
v4x = 0.175826619093915 m/sek;
v4y = 0 m/sek;
v4z = 1.02361310165051 m/sek;
 
(* Konstanten *)
 
G  = 667384/10^16 m^3/kg/sek^2;
ε0 = 8.854187817*^-12 Amp^2 sek^4/kg/m^3;
Au = 149597870690 m;
m0 = 1988*10^27 kg;
yr = 36525*24*36 sek;
 
(* Differentialgleichung *)
 
nds=NDSolve[{
 
x1'[t] == vx1[t], y1'[t] == vy1[t], z1'[t] == vz1[t],
x2'[t] == vx2[t], y2'[t] == vy2[t], z2'[t] == vz2[t],
x3'[t] == vx3[t], y3'[t] == vy3[t], z3'[t] == vz3[t],
x4'[t] == vx4[t], y4'[t] == vy4[t], z4'[t] == vz4[t], 
 
vx1'[t] ==
 
(G m2 (x2[t]-x1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(G m3 (x3[t]-x1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3]+
(G m4 (x4[t]-x1[t]))/Sqrt[((x4[t]-x1[t])^2+(y4[t]-y1[t])^2+(z4[t]-z1[t])^2)^3]+
 
(-q1*q2/(4Pi ε0 )/m1 (x2[t]-x1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(-q1*q3/(4Pi ε0 )/m1 (x3[t]-x1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3]+
(-q1*q4/(4Pi ε0 )/m1 (x4[t]-x1[t]))/Sqrt[((x4[t]-x1[t])^2+(y4[t]-y1[t])^2+(z4[t]-z1[t])^2)^3],
 
 vy1'[t] ==
 
(G m2 (y2[t]-y1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(G m3 (y3[t]-y1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3]+
(G m4 (y4[t]-y1[t]))/Sqrt[((x4[t]-x1[t])^2+(y4[t]-y1[t])^2+(z4[t]-z1[t])^2)^3]+
 
(-q1*q2/(4Pi ε0 )/m1 (y2[t]-y1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(-q1*q3/(4Pi ε0 )/m1 (y3[t]-y1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3]+
(-q1*q4/(4Pi ε0 )/m1 (y4[t]-y1[t]))/Sqrt[((x4[t]-x1[t])^2+(y4[t]-y1[t])^2+(z4[t]-z1[t])^2)^3],
 
vz1'[t] ==
 
(G m2 (z2[t]-z1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(G m3 (z3[t]-z1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3]+
(G m4 (z4[t]-z1[t]))/Sqrt[((x4[t]-x1[t])^2+(y4[t]-y1[t])^2+(z4[t]-z1[t])^2)^3]+
 
(-q1*q2/(4Pi ε0 )/m1 (z2[t]-z1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(-q1*q3/(4Pi ε0 )/m1 (z3[t]-z1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3]+
(-q1*q4/(4Pi ε0 )/m1 (z4[t]-z1[t]))/Sqrt[((x4[t]-x1[t])^2+(y4[t]-y1[t])^2+(z4[t]-z1[t])^2)^3],
 
vx2'[t] ==
 
(G m1 (x1[t]-x2[t]))/Sqrt[((x1[t]-x2[t])^2+(y1[t]-y2[t])^2+(z1[t]-z2[t])^2)^3]+
(G m3 (x3[t]-x2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3]+
(G m4 (x4[t]-x2[t]))/Sqrt[((x4[t]-x2[t])^2+(y4[t]-y2[t])^2+(z4[t]-z2[t])^2)^3]+
 
(-q2*q1/(4Pi ε0 )/m2 (x1[t]-x2[t]))/Sqrt[((x1[t]-x2[t])^2+(y1[t]-y2[t])^2+(z1[t]-z2[t])^2)^3]+
(-q2*q3/(4Pi ε0 )/m2 (x3[t]-x2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3]+
(-q2*q4/(4Pi ε0 )/m2 (x4[t]-x2[t]))/Sqrt[((x4[t]-x2[t])^2+(y4[t]-y2[t])^2+(z4[t]-z2[t])^2)^3],
 
 vy2'[t] ==
 
(G m1 (y1[t]-y2[t]))/Sqrt[((x1[t]-x2[t])^2+(y1[t]-y2[t])^2+(z1[t]-z2[t])^2)^3]+
(G m3 (y3[t]-y2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3]+
(G m4 (y4[t]-y2[t]))/Sqrt[((x4[t]-x2[t])^2+(y4[t]-y2[t])^2+(z4[t]-z2[t])^2)^3]+
 
(-q2*q1/(4Pi ε0 )/m2 (y1[t]-y2[t]))/Sqrt[((x1[t]-x2[t])^2+(y1[t]-y2[t])^2+(z1[t]-z2[t])^2)^3]+
(-q2*q3/(4Pi ε0 )/m2 (y3[t]-y2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3]+
(-q2*q4/(4Pi ε0 )/m2 (y4[t]-y2[t]))/Sqrt[((x4[t]-x2[t])^2+(y4[t]-y2[t])^2+(z4[t]-z2[t])^2)^3],
 
vz2'[t] ==
 
(G m1 (z1[t]-z2[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(G m3 (z3[t]-z2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3]+
(G m4 (z4[t]-z2[t]))/Sqrt[((x4[t]-x2[t])^2+(y4[t]-y2[t])^2+(z4[t]-z2[t])^2)^3]+
 
(-q2*q1/(4Pi ε0 )/m2 (z1[t]-z2[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(-q2*q3/(4Pi ε0 )/m2 (z3[t]-z2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3]+
(-q2*q4/(4Pi ε0 )/m2 (z4[t]-z2[t]))/Sqrt[((x4[t]-x2[t])^2+(y4[t]-y2[t])^2+(z4[t]-z2[t])^2)^3],
 
vx3'[t] ==
 
(G m1 (x1[t]-x3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(G m2 (x2[t]-x3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3]+
(G m4 (x4[t]-x3[t]))/Sqrt[((x4[t]-x3[t])^2+(y4[t]-y3[t])^2+(z4[t]-z3[t])^2)^3]+
 
(-q3*q1/(4Pi ε0 )/m3 (x1[t]-x3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(-q3*q2/(4Pi ε0 )/m3 (x2[t]-x3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3]+
(-q3*q4/(4Pi ε0 )/m3 (x4[t]-x3[t]))/Sqrt[((x4[t]-x3[t])^2+(y4[t]-y3[t])^2+(z4[t]-z3[t])^2)^3],
 
vy3'[t] ==
 
(G m1 (y1[t]-y3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(G m2 (y2[t]-y3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3]+
(G m4 (y4[t]-y3[t]))/Sqrt[((x4[t]-x3[t])^2+(y4[t]-y3[t])^2+(z4[t]-z3[t])^2)^3]+
 
(-q3*q1/(4Pi ε0 )/m3 (y1[t]-y3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(-q3*q2/(4Pi ε0 )/m3 (y2[t]-y3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3]+
(-q3*q4/(4Pi ε0 )/m3 (y4[t]-y3[t]))/Sqrt[((x4[t]-x3[t])^2+(y4[t]-y3[t])^2+(z4[t]-z3[t])^2)^3],
 
vz3'[t] ==
 
(G m1 (z1[t]-z3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(G m2 (z2[t]-z3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3]+
(G m4 (z4[t]-z3[t]))/Sqrt[((x4[t]-x3[t])^2+(y4[t]-y3[t])^2+(z4[t]-z3[t])^2)^3]+
 
(-q3*q1/(4Pi ε0 )/m3 (z1[t]-z3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(-q3*q2/(4Pi ε0 )/m3 (z2[t]-z3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3]+
(-q3*q4/(4Pi ε0 )/m3 (z4[t]-z3[t]))/Sqrt[((x4[t]-x3[t])^2+(y4[t]-y3[t])^2+(z4[t]-z3[t])^2)^3],
 
vx4'[t] ==
 
(G m1 (x1[t]-x4[t]))/Sqrt[((x1[t]-x4[t])^2+(y1[t]-y4[t])^2+(z1[t]-z4[t])^2)^3]+
(G m2 (x2[t]-x4[t]))/Sqrt[((x2[t]-x4[t])^2+(y2[t]-y4[t])^2+(z2[t]-z4[t])^2)^3]+
(G m3 (x3[t]-x4[t]))/Sqrt[((x3[t]-x4[t])^2+(y3[t]-y4[t])^2+(z3[t]-z4[t])^2)^3]+
 
(-q4*q1/(4Pi ε0 )/m4 (x1[t]-x4[t]))/Sqrt[((x1[t]-x4[t])^2+(y1[t]-y4[t])^2+(z1[t]-z4[t])^2)^3]+
(-q4*q2/(4Pi ε0 )/m4 (x2[t]-x4[t]))/Sqrt[((x2[t]-x4[t])^2+(y2[t]-y4[t])^2+(z2[t]-z4[t])^2)^3]+
(-q4*q3/(4Pi ε0 )/m4 (x3[t]-x4[t]))/Sqrt[((x3[t]-x4[t])^2+(y3[t]-y4[t])^2+(z3[t]-z4[t])^2)^3],
 
vy4'[t] ==
 
(G m1 (y1[t]-y4[t]))/Sqrt[((x1[t]-x4[t])^2+(y1[t]-y4[t])^2+(z1[t]-z4[t])^2)^3]+
(G m2 (y2[t]-y4[t]))/Sqrt[((x2[t]-x4[t])^2+(y2[t]-y4[t])^2+(z2[t]-z4[t])^2)^3]+
(G m3 (y3[t]-y4[t]))/Sqrt[((x3[t]-x4[t])^2+(y3[t]-y4[t])^2+(z3[t]-z4[t])^2)^3]+
 
(-q4*q1/(4Pi ε0 )/m4 (y1[t]-y4[t]))/Sqrt[((x1[t]-x4[t])^2+(y1[t]-y4[t])^2+(z1[t]-z4[t])^2)^3]+
(-q4*q2/(4Pi ε0 )/m4 (y2[t]-y4[t]))/Sqrt[((x2[t]-x4[t])^2+(y2[t]-y4[t])^2+(z2[t]-z4[t])^2)^3]+
(-q4*q3/(4Pi ε0 )/m4 (y3[t]-y4[t]))/Sqrt[((x3[t]-x4[t])^2+(y3[t]-y4[t])^2+(z3[t]-z4[t])^2)^3],
 
vz4'[t] ==
 
(G m1 (z1[t]-z4[t]))/Sqrt[((x1[t]-x4[t])^2+(y1[t]-y4[t])^2+(z1[t]-z4[t])^2)^3]+
(G m2 (z2[t]-z4[t]))/Sqrt[((x2[t]-x4[t])^2+(y2[t]-y4[t])^2+(z2[t]-z4[t])^2)^3]+
(G m3 (z3[t]-z4[t]))/Sqrt[((x3[t]-x4[t])^2+(y3[t]-y4[t])^2+(z3[t]-z4[t])^2)^3]+
 
(-q4*q1/(4Pi ε0 )/m4 (z1[t]-z4[t]))/Sqrt[((x1[t]-x4[t])^2+(y1[t]-y4[t])^2+(z1[t]-z4[t])^2)^3]+
(-q4*q2/(4Pi ε0 )/m4 (z2[t]-z4[t]))/Sqrt[((x2[t]-x4[t])^2+(y2[t]-y4[t])^2+(z2[t]-z4[t])^2)^3]+
(-q4*q3/(4Pi ε0 )/m4 (z3[t]-z4[t]))/Sqrt[((x3[t]-x4[t])^2+(y3[t]-y4[t])^2+(z3[t]-z4[t])^2)^3],
 
x1[0] == x1x, y1[0] == y1y, z1[0] == z1z,
x2[0] == x2x, y2[0] == y2y, z2[0] == z2z,
x3[0] == x3x, y3[0] == y3y, z3[0] == z3z,
x4[0] == x4x, y4[0] == y4y, z4[0] == z4z,
 
vx1[0] == v1x, vy1[0] == v1y, vz1[0] == v1z,
vx2[0] == v2x, vy2[0] == v2y, vz2[0] == v2z,
vx3[0] == v3x, vy3[0] == v3y, vz3[0] == v3z,
vx4[0] == v4x, vy4[0] == v4y, vz4[0] == v4z},
 
{x1, x2, x3, x4, y1, y2, y3, y4, z1, z2, z3, z4,
vx1, vx2, vx3, vx4, vy1, vy2, vy3, vy4, vz1, vz2, vz3, vz4},
 
{t, 0, Tmax},

WorkingPrecision-> wp,
MaxSteps-> Infinity,
Method-> mta,
InterpolationOrder-> All];
 
(* Position, Geschwindigkeit *)
 
f2p={{x1[t], y1[t], z1[t]}, {x2[t], y2[t], z2[t]}, {x3[t], y3[t], z3[t]}, {x4[t], y4[t], z4[t]}}/.nds[[1]];
f2v={{vx1[t], vy1[t], vz1[t]}, {vx2[t], vy2[t], vz2[t]}, {vx3[t], vy3[t], vz3[t]}, {vx4[t], vy4[t], vz4[t]}}/.nds[[1]];
 
(* Formatierung *)
 
s[text_]:=Style[text, FontSize->11];
sw[text_]:=Style[text, White, FontSize->11];
colorfunc[n_]:=Function[{x, y, z, t},
Hue[0, n, 0.5,
If[Tmax<0, Max[Min[(+T+(-t+trail))/trail, 1], 0],
Max[Min[(-T+(t+trail))/d1, 1], 0]]]];
 
(* Animation *)
 
Do[Print[Rasterize[
Grid[{{
Show[

If[T == 0, {},

ParametricPlot3D[Evaluate[f2p],
{t, Max[0, T-trail], T},

PlotStyle->{
{Thickness[thk], Red},
{Thickness[thk], Blue},
{Thickness[thk], Green},
{Thickness[thk], Magenta}},

PlotRange->plotrange, AspectRatio->1, MaxRecursion->15, Axes->True, ImageSize->imagesize]],
 
Graphics3D[{
 
{PointSize[2point/3], Lighter[Red],     Point[{x1x, y1y, z1z}]},
{PointSize[2point/3], Lighter[Blue],    Point[{x2x, y2y, z2z}]},
{PointSize[2point/3], Lighter[Green],   Point[{x3x, y3y, z3z}]},
{PointSize[2point/3], Lighter[Magenta], Point[{x4x, y4y, z4z}]}},

PlotRange->plotrange, AspectRatio->1, Axes->True, ImageSize->imagesize],
 
Graphics3D[{PointSize[point], Red,      Point[With[{t=T}, Evaluate[f2p]][[1]]]}],
Graphics3D[{PointSize[point], Blue,     Point[With[{t=T}, Evaluate[f2p]][[2]]]}],
Graphics3D[{PointSize[point], Green,    Point[With[{t=T}, Evaluate[f2p]][[3]]]}],
Graphics3D[{PointSize[point], Magenta,  Point[With[{t=T}, Evaluate[f2p]][[4]]]}],
 
ViewPoint->{0, Infinity, 0}]},
 
{ },
{s["t"->N[T]]},
{ },
{s["p1{x,y,z}"->With[{t=T}, Evaluate[f2p[[1]]]]], sw[1/2]},
{s["v1{x,y,z}"->With[{t=T}, Evaluate[f2v[[1]]]]], sw[1/2]},
{s["v1{total}"->{With[{t=T}, Evaluate[Norm[f2v[[1]]]]]}], sw[1/2]},
{ },
{s["p2{x,y,z}"->With[{t=T}, Evaluate[f2p[[2]]]]], sw[1/2]},
{s["v2{x,y,z}"->With[{t=T}, Evaluate[f2v[[2]]]]], sw[1/2]},
{s["v2{total}"->{With[{t=T}, Evaluate[Norm[f2v[[2]]]]]}], sw[1/2]},
{ },
{s["p3{x,y,z}"->With[{t=T}, Evaluate[f2p[[3]]]]], sw[1/2]},
{s["v3{x,y,z}"->With[{t=T}, Evaluate[f2v[[3]]]]], sw[1/2]},
{s["v3{total}"->{With[{t=T}, Evaluate[Norm[f2v[[3]]]]]}], sw[1/2]},
{ },
{s["p4{x,y,z}"->With[{t=T}, Evaluate[f2p[[4]]]]], sw[1/2]},
{s["v4{x,y,z}"->With[{t=T}, Evaluate[f2v[[4]]]]], sw[1/2]},
{s["v4{total}"->{With[{t=T}, Evaluate[Norm[f2v[[4]]]]]}], sw[1/2]}
}, Alignment->Left]]],
 
(* Zeitregler *)
 
{T, 0, Tmax, Tmax/5}]

(* Export als HTML Dokument *)
(* Export["dateiname.html", EvaluationNotebook[], "GraphicsOutput" -> "PNG"] *)
(* Export direkt als Bildsequenz *)
(* Do[Export["dateiname" <> ToString[T] <> ".png", Rasterize[...] ], {T, 0, 10, 5}] *)










5 Körper:

Code: Alles auswählen

(* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)
(* ||| Mathematica Syntax || yukterez.net || 5 Body Newtonian Mass & Charge Simulator ||| *)
(* |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| *)

ClearAll["Global`*"]
Amp = 1; kg = 1; m = 1; sek = 1; km = 1000 m; (* SI Einheiten *)
 
mt1 = {"StiffnessSwitching", Method-> {"ExplicitRungeKutta", Automatic}};
mt2 = {"ImplicitRungeKutta", "DifferenceOrder"-> 20};
mt3 = {"EquationSimplification"-> "Residual"};
mt0 = Automatic;
mta = mt2;
wp  = MachinePrecision;
 
(* Plot Optionen *)
 
Tmax      = 2π sek;
trail     = π/3 sek;
point     = 0.015;
thk       = 0.004;
plotrange = {{-1.2 m, +1.2 m}, {-1.2 m, +1.2 m}, {-1.2 m, +1.2 m}};
imagesize = 430;
 
(* Masse 1 *)
 
m1  = 1 m^3/sek^2/G;
q1  = 0 Amp sek;
 
x1x = 0.3673027525587564 m;
y1y = 0 m;
z1z = 0 m;
 
v1x = 0 m/sek;
v1y = 0 m/sek;
v1z = 2.134808238913507 m/sek;
 
(* Masse 2 *)
 
m2  = 1 m^3/sek^2/G;
q2  = 0 Amp sek;
 
x2x = 0.6760231721222919 m;
y2y = 0 m;
z2z = -0.1852023639618821 m;
 
v2x = -2.138115663581648 m/sek;
v2y = 0 m/sek;
v2z = -0.7311482987140561 m/sek;
 
(* Masse 3 *)
 
m3  = 1 m^3/sek^2/G;
q3  = 0 Amp sek;
 
x3x = -0.859674548416703 m;
y3y = 0 m;
z3z = -0.3897351750882059 m;
 
v3x = 0.5365854530502069 m/sek;
v3y = 0 m/sek;
v3z = -0.3362558207426968 m/sek;
 
(* Masse 4 *)
 
m4  = 1 m^3/sek^2/G;
q4  = 0 Amp sek;
 
x4x = -0.8596745484016705 m;
y4y = 0 m;
z4z = 0.3897351750882059 m;
 
v4x = -0.5365854530502073 m/sek;
v4y = 0 m/sek;
v4z = -0.3362558207426963 m/sek;
 
(* Masse 5 *)
 
m5  = 1 m^3/sek^2/G;
q5  = 0 Amp sek;
 
x5x = 0.6760231721222916 m;
y5y = 0 m;
z5z = 0.1852023639618821 m;
 
v5x = 2.138115663581649 m/sek;
v5y = 0 m/sek;
v5z = -0.7311482987140556 m/sek;
 
(* Konstanten *)
 
G  = 667384/10^16 m^3/kg/sek^2;
ε0 = 8.854187817*^-12 Amp^2 sek^4/kg/m^3;
Au = 149597870690 m;
m0 = 1988*10^27 kg;
yr = 36525*24*36 sek;
 
(* Differentialgleichung *)
 
nds=NDSolve[{
 
x1'[t] == vx1[t], y1'[t] == vy1[t], z1'[t] == vz1[t],
x2'[t] == vx2[t], y2'[t] == vy2[t], z2'[t] == vz2[t],
x3'[t] == vx3[t], y3'[t] == vy3[t], z3'[t] == vz3[t],
x4'[t] == vx4[t], y4'[t] == vy4[t], z4'[t] == vz4[t],
x5'[t] == vx5[t], y5'[t] == vy5[t], z5'[t] == vz5[t],
 
vx1'[t] ==
 
(G m2 (x2[t]-x1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(G m3 (x3[t]-x1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3]+
(G m4 (x4[t]-x1[t]))/Sqrt[((x4[t]-x1[t])^2+(y4[t]-y1[t])^2+(z4[t]-z1[t])^2)^3]+
(G m5 (x5[t]-x1[t]))/Sqrt[((x5[t]-x1[t])^2+(y5[t]-y1[t])^2+(z5[t]-z1[t])^2)^3]+
 
(-q1*q2/(4Pi ε0 )/m1 (x2[t]-x1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(-q1*q3/(4Pi ε0 )/m1 (x3[t]-x1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3]+
(-q1*q4/(4Pi ε0 )/m1 (x4[t]-x1[t]))/Sqrt[((x4[t]-x1[t])^2+(y4[t]-y1[t])^2+(z4[t]-z1[t])^2)^3]+
(-q1*q5/(4Pi ε0 )/m1 (x5[t]-x1[t]))/Sqrt[((x5[t]-x1[t])^2+(y5[t]-y1[t])^2+(z5[t]-z1[t])^2)^3],
 
 vy1'[t] ==
 
(G m2 (y2[t]-y1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(G m3 (y3[t]-y1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3]+
(G m4 (y4[t]-y1[t]))/Sqrt[((x4[t]-x1[t])^2+(y4[t]-y1[t])^2+(z4[t]-z1[t])^2)^3]+
(G m5 (y5[t]-y1[t]))/Sqrt[((x5[t]-x1[t])^2+(y5[t]-y1[t])^2+(z5[t]-z1[t])^2)^3]+
 
(-q1*q2/(4Pi ε0 )/m1 (y2[t]-y1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(-q1*q3/(4Pi ε0 )/m1 (y3[t]-y1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3]+
(-q1*q4/(4Pi ε0 )/m1 (y4[t]-y1[t]))/Sqrt[((x4[t]-x1[t])^2+(y4[t]-y1[t])^2+(z4[t]-z1[t])^2)^3]+
(-q1*q5/(4Pi ε0 )/m1 (y5[t]-y1[t]))/Sqrt[((x5[t]-x1[t])^2+(y5[t]-y1[t])^2+(z5[t]-z1[t])^2)^3],
 
vz1'[t] ==
 
(G m2 (z2[t]-z1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(G m3 (z3[t]-z1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3]+
(G m4 (z4[t]-z1[t]))/Sqrt[((x4[t]-x1[t])^2+(y4[t]-y1[t])^2+(z4[t]-z1[t])^2)^3]+
(G m5 (z5[t]-z1[t]))/Sqrt[((x5[t]-x1[t])^2+(y5[t]-y1[t])^2+(z5[t]-z1[t])^2)^3]+
 
(-q1*q2/(4Pi ε0 )/m1 (z2[t]-z1[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(-q1*q3/(4Pi ε0 )/m1 (z3[t]-z1[t]))/Sqrt[((x3[t]-x1[t])^2+(y3[t]-y1[t])^2+(z3[t]-z1[t])^2)^3]+
(-q1*q4/(4Pi ε0 )/m1 (z4[t]-z1[t]))/Sqrt[((x4[t]-x1[t])^2+(y4[t]-y1[t])^2+(z4[t]-z1[t])^2)^3]+
(-q1*q5/(4Pi ε0 )/m1 (z5[t]-z1[t]))/Sqrt[((x5[t]-x1[t])^2+(y5[t]-y1[t])^2+(z5[t]-z1[t])^2)^3],
 
vx2'[t] ==
 
(G m1 (x1[t]-x2[t]))/Sqrt[((x1[t]-x2[t])^2+(y1[t]-y2[t])^2+(z1[t]-z2[t])^2)^3]+
(G m3 (x3[t]-x2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3]+
(G m4 (x4[t]-x2[t]))/Sqrt[((x4[t]-x2[t])^2+(y4[t]-y2[t])^2+(z4[t]-z2[t])^2)^3]+
(G m5 (x5[t]-x2[t]))/Sqrt[((x5[t]-x2[t])^2+(y5[t]-y2[t])^2+(z5[t]-z2[t])^2)^3]+
 
(-q2*q1/(4Pi ε0 )/m2 (x1[t]-x2[t]))/Sqrt[((x1[t]-x2[t])^2+(y1[t]-y2[t])^2+(z1[t]-z2[t])^2)^3]+
(-q2*q3/(4Pi ε0 )/m2 (x3[t]-x2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3]+
(-q2*q4/(4Pi ε0 )/m2 (x4[t]-x2[t]))/Sqrt[((x4[t]-x2[t])^2+(y4[t]-y2[t])^2+(z4[t]-z2[t])^2)^3]+
(-q2*q5/(4Pi ε0 )/m2 (x5[t]-x2[t]))/Sqrt[((x5[t]-x2[t])^2+(y5[t]-y2[t])^2+(z5[t]-z2[t])^2)^3],
 
 vy2'[t] ==
 
(G m1 (y1[t]-y2[t]))/Sqrt[((x1[t]-x2[t])^2+(y1[t]-y2[t])^2+(z1[t]-z2[t])^2)^3]+
(G m3 (y3[t]-y2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3]+
(G m4 (y4[t]-y2[t]))/Sqrt[((x4[t]-x2[t])^2+(y4[t]-y2[t])^2+(z4[t]-z2[t])^2)^3]+
(G m5 (y5[t]-y2[t]))/Sqrt[((x5[t]-x2[t])^2+(y5[t]-y2[t])^2+(z5[t]-z2[t])^2)^3]+
 
(-q2*q1/(4Pi ε0 )/m2 (y1[t]-y2[t]))/Sqrt[((x1[t]-x2[t])^2+(y1[t]-y2[t])^2+(z1[t]-z2[t])^2)^3]+
(-q2*q3/(4Pi ε0 )/m2 (y3[t]-y2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3]+
(-q2*q4/(4Pi ε0 )/m2 (y4[t]-y2[t]))/Sqrt[((x4[t]-x2[t])^2+(y4[t]-y2[t])^2+(z4[t]-z2[t])^2)^3]+
(-q2*q5/(4Pi ε0 )/m2 (y5[t]-y2[t]))/Sqrt[((x5[t]-x2[t])^2+(y5[t]-y2[t])^2+(z5[t]-z2[t])^2)^3],
 
vz2'[t] ==
 
(G m1 (z1[t]-z2[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(G m3 (z3[t]-z2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3]+
(G m4 (z4[t]-z2[t]))/Sqrt[((x4[t]-x2[t])^2+(y4[t]-y2[t])^2+(z4[t]-z2[t])^2)^3]+
(G m5 (z5[t]-z2[t]))/Sqrt[((x5[t]-x2[t])^2+(y5[t]-y2[t])^2+(z5[t]-z2[t])^2)^3]+
 
(-q2*q1/(4Pi ε0 )/m2 (z1[t]-z2[t]))/Sqrt[((x2[t]-x1[t])^2+(y2[t]-y1[t])^2+(z2[t]-z1[t])^2)^3]+
(-q2*q3/(4Pi ε0 )/m2 (z3[t]-z2[t]))/Sqrt[((x3[t]-x2[t])^2+(y3[t]-y2[t])^2+(z3[t]-z2[t])^2)^3]+
(-q2*q4/(4Pi ε0 )/m2 (z4[t]-z2[t]))/Sqrt[((x4[t]-x2[t])^2+(y4[t]-y2[t])^2+(z4[t]-z2[t])^2)^3]+
(-q2*q5/(4Pi ε0 )/m2 (z5[t]-z2[t]))/Sqrt[((x5[t]-x2[t])^2+(y5[t]-y2[t])^2+(z5[t]-z2[t])^2)^3],
 
vx3'[t] ==
 
(G m1 (x1[t]-x3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(G m2 (x2[t]-x3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3]+
(G m4 (x4[t]-x3[t]))/Sqrt[((x4[t]-x3[t])^2+(y4[t]-y3[t])^2+(z4[t]-z3[t])^2)^3]+
(G m5 (x5[t]-x3[t]))/Sqrt[((x5[t]-x3[t])^2+(y5[t]-y3[t])^2+(z5[t]-z3[t])^2)^3]+
 
(-q3*q1/(4Pi ε0 )/m3 (x1[t]-x3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(-q3*q2/(4Pi ε0 )/m3 (x2[t]-x3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3]+
(-q3*q4/(4Pi ε0 )/m3 (x4[t]-x3[t]))/Sqrt[((x4[t]-x3[t])^2+(y4[t]-y3[t])^2+(z4[t]-z3[t])^2)^3]+
(-q3*q5/(4Pi ε0 )/m3 (x5[t]-x3[t]))/Sqrt[((x5[t]-x3[t])^2+(y5[t]-y3[t])^2+(z5[t]-z3[t])^2)^3],
 
vy3'[t] ==
 
(G m1 (y1[t]-y3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(G m2 (y2[t]-y3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3]+
(G m4 (y4[t]-y3[t]))/Sqrt[((x4[t]-x3[t])^2+(y4[t]-y3[t])^2+(z4[t]-z3[t])^2)^3]+
(G m5 (y5[t]-y3[t]))/Sqrt[((x5[t]-x3[t])^2+(y5[t]-y3[t])^2+(z5[t]-z3[t])^2)^3]+
 
(-q3*q1/(4Pi ε0 )/m3 (y1[t]-y3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(-q3*q2/(4Pi ε0 )/m3 (y2[t]-y3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3]+
(-q3*q4/(4Pi ε0 )/m3 (y4[t]-y3[t]))/Sqrt[((x4[t]-x3[t])^2+(y4[t]-y3[t])^2+(z4[t]-z3[t])^2)^3]+
(-q3*q5/(4Pi ε0 )/m3 (y5[t]-y3[t]))/Sqrt[((x5[t]-x3[t])^2+(y5[t]-y3[t])^2+(z5[t]-z3[t])^2)^3],
 
vz3'[t] ==
 
(G m1 (z1[t]-z3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(G m2 (z2[t]-z3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3]+
(G m4 (z4[t]-z3[t]))/Sqrt[((x4[t]-x3[t])^2+(y4[t]-y3[t])^2+(z4[t]-z3[t])^2)^3]+
(G m5 (z5[t]-z3[t]))/Sqrt[((x5[t]-x3[t])^2+(y5[t]-y3[t])^2+(z5[t]-z3[t])^2)^3]+
 
(-q3*q1/(4Pi ε0 )/m3 (z1[t]-z3[t]))/Sqrt[((x1[t]-x3[t])^2+(y1[t]-y3[t])^2+(z1[t]-z3[t])^2)^3]+
(-q3*q2/(4Pi ε0 )/m3 (z2[t]-z3[t]))/Sqrt[((x2[t]-x3[t])^2+(y2[t]-y3[t])^2+(z2[t]-z3[t])^2)^3]+
(-q3*q4/(4Pi ε0 )/m3 (z4[t]-z3[t]))/Sqrt[((x4[t]-x3[t])^2+(y4[t]-y3[t])^2+(z4[t]-z3[t])^2)^3]+
(-q3*q5/(4Pi ε0 )/m3 (z5[t]-z3[t]))/Sqrt[((x5[t]-x3[t])^2+(y5[t]-y3[t])^2+(z5[t]-z3[t])^2)^3],
 
vx4'[t] ==
 
(G m1 (x1[t]-x4[t]))/Sqrt[((x1[t]-x4[t])^2+(y1[t]-y4[t])^2+(z1[t]-z4[t])^2)^3]+
(G m2 (x2[t]-x4[t]))/Sqrt[((x2[t]-x4[t])^2+(y2[t]-y4[t])^2+(z2[t]-z4[t])^2)^3]+
(G m3 (x3[t]-x4[t]))/Sqrt[((x3[t]-x4[t])^2+(y3[t]-y4[t])^2+(z3[t]-z4[t])^2)^3]+
(G m5 (x5[t]-x4[t]))/Sqrt[((x5[t]-x4[t])^2+(y5[t]-y4[t])^2+(z5[t]-z4[t])^2)^3]+
 
(-q4*q1/(4Pi ε0 )/m4 (x1[t]-x4[t]))/Sqrt[((x1[t]-x4[t])^2+(y1[t]-y4[t])^2+(z1[t]-z4[t])^2)^3]+
(-q4*q2/(4Pi ε0 )/m4 (x2[t]-x4[t]))/Sqrt[((x2[t]-x4[t])^2+(y2[t]-y4[t])^2+(z2[t]-z4[t])^2)^3]+
(-q4*q3/(4Pi ε0 )/m4 (x3[t]-x4[t]))/Sqrt[((x3[t]-x4[t])^2+(y3[t]-y4[t])^2+(z3[t]-z4[t])^2)^3]+
(-q4*q5/(4Pi ε0 )/m4 (x5[t]-x4[t]))/Sqrt[((x5[t]-x4[t])^2+(y5[t]-y4[t])^2+(z5[t]-z4[t])^2)^3],
 
vy4'[t] ==
 
(G m1 (y1[t]-y4[t]))/Sqrt[((x1[t]-x4[t])^2+(y1[t]-y4[t])^2+(z1[t]-z4[t])^2)^3]+
(G m2 (y2[t]-y4[t]))/Sqrt[((x2[t]-x4[t])^2+(y2[t]-y4[t])^2+(z2[t]-z4[t])^2)^3]+
(G m3 (y3[t]-y4[t]))/Sqrt[((x3[t]-x4[t])^2+(y3[t]-y4[t])^2+(z3[t]-z4[t])^2)^3]+
(G m5 (y5[t]-y4[t]))/Sqrt[((x5[t]-x4[t])^2+(y5[t]-y4[t])^2+(z5[t]-z4[t])^2)^3]+
 
(-q4*q1/(4Pi ε0 )/m4 (y1[t]-y4[t]))/Sqrt[((x1[t]-x4[t])^2+(y1[t]-y4[t])^2+(z1[t]-z4[t])^2)^3]+
(-q4*q2/(4Pi ε0 )/m4 (y2[t]-y4[t]))/Sqrt[((x2[t]-x4[t])^2+(y2[t]-y4[t])^2+(z2[t]-z4[t])^2)^3]+
(-q4*q3/(4Pi ε0 )/m4 (y3[t]-y4[t]))/Sqrt[((x3[t]-x4[t])^2+(y3[t]-y4[t])^2+(z3[t]-z4[t])^2)^3]+
(-q4*q5/(4Pi ε0 )/m4 (y5[t]-y4[t]))/Sqrt[((x5[t]-x4[t])^2+(y5[t]-y4[t])^2+(z5[t]-z4[t])^2)^3],
 
vz4'[t] ==
 
(G m1 (z1[t]-z4[t]))/Sqrt[((x1[t]-x4[t])^2+(y1[t]-y4[t])^2+(z1[t]-z4[t])^2)^3]+
(G m2 (z2[t]-z4[t]))/Sqrt[((x2[t]-x4[t])^2+(y2[t]-y4[t])^2+(z2[t]-z4[t])^2)^3]+
(G m3 (z3[t]-z4[t]))/Sqrt[((x3[t]-x4[t])^2+(y3[t]-y4[t])^2+(z3[t]-z4[t])^2)^3]+
(G m5 (z5[t]-z4[t]))/Sqrt[((x5[t]-x4[t])^2+(y5[t]-y4[t])^2+(z5[t]-z4[t])^2)^3]+
 
(-q4*q1/(4Pi ε0 )/m4 (z1[t]-z4[t]))/Sqrt[((x1[t]-x4[t])^2+(y1[t]-y4[t])^2+(z1[t]-z4[t])^2)^3]+
(-q4*q2/(4Pi ε0 )/m4 (z2[t]-z4[t]))/Sqrt[((x2[t]-x4[t])^2+(y2[t]-y4[t])^2+(z2[t]-z4[t])^2)^3]+
(-q4*q3/(4Pi ε0 )/m4 (z3[t]-z4[t]))/Sqrt[((x3[t]-x4[t])^2+(y3[t]-y4[t])^2+(z3[t]-z4[t])^2)^3]+
(-q4*q5/(4Pi ε0 )/m4 (z5[t]-z4[t]))/Sqrt[((x5[t]-x4[t])^2+(y5[t]-y4[t])^2+(z5[t]-z4[t])^2)^3],
 
vx5'[t] ==
 
(G m1 (x1[t]-x5[t]))/Sqrt[((x1[t]-x5[t])^2+(y1[t]-y5[t])^2+(z1[t]-z5[t])^2)^3]+
(G m2 (x2[t]-x5[t]))/Sqrt[((x2[t]-x5[t])^2+(y2[t]-y5[t])^2+(z2[t]-z5[t])^2)^3]+
(G m3 (x3[t]-x5[t]))/Sqrt[((x3[t]-x5[t])^2+(y3[t]-y5[t])^2+(z3[t]-z5[t])^2)^3]+
(G m4 (x4[t]-x5[t]))/Sqrt[((x4[t]-x5[t])^2+(y4[t]-y5[t])^2+(z4[t]-z5[t])^2)^3]+
 
(-q5*q1/(5Pi ε0 )/m5 (x1[t]-x5[t]))/Sqrt[((x1[t]-x5[t])^2+(y1[t]-y5[t])^2+(z1[t]-z5[t])^2)^3]+
(-q5*q2/(5Pi ε0 )/m5 (x2[t]-x5[t]))/Sqrt[((x2[t]-x5[t])^2+(y2[t]-y5[t])^2+(z2[t]-z5[t])^2)^3]+
(-q5*q3/(5Pi ε0 )/m5 (x3[t]-x5[t]))/Sqrt[((x3[t]-x5[t])^2+(y3[t]-y5[t])^2+(z3[t]-z5[t])^2)^3]+
(-q5*q4/(5Pi ε0 )/m5 (x4[t]-x5[t]))/Sqrt[((x4[t]-x5[t])^2+(y4[t]-y5[t])^2+(z4[t]-z5[t])^2)^3],
 
vy5'[t] ==
 
(G m1 (y1[t]-y5[t]))/Sqrt[((x1[t]-x5[t])^2+(y1[t]-y5[t])^2+(z1[t]-z5[t])^2)^3]+
(G m2 (y2[t]-y5[t]))/Sqrt[((x2[t]-x5[t])^2+(y2[t]-y5[t])^2+(z2[t]-z5[t])^2)^3]+
(G m3 (y3[t]-y5[t]))/Sqrt[((x3[t]-x5[t])^2+(y3[t]-y5[t])^2+(z3[t]-z5[t])^2)^3]+
(G m4 (y4[t]-y5[t]))/Sqrt[((x4[t]-x5[t])^2+(y4[t]-y5[t])^2+(z4[t]-z5[t])^2)^3]+
 
(-q5*q1/(5Pi ε0 )/m5 (y1[t]-y5[t]))/Sqrt[((x1[t]-x5[t])^2+(y1[t]-y5[t])^2+(z1[t]-z5[t])^2)^3]+
(-q5*q2/(5Pi ε0 )/m5 (y2[t]-y5[t]))/Sqrt[((x2[t]-x5[t])^2+(y2[t]-y5[t])^2+(z2[t]-z5[t])^2)^3]+
(-q5*q3/(5Pi ε0 )/m5 (y3[t]-y5[t]))/Sqrt[((x3[t]-x5[t])^2+(y3[t]-y5[t])^2+(z3[t]-z5[t])^2)^3]+
(-q5*q4/(5Pi ε0 )/m5 (y4[t]-y5[t]))/Sqrt[((x4[t]-x5[t])^2+(y4[t]-y5[t])^2+(z4[t]-z5[t])^2)^3],
 
vz5'[t] ==
 
(G m1 (z1[t]-z5[t]))/Sqrt[((x1[t]-x5[t])^2+(y1[t]-y5[t])^2+(z1[t]-z5[t])^2)^3]+
(G m2 (z2[t]-z5[t]))/Sqrt[((x2[t]-x5[t])^2+(y2[t]-y5[t])^2+(z2[t]-z5[t])^2)^3]+
(G m3 (z3[t]-z5[t]))/Sqrt[((x3[t]-x5[t])^2+(y3[t]-y5[t])^2+(z3[t]-z5[t])^2)^3]+
(G m4 (z4[t]-z5[t]))/Sqrt[((x4[t]-x5[t])^2+(y4[t]-y5[t])^2+(z4[t]-z5[t])^2)^3]+
 
(-q5*q1/(5Pi ε0 )/m5 (z1[t]-z5[t]))/Sqrt[((x1[t]-x5[t])^2+(y1[t]-y5[t])^2+(z1[t]-z5[t])^2)^3]+
(-q5*q2/(5Pi ε0 )/m5 (z2[t]-z5[t]))/Sqrt[((x2[t]-x5[t])^2+(y2[t]-y5[t])^2+(z2[t]-z5[t])^2)^3]+
(-q5*q3/(5Pi ε0 )/m5 (z3[t]-z5[t]))/Sqrt[((x3[t]-x5[t])^2+(y3[t]-y5[t])^2+(z3[t]-z5[t])^2)^3]+
(-q5*q4/(5Pi ε0 )/m5 (z4[t]-z5[t]))/Sqrt[((x4[t]-x5[t])^2+(y4[t]-y5[t])^2+(z4[t]-z5[t])^2)^3],
 
x1[0] == x1x, y1[0] == y1y, z1[0] == z1z,
x2[0] == x2x, y2[0] == y2y, z2[0] == z2z,
x3[0] == x3x, y3[0] == y3y, z3[0] == z3z,
x4[0] == x4x, y4[0] == y4y, z4[0] == z4z,
x5[0] == x5x, y5[0] == y5y, z5[0] == z5z,
 
vx1[0] == v1x, vy1[0] == v1y, vz1[0] == v1z,
vx2[0] == v2x, vy2[0] == v2y, vz2[0] == v2z,
vx3[0] == v3x, vy3[0] == v3y, vz3[0] == v3z,
vx4[0] == v4x, vy4[0] == v4y, vz4[0] == v4z,
vx5[0] == v5x, vy5[0] == v5y, vz5[0] == v5z},
 
{x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5,
vx1, vx2, vx3, vx4, vx5, vy1, vy2, vy3, vy4, vy5, vz1, vz2, vz3, vz4, vz5},
 
{t, 0, Tmax},

WorkingPrecision-> wp,
MaxSteps-> Infinity,
Method-> mta,
InterpolationOrder-> All];
 
(* Position, Geschwindigkeit *)
 
f2p[t_]:={{x1[t], y1[t], z1[t]}, {x2[t], y2[t], z2[t]}, {x3[t], y3[t], z3[t]}, {x4[t], y4[t], z4[t]}, {x5[t], y5[t], z5[t]}}/.nds[[1]];
f2v[t_]:={{vx1[t], vy1[t], vz1[t]}, {vx2[t], vy2[t], vz2[t]}, {vx3[t], vy3[t], vz3[t]}, {vx4[t], vy4[t], vz4[t]}, {vx5[t], vy5[t], vz5[t]}}/.nds[[1]];
sp[T_] :=(m1 Evaluate[f2p[T][[1]]]+m2 Evaluate[f2p[T][[2]]]+m3 Evaluate[f2p[T][[3]]]+m4 Evaluate[f2p[T][[4]]]+m5 Evaluate[f2p[T][[5]]])/(m1+m2+m3+m4+m5)
 
(* Formatierung *)
 
s[text_]:=Style[text, FontSize->11];
sw[text_]:=Style[text, White, FontSize->11];
colorfunc[n_]:=Function[{x, y, z, t},
Hue[0, n, 0.5,
If[Tmax<0, Max[Min[(+T+(-t+trail))/trail, 1], 0],
Max[Min[(-T+(t+trail))/d1, 1], 0]]]];
 
(* Animation *)
 
Do[Print[Rasterize[
Grid[{{
Show[

If[T == 0, {},

ParametricPlot3D[Evaluate[f2p[T]],
{t, Max[0, T-trail], T},

PlotStyle->{
{Thickness[thk], Red},
{Thickness[thk], Blue},
{Thickness[thk], Green},
{Thickness[thk], Magenta},
{Thickness[thk], Cyan}},

PlotRange->plotrange, AspectRatio->1, MaxRecursion->15, Axes->True, ImageSize->imagesize]],
 
Graphics3D[{
 
{PointSize[2point/3], Lighter[Red],     Point[{x1x, y1y, z1z}]},
{PointSize[2point/3], Lighter[Blue],    Point[{x2x, y2y, z2z}]},
{PointSize[2point/3], Lighter[Green],   Point[{x3x, y3y, z3z}]},
{PointSize[2point/3], Lighter[Magenta], Point[{x4x, y4y, z4z}]},
{PointSize[2point/3], Lighter[Cyan],    Point[{x5x, y5y, z5z}]}},

PlotRange->plotrange, AspectRatio->1, Axes->True, ImageSize->imagesize],
 
Graphics3D[{PointSize[point], Red,      Point[With[{t=T}, Evaluate[f2p[T]]][[1]]]}],
Graphics3D[{PointSize[point], Blue,     Point[With[{t=T}, Evaluate[f2p[T]]][[2]]]}],
Graphics3D[{PointSize[point], Green,    Point[With[{t=T}, Evaluate[f2p[T]]][[3]]]}],
Graphics3D[{PointSize[point], Magenta,  Point[With[{t=T}, Evaluate[f2p[T]]][[4]]]}],
Graphics3D[{PointSize[point], Cyan,     Point[With[{t=T}, Evaluate[f2p[T]]][[5]]]}],
 
ViewPoint->{0, Infinity, 0}]},
 
{ },
{s["t"->N[T]]},
{ },
{s["p1{x,y,z}"-> With[{t=T}, Evaluate[f2p[T][[1]]]]],        sw[1/2]},
{s["v1{x,y,z}"-> With[{t=T}, Evaluate[f2v[T][[1]]]]],        sw[1/2]},
{s["v1{total}"->{With[{t=T}, Evaluate[Norm[f2v[T][[1]]]]]}], sw[1/2]},
{ },
{s["p2{x,y,z}"-> With[{t=T}, Evaluate[f2p[T][[2]]]]],        sw[1/2]},
{s["v2{x,y,z}"-> With[{t=T}, Evaluate[f2v[T][[2]]]]],        sw[1/2]},
{s["v2{total}"->{With[{t=T}, Evaluate[Norm[f2v[T][[2]]]]]}], sw[1/2]},
{ },
{s["p3{x,y,z}"-> With[{t=T}, Evaluate[f2p[T][[3]]]]],        sw[1/2]},
{s["v3{x,y,z}"-> With[{t=T}, Evaluate[f2v[T][[3]]]]],        sw[1/2]},
{s["v3{total}"->{With[{t=T}, Evaluate[Norm[f2v[T][[3]]]]]}], sw[1/2]},
{ },
{s["p4{x,y,z}"-> With[{t=T}, Evaluate[f2p[T][[4]]]]],        sw[1/2]},
{s["v4{x,y,z}"-> With[{t=T}, Evaluate[f2v[T][[4]]]]],        sw[1/2]},
{s["v4{total}"->{With[{t=T}, Evaluate[Norm[f2v[T][[4]]]]]}], sw[1/2]},
{ },
{s["p5{x,y,z}"-> With[{t=T}, Evaluate[f2p[T][[5]]]]],        sw[1/2]},
{s["v5{x,y,z}"-> With[{t=T}, Evaluate[f2v[T][[5]]]]],        sw[1/2]},
{s["v5{total}"->{With[{t=T}, Evaluate[Norm[f2v[T][[5]]]]]}], sw[1/2]},
{ },
{s["sp{x,y,z}"-> With[{t=T}, sp[T]]],                        sw[1/2]},
{s["sv{x,y,z}"-> With[{t=T}, sp'[T]]],                       sw[1/2]},
{s["sv{total}"->{With[{t=T}, Norm[sp'[T]]]}],                sw[1/2]}
}, Alignment->Left]]],
 
(* Zeitregler *)
 
{T, 0, Tmax, Tmax/5}]

(* Export als HTML Dokument *)
(* Export["dateiname.html", EvaluationNotebook[], "GraphicsOutput" -> "PNG"] *)
(* Export direkt als Bildsequenz *)
(* Do[Export["dateiname" <> ToString[T] <> ".png", Rasterize[...] ], {T, 0, 10, 5}] *)









Bild
Simon Tyran aka Симон Тыран @ vk || wikipedia || stackexchange || wolfram


Zurück zu „Yukterez Notizblock“

Wer ist online?

Mitglieder in diesem Forum: 0 Mitglieder und 3 Gäste