kleiner planetarer Simulator

Kibo

Registriertes Mitglied
Wie im anderen Thread besprochen, eröffne ich einen neuen Thread um keinen zu stören.
Also es gibt eine neue Version (1.9.1)

Hier die Neuerungen:

Die Buttons wurden durch Kreise ersetzt was viele Vorteile hat.
Die Initialisierung mit dynamischen Objekten ist wesentlich einfacher, das findet sich im Code noch nicht wieder.
Neben dem ästhetischen Faktor wird über den Kugelradius auch die Masse der Planeten wiedergegeben. Das läuft über die nicht ganz realitätsgetreue Formel.:
r :=sqrt(m)+5

Das Programm läuft bis auf Weiteres noch in 2D und mit GDI+. Der Downloadlink ist immer noch der Gleiche über Skydrive


mfg
 
Zuletzt bearbeitet:

Bernhard

Registriertes Mitglied
Also es gibt eine neue Version (1.9.1)
Hi Kibo,

nettes Programm. Man sieht sogar eine Perihelbewegung des innersten Planeten. Die rote Linie gehört vermutlich also zu Merkur? Hast Du den BS-Integrator verwendet?

Das läuft über die nicht ganz realitätsgetreue Formel.: r :=sqrt(m)+5
Die Sonne würde ich etwas kleiner darstellen, weil die Merkurbahn schon fast die "Sonnenoberfläche" streift.
Gruß
 

Kibo

Registriertes Mitglied
Hallo Bernhard,

Der Code für die physikalischen Berechnungen ist komplett selbst geschrieben, ich wollte einfach meine Kenntnisse auf dem Gebiet auffrischen und habe mir daher meine eigenen Gedanken gemacht weil man dabei mehr lernt.
Wäre ich auf größere Hürden gestoßen, hätte ich mir natürlich euren Quellcode genauer angeschaut. Wenn du möchtest, schicke ich dir mal meinen Code zu.
Was die Größenverhältnisse angeht, steht da noch nicht fest wie ich es am Ende haben möchte. Die Körper sollen halt alle gut sichtbar sein aber eben auch die Massenverhältnisse wiederspiegeln.

Das Beispielsystem soll nicht unser Sonnensystem darstellen, ich habe einfach ein bisschen mit Zahlen herumgespielt bis ein halbwegs stabiles System da stand.
Die Werte für die Planeten können in der Datei Planeten.txt nach belieben angepasst werden. Man kann auch noch mehr Objekte hinzufügen (wobei das noch nicht ganz nach meiner Zufriedenheit funktioniert).

mfg
 
Zuletzt bearbeitet:

Kibo

Registriertes Mitglied
Nun ich habe 2 Schleifen, dort nehme ich einen Planeten, berechne Entfernung, Masse und Richtung zu einem anderen Planeten und addiere diesen Vektor mit dem schon vorhandenen Impuls. Das mache ich so lange bis jeder Körper zu jedem anderen in Beziehung steht. Das ist so in etwa der Kern.
 

Kibo

Registriertes Mitglied
Es gibt ein etwas größeres Update auf Version 1.10!

Es wurde ein kleines Formular zur direkten Manipulation der Himmelskörper integriert, ausserdem kann man jetzt mit den hoch und runter Tasten zoomen.

Die Änderungen gehen mit einem erheblichen Gewinn an Benutzerfreundlichkeit (und Spaß ) einher, von daher kann ich den Download nur wärmstens empfehlen!

Downloadlink

mfg
 

Kibo

Registriertes Mitglied
Guten Abend,

Es gibt mal wieder ein Update, Ich bin mitlerweile bei Version 1.12 angelangt.
Diese ist nun endlich in der Lage mit einer beliebigen Anzahl an Testkörpern umzugehen. Ich habe den Physikteil erheblich verbessert, sodass sich keine Ungenauigkeiten aufgrund von veränderten Planetenanzahlen ergeben sollten. Ich habe soweit Ich das sehe nun alle Bugs behoben, damit ist das jetzt der stabile Release von Version 1.
Das ist auch der letzte Release der GDI verwenden wird für Version 2 in 3d steige ich dann auf OpenGl um.

Falsl jemand Lust hat, würde Ich gerne den Integrator besprechen, Ich wäre interessiert mehr zu dem Thema zu erfahren, mein jetziger Integrator ist ja komplett von mir selbst zusammen geschustert ohne das ich dazu auch nur eine Seite Theorie zu gelesen habe, daher weiß ich nicht einmal auf welchem Niveau an Genauigkeit und Schnelligkeit der nun arbeitet.
 

Bernhard

Registriertes Mitglied
Falsl jemand Lust hat, würde Ich gerne den Integrator besprechen, Ich wäre interessiert mehr zu dem Thema zu erfahren, mein jetziger Integrator ist ja komplett von mir selbst zusammen geschustert ohne das ich dazu auch nur eine Seite Theorie zu gelesen habe, daher weiß ich nicht einmal auf welchem Niveau an Genauigkeit und Schnelligkeit der nun arbeitet.
Hi Kibo,

ohne explizite Vorkenntnisse hast Du ziemlich sicher das Eulersche Polygonzugverfahren implementiert. Sobald man weiß, was die Ableitung einer Funktion einer Veränderlichen ist, kann man das mehr oder weniger intuitiv verstehen. Etwas ausgefeilter ist das Runge-Kutta-Verfahren. Versuche also erst mal Dein aktuell eingesetztes Verfahren anhand dieser Links zu verbessern.
 

Kibo

Registriertes Mitglied
Hallo Bernhard,

Ich weiß nicht so genau ob das Eulersche Polygonzugverfahren meinem entspricht. Ich fasse einfach mal mein Verfahren zusammen, sonst kann man ja schlecht drüber reden.

Also ich habe 2 Schleifen, die dienen zum bestimmen welcher Planet gerade mit der Berechnung dranne ist. Zum beispiel wird gerade der Planet a von b beinflusst.
Ich speichere zuerst den jetzigen Impuls von a (liegt in a.ix und iy) in die Hilfsvariablen a.hx und a.hy, wichtig später für die Trägheit. Dan kommt die Gravitationsberechnung die b auf a ausübt, die läuft über eine Methode die in der Planetenklasse steht. Ich bestimme den Vektor zwischen a und b mit
this.xd:=gegner.x-this.x
this.yd:=gegner.y-this.y

Daraus berechne ich den Abstand zwischen a und b.
abstand:=sqrt(this.xd*this.xd+this.yd*this.yd)

aus dem Vektor und dem Abstand kann na dann leicht den Einheitsvektor bestimmen.
this.xn:=this.xd/abstand
this.yn:=this.yd/abstand

Ich berechne dann die Stärke der Gravitation aus dem Abstand und der Masse der Planeten
this.s:=0.667384*gegner.m*this.m/(abstand*abstand)

Als Nächste schreibe ich die Gravitation noch in einen Vektor...
this.gx:=this.s*this.xn
this.gy:=this.s*this.yn

und addiere den neuen Impuls auf den alten auf.
%Planeten%.ix := %Planeten%.hx + %Planeten%.gx /(genauigkeit * %Planeten%.m)
%Planeten%.iy := %Planeten%.hy + %Planeten%.gy /(genauigkeit * %Planeten%.m)
Genauigkeit soll hier die Schrittweite bestimmen. %Planeten% ist der Dynamische Aufruf für den jetzigen Planeten also a (innerhalb der Klasse war das dann this für a und Gegner für b)

das war der inhalt der inneren Schleife, die wird nun für jeden anderen Planeten noch für a ausgeführt und am Schluss, wenn man alle Gravitationseinflüsse für a berechnet hat, werden die neuen Koordinaten für a bestimmt
%Planeten%.x:=%Planeten%.x+ %Planeten%.ix/genauigkeit
%Planeten%.y:=%Planeten%.y+ %Planeten%.iy/genauigkeit

das Ganze wird durch die äußere Schleife noch für alle anderen in die Planeten.txt eingetragenen Planeten wiederholt.

Wäre nett wenn du mal drüber schauen würdest Bernhard:)
 

Bernhard

Registriertes Mitglied
aus dem Vektor und dem Abstand kann na dann leicht den Einheitsvektor bestimmen.
this.xn:=this.xd/abstand
this.yn:=this.yd/abstand
Hallo Kibo,

den Einheitsvektor braucht man eigentlich gar nicht. Besser man speichert die Differenz zwischen dem Ort des aktuellen Planeten und dem Ort des "Gegners" in einem Hilfsvektor h. Dann berechnet man über das Skalarprodukt dieses Vektors mit sich selbst das Quadrat des Abstandes sqr_abstand und zieht dann die Wurzel aus diesem Wert, um zuletzt den Abstand zu bekommen.

Ich berechne dann die Stärke der Gravitation aus dem Abstand und der Masse der Planeten
this.s:=0.667384*gegner.m*this.m/(abstand*abstand)
Mit den eben erwähnten drei Variablen sieht das dann folgendermaßen aus:

Beschleunigung(Planet) = 0.667384 * gegner.m * h / (sqr_abstand * abstand)

Man bekommt dann die Beschleunigung gleich als Vektor und vermeidet auch, dass man aus einer Größe zuerst die Wurzel zieht und diese dann wieder quadriert.

ZUsätzlich kann man anstelle von der Masse eines Planeten mit der kombinierten Variablen 0.667384 * gegner.m rechnen. Die Gravitationskonstante wird also einfach in die Masse fest rein multipliziert. Das spart haufenweise überflüssige Multiplikationen. Die Masse alleine wird im gesamten Programm nämlich gar nicht benötigt.

Als Nächste schreibe ich die Gravitation noch in einen Vektor...
this.gx:=this.s*this.xn
this.gy:=this.s*this.yn
Dieser Schritt entfällt dann.

und addiere den neuen Impuls auf den alten auf.
%Planeten%.ix := %Planeten%.hx + %Planeten%.gx /(genauigkeit * %Planeten%.m)
%Planeten%.iy := %Planeten%.hy + %Planeten%.gy /(genauigkeit * %Planeten%.m)
Genauigkeit soll hier die Schrittweite bestimmen. %Planeten% ist der Dynamische Aufruf für den jetzigen Planeten also a (innerhalb der Klasse war das dann this für a und Gegner für b)
Die Genauigkeit ist offensichtlich das Inverse der Schrittweite dt. Hier würde ich besser direkt mit der Schrittweite dt im Zähler rechnen. Das macht den Code übersichtlicher und vermeidet erneut unnötiges Ausrechnen. Wenn Du dem Anwender als Parameter lieber das Inverse der Schrittweite geben willst, musst Du eben beim Start des Programmes einmal den Kehrwert berechnen, um dt zu kennen.

das war der inhalt der inneren Schleife, die wird nun für jeden anderen Planeten noch für a ausgeführt und am Schluss, wenn man alle Gravitationseinflüsse für a berechnet hat, werden die neuen Koordinaten für a bestimmt
%Planeten%.x:=%Planeten%.x+ %Planeten%.ix/genauigkeit
%Planeten%.y:=%Planeten%.y+ %Planeten%.iy/genauigkeit
Das ist das genau das Euler-Verfahren, falls genauigkeit = 1.0 / dt.

das Ganze wird durch die äußere Schleife noch für alle anderen in die Planeten.txt eingetragenen Planeten wiederholt.
Hier fehlt noch das dritte newtonsche Gesetz (actio gegengleich reactio). Du berechnest in jedem Schritt sowohl actio, als auch reactio getrennt. Durch geschicktes Ordnen der Schleife berechnet man einmal actio (oder reactio) und verwendet diesen Wert dann sowohl mit positivem als auch mit negativem Vorzeichen. Man braucht dazu eigentlich nur für jeden Planeten eine vektorielle Variable für die Beschleunigung und läßt die äußere Schleife immer ab dem aktuellen Planeten starten. Die äußere Schleife wird dann mit fortschreitender innerer Schleife immer kleiner und man berechnet keine Zwischenergebnisse mehr doppelt.

Sämtliche Tips (mit unterschiedlichem Schwierigkeitsgrad) kann man in den genaueren Verfahren (Runge-Kutte, Bulirsch-Stoer, Symplektischer Integrator) beibehalten.
Gruß
 

Kibo

Registriertes Mitglied
Hallo Bernhard,

Erst einmal danke für die Tipps. Die Masse außerhalb der Schleife mit der Gravitationskonstante zu verrechnen ist kein Problem, auch die kleine Anpassung mit der Schrittweite stellt mich vor keinem Hindernis, das mit der Beschleunigung direkt als Vektor ist allerdings neu für mich. Wie sieht denn so ein Beschleunigungsvektor aus, ist das ein Array?
Hab ich deinen letzten Tipp richtig verstanden? Du berechnest quasi gleichzeitig die Beeinflussung von Planet a durch b und b durch a ja? Das leuchtet mir ja ein, aber wozu ist da zwingend eine vektorielle Beschleunigungsvariable nötig.

Ich setze erstmal deine Tipps um, so wie ich sie verstanden habe, parallel dann zu OpenGl und Kollisionsdetektion.

mfg
 
Zuletzt bearbeitet:

Bernhard

Registriertes Mitglied
Wie sieht denn so ein Beschleunigungsvektor aus, ist das ein Array?
Von der Mathematik her sind das nur drei reelle Zahlen für die drei Raumdimensionen. SRMeister hatte eine eigene Klasse PVektor dazu geschrieben. Wie Du das datentechnisch im Detail umsetzt, bleibt eigentlich Dir überlassen.

Hab ich deinen letzten Tip richtig verstanden? Du berechnest quasi gleichzeitig die Beeinflussung von Planet a durch b und b durch a ja?
Genau. Es geht allerdings hauptsächlich darum Zwischenergebnisse nicht doppelt zu berechnen, da sonst das Programm unnötig langsam wird.

Das leuchtet mir ja ein, aber wozu ist da zwingend eine vektorielleBeschleunigungsvariable nötig.
Der Beschleunigungsvektor dient nur der Übersichtlichkeit. Du musst das nicht so machen.
Gruß
 

SRMeister

Registriertes Mitglied
Das Programm gefällt mir auch.
So wie du, hatte ich auch zu Beginn den intuitiven Ansatz gewählt und mich Schritt für Schritt weiter vertieft. Auch ich/wir haben nicht alles an einem Tag geschafft, also, es ist ganz normal, wenn man nicht alles sofort umsetzen kann. (Bernhard hat dir immerhin die komprimierte Schlussfolgerung aus mehreren Monaten gegeben)
Eine Anmerkung möchte ich noch zum Eulerverfahren machen: Prinzipbedingt summieren sich die Ungenauigkeiten "in eine Richtung" - das Bedeutet, die Fehler sind immer positiv und addieren sich ständig. Besser ist ein Verfahren welches zufällig Fehler verursacht, die mal positiv und mal negativ sind, damit sie sich in Summe eher aufheben.

So ein Verfahren ist das Midpoint-Verfahren und das sollte auch dein nächster Schritt sein. Hast du die Verallgemeinerung von Euler zu Midpoint verstanden/umgesetzt, sind die anderen Integratoren, wie Runge-Kutta und Bulirsch-Stoer nicht mehr weit entfernt.

Viel Spaß
 

Kibo

Registriertes Mitglied
Danke SRMeister,

das wäre eh mein nächster Schritt gewesen, das Verfahren zu wechseln. Nun weiß ich dann auch welches Verfahren ich als Nächstes nehme.

mfg
 

Bernhard

Registriertes Mitglied
Besser ist ein Verfahren welches zufällig Fehler verursacht, die mal positiv und mal negativ sind, damit sie sich in Summe eher aufheben.
Die Gesamtenergie erfährt beim Mittelpunktsverfahren angeblich ebenfalls eine Fehlerdrift, weswegen ein symplektischer Integrator zusätzliche Vorteile bringt. Diese Integratoren besitzen, so wie die anderen Verfahren auch, eine bestimmte Fehlerordnung, weswegen man auch beim symplektischen Integrator ein Extrapolationspolynom verwenden kann (s.a. Richardson-Extrapolation).

EDIT: Um "Berührungsängste" mit dem symplektischen Integrator zu vermeiden kann man sich vergegenwärtigen, dass es sich dabei auch nur um sehr geschickte Differenzenverfahren handelt.
 
Zuletzt bearbeitet:

Kibo

Registriertes Mitglied
Guten Abend,

Ich habe mich dieses Wochenende mal wieder kurz ran setzen können.
Das ist dabei raus gekommen. Wie man sieht, jetzt mit Open-Gl in 3 Dimensionen. Die Texturierung funktioniert noch nicht wie sie soll und bei der Platzierung der Sphären scheint es auch noch einen Fehler zu geben. Des weiteren gab es ein paar kleinere Probleme mit dem klassischen Runge-Kutta, weswegen in der Version das gute alte stabile Euler-verfahren werkelt, auf Grundlage von Bernhards Tips angepasst. Die Kamerasteuerung war für mich daher noch zweitrangig weswegen man mit WSAD nur zoomen und drehen kann (eher weniger als mehr bugfrei). Wie man sieht, hab ich im Moment nicht so viel Zeit, sobald mich die Muse küsst werden die aufgezählten Mängel demnächst behoben.

mfg
 

Kibo

Registriertes Mitglied
Realitätstreue wäre programmiertechnisch kein Problem, blos dann erscheinen die Planeten im Vergleich zu den unendlichen Weiten des Weltalls doch wirklich sehr klein, nämlich um die Pixelgröße. Vielleicht bau ich noch einen Regler im GUI für ein. Wie gesagt, ich programmier das Teil nur so zum Spaß ohne ein bestimmtes Motiv dahinter.
Ich komm im Moment mit Quadric Texturen nicht wirklich klar, hab halt keine Ahnung von OpenGl ^^''
Habe im Moment relativ wenig Freizeit, so dass fürs Programmieren eher kaum was am Tag übrig bleibt. Ausserdem arbeite ich gerade noch an einer Client-Server Geschichte, die ist Momentan interessanter, von daher ruht der Simu jetzt etwas. Guter Hefeteig braucht ja auch viel Zeit zum aufgehen.:)

mfg
 

Kibo

Registriertes Mitglied
Hallo,

Ich bin mal wieder dabei meinen Simulator ein bisschen aufzupolieren und rechne gerade am Beispielsystem Erde-Mond.

Bei mir hakt es im Moment an der Formel der Bahngeschwindigkeit.

wurzel((GM+m)(2/r-1/a))
G=Gravitationskonstante
M=Masse Erde
m= Masse Mond
r=a= große Halbachse

ich setze ein: G*(5,974*10^24+7,349*10^22)*(1/384400)= 32403

Richtige Lösung wäre meiner Meinung nach 1,02 km/s bzw 3672 km/h. Sieht wer meinen Fehler?

mfg
 
Oben