Was ist eigentlich Gravitation?

[Yukterez]

Solange die App einvernehmlich funzt, muss ausser ich den Code niemand verstehen.

Naja noch arbeitest du ja daran sie zu vervollkommnen, da könnte es schon nutzen wenn auch andere den Code nachvollziehen können.

Die 1 im 1-a auch nicht verstehend,

[Yukterez]

Also innerhalb des Codes einer Methode wüsste ich nicht, wie ich da die Einheiten mitschleifen soll, schließlich soll man ja nicht jede einzelne Zeile kommentieren.

Gannz einfach:

Code: Alles auswählen

(* Yukterez Code, Gravity Simulator *)

(* Einstellungen *)

M1 = mSol; M2 = 2 mSol; M3 = 3 mSol; M4 = 4 mSol; (* Massen *)

P1x = Au; P2x = Au; P3x = -Au; P4x = -Au;  (* Anfangspositionen *)
P1y = Au; P2y = -Au; P3y = Au; P4y = -Au;
P1z = 0; P2z = 0; P3z = 0; P4z = 0;

v1x = 0; v2x = 0; v3x = 0; v4x = 0; (* Anfangsgeschwindigkeiten *)
v1y = 0; v2y = 0; v3y = 0; v4y = 0;
v1z = 0; v2z = 0; v3z = 0; v4z = 0;

(* Eckdaten *)

G = 667384/10^16 m³/kg/sek²; (* Gravitationskonstante *)
Au = 149597870690 m; (* Astronomische Einheit *)
mSol = 1988*10^27 kg; (* Sonnenmasse *)

(* Dimensionen *)

kg=1; m=1; sek=1;

(* Container *)

Funktion[{{x1x_, y1y_, z1z_}, {x2x_, y2y_, z2z_}, {x3x_, y3y_,
    z3z_}, {x4x_, y4y_, z4z_}}, {{vx1x_, vy1y_, vz1z_}, {vx2x_, vy2y_,
     vz2z_}, {vx3x_, vy3y_, vz3z_}, {vx4x_, vy4y_, vz4z_}}, {m1_, m2_,
    m3_, m4_}, T_, plotType : ("x" | "v"), plotOptions___] :=
Module[{nds, Tmax, funcToPlot},
 
(* DGL *)

   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],
       
       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],
         
       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],
     
      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],
         
       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],
     
       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],
     
      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],
       
       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],
     
       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],   
     
       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],   
     
       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],     
     
       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],
     
      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] == vx1x, vy1[0] == vy1y, vz1[0] == vz1z,
      vx2[0] == vx2x, vy2[0] == vy2y, vz2[0] == vz2z,
      vx3[0] == vx3x, vy3[0] == vy3y, vz3[0] == vz3z,
      vx4[0] == vx4x, vy4[0] == vy4y, vz4[0] == vz4z},
     
     {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, T}];

(* Solve *)
   
   If[Head[nds] =!= NDSolve, Tmax = nds[[1, 1, 2, 1, 1, 2]];
   
    funcToPlot =
     If[plotType ===
        "x", {{x1[t], y1[t], z1[t]}, {x2[t], y2[t], z2[t]}, {x3[t],
         y3[t], z3[t]}, {x4[t], y4[t], z4[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]}}] /. nds[[1]];
   
(* Plot *)
     
    ParametricPlot3D[Evaluate[funcToPlot], {t, 0, Tmax},
     PlotStyle -> {{Red, Thick}, {Blue, Thick}, {Green,
        Thick}, {Orange, Thick}},
     
(* Plot Range *)
   
     PlotRange -> {{-2 Au, 2 Au}, {-2 Au, 2 Au}, {-2 Au, 2 Au}},
     AspectRatio -> 1, MaxRecursion -> ControlActive[3, 100],
     Axes -> True, plotOptions], Text["Yukterez Mod."]]] // Quiet

(* Manipulate *)
   
    plot = Manipulate[Show[Funktion[

   {{P1x, P1y, P1z}, {P2x, P2y, P2z}, {P3x, P3y, P3z}, {P4x, P4y, P4z}}, 
   {{v1x, v1y, v1z}, {v2x, v2y, v2z}, {v3x, v3y, v3z}, {v4x, v4y, v4z}},
   {M1, M2, M3, M4},
    T, xv,
   
    ImageSize -> {440, 440}],

(* Startpositionen *)
   
  Graphics3D[{
  {Red, Point[{P1x, P1y, P1z}]},
  {Blue, Point[{P2x, P2y, P2z}]},
  {Green, Point[{P3x, P3y, P3z}]},
  {Orange, Point[{P4x, P4y, P4z}]}}],

(* Perspektive *)

  ViewPoint -> {0, 0, 2}],
 
{{xv, "x", "Position"}, {"x" -> "Position"}},

(* Zeitregler *)
 
{{T, 1*^4, "Time"}, 1*^4, 1*10^7},

ControlPlacement -> {Bottom}]

Export["plot.avi", plot, "FrameRate" -> 75]

Zumindest hoffend dass man hier trotz relativ unleserlichem .txt Format die Einheiten alle nachvollziehen kann,

[contravariant]

Deswegen schrieb ich ja, dass die Anwendung einvernehmlich laufen muss. Wenn nur ich damit einverstanden wär, könnte ich sie mir auch an die Backe tackern. Deswegen gibt es inzwischen ja auch schon zwei davon und ihr (jedenfalls die, die sie zum Laufen bekommen) sollt mir sagen, welche davon richtiger rechnet. Das es noch nicht komplett richtig ist hatte ich ja bereits erwähnt, zumindest kam es schon mal bei Yukterez so an.

Beide rechnen falsch, da in beiden die Berechnung der neuen Positionen der Massen denselben Fehler enthalten.

[contravariant]

Diese Aussage ist unbegründet, solange du dich von der Rechnung nicht überzeugt, also die Anwendung laufen lassen hast. Das du den Code nicht nachvollziehen kannst ist nach wie vor nicht mein Problem. Mein letztes Wort zu diesem Thema.

Wenn die Physik, die einem solchen Code zugrunde liegen sollte, dich nicht interessiert, bitte.

[McDaniel-77]

Endlich stimmt's Spacerat,

bis auf die rote Masse. Diese bewegt sich kaum, auch wenn sie genauso schwer wie die blaue Masse ist.

Die großen Massen sind schneller als die kleine Masse. Die kleinste Masse fällt am langsamsten und ist deshalb meistens zu spät.

McDaniel-77

[Yukterez]

Das ist jetzt aber ein bißchen sehr widersprüchlich: wenn es stimmt, kann das von dir beschriebene Verhalten nicht eintreten, und wenn es so ist wie du beschreibst, dann kann es nicht stimmen.

Doch nicht beides zugleich haben könnend,

[Yukterez]

Der Fehler hat aber nichts mit Einheiten oderso zu tun... Ich habe mich schlicht beim Berechnen der gemeinsamen Massezentren verhauen.

Berechne doch einfach die Richtungskomponenten für x, y und z getrennt! Das zeitabhängige Differential



aus meiner Formel ist in abstandsabhängiger Notation nichts anderes als

,

wo links nur die x,y,z-Koordinatenversion von G·M/r² steht, und rechts die Richtungskomponente nach Pythagoras (hier für den Körper 1 auf der Achse x). Das gleiche machst du dann noch für die y und z-Achse, dann steht rechts einfach y oder z statt x im Zähler.

Ich bin mir nicht sicher, ob sich das Drei- und Mehrkörperproblem so ganz ohne Pythagoras überhaupt lösen lässt, deinen Ansatz mit den Vektoren kann ich leider nicht nachvollziehen (aber da bist du z.T. auch dran schuld weil nicht in jeder Zeile deines Codes ein Kommentar steht :)

,

[Yukterez]

Was also noch fehlt ist die Masse, die an diesem Punkt wirkt (Wirkmasse). Das kann doch schlecht m1+m2 sein oder doch?

Für den ersten Körper mit der Masse m1 ist die Wirkmasse m2 und m3, also rechnest du mit G·m2/{r_(m1→m2)}²×Vektor_(m1→m2) und G·m3/{r_(m1→m3)}²×Vektor_(m1→m3). Dann hast du für die Masse m1 zwei verschieden starke Vektoren, die du dann noch zusammenrechnen musst. Wie du das ohne x,y,z-Komponenten als Basis machst musst du allerdings selber wissen, wie man das in Java programmiert habe ich leider keine Ahnung.

Hoffend, mit "Wirkmasse" nicht was anderes zu meinen als du,

[Yukterez]

Nee, meinst du nicht... Habe es inzwischen schon selbst gesehen. Hatte wohl einen Gemüsegarten auf den Augen.

Dafür hatte ich einen solchen anscheinend im Gehirn: statt

Yukterez hat geschrieben:
dann steht rechts einfach y oder z statt x im Nenner

meinte ich natürlich

Yukterez wollte damit sagen:
dann steht rechts einfach y oder z statt x im Zähler

Annehmend, dass das aber sowieso aus dem Kontext heraus ersichtlich war,

[Mikesch]

Jede einzelne Zeile kommentieren ist dagegen sehr schlechter Stil. Aber das sind CodingConventions und sollen hier nicht interessieren. Mal sehen, wenn ich Zeit und Lust habe oder gar nicht weiter komme, werde ich es wohl doch noch kommentieren.

Das zeigt doch sehr schön den Unterschied zwischen Profis und Laien. Einer der ersten Qualitäts- und Testkriterien ist das Vorliegen einer anständigen Dokumentation und notfalls erfolgt dies eben auch zeilenweise. Die Tests werden durch unabhängige Tester durchgezogen, die müssen anhand des Quellcodes und der Kommentierung und weiterer Dokumentationen (UML bspw), auch ohne das der Programmierer daneben sitzt, grundsätzlich erkennen, was geplant war und was tatsächlich umgesetzt wurde. Ein Try-and-Error ist ganz schlechter Entwicklungsstil. Sie hätten noch nicht einmal die automatischen Tests bestanden, um ihren Code in eine Testumgebung schieben zu dürfen. Eine Methodik, wie z.B. Pair-Programming, könnte bei ihrem Stil (wenn man das als Stil bezeichnen möchte) nicht angewendet werden.
Wenn Sie meinen, das Sie eine vernünftige Kommentierung nicht nötig haben, und das der Leser gefälligst selber darauf kommen soll, was Sie gemeint haben könnten und dass der Code selbstverständlich fehlerfrei und eine Transparenz daher völlig unnötig ist, dann zeugt das von einer Arroganz dem Gesprächspartner gegenüber und das auch noch gepaart bei völliger Selbstüberschätzung Ihrer Fähigkeiten.

Mike

PS:
Zum Thema Ausführungsgeschwindigkeit von diversen Programmiersprachen, auch eines ihrer Null-Ahnung-Themen, schauen Sie nach bei http://benchmarksgame.alioth.debian.org/. Dort sehen Sie, wie man dokumentiert. Nehmen Sie sich das als Vorbild.