Als Himmelsmechaniker der besonders an der Bewegung von Planeten und Asteroiden interessiert ist, bin ich auf Simulationen angewiesen. Wir können (zumindest in unserem Sonnensystem) die Himmelskörper zwar beobachten – wenn auch nicht alle – aber ihre Bewegung ist zu langsam, um ihre Bewegung direkt und rein aus Beobachtungsdaten vernünftig analysieren zu können.
Der äußerste Planet unseres Sonnensystems, Neptun, braucht für einen Umlauf beispielsweise 165 Jahre. Er wurde aber erst 1846 entdeckt. Und obwohl das schon 164 Jahre her ist, ist es doch noch nicht lang genug um einen kompletten Neptunumlauf beobachten zu können. Glücklicherweise müssen wir das aber auch nicht. Wir kennen ja die Kräfte, die für die Bewegung der Himmelskörper verantwortlich sind. Diese Kräfte können wir mathematisch beschreiben und damit auch Modelle eines Planetensystems basteln das wir dann am Computer betrachten. Hier können wir die Zeit schneller laufen lassen und zusehen, wie sich die Planeten und Asteroiden bewegen würden, wenn wir lange genug beobachten könnten.
Kleine Programme, die solche Simulationen durchführen gibt es haufenweise (hier habe ich mal eines vorgestellt). Die sind aber eher nette wissenschaftliche Spielereien und für die echte Forschung nicht wirklich brauchbar. Dabei gibt es auch professionelle Programme frei verfügbar – und die sind oft nichtmal kompliziert.
Das Programm, das ich hier vorstellen möchte heisst „Mercury„. Entwickelt wurde es von John Chambers auf dessen Homepage man das Paket auch runterladen kann.
Was macht Mercury?
Mit Mercury lässt sich die Bewegung eines beliebigen Planetensystems simulieren. Man hat einen Zentralkörper (den Stern) und dazu große und kleine Himmelskörper wie Planeten oder Asteroiden. Man kann sich mit Mercury nun ansehen, wie deren Bewegung unter ihrem gegenseitigen gravitativen Einfluss aussieht.
Ich habe den ganzen Prozess einer himmelsmechanischen Simulation früher schonmal in einer kleinen Serie (Teil 1, Teil 2, Teil 3, Teil 4) beschrieben und möchte mich hier nicht wiederholen.
Ich möchte auch nicht auf die genaue mathematische Methode eingehen, mit der Mercury die der Planetenbewegung zugrunde liegenden Differentialgleichungen löst (es handelt sich dabei nicht um die Methode der Lie-Reihen, die ich hier schonmal vorgestellt habe). Wer die Details wissen möchte, der kann das in Chambers‘ Artikel nachlesen („A hybrid symplectic integrator that permits close encounters between massive bodies“, MNRAS 304, 793 (1999))
Ich möchte hier erklären, wie jeder dieses Programm benutzen kann, um selbst die Bewegung der Himmelskörper auf professionellen Niveau zu untersuchen
Installation
Bevor man Mercury verwenden kann, muss man es erstmal installieren. Das ist nicht so einfach wie bei „normalen“ Programmen. Es gibt keine schicken grafischen Installationsassistenten oder ähnliches – hier ist Handarbeit gefragt 😉
Zuerst einmal lädt man sich das Archiv mit den Dateien runter und entpackt es. Linux-Usern muss ich wohl nicht erklären, wie man ein tar-file entpackt und die gängigen Windowsprogramme sollten das von alleine schaffen.
Ihr habt nun einige Dateien vor euch. Diejenigen, deren Endung „.for“ lautet, stellen die Komponenten des eigentlichen Programs dar. Es ist der Fortran-Quellcode von Mercury6 und wer nicht wirklich weiß was er tut, sollte diese Dateien lieber nicht bearbeiten. Was man aber tun muss, ist die Dateien zu kompilieren. Man muss aus dem Quellcode ein ausführbares Programm machen. Unter Linux sollte sich leicht ein passender Kompiler finden lassen Ich verwende meistens „g77“, der mittlerweile aber gfortran heisst soweit ich das mitbekommen habe. Mit dem Befehl
g77 -O -o mercury6 mercury6_2.for
wird aus dem Quellcode („mercury6_2.for“) eine ausführbare Datei („mercury6“) erzeugt. Auf diese Weise kompiliert man am besten gleich auch noch „element6.for“ und „close6.for“. Natürlich kann man die Fortran-Dateien auch unter Windows kompilieren. Ich weiß gerade nicht, welche Programme hier momentan empfehlenswert sind – ich glaube, ich habe Mercury mal mit dem Watcom Fortran Compiler kompiliert.
Ok – gehen wir davon aus, dass wir nun ein lauffähige Versionen von Mercury haben. Nun müssen wir dem Program natürlich noch sagen, was wir genau simulieren wollen. Dafür gibt es eine Reihe von Inputdateien, die wir editieren müssen.
Was soll Mercury machen
Wichtig sind hier vor allem die Dateien mit den Namen „big.in“ und „small.in“. So wie die meisten numerischen Simulationsprogramme unterscheidet auch Mercury zwischen „großen“ und „kleinen“ Himmelskörper. Große Objekte sind massiv genug, um alle anderen Himmelskörper durch ihre Gravitationskraft zu beeinflussen. Kleine Objekte spüren nur die Gravitationskraft aller anderen Körper; sind selbst aber zu winzig um einen merkbaren Gravitationseinfluss auszuüben; dieser wird vernachlässigt.
Typischerweise trägt man bei „big.in“ die Planeten ein, die man simulieren will; bei „small.in“ finden sich dann Asteroiden, Kometen oder anderer „Kleinkram“. Man sollte sich hier gut überlegen, was man alles inkludieren will. Während man von den „kleinen“ Objekten durchaus sehr viele gleichzeitig simulieren kann (sie beeinflussen einander ja nicht) steigt der Rechenaufwand enorm mit der Zahl der „großen“ Körper – denn hier muss ja immer ausgerechnet werden wie jedes einzelne Objekt jedes andere beeinflusst.
Mercury kommt gleich mit voll funktionsfähigen Inputdateien. In der mitgelieferten „big.in“-Datei sind die passenden Werte für alle Planeten (plus Pluto) angegeben. Das sieht so aus:
Die ersten beiden wichtigen Angaben, die man hier machen muss, sind die für „style“ und „epoch“. Unter „style“ spezifiziert man, wie genau man die Parameter der Himmelskörper angegeben hat. Im Beispiel ist „Cartesian“ eingestellt – das bedeutet die Planeten werden durch die Angabe ihrer heliozentrischen x,y und z-Koordinaten und durch die Angabe ihre x, y und z-Geschwindigkeitskomponenten beschrieben (zu den anderen Optionen sage ich später etwas). Der Punkt „epoch“ ist der Zeitpunkt – aufgeschrieben im julianischen Datum, für den die Angaben gültig sind. Diese Information wird dann später noch nützlich sein…
Dann kommen die Daten der Planeten selbst. In der ersten Zeile steht jeweils der Name des Objekts, seine Masse (in Einheiten der Sonnenmasse), ein Wert für „r“ der uns vorerst nicht weiter zu interessieren braucht und die Dichte des Planeten (die eigentlich auch nicht relevant für die allermeisten Berechnungen ist). Danach folgen die schon oben genannten Werte für Position und Geschwindigkeit und schließlich drei weitere Zahlen die im Beispiel alle den Wert „0“ haben (die haben mit dem Drehmoment des Körpers zu tun und sind fürs erste auch uninteressant).
Die Datei „small.in“ sieht ähnlich aus:
Auch hier wird am Anfang spezifiziert, in welchen Format die Daten angegeben sind. In diesem Fall ist es „asteroidal“ – das bedeutet, dass man nicht x,y und z-Koordinaten angibt sondern dessen Bahnelemente.
Der Zeitpunkt für den diese Elemente gültig sind („epoch“) wird diesmal für jeden kleinen Körper einzeln angegeben. Das ist der Wert der bei „ep“ gleich neben dem Namen des Objekts steht. Diese Zeitpunkte können identisch sein – müssen aber nicht. Sind die Zeiten verschiedenen, dann rechnet der Integrator erstmal alle Parameter auf den gleichen Zeitpunkt um bevor die eigentliche Simulation startet.
Jetzt hat man eigentlich schon das wichtigste geschafft – es fehlen allerdings noch einige allgemeine Angaben zur Simulation. Die gibt man in der Datei „param.in“ an:
Wichtig sind hier vor allem die ersten 6 Zeilen. Mit „algorithm“ wählt man die numerische Methode aus, die für die Simulation verwendet werden soll. Hier ist für die meisten Anwendungen „HYB“ die beste Wahl. Dieser „Hybrid-Integrator“ wählt je nach den aktuellen Bedingungen automatisch die beste Methode aus. Natürlich muss das Programm auch wissen, wie lange es simulieren soll. Dafür gibt man Start- und Stopzeitpunkt (wieder im Julianischen Datum) an. Das „Output-Interval“ ist der Zeitraum, für den man seine Datenpunkte bekommen möchte. Wenn ich beispielsweise die Bewegung der Himmelskörper eine Million Jahre lang berechne, dann möchte ich vielleicht nicht für jeden einzelnen Tag alle Positionen aller Objekte erfahren (die Datenmenge wäre enorm groß) sondern z.B. nur alle 10000 Jahre. Dieser Zeitraum kann hier spezifiziert werden. Sehr wichtig ist der Wert bei „timestep“. Das ist der Zeitschritt der Integration und diesen Wert auszuwählen ist knifflig. Ich habe das hier schonmal genauer beschrieben. Als Richtwert (der aber wirklich nur als solcher zu betrachten ist) kann man hier etwa 10% der kürzesten Umlaufzeit eines in der Simulation auftretenden Himmelskörpers nehmen.
Den „accuracy parameter“ kann man fürs erste ignorieren – genauso wie die Einstellungen, die man in den nächsten Zeilen machen kann. Wichtig sind nur zwei Parameter aus dem letzten Absatz. Bei „radius of central body“ und „central mass“ muss man die entsprechenden Werte für den Zentralstern des Planetensystems eingeben.
Rechnen und Auswerten
Hat man das alles eingestellt (bzw. nur angesehen – denn die Beispieldateien sind ja schon korrekt ausgefüllt), dann kann es eigentlich losgehen. Will man die Namen der ganzen Inputdateien ändern, dann muss man das noch in „files.in“ angeben und danach kann man „mercury6“ laufen lassen. Macht man das mit den Beispieldateien, sollte das Ergebnis so aussehen:
Die ganzen Zahlen die hier stehen, braucht man auch erstmal nicht genau zu verstehen – interessant ist vielleicht das, was bei „dE/E“ und „dL/L“ steht. Hier wird angegeben wie sich Gesamtenergie und -drehimpuls des Systems im Laufe der Zeit geändert haben. Da die Energie- bzw. Drehimpulserhaltung gilt, sollten diese Werte natürlich sehr klein sein!
Aber was hat das Programm denn nun simuliert? Die Ergebnisse finden sich in einer neuen Datei namens „xv.out“ die man erst bearbeiten muss, um sie betrachten zu können. Dazu braucht man das Programm „element6“, das aus dieser Datei die richtigen Ergebnisse extrahiert. Man editiert dazu die Datei „element.in“ – hier muss man angeben, welche Datei bearbeite werden soll („xv.out“), in welchen Format die Daten ausgegeben werden sollen, welcher Zeitraum zwischen den Ausgabeschritten liegen soll und natürlich welche Daten ausgegeben werden sollen! Hier kann man aus einer langen Liste von Möglichkeiten wählen (die einzelnen Parameter sind alle in der Bedinungsanleitung „mercury6.man“ aufgeführt). Im Beispiel findet man hier folgende Zeile:
a8.5 e8.6 i8.4 g8.4 n8.4 l8.4 m13e
Sieht kryptisch aus, ist aber simpel. „a“ steht für die große Halbachse des Objekts – „a8.5“ heisst, dass der Wert auf 8 Stellen genau, 5 davon hinter dem Komma ausgegeben werden soll. „e“ steht für die Exzentrizität der Bahn; „i“ für die Inklination, etc.
Lässt man nun das Programm „element6“ laufen erhält man einen Haufen neue Dateien mit der Endung „.aei“. Das sind nun die eigentlichen Resultate der Simulation. So sehen sie beispielsweise für den Asteroiden APOLLO (einer aus der Datei „small.in“) aus:
Man sieht hier genau, was sich im Simulationszeitraum von 30 Jahren getan hat. Die große Halbachse der Bahn des Asteroiden ist beispielsweise von 1.47105 astronomischen Einheiten auf 1.47075 gewachsen während die Bahnexzentrizität von 0.560011 auf 0.559899 gesunken ist.
Gut – das alles ist noch nicht wirklich schön bunt und spektakulär. Man kann diesen Output mit den Rohdaten vergleichen, die ein beobachtender Astronom von seinem Teleskop bekommt. Da muss man nun noch analysieren, reduzieren und ein paar nette grafische Programme benutzen um schließlich schöne Plots zu bekommen. Aber dann sieht man zum Beispiel so wie im nächsten Bild sehr schön wie sich die Exzentrizitäten der Erd- und Venusbahn gekoppelt und gegengleich im Laufe der Jahrtausende ändern:
Soweit zum Einsteigercrashkurs „Sonnensystemsimulation“ 😉 Mit mercury6 kann man noch viel mehr machen – ich habe z.B. noch nicht erklärt, wozu das Program „close6“ gut ist und wie das mit den nahen Begegnungen und den Kollisionen ist. Aber das muss man auch nicht wissen um das Program erstmal kennenzulernen – und bei Interesse erkläre ich das gern in einem späteren Beitrag. In der Bedienungsanleitung („mercury6.man“) werden alle Aspekte und Details nochmal genau erklärt – und wenn noch Fragen offen sind, dann antworte ich natürlich auch gerne! Wer also bisher durchgehalten und immer noch Lust hat, solche Simulationen zu rechnen, dem wünsche ich viel Erfolg!
Und wer dabei auf den Geschmack gekommen ist, der soll mir Bescheid sagen! Ich kann immer Leute (und vor allem Computer 😉 ) brauchen die mir bei meinen Simulationen helfen (die dann doch etwas umfangreicher sind als die paar Asteroiden und 30 Jahre im Beispiel 😉 ).
Hui, das erinnert mich daran, wie ich vor etlichen Jahren einmal ein C-Programm schnell zusammengehackt habe, um einen Eindruck von der Newtonschen Mechanik zu bekommen — Simulation eines n-Körper-Problemes. (Differentialgleichungen mit dem klassischen Runge-Kutta-Verfahren „gelöst“, weil ich für einen schnellen Hack keine Lust hatte, ein brauchbareres Verfahren zu coden. Und „natürlich“ war das Modell strikt klassisch und hat sich nicht um relativistische Effekte gekümmert, wie sie sich etwa schon beim Merkur niederschlagen.)
Und tatsächlich, ich habe etwas dabei gelernt: Nämlich, was sich aufschaukelnde Rundungsfehler sind und wie sie das Ergebnis eines fröhlichen Tages mit emacs und gcc völlig wertlos machen können. 😉 Und dabei kam ich mir schon gewitzt vor, weil ich gleich mit Runge-Kutta anfing und nicht die banalste Form der Näherung nahm…
Danke für den Tipp mit Mercury, ich glaube, ich habe jetzt einen Grund, mir mal einen Fortran-Compiler zu installieren.
„hybrid symplectic integrator“ ist aber auch ein Ungetüm von einem Namen ;-). Die Methode ist mal eine hübsche Anwendung der Poisson-Klammer, schick!
Hm… ich überlege gerade, ob man mit dem Programm wohl auch das fiktive „Cyrannus Sternensystem“ simulieren könnte?
https://de.battlestarwiki.org/wiki/Cyrannus_Sternensystem
Aber leider fehlen auf der oben verlinketen Seite eine Menge Daten dazu. Und ob Kevin Grazier, „science advisor“ für die Serie Battlestar Galactica, „who works at NASA’s Jet Propulsion Laboratory“ sich in seinem Buch „The Science of Battlestar Galactica“, näher dazu ausgelassen hat, weis ich leider nicht, weil ich das Buch (noch) nicht habe. 🙁
Anyway: Cool, dass es sowas inzwischen gibt.
@Hans: „Anyway: Cool, dass es sowas inzwischen gibt.“
Naja, Programme wie das hier gibt es seit den 70er Jahren. Mercury hat auch schon bald 15 Jahre auf dem Buckel. Sooo extrem neu ist das jetzt nicht.
Okay, auch wenn es jetzt nicht mehr so ganz neu ist, so wusste ich bisher trotzdem nichts von deren Existens. Zumindest nicht für den Gebrauch auf PCs und das meinte ich. Das es die für Main Frames schon länger gibt, ist mir klar. Aber damit wären wir dann auch bei der Frage nach der benötigten Rechenleistung. Und die ist ja erst in den letzten 15 Jahren in Bereiche vorgerückt, wo es auch für „ambitionierte“ bzw. „semiprofessionelle Amateure“ interessant wird, sich an Numbercrunching Aufgaben im grösseren Stil zu wagen. Und ich denke, ein N-Körper Problem wie es sich aus dem oben genannten Battlestar-System ergibt, ist ’ne Numbercrunchingaufgabe im grossen Stil. Vor allem, wenn man wissen will, ob das System auch über eine Millionen Jahre hinweg stabil ist oder warum es sehr wahrscheinlich nicht funtionieren würde.
Dahinter steckt die Frage, wieviel reale Science in der Fiktion drin ist? Freilich wird dabei auch getrickst, weil es die meissten Zuschauer eher weniger interessiert solange es vernünftig klingt. Auch wenn es aus wissenschaftlicher Sicht eigentlich grosser Unfug ist. Aber das ist ’ne andere Baustelle.
Hallo,
Ich habe versucht obliquity auszugeben und erhalten leider nur die gleihen Daten wie bei inclination.
Eventuell ist die Kompilierung fehlgeschlagen.
Ich erhalte dort Warnings.
Das sind meine outputparameter. Der letzte ist die obliquity
a8.5 e8.6 i8.4 g8.4 n8.4 l8.4 m13e o8.4
Im Code steht – vereinfacht:
obliquity = acos (cos(inclination)*cos(is(k))
+ sin(inclination)*sin(is(k))*cos(ns(k) – longitude of ascending node))
is und ns sind 0, da sqrt(s(1)*s(1) + s(2)*s(2) + s(3)*s(3)) kleiner gleich 0 ist. Die Werte sind die 3 Drehmoments-werte aus big.in und die habe ich auf 0 belassen.
Das sind meine param.in werte:
algorithm (MVS, BS, BS2, RADAU, HYBRID etc) = hyb
start time (days)= 2451179.5
stop time (days)= 2564409.5
output interval (days) = 3652.50d0
timestep (days) = 80
accuracy parameter=1.d-10
Vll sieht jemand den Fehler
@Andreas: Hast du denn im Manual nachgesehen, ob die Obliquity überhaupt schon komplett im Code implementiert ist? Ich bin mir nicht sicher, ob das mit mercury6 geht.
Ich habe den Wert aus dem Manual.
Von nicht komplett stand dort jetzt nichts.
Nach etwas debuggen und probieren, sind die letzten 3 Nuller-werte in – für Erde – Big.in verantwortlich…. aus dem Manual:
„After these numbers you should give the 3 components
of spin angular momentum for the body, in units of solar masses AU^2 per
day (if in doubt enter these as 0)“
Kann man die für die Erde irgendwo nachschlagen?
tschau
@Andreas: Diese Funktion hab ich nie benutzt. Vor allem weil die Werte für die meisten Planeten eh ncht bzw. nur unvollständig bekannt sind. MMn muss man sowieso andere Methoden verwenden, wenn man die Planeten nicht als Massepunkte rechnet sondern als ausgedehnte Körper. Da brauchts dann andere Diffgleichungen…
Hallo,
Ich habe Daten vom JPL Horizons und die Aussage das der „spin angular momentum“ nicht zur Verfügung steht. Der Rest ist schön zu übertragen.
Kennen Sie eine Formel um den Wert zu errechnen?
https://ssd.jpl.nasa.gov/horizons.cgi
@Andreas
Klingt wie der Rotations-Drehimpuls des Planeten (es gibt auch einen Bahn-Drehimpuls, aber „Spin“ bedeutet ja Eigendrehung). Der ist nicht schwer zu berechnen:
L = ϴ*ω
mit dem Trägheitsmoment ϴ und der Winkelgeschwindigkeit ω. Da Planeten mit guter Näherung massive Kugeln sind, ist ϴ ≈ 2/5 m/r² (wenn es genauer sein muss, dann musst Du allerdings die Formel für ein Rotationsellipsoid mit entsprechender Abplattung suchen, das um seine kurze Achse rotiert; wenn Du es noch genauer haben willst, musst Du auch noch berücksichtigen, dass die Planeten innen dichter als außen sind, da braucht es dann allerdings wohl an den Planeten gemessene Werte, wenn es solche gibt). Die Massen m und Radien r für alle Planeten des Sonnensystems finden sich in den zugehörigen Wikpedia-Einträgen.
Wenn Du die wie in #8 zitiert in merkwürdigen Einheiten brauchst, musst Du sie vorher entsprechend durch diese dividieren. Beispiel Erde:
Erdmasse = 5,974E+24 kg, Sonnenmasse = 1,989E+30 kg; Erdmasse in Sonnenmassen = 5,974E+24 kg / 1,989E+30 kg =3,0035E-6
Erdradius = 6378 km, AU = 149 597 871 km; Erdradius in AU = 4,26343E-5.
Für ω ist die siderische Rotationsperiode anzusetzen, also eine Drehung bzgl. des Sternenhimmel (die auf die Sonne bezogene bürgerliche Tageslänge von 24 h ist 3 min 56 s länger), das sind 86164 s oder in Tagen 86164/86400 = 0,997269 Tage.
So, jetzt brauchst Du die Werte nur oben einzusetzen (bei r das Quadrieren nicht vergessen).
@Andreas
Ich hab‘ noch vergessen zu erwähnen, dass ω = 2π/T ist, und T die oben genannte Rotationsdauer. Die endgültige Formel lautet also
L = 2/5 m/r² * 2π/T
mit m, r und T wie oben angegeben.
Das Problem ist leider, dass man 3 Werte für einen Vektor in 3D benötigt.
@Andreas
Du meinst für $latex \vec{\omega}$?
Müsste folgendermaßen gehen für Winkelgeschwindigkeit und Achsenrichtung:
Der Betrag ist die Winkelgeschwindigkeit, wie Alderamin sie angegeben hat, die Richtung parrallel zur Drehachse und so, das es eine Rechtsdrehung ist. Außerdem gilt $latex \vec{v}=\vec{\omega}x\vec{r}$, dabei ist x das Vektorprodukt.
Also wenn das Ding parrallel zum Daumen an der rechten Hand ist, dann dreht sich der Planet in die Richtung, in die deine Finger zeigen.
@Andreas
Die Richtung des Drehmoments ist diejenige der Drehachse der Erde, der Betrag entspricht dem der obengenannten Formel. Du kannst also einfach einen Einheitsvektor in der Richtung der Erdachse mit dem Betrag des Drehimpulses multiplizieren (also komponentenweise L = (x*|L|, y*|L|, z*|L|), wenn (x,y,z) ein Einheitsvektor ist).
Fragt sich noch, wie der Einheitsvektor der Erdachse ausschaut. Das hängt davon ab, welches Koordinatensystem verwendet wird, was die x-, y- und z-Achsen sind. Sind das ekliptikale Koordinaten, wo die x-y-Ebene die Ebene der Erdbahn ist und z die Achse senkrecht dazu? Zeigt eine die x-Achsen zum Frühlingspunkt (engl. vernal equinoix point)? Dann kann wäre x=0, y=sin i und z=cos i mit der Schiefe der Erdachse i = 23,44°.
Ggf. könnte auch eine andere Anordnung der Koordinaten gelten, y nach oben und eine x-z-Ebene, oder y zeigt zum Frühlingspunkt, dann wären die Werte entsprechend zu vertauschen. Oder es werden äquatoriale Koordinaten verwendet, dann wäre einfach x=0, y=0 und z=1 (wenn z die Nordrichtung wäre).
@Andreas
Da sind genug Fehler drin, das nochmals zu schreiben:
Die Richtung des Drehimpulses ist diejenige der Drehachse der Erde, der Betrag entspricht dem der obengenannten Formel. Du kannst also einfach einen Einheitsvektor in der Richtung der Erdachse mit dem Betrag des Drehimpulses multiplizieren (also komponentenweise L = (x*|L|, y*|L|, z*|L|), wenn (x,y,z) ein Einheitsvektor ist).
Fragt sich noch, wie der Einheitsvektor der Erdachse ausschaut. Das hängt davon ab, welches Koordinatensystem verwendet wird, was die x-, y- und z-Achsen sind. Sind das ekliptikale Koordinaten, wo die x-y-Ebene die Ebene der Erdbahn ist und z die Achse senkrecht dazu? Zeigt die x-Achse zum Frühlingspunkt (engl. vernal equinox point)? Dann kann wäre x=0, y=sin i und z=cos i mit der Schiefe der Erdachse i = 23,44°.
Ggf. könnte auch eine andere Anordnung der Koordinaten gelten, y nach oben und eine x-z-Ebene, oder y zeigt zum Frühlingspunkt, dann wären die Werte entsprechend zu vertauschen. Oder es werden äquatoriale Koordinaten verwendet, dann wäre einfach x=0, y=0 und z=1 (wenn z die Nordrichtung wäre).
Vielen Dank euch.
Jetzt macht’s mehr und mehr Sinn. Ich denke so muss es gehen.
Ephemeris said:
Coordinate systm: „Ecliptic and Mean Equinox of Reference Epoch“
Klingt für mich nach ekliptikalen Koordniatensystem.
/off
Nach einer befreienden Diagnose, werde ich mich mal hoffmännisch auf’s Rad schwingen und schauen wohin michs so treibt.
//
Später werde ich versuchen den gleichen obliquity-wert zu berechnen, den ich hier als i mit in die Berechnung einfliessen lasse. Aber dann kann man wohl hoffen,
dass das Programm funktioniert und ich mich nicht verechnet habe :).
Adios
Hallo,
Sie meinen 2/5 m*r² * 2π/T?!
Dann ist 2/5 * 3,0035193565E-6 * 4.26342965803E-5² * 2π/0,997269
=
0.0000000000000137586612968548701011352894601879605410688916810119
=
1.375866129685487E-14
Dann erhalte ich für einen Vektor von
(0, 5.473037341070530E-15, 1.262325722406288E-14)
am 2012-Dec-12
eine obliquity von 18.4367
Ich habe auf 23.44 gehöfft.
Das die Berechnung funktioniert und ich den Code richtig gelesen habe, hat mit John Chambers per Mail bestätigt.
1.) Rechenfehler
2.) Fehldeutung vom Koordniatensystem (Da werde ich jetzt nochmal nachschauen)
3.) Programm, bzw. Kompilierung fehlerhaft.
zu 2) Ich vermute ein Koordinatensystem, was eher „heliozentrisch“ und nicht „geozentrisches“ ist. Grob gedacht müsste die inclination von 7.155° noch reinkommen.
Alles Gute
@Andreas
Heliozentrische Koordinaten sind auch nur ekliptikale, nur mit anderem Ursprung. Für den Winkel der Erdachse sollte das egal sein. Nennt man die Achsneigung nicht „inclination“? Ist mit „obliquity“ vielleicht was anderes gemeint? Die Achsneigung darf sich über ein paar Jahre nicht ändern.
@Andreas
Ha, die Formel ist L = 2/5 m / r² * 2π/T
m / r², nicht m * r². Das macht dann am Ende 4,164278E-7.
Schon mal ein Rechenfehler.
echt /r²
okay probiers nochmal
mal schaun was ich da raus krieg
komisch die Einheit ist doch auch solar mass*AU² per day
Nein ich denke die „inclination“ ist die Neigung des Orbits um die Sonne und „obliquity“ ist die Neigung der Rotation um sich selbst.
https://en.wikipedia.org/wiki/Orbital_inclination.
https://en.wikipedia.org/wiki/Obliquity
Laut wiki wird der Begriff inclination wohl für beides benutzt.
Hier wird das Koordinatensystem von Horizons erklärt.
https://naif.jpl.nasa.gov/pub/naif/toolkit_docs/FORTRAN/req/frames.html#Mean%20Ecliptic%20and%20Equinox%20of%20Date%20Frames
Das mit dem heliozentrisch… war bla bla .. sry… war vorschnell. Die Positionsdaten würden dann wenig Sinn machen.
Zum Koordinatensystem:
Ecliptic and Mean Equinox of Reference Epoch
Reference epoch: J2000.0
xy-plane: plane of the Earth’s orbit at the reference epoch
x-axis : out along ascending node of instantaneous plane of the Earth’s
orbit and the Earth’s mean equator at the reference epoch
z-axis : perpendicular to the xy-plane in the directional (+ or -) sense
of Earth’s north pole at the reference epoch.
Der Begriff „J2000“ wird dann so erklärt:
The x-axis is aligned with the mean equinox. The z-axis is aligned with the Earth’s spin axis or celestial North Pole. The y-axis is rotated by 90° East about the celestial equator.
Die Abbildung könnte sein:
1 0 0
R= 0 cos -sin
0 sin cos
(Quelle:https://gmat.gsfc.nasa.gov/docs/GMATMathSpec.pdf)
Ich dachte an einen Schreibfehler bei der Formel,
weil
https://de.wikipedia.org/wiki/Tr%C3%A4gheitsmoment#cite_note-dem149-6
und wg der Einheiten.
tschau
@Andreas
Ja Du hast Recht, ich hab‘ nochmal bei Wikipedia nachgeschaut, es ist tatsächlich L = 2/5 m * r² * 2π/T. Mein Fehler.
Ja, das ist umständlich formuliert aber beschreibt den Frühlingspunkt (da steigt der Äquator über die Ebene der Erdbahn, nennt sich deswegen „aufsteigender Knoten“, ascending node). Also die x-Achse geht erwartungsgemäß zum Frühlingspunkt, die y-Achse mathematisch positiv im 90°-Winkel dazu in der Ebene der Erdbahn, die z-Achse senkrecht dazu in Nordrichtung, gerade so, wie ich’s ursprünglich vermutete. Dann müsste x=0, y=sin i und z=cos i richtig sein.
Der Winkel i wird gemessen zwischen z-Achse und Erdachse, es ist der selbe, der auch zwischen Äquator und Bahnebene (x-y) besteht. Wäre er 0 (Erdachse parallel zur z-Achse, senkrecht auf der Bahnebene), dann wären x = y = 0 und z = 1. x ist deswegen immer 0, weil der Frühlingspunkt gerade ein Ende der Achse ist, um die der Äquator gegen die Erdbahn verkippt ist – so ist er definiert (der Äquator schneidet dort die Bahnebene). Also neigt sich die Erdachse senkrecht zur x-Achse nur in der y-z-Ebene.
Diese Transformationsmatrix angewendet auf die äquatorialen Koordinaten des Himmelspols (0,0,1) ergibt (0, -sin i, cos i). Vielleicht hilft das Minuszeichen, aber ich bezweifle es. Dann weiß ich auch nicht weiter. Oder spielt doch die Abplattung und Dichtevariation der Erde eine Rolle, die das Trägheitsmoment beeinflusst? Die Formel gilt ja nur für eine homogene Kugel. Allerdings schreibt Wikipedia auch, dass die Abweichung nur 0,3 % sei.
???
okay – es ist also die transformation vom Koordinatensystem sozusagen. Dann verstehe ich es besser. Das Minus macht einen gr. Unterschied.
2 Sachen.
1.) So schlecht schätze ich die Werte nicht ein, weil ansonsten gar nicht mehr neu berechnet wird. Da ist das Programm sehr benutzerunfreundlich.
2.)….. ääää(*) … okay ich sprech mal weiter….
Sind die Daten vll falsch, da ich einen schlechten inclinations-wert habe.
Ich hatte John Chambers nochmal gebeten die Daten zu prüfen oder mir geprüfte Testdaten zu schicken. Ich denke aber er hat besseres zu tun.
(*) Deshalb habe ich die Testdaten vom Programm wieder eingetragen und nicht die Daten vom Ephemeris-programm und ich erhalte … tata … 23.4417 und in 580 Jahren 23.5167.
Heureka 🙂
Allerdings habe ich die Zeit noh nicht angepasst. Vorerst bin ich zufrieden.
Bis bald
Hallo,
Das davor keine ordentliche Werte rausgekommen sind, lag daran dass ich nur die Erddaten mit denen von Ephemeris ersetzt habe und die anderen auf den Beispieldaten belassen habe. Na ja.
Ich denke das hat mit dem Uranus auch ganu gut geklappt.
Ich habe beide „obliquities“ gedruckt und abgelegt.
Achsenneigung Erde:
https://www.klingbim.de/images/earth_obliquity.jpg
Achsenneigung Uranus:
https://www.klingbim.de/images/uranus_obliquity.jpg
tschau
p.s.: es ist plus sinus.
Zumindest denke ich das nach einem Kurvenvergleich.
Kommen die Daten aus Ephemeris müsste es -sin sein. Die Transformation scheint mir hier richtig erklärt.
Das PDF stützt sich auf das Programm.
Ich vermute das die Beispieldaten aus einem anderen Koordinatensysten kommen.
Das reicht jetzt aber.
Tschau und danke nochmal.
@Andreas
Na, herzlichen Glückwunsch. Scheint ja ein dolles Programm zu sein.
Und Respekt dafür, dass Du Dich so in die Materie einarbeitest, das hat man sonst eigentlich nie. Genau so muss das sein. Alles selbst nachprüfen. Super.
https://universesandbox.com/
Das macht auch ziemlich „Spass“. Weiß leider nur nicht, wie akkurat das ganze ist. Aber es gibt fast alle Planeten/Monde/Kleinkörper/etc.