kleiner planetarer Simulator

[Kibo]

Stimmt, danke

[Kibo]

Hmm erstaunlicherweise werden beim Runge Kutta die Bahnen ganz schön instabil im vergleich zum Euler.

Wer mal testen will hier die 2 verschiedenen Versionen

Mirrors

Rapidshare
Runge Kutta
Euler
Planeten Datei für Erde-Mond


Skydrive
Alle 3 Dateien

Anleitung: Jeweilige exe und die Planeten.txt im selben Ordner herunterladen und exe starten. Wenn die Planeten.txt für Erde-Mond nicht bei der Exe Datei mit zu liegt, wird ein Standartsystem erstellt. Die Szenerie wird mit WASD QY und den Pfeiltasten gesteuert.

mfg

[Bernhard]

Hallo Kibo,

Hmm erstaunlicherweise werden beim Runge Kutta die Bahnen ganz schön instabil im vergleich zum Euler.

ich vermute da eher einen Programmierfehler. Poste hier, bei Bedarf, die entscheidenden Stellen deines Quelltextes. Dann kann man sich das genauer ansehen.

Runge Kutta
Euler
Planeten Datei für Erde-Mond

Diese drei Links funktionieren bei mir nicht.
MfG

[Kibo]

Ok, Ist gefixt. Da ich die Links nicht mehr ändern kann muss ich sie halt noch mal neu rein stellen:

Rapidshare

Euler

Runge Kutta

Planeten.txt

mfg

[Bernhard]

Ich komme mit der Bedienung des Programmes überhaupt nicht klar.
In welchen physikalischen Einheiten arbeitet das Programm?
Wie kann man die Startbedingungen vorgeben?
Was bedeuten Zeitangaben in Millisekunden (!)?
MfG

EDIT: Hier noch ein Link auf eine vorbildliche Lösung deiner Aufgabenstellung: http://phet.colorado.edu/sims/my-sol...system_en.html .

[Kibo]

Hallo Bernhard,

Die Startbedingungen gibt man über die Datei Planeten.txt vor:

Koordinaten in x y z Impuls in x y z Masse s v Dichte des Planeten

Die Koordinaten sind in Kilometern anzugeben, der Impuls gibt quasi Startgeschwindigkeit und Richtung eines Planeten vor und bezieht sich im Moment noch auf das Kilometer/Berechnungsschritt. Die Masse sollte in Kilogramm angeben werden und ist im Moment noch der Übersichtlichkeit halber um 20 Nullen gekürzt, man kann auch solche Werte wie 5.7e11 eingeben, daher wird das in der nächsten Version angepasst. Die werte s v bitte auf 0 lassen, das funktioniert noch nicht so richtig. Die Dichte wird in g/cm³ angegeben.

Das Formular im Programm ist analog zu der Textdatei zu betrachten, und ist gedacht um die aktuell berechneten Werte abzulesen und anzupassen.

Das Feature Maus funktioniert nicht mehr, bitte ignorieren.
Der Slider Genauigkeit ermöglicht das Anpassen der Schrittweite, Der Slider geht von 1-100 und wird im Programm durch 5 gerechnet und dann als Kehrwert für die Schrittweite genommen. Das bedeutet Geneauigkeit 1 entspricht einer Schrittweite von 5 Einheiten (Kilometern?) und 100 entspricht der Schrittweite 1/ (100/5) 0.05.

Für die schlechte Beschriftung möchte ich mich entschuldigen, das Programm ist ja noch in Arbeit.

mfg

[Bernhard]

Hallo Kibo,

Die Startbedingungen gibt man über die Datei Planeten.txt vor:

die hatte mir gefehlt und damit wird jetzt alles gleich viel verständlicher . Nett wäre, wenn die Programme genau diese Datei selbst generieren, falls sie fehlt. Die aktuelle Default-Datei liefert ein Szenario, das eher an die Ergebnisse eines Teilchenbeschleunigers erinnert.

Die Masse sollte in Kilogramm angeben werden und ist im Moment noch der Übersichtlichkeit halber um 20 Nullen gekürzt, man kann auch solche Werte wie 5.7e11 eingeben, daher wird das in der nächsten Version angepasst.

Solche Spezialitäten solltest Du unbedingt vermeiden, da man die bei längeren Pausen sehr gerne vergisst. Es sei denn, Du machst Dir sehr genaue Notizen. Mir ist die letzten Tage selbst aufgefallen, wie wichtig Kommentare gerade in privaten Quelltexten sind.

Die werte s v bitte auf 0 lassen, das funktioniert noch nicht so richtig. Die Dichte wird in g/cm³ angegeben.

Das Formular im Programm ist analog zu der Textdatei zu betrachten, und ist gedacht um die aktuell berechneten Werte abzulesen und anzupassen.

Das Feature Maus funktioniert nicht mehr, bitte ignorieren.
Der Slider Genauigkeit ermöglicht das Anpassen der Schrittweite, Der Slider geht von 1-100 und wird im Programm durch 5 gerechnet und dann als Kehrwert für die Schrittweite genommen. Das bedeutet Geneauigkeit 1 entspricht einer Schrittweite von 5 Einheiten (Kilometern?) und 100 entspricht der Schrittweite 1/ (100/5) 0.05.

Für die schlechte Beschriftung möchte ich mich entschuldigen, das Programm ist ja noch in Arbeit.

OK, Rest ist erst mal Baustelle.

Jetzt zum Programm:
Bei dieser Auflösung (=visuell) und einem Umlauf sollte das Euler-Verfahren und das Runge-Kutta-Verfahren genau das gleiche Ergebnis berechnen. Du solltest den Quelltext des Runge-Kutta-Verfahrens nach kleinen aber wichtigen Fehlern durchsuchen. Wenn Du nicht weiterkommst, können wir den ersten Zeitschritt (ein Delta t) hier diskutieren, denn den kann man auch mit dem Taschenrechner nachrechnen.
MfG

[Bernhard]

Hallo Kibo,

um den Fehler in deinem Programm zu finden schreibe ich erst mal die DGL an (leider ohne LaTeX):

\dot{\vec{x}_1} = \vec{v}_1
\dot{\vec{v}_1} = G m_2 \frac{\vec{x_2}-\vec{x_1}}{\| \vec{x_2}-\vec{x_1} \|^3}

\dot{\vec{x}_2} = \vec{v}_2
\dot{\vec{v}_2} = G m_1 \frac{\vec{x_1}-\vec{x_2}}{\| \vec{x_1}-\vec{x_2} \|^3}

Es ist eine gewöhnliche DGL 1. Ordnung mit den 12 Variablen x_1, v_1, x_2 und v_2. Der erste Schritt berechnet den Wert dieser Variablen nach der Zeit dt. Bei einem ersten Test, um den Fehler zu finden, würde ich dt vergleichsweise groß wählen, z.B. 1 Tag = 86400 Sekunden.

Man kann dann die vier Hilfsvariablen k_1, k_2, k_3 und k_4 des Runge-Kutta-Verfahrens berechnen und testen, ob das Programm diese Variablen korrekt berechnet.

Beim ersten Überfliegen der Gleichungen sollte man bemerken, dass in allen vier Teilgleichungen der DGL die Gravitationskonstante immer als Produkt mit der Masse der Erde oder der Masse des Mondes vorkommt. Man kann die Gravitationskonstante also mit diesen zwei Massen schon vor der eigentlichen Berechnung multiplizieren und das Ergebnis als Masse verwenden. Die Gravitationskonstante taucht im Programm dann vorerst gar nicht mehr auf, was die Fehlersuche etwas vereinfacht.

Dann sollte man noch bemerken, dass man direkt mit Orts- und Geschwindigkeitsvektoren rechnen kann. Der Umweg über eine Impulsvariable ist momentan eher hinderlich. Die Werte aus der Datei planeten.txt kann ich zudem nicht nachvollziehen. Bei einer mittleren Orbitalgeschwindigkeit des Mondes von 1 km/s und einer Masse von 7.xeXX würde ich einen Wert von 7.xeXX oder etwas vergleichbares als Impuls erwarten.
MfG

[Kibo]

Hallo Bernhard,

Also deine Differentialgleichung kommt mir irgendwie spanisch vor Ich poste einfach mal den unterschiedlichen code. Ich hab ihn mal vonn Loops und Klassen befreit so dass man nicht viel umdenken muss, das ändert nichts am Ergebnis. Kommentare sind mit ; abgetrennt

Euler:

%Planeten%.hx:= %Planeten%.ix ; speichern des alten Impulses für trägheit
%Planeten%.hy:= %Planeten%.iy
%Planeten%.hz:= %Planeten%.iz

;steht in klasse this = planet um den es geht, gegner= beeinflussender planet

this.xd:=gegner.x-this.x ;Berechnung des Abstands zwischen den Planeten
this.yd:=gegner.y-this.y
this.zd:=gegner.z-this.z

abstand:=sqrt(this.xd*this.xd+this.yd*this.yd+this .zd*this.zd) ;bis hier

this.xn:=this.xd/abstand ; Einheitsvektor
this.yn:=this.yd/abstand
this.zn:=this.zd/abstand

this.s:=gegner.gm/(abstand*abstand) ; Berechnen der einwirkinden Gravitationskräfte

this.gx:=this.s*this.xn ; Verechnung der Grav mit der Richtung
this.gy:=this.s*this.yn
this.gz:=this.s*this.zn

;Ende Klassenteil

if (plzB=1)
%Planeten%.abstand:=abstand ; Speichert den Abstand zum Hauptkörper (idR Stern)


;aufaddieren des berechneten kraftvektors auf den vorhandenen impuls
%Planeten%.ix := %Planeten%.hx + %Planeten%.gx/genauigkeit
%Planeten%.iy := %Planeten%.hy + %Planeten%.gy/genauigkeit
%Planeten%.iz := %Planeten%.hz + %Planeten%.gz/genauigkeit

;Oben = Innerer schleifenteil, wird für jeden anderen Körper im System durchgeführt

%Planeten%.x:=%Planeten%.x+ %Planeten%.ix/genauigkeit ; Berechnung der neuen Koordinaten aus den summierten Kräften
%Planeten%.y:=%Planeten%.y+ %Planeten%.iy/genauigkeit
%Planeten%.z:=%Planeten%.z+ %Planeten%.iz/genauigkeit

Teil 2 folgt

[Kibo]

Nächste Ladung, klassisches Runge Kutta:

%Planeten%.hx:= %Planeten%.ix ; speichern des alten Impulses für trägheit
%Planeten%.hy:= %Planeten%.iy
%Planeten%.hz:= %Planeten%.iz

;Klassenteil ist der selbe daher siehe Teil 1

if (plzB=1)
%Planeten%.abstand:=abstand


%Planeten%.P1ix := %Planeten%.hx + %Planeten%.gx/(genauigkeit*2) ; %Planeten%.g/(genauigkeit*2 )=m1 prognose punkt 1 impuls (halbe schrittweite)
%Planeten%.P1iy := %Planeten%.hy + %Planeten%.gy/(genauigkeit*2)
%Planeten%.P1iz := %Planeten%.hz + %Planeten%.gz/(genauigkeit*2)

%Planeten%.P1x:=%Planeten%.x+ %Planeten%.P1ix/(genauigkeit*2) ;prognose punkt 1 koordinaten
%Planeten%.P1y:=%Planeten%.y+ %Planeten%.P1iy/(genauigkeit*2)
%Planeten%.P1z:=%Planeten%.z+ %Planeten%.P1iz/(genauigkeit*2)

%Planeten%.xd:=%PlanetenB%.x-%Planeten%.P1x ; abstand
%Planeten%.yd:=%PlanetenB%.y-%Planeten%.P1y
%Planeten%.zd:=%PlanetenB%.z-%Planeten%.P1z

abstand:=sqrt(%Planeten%.xd*%Planeten%.xd+%Planete n%.yd*%Planeten%.yd+%Planeten%.zd*%Planeten%.zd) ; immernoch abstand
%Planeten%.xn:=%Planeten%.xd/abstand
%Planeten%.yn:=%Planeten%.yd/abstand
%Planeten%.zn:=%Planeten%.zd/abstand

%Planeten%.s:=%PlanetenB%.gm/(abstand*abstand) ; gravitation

%Planeten%.g2x:=%Planeten%.s*%Planeten%.xn ;%Planeten%.g/(genauigkeit*2 )=m2 ; gravitationsvektor
%Planeten%.g2y:=%Planeten%.s*%Planeten%.yn
%Planeten%.g2z:=%Planeten%.s*%Planeten%.zn

%Planeten%.P2ix := %Planeten%.hx + %Planeten%.g2x/(genauigkeit*2) ; %Planeten%.g/(genauigkeit*2 )=m2 ;Prognosepunkt 2 impuls
%Planeten%.P2iy := %Planeten%.hy + %Planeten%.g2y/(genauigkeit*2)
%Planeten%.P2iz := %Planeten%.hz + %Planeten%.g2z/(genauigkeit*2)

%Planeten%.P2x:=%Planeten%.x+ %Planeten%.P2ix/(genauigkeit*2) ;teil 2 ; Prognosepunkt 2 koordninaten
%Planeten%.P2y:=%Planeten%.y+ %Planeten%.P2iy/(genauigkeit*2)
%Planeten%.P2z:=%Planeten%.z+ %Planeten%.P2iz/(genauigkeit*2)

%Planeten%.xd:=%PlanetenB%.x-%Planeten%.P2x ; siehe oben usw und so fort
%Planeten%.yd:=%PlanetenB%.y-%Planeten%.P2y
%Planeten%.zd:=%PlanetenB%.z-%Planeten%.P2z

abstand:=sqrt(%Planeten%.xd*%Planeten%.xd+%Planete n%.yd*%Planeten%.yd+%Planeten%.zd*%Planeten%.zd)
%Planeten%.xn:=%Planeten%.xd/abstand
%Planeten%.yn:=%Planeten%.yd/abstand
%Planeten%.zn:=%Planeten%.zd/abstand

%Planeten%.s:=%PlanetenB%.gm/(abstand*abstand)

%Planeten%.g3x:=%Planeten%.s*%Planeten%.xn ;%Planeten%.g/(genauigkeit*2 )=m3 teil 1
%Planeten%.g3y:=%Planeten%.s*%Planeten%.yn
%Planeten%.g3z:=%Planeten%.s*%Planeten%.zn

%Planeten%.P3ix := %Planeten%.hx + %Planeten%.g3x/genauigkeit ; %Planeten%.g/(genauigkeit*2 )=m3
%Planeten%.P3iy := %Planeten%.hy + %Planeten%.g3y/genauigkeit
%Planeten%.P3iz := %Planeten%.hz + %Planeten%.g3z/genauigkeit

%Planeten%.P3x:=%Planeten%.x+ %Planeten%.P3ix/genauigkeit ;teil 2
%Planeten%.P3y:=%Planeten%.y+ %Planeten%.P3iy/genauigkeit
%Planeten%.P3z:=%Planeten%.z+ %Planeten%.P3iz/genauigkeit

%Planeten%.xd:=%PlanetenB%.x-%Planeten%.P3x
%Planeten%.yd:=%PlanetenB%.y-%Planeten%.P3y
%Planeten%.zd:=%PlanetenB%.z-%Planeten%.P3z

abstand:=sqrt(%Planeten%.xd*%Planeten%.xd+%Planete n%.yd*%Planeten%.yd+%Planeten%.zd*%Planeten%.zd)
%Planeten%.xn:=%Planeten%.xd/abstand
%Planeten%.yn:=%Planeten%.yd/abstand
%Planeten%.zn:=%Planeten%.zd/abstand

%Planeten%.s:=%PlanetenB%.gm/(abstand*abstand)

%Planeten%.g4x:=%Planeten%.s*%Planeten%.xn ;%Planeten%.g/(genauigkeit*2 )=m4 teil 1
%Planeten%.g4y:=%Planeten%.s*%Planeten%.yn
%Planeten%.g4z:=%Planeten%.s*%Planeten%.zn

%Planeten%.P4ix := %Planeten%.hx + %Planeten%.g4x/genauigkeit ; %Planeten%.g/(genauigkeit*2 )=m4
%Planeten%.P4iy := %Planeten%.hy + %Planeten%.g4y/genauigkeit
%Planeten%.P4iz := %Planeten%.hz + %Planeten%.g4z/genauigkeit

%Planeten%.ix:=1/6*(%Planeten%.P1ix + 2 * %Planeten%.P2ix + 2*%Planeten%.P3ix + %Planeten%.P4ix) Verrechnung der 3 Prognosepunkte zum dann verwendeten Impuls
%Planeten%.iy:=1/6*(%Planeten%.P1iy + 2 * %Planeten%.P2iy + 2*%Planeten%.P3iy + %Planeten%.P4iy)
%Planeten%.iz:=1/6*(%Planeten%.P1iz + 2 * %Planeten%.P2iz + 2*%Planeten%.P3iz + %Planeten%.P4iz)

; ende innere schleife, ausgeführt für jeden gegnerischen Körper

%Planeten%.x:=%Planeten%.x+ %Planeten%.ix/genauigkeit ; Verrechnen des summierten neuen Impuls mit den Koordinaten
%Planeten%.y:=%Planeten%.y+ %Planeten%.iy/genauigkeit
%Planeten%.z:=%Planeten%.z+ %Planeten%.iz/genauigkeit

Ich hoffe es ist einigermaßen verständlich. Falls fragen sein sollten immer her damit.