Das Mehrkörperproblem in der Astronomie/ Allgemeine Lösungsmethoden/ Tests

Im Folgenden wird die bislang nur theoretisch skizzierte Güte der hier vorgestellten Lösungsmethoden einer praktischen Prüfung unterzogen. Das naheliegendste Testverfahren besteht darin, die simulierten Bewegungen mit den vom Newtonschen Gravitationsgesetz vorhergesagten genauen Bahnen zu vergleichen. Leider ist dies wie bereits erwähnt im Allgemeinen nur für sehr einfache, lediglich aus zwei Massenpunkten bestehende Systeme möglich. Zur Prüfung der Simulation eines realistischen Mehrkörperensembles muss auf allgemeingültige Eigenschaften wie die Erhaltung der Gesamtenergie zurückgegriffen werden. In der Praxis betrachtet man weitere Größen, z.B. die räumliche Verteilung der Sterne in einem Haufen oder einer Galaxie sowie deren Geschwindigkeiten und vergleicht diese ebenfalls mit den Vorhersagen, was im letzten Kapitel dieses Buches ausführlich behandelt wird.

Bevor die Ergebnisse einfacher Tests vorgestellt werden, seien auch noch einige Worte zur Wahl der Maßeinheiten angebracht. Alltäglichen Einheiten wie Kilogramm, Sekunde und Meter sind für die Behandlung astronomischer Systeme selbstverständlich äußerst unhandlich. Viel praktischer ist es, z.B. Massen in Sonnenmassen (1.989 10³⁰ kg), Zeiten in Jahren (31558150 s) und Längen in Astronomischen Einheiten (149597871 km, die mittlere Entfernung der Erde von der Sonne) anzugeben. Die Gravitationskonstante hat dann den Wert 39.4899. Geschwindigkeiten werden in der Astronomie üblicherweise in km/s angegeben, es gilt 1 Astronomische Einheit/Jahr = 4.740388 km/s.

Für Längen wird oft auch die Einheit Parallaxensekunde Parsec (die Entfernung, unter welcher die Strecke Erde-Sonne unter einem Winkel von 1 Bogensekunde erscheint) benutzt, es gilt 1 Parsec = 206265 Astronomische Einheiten (= 3.26 Lichtjahre). Drückt man Massen weiterhin in Sonnenmassen und Zeiten weiterhin in Jahren aus, so nimmt die Gravitationskonstante nun den Wert 4.49814 10^-15 an. Ein Parsec/Jahr entspricht 977813 km/s.

Euler-Verfahren

Die in der Theorie schon ausführlich diskutierte Erdbahn sei auch für die Praxistests herangezogen. Die Startposition der Erde sei (1 AE / 0), ihre Anfangsgeschwindigkeit (0 / 30 km/s). Die Sonne steht zu Beginn mit der Geschwindigkeit (0 / 0) im Ursprung. Man erwartet einen stabilen, fast kreisförmigen Orbit.

Simulation der Erdbahn mit dem Euler-Verfahren über einen Zeitraum von 10 Jahren mit Zeitschritten von 1 Tag (blaue Kurve) und 1/10 Tag (rote Kurve)

Die Simulation liefert jedoch ein ganz anderes Bild. Wie befürchtet bewegt sich die Erde auf einer Spiralbahn rasch von der Sonne weg. Verwendet man Zeitschritte von 1 Tag, verdoppelt sich ihre Entfernung innerhalb von nur 10 Jahren, was einem jährlichen relativen Fehler von ungefähr 10% entspricht. Die Euler-Methode schneidet in diesem Beispiel somit noch schlechter ab als vorhergesagt. Kleinere Schritte von 1/10 Tag bringen wie erwartet keine entscheidende Verbesserung. Die Erde entfernt sich nach 10 Jahren bis auf 1.2 AE, der jährliche relative Fehler ist mit etwa 2% noch immer viel zu hoch.

Verletzung der Energieerhaltung durch eine Simulation der Erdbahn mit dem Euler-Verfahren über einen Zeitraum von 10 Jahren mit Zeitschritten von 1 Tag (blaue Kurve) und 1/10 Tag (rote Kurve). Die zeitliche Entwicklung der Gesamtenergie ist relativ zu ihrem Anfangswert angegeben.

Angesichts der nach "außen" driftenden Erde ist es nicht verwunderlich, dass die Simulation auch keine konstante Gesamtenergie wiederzugeben vermag. Diese nimmt vielmehr mit der Zeit stetig zu, da die kinetische Energie zu hoch und die potentielle Energie betragsmäßig zu niedrig angesetzt wird. Mit einer Schrittweite von 1 Tag verändert sich die Energie binnen 10 Jahren etwa um die Hälfte ihres ursprünglichen Betrages, also um ungefähr 5% pro Jahr. Mit Schritten von 1/10 Tag stellt sich eine gesamte Änderung von circa 1/6 des anfänglichen Betrages bzw. von etwa 1.5% pro Jahr ein.

Zwischenschritt-Verfahren

Leapfrog-Verfahren

Die Leapfrog-Methode erfüllt die in ihr gesetzten Erwartungen. Eine Simulation der Erdbahn mit Zeitschritten von 1 Tag bleibt selbst über 100000 Jahre weitgehend stabil, der Orbit weitet sich nur um 0.03 AE. Der jährliche relative Fehler beläuft sich auf nur 3 Zehnmillionstel und ist damit noch deutlich geringer als theoretisch abgeschätzt. Im Gegensatz zur Euler-Methode neigen die individuellen Fehler der Einzelschritte nicht dazu, sich komplett aufzuaddieren, sondern heben sich gegenseitig zum Teil wieder auf. Überraschenderweise zeigt sich hinsichtlich der prognostizierten Bahn nur eine geringfügige Verschlechterung des Ergebnisses, wenn man die Schrittweite $\Delta t$ auf 10 Tage erhöht. Angesichts der starken nominellen Abhängigkeit des Fehlers eines Schrittes von $\Delta t$ - es gilt ja $\Delta R/R\sim (\Delta t)^{4}$ - sollte die Simulation sehr viel instabiler werden. Tatsächlich aber stellt sich nach 100000 Jahren mit 0.04 AE nur eine geringfügig größere Abweichung ein. Die viel zahlreicheren Rundungsfehler mit endlich genauen Gleitkommazahlen bei Verwendung einer kleineren Schrittweite kommen hier deutlich zur Geltung. Die gegenseitige Kompensation der Einzelfehler spielt für das gute Ergebnis auch mit grober Schrittweite ebenfalls eine wichtige Rolle.

Simulation der Erdbahn mit dem Leapfrog-Verfahren über einen Zeitraum von 100000 Jahren mit Zeitschritten von 1 Tag (blaue Kurve) und 10 Tagen (rote Kurve)

Erwartungsgemäß zeigt sich die Leapfrog-Methode auch hinsichtlich der Energieerhaltung dem Euler-Verfahren als klar überlegen. Zwar liefert die Simulation auch jetzt keine konstante Gesamtenergie, doch driftet diese hier nicht zu immer höheren Werten hin ab, sondern pendelt jeweils im Verlauf eines Umlaufs der Erde um die Sonne zwischen stabilen Grenzen hin und her. Vergleicht man die beiden Berechnungen erneut, zeigt diejenige mit einer Schrittweite von 10 Tagen nun doch eine klar schlechtere Qualität. Bei einer Schrittweite von 1 Tag schwankt die Gesamtenergie um etwa 1 / 2000, mit den vergröberten Schritten von jeweils 10 Tagen jedoch schon um etwa 1 / 250 ihres Absolutbetrages. Die regelmäßigen Oszillationen der Gesamtenergie bei der Simulation eines Systems aus zwei Massenpunkten sind für das Leapfrog-Verfahren typisch und von wissenschaftlichen Arbeiten wiederholt bestätigt worden.

Verletzung der Energieerhaltung durch eine Simulation der Erdbahn mit dem Leapfrog-Verfahren über einen Zeitraum von 100000 Jahren mit Zeitschritten von 1 Tag (blaue Kurve) und 10 Tagen (rote Kurve). Die zeitliche Entwicklung der Gesamtenergie ist relativ zu ihrem Anfangswert angegeben. Die Schwankungen bleiben über die gesamte Simulation hinweg stabil, so dass es genügt, nur die ersten Jahre darzustellen

Runge-Kutta-Verfahren

Bei geeigneter Wahl der Zeitschritte liefert das Runge-Kutta-Verfahren eine weitere enorme Steigerung der Vorhersagegenauigkeit. Setzt man für die Simulation der Erdbahn $\Delta t$ wiederum auf 1 Tag, so ist auch nach 100000 Jahren praktisch kein Fehler zu erkennen. Bei einer Vergröberung der Zeitschritte bricht die Runge-Kutta-Methode aber im Gegensatz zum Leapfrog-Verfahren rasch zusammen. Mit einer Schrittweite von 10 Tagen wird die Simulation der Erdbahn nach einigen Jahrtausenden völlig instabil und damit unbrauchbar. Die noch stärkere nominelle Abhängigkeit des Fehlers eines Schritts von $\Delta t$ gemäß $\Delta R/R\sim (\Delta t)^{5}$ behält hier gegenüber der geringeren Zahl von Rechenoperationen (und damit Rundungsfehlern) sowie einer gegenseitigen Kompensation von Einzelfehlern die Oberhand.

Simulation der Erdbahn mit dem Runge-Kutta-Verfahren über einen Zeitraum von 100000 Jahren mit Zeitschritten von 1 Tag.

Dass die Simulation mit einer Schrittweite von 1 Tag überhaupt einen Fehler aufweist, zeigt erst die Betrachtung der Gesamtenergie. Sie pendelt nicht wie bei dem Leapfrog-Verfahren zwischen stabilen Grenzen hin und her, sondern weist eine Verschiebung hin zu kleineren Werten auf. Diese ist allerdings außerordentlich gering - in 100000 Jahren verändert sich die Gesamtenergie nur um etwa 2 Hunderttausendstel bzw pro Jahr um 2 Zehnmilliardstel ihres Absolutwertes. Die langsame Drift der Gesamtenergie ist für die Runge-Kutta-Methode charakteristisch und wie die vom Leapfrog-Verfahren erzeugten Oszillationen von Forschungsarbeiten ebenfalls bestätigt worden.

Verletzung der Energieerhaltung durch eine Simulation der Erdbahn mit dem Runge-Kutta-Verfahren über einen Zeitraum von 100000 Jahren mit Zeitschritten von 1 Tag. Die zeitliche Entwicklung der Gesamtenergie ist relativ zu ihrem Anfangswert angegeben.

Prädiktor-Korrektor-Verfahren

Leapfrog-Verfahren

Mit dem Umbau zu einem Prädiktor-Korrektor-Verfahren reicht die Leapfrog-Methode hinsichtlich ihrer Genauigkeit fast an das Runge-Kutta-Verfahren heran, sofern die Korrektur genügend oft als Iteration durchgeführt wird. Wie bei der klassischen Formulierung oszilliert die Gesamtenergie während eines Umlaufs der Erde um die Sonne zwischen gewissen Grenzen hin und her. Diese Schwankungen fallen jetzt aber (erneut ist $\Delta t$ auf 1 Tag gesetzt) mit nur etwa 1 / 100000 (statt 1 / 2000) des Absolutwertes sehr viel schwächer aus. Andererseits sind die Extremwerte nicht mehr stabil, sondern wandern langsam hin zu höheren Werten. Der dadurch bewirkte relative Fehler beläuft sich nach 100000 Jahren jedoch nur auf 5 Hunderttausendstel der Gesamtenergie bzw. 5 Zehnmilliardstel pro Jahr.

Verletzung der Energieerhaltung durch eine Simulation der Erdbahn mit dem Leapfrog-Verfahren über einen Zeitraum von 100000 Jahren mit Zeitschritten von 1 Tag. Die zeitliche Entwicklung der Gesamtenergie ist relativ zu ihrem Anfangswert angegeben.

Mit dem Prädiktor-Korrektor-Algorithmus reagiert das Leapfrog-Verfahren wesentlich empfindlicher auf eine größere Schrittweite als gemäß der Zwischenschritt-Form. Mit Schritten von 10 Tagen bricht die Simulation zwar nicht völlig zusammen, doch variiert die Gesamtenergie um etwa 1/1000 pro Umlauf. Zudem steigt sie im Verlauf der Berechnungen um mehr als 1/3 an, d.h. um circa 3-4 Hunderttausendstel pro Jahr.

Hermite-Polynome- und Mehrschrittverfahren

Abermals eine optimale zeitliche Schrittweite vorausgesetzt, gewinnt man selbst mit einem nur einmaligen Durchlauf der Korrekturformeln gegenüber den bereits sehr guten Runge-Kutta- und Leapfrog-Verfahren noch an Stabilität. Eine Simulation der Erdbahn auf Grundlage der Hermite-Polynome mit einer Schrittweite von 1 Tag lässt erneut eine allmähliche Zunahme der Gesamtenergie erkennen, doch liegt diese nach 100000 Jahren lediglich bei etwa 4 Millionstel ihres Absolutwertes. Damit beträgt der jährliche relative Fehler gerade einmal 4 Hundertmilliardstel. Noch besser schneidet das Mehrschrittverfahren ab, das nach 100000 Jahren einen relativem Fehler der Gesamtenergie von nur 3 Zehnmillionstel liefert - entsprechend 3 Billionstel pro Jahr!

Auf eine Vergröberung der Zeitschritte reagieren beide Verfahren allerdings ebenso sensibel wie die Runge-Kutta-Methode. Mit $\Delta t$ = 10 Tage wird die Simulation nach mehreren Jahrtausenden erneut völlig instabil.

Beispiel eines Drei-Körper-Systems

Am Ende dieses Abschnitts seien die Leapfrog-Methode (als Prädiktor-Korrektor-Algorithmus) und das Hermite-Polynome-Verfahren anhand eines aus drei Körpern von je 1 Sonnenmasse bestehenden Systems vorgeführt. Runge-Kutta- und Mehrschritt-Verfahren liefern fast exakt die gleichen Resultate wie die Hermite-Polynome-Methode, so dass darauf nicht weiter eingegangen werden muss. Die Anfangspositionen (vor der Umrechnung ins Schwerpunktsystem) sind (-10 AE / 0), (10 AE / 0) und (0 / 10 AE), die Geschwindigkeiten zu Beginn (0 / -5 km/s), (0 / 5 km/s) und (-5 km/s / 0). Die zeitliche Schrittweite ist wiederum auf 1 Tag gesetzt.

Simulation eines Drei-Körper-Systems mit dem Leapfrog- (blaue Kurve) und Hermite-Polynome-Verfahren (rote Kurve) über einen Zeitraum von 50 Jahren mit festen Zeitschritten von 1 Tag.

Schon dieses scheinbar so einfache System bietet dem Betrachter hinsichtlich der Bahnen seiner Mitglieder ein geradezu sinnverwirrendes Bild. Erschwerend kommt hinzu, dass die Leapfrog-Methode schon nach kurzer Zeit völlig andere Ergebnisse liefert als das Hermite-Polynome-Verfahren. Wirft man einen Blick auf die Energieerhaltung, so mag man Letzterem mehr Vertrauen schenken. Wie für die Erdbahn liefern dieses eine sehr stabile Gesamtenergie, während das Leapfrog-Verfahren bald unter heftigen Ausschlägen leidet. Nach 7 Jahren zeigt dieses kurzzeitig einen relativen Fehler von 1 / 50 des Absolutwertes. Nach 33 Jahren steigt dieser abrupt sogar auf ein Vielfaches der Gesamtenergie.

Eine genaue Betrachtung zeigt, dass solche Fehler stets dann auftreten, wenn zwei Massenpunkte sich sehr nahe kommen. Gerade dann aber nimmt der relative Fehler eines Simulationsschritts dramatisch zu. Das Hermite-Polynome-Verfahren scheint in vorliegendem Beispiel mit solch engen Begegnungen noch fertig werden zu können, die Leapfrog-Methode trotz ihrer nur etwas schlechteren Stabilität dagegen nicht mehr. Der erste Beinahezusammenstoß sorgt bereits dafür, dass die Simulationen sich mit fortschreitender Zeit immer mehr auseinanderentwickeln. Beim zweiten derartigen Ereignis überschätzt das Leapfrog-Verfahren die kinetische Energie der Stoßpartner dramatisch, wodurch sich das Ensemble schlagartig auflöst.

Verletzung der Energieerhaltung durch eine Simulation eines Drei-Körper-Systems mit dem Leapfrog- (blaue Kurve) und Hermite-Polynome-Verfahren (rote Kurve) über einen Zeitraum von 50 Jahren mit festen Zeitschritten von 1 Tag.

Doch auch die zunächst zuverlässig erscheinenden Hermite-Polynome stoßen rasch an Grenzen. Lässt man obige Simulation weiterlaufen, tritt nach 70 Jahren eine so dichte Begegnung zweier Körper auf, dass diese mit der festen Schrittweite von 1 Tag auch mit der genaueren Methode nicht mehr behandelt werden kann. Wieder fliegen die davon betroffenen Massenpunkte unvermittelt auseinander. Da dadurch das dritte Mitglied praktisch keine Anziehungskraft mehr durch die anderen beiden Körper verspürt, bewegt sich auch dieses schließlich ungehindert fort.

Simulation eines Drei-Körper-Systems mit dem Hermite-Polynome-Verfahren über einen Zeitraum von 100 Jahren mit festen Zeitschritten von 1 Tag. Nach einem Beinahezusammenstoß zweier Massenpunkte nach 70 Jahren löst sich das System schlagartig auf, weil die sehr enge Begegnung mit konstanter Schrittweite nicht mehr korrekt wiedergegeben werden kann.

Schaut man sich die Beinahekollision genau an, so erkennt man das mit einer festen Schrittweite einhergehende Problem. Kommen zwei Körper sich sehr nahe, so üben sie starke Anziehungskräfte aufeinander aus und werden dementsprechend auf hohe Geschwindigkeiten beschleunigt. Dies aber führt dazu, dass sie nach der engsten Passage mit eben diesen Geschwindigkeiten weit auseinander gesetzt werden. Mit diesem größeren Abstand jedoch müssten sie wieder weit geringere Geschwindigkeiten aufweisen. Die kinetische Energie der Massenpunkte wird nach deren Vorübergang also maßlos überschätzt, wodurch sie sich gegenseitig nicht mehr festhalten können und dem System irrtümlich entkommen. Wie diese Instabilität verhindert werden kann, ist Gegenstand des nächsten Kapitels.

Simulation eines Drei-Körper-Systems mit dem Hermite-Polynome-Verfahren über einen Zeitraum von 100 Jahren mit festen Zeitschritten von 1 Tag. Der Zoom zeigt, wie bei einer Beinahekollision die beteiligten Massenpunkte in demselben Schritt stark beschleunigt und mit hoher Geschwindigkeit unrealistisch weit auseinandergesetzt werden.