Warum funktioniert mein Zweikörpermodell nicht? [geschlossen]

In Situationen wie dieser ist es eine gute Idee, den Zeitschritt basierend auf dem Gradienten der Kraft anzupassen - denn das gesamte Konzept der numerischen Integration besagt, dass "sich die Dinge von jetzt an bis zum nächsten Zeitschritt nicht zu sehr ändern" und Diese Annahme wird verletzt, wenn Sie sich schnell durch eine Region mit sich schnell ändernder Kraft bewegen.

Dies hat einen riskanten Nebeneffekt: Wenn Sie den Zeitschritt "klein genug" machen, um mit der schnellen Variation fertig zu werden, und die Variation "unendlich schnell" wird, kann Ihr Algorithmus "stecken bleiben" - Sie könnten auf Xenos Paradoxon von stoßen Achilles und die Schildkröte.

Um nicht stecken zu bleiben, sind einige numerische Integrationen (wie Runge-Kutta) besser darin, die Krümmung höherer Ordnung zu berücksichtigen - was einen größeren Schritt ohne Genauigkeitsverlust ermöglicht.

In Ihrem Fall wird die Beschleunigung für den gesamten nächsten Zeitschritt jedoch durch die aktuelle Position bestimmt. Wenn Ihre Partikel also zu Beginn des Schritts nahe beieinander liegen, werden sie weit weggeschleudert, und da ihre Anziehungskraft dies tut dann stark verringert werden, werden sie es schwer haben, zurückzukommen.

Ein paar Kritikpunkte zu deinem Code:

1. Sie definieren particle1[4], verwenden aber nur zwei Komponenten
2. Sie codieren den Zeitschritt hart, anstatt eine Variable wie zu verwendendt
3. Sie verwenden die grundlegendste Integrationsmethode ... bitte informieren Sie sich über andere
4. Sie definieren implizit Masse = 1 und Kraftkonstante = 1; erwägen Sie, diese Variablen zu erstellen (auch wenn Sie sie auf 1 setzen)

In Bezug auf den ersten Punkt würde ich a structfür meine Partikel verwenden:

typedef struct{
  double x;
  double v;
} PARTICLE;

und dann können Sie mit auf ihre Eigenschaften zugreifen

PARTICLE p1, p2;
p1.v = 0.;
p1.x = 0.;
p2.v = 0.;
p2.x = 1.;

Was sofort besser lesbar ist ... Wenn Sie in zwei Dimensionen arbeiten möchten, können Sie sie entweder explizit als Teil Ihres struct: ( double x; double y; double vx; double vy;) schreiben oder sie zu Arrays machen.

Aber wirklich - variable Zeitschritte und Interpolationen höherer Ordnung würden zur Stabilität Ihrer Lösung beitragen.

Sie erhalten Partikel, die 0 Abstand voneinander haben und daher eine unendliche Beschleunigung haben. Der Computer kann dies nicht auflösen und explodiert daher, je nachdem, wie nahe der Schritt an Null herankommt.

Ich bin mir nicht sicher, was "physikalisches" Verhalten wäre, da Punktmassen in Nullentfernung schwer zu untersuchen sind. Wenn Sie jedoch eine y-Koordinate hinzufügen und ihr eine Anfangsgeschwindigkeit geben, erhalten Sie möglicherweise bessere Ergebnisse (im Grunde eine längliche Ellipse anstelle einer oszillierenden Linie).

Sie könnten die anfängliche y-Geschwindigkeit sehr klein machen, sodass Sie sich der X-Bewegung annähern. Oder fügen Sie Ihrer Division einfach einen kleinen Offset hinzu.