Welche Variablen werden benötigt, um einfache Hufeisenumlaufzeiten zu berechnen?

Das hängt von der Genauigkeit ab, mit der Sie arbeiten möchten. In nullter Ordnung, wie in Murray & Dermott, "Solar System dynamics", Kap. 3. beschrieben, können Sie Folgendes tun:

• Die Nullgeschwindigkeitskonturen, die in Ihrem Bild gezeichnet sind, werden nicht mit Teilchenumlaufbahnen mit unendlicher Genauigkeit übereinstimmen, aber sie sind eine gute Annäherung nullter Ordnung für Objekte mit geringer Exzentrizität in Bezug auf den Stern ($v_{r}/v_{\theta} \ll 1$)
• Ein Teilchen auf diesen Bahnen befindet sich auf einer regulären Keplerbahn mit Radius$r_{\rm H}$, außerhalb der Reichweite des Gravitationseinflusses des Planeten bei$r_{\rm P}$. Um also den größten "Bruch" der Reisezeit eines Teils der Hufeisenbahn zu erhalten, entweder des inneren oder des äußeren Teils, können Sie mit relativen Geschwindigkeiten und der Annahme von keplerschen Geschwindigkeiten arbeiten.
• Vorsicht ist geboten, an welcher Umlaufzeit Sie interessiert sind: Wenn$a$ist die Entfernung der Erde von der Sonne, und das Hufeisen ist$d$weg von einer perfekt kreisförmigen Umlaufbahn, wodurch zwischen den Entfernungen gewechselt wird$a\pm d$, dann ist die Umlaufzeit im Ruhesystem der Sonne$v_{-}=\sqrt{\frac{GM}{a-d}}$, für eine große Anzahl von Umlaufbahnen bis zur engen Begegnung und danach$v_{+}=\sqrt{\frac{GM}{a+d}}$.
• Daher gilt für die meiste Zeit der vereinfachte Fall, dass sich der Satellit auf einer eigenen Keplerbahn befindet. Mit diesem Wissen bewaffnet, können wir die planetozentrische Wiederkehrzeit annähern als$t_{\rm rec}=\frac{2\pi a}{v_{K} - v_{-}}$als einfache Aufholzeit zwischen Objekten auf unterschiedlichen Umlaufbahnen.
• Die Relativgeschwindigkeit$v_{\rm rel} = v_{K} - v_{-} = \sqrt{\frac{GM}{a}}-\sqrt{\frac{GM}{a-d}}$kann in der Grenze von erweitert werden$d/a\ll1$hinein$v_{\rm rel} \approx \frac{1}{2} \sqrt\frac{GM}{a} \frac{d}{a}$und so leite ich daraus ab$t_{rec}\approx \frac{4\pi}{\sqrt{GM}} \frac{a^{5/2}}{d}$. Wie erwartet divergiert die Zykluszeit für$d=0$, da im Fall der gemeinsamen Umlaufbahn mit der Erde diese Zeit unendlich sein muss. Eine Obergrenze an$d$daraus nicht abgeleitet werden kann, müsste man sich dafür der vollständigen Lösung zuwenden.

Nur aus Neugier habe ich ein paar Werte in diese Formel eingefügt und schnell etwas in Python geschrieben:

import numpy as np

#Basic physics quantities
G      = 6.678e-8 #cgs units
pi     = 3.141592
navo   = 6e23 # particles per mole
sigma  = 5.67e-5   #erg cm-2 s-1 K-4
kb     = 1.38e-16  #erg/K
km     = 1e5 #kilometers in cm
mearth = 5.98e27  #g
msun   = 2.0e33   #g
au     = 1.49e13  #cm
yr   = 365*24*3600
rearth   = 6370e5
rjupiter = 74000*km

#
# Returns the approximate horseshoe-cycle time in the reduced 3body problem
# Masses of bodies: m0>>m1>>m2
# Semimajor axis distance is from m0 to m1, radial distance is a(m0->m1)-a(m0->m2)
#
def hs_cycle(mcentral, semimajor_axis, radial_distance):
    return 4*pi/np.sqrt(G*mcentral)*semimajor_axis**(5./2.)/radial_distance/yr

#
# https://en.wikipedia.org/wiki/(419624)_2010_SO16 around the Sun
#
# Quoted cycle time ~350 years, with d=0.004 AU
#
print("Predicted 2010_SO16 cycletime [years] = ", hs_cycle(msun, 1*au,0.004*au), " predicted = 350 yrs")

#
# Janus/Epimetheus around Saturn
#
# a = 151410 km, d = 25 km, as stated in https://en.wikipedia.org/wiki/Epimetheus_(moon)
# Quoted cycle time = 8 years (from comments)
#
print("Predicted Janus/Epimetheus cycletime [years] = ", hs_cycle(95*mearth, 151410*km,50*km), " predicted = 4 yrs")


#
# 3753 Cruithne
#
# a = 1 AU and semimajor axis difference from https://en.wikipedia.org/wiki/3753_Cruithne
# Quoted cycle time = 770 years
#
print("Predicted 3753 Cruithne cycletime [years] = ",hs_cycle(msun, 1*au, (1.0-0.99774)*au), " predicted = 770 yrs")

und die Ergebnisse, die ich bekomme, sind

Predicted 2010_SO16 cycletime [years] =  495.7747141830971  predicted = 350 yrs
Predicted Janus/Epimetheus cycletime [years] =  11.542076781209305  predicted = 8 yrs
Predicted 3753 Cruithne cycletime [years] =  877.4773702355546  predicted = 770 yrs

Die Formel kann also bis zu einem Faktor von ~2 abweichen. Dies liegt sicherlich einfach daran, dass die Realität komplexer ist als eine einfache Kreisbahnnäherung, aber auch an der Qualität der verwendeten Werte. Wikipedia ist nicht dafür bekannt, bestimmte Werte gut zu recherchieren. Ich nahm die, die ich dort fand. Für SO16 war die Auswahl besonders verwirrend, also nahm ich die beiden, die in derselben Textzeile erwähnt wurden, in der Hoffnung, dass sie aus derselben Quelle stammen würden.

Wer konsistentere Werte findet, kann gerne kommentieren.

Ernest W. Browns On a New Family of Periodic Orbits in the Problem of Three Bodies: (Tafeln 6, 7.) in MNRAS, 71, (5), S. 438–454, veröffentlicht am 10. März 1911, scheint dort zu sein, wo Hufeisenbahnen waren zuerst vorgeschlagen. ( Auch hier erhältlich). Es beginnt:

Es gibt vier bekannte Asteroiden, die um den einen oder anderen Eckpunkt der beiden gleichseitigen Dreiecke zu oszillieren scheinen, deren Basis die Verbindungslinie zwischen Jupiter und Sonne ist. Diese Scheitelpunkte sind die wohlbekannten Positionen des relativen Gleichgewichts. Der heliozentrische Vektor eines dieser Asteroiden kann sich anscheinend bis zu 17° von seiner Gleichgewichtslage entfernen.* Die Schwingungen können daher nicht als sehr klein angesehen werden. Man fragt sich natürlich, ob Schwingungen dieser Art in noch grösseren Bögen möglich sind; und wenn ja, auf welche Weise kann er die Umlaufbahnen am bequemsten erhalten?

*LJ Linders, Arhivfor Mat., Ast. och Fys., So. Tierarzt. Ak. Ich Stockholm, Bd. 4, Nr. 20.

---

---

Ich werde einige Hufeisenbahnen im Circular Restricted Three-Body Problem-Formalismus erstellen und sie in Python darstellen und dann mit der Schätzung der synodischen Periode vergleichen, die in der Antwort von @AtmosphericPrisonEscape beschrieben wird .

tl;dr: Gute qualitative Übereinstimmung, keine Überraschungen!

---

Eine kurze Zusammenfassung der CR3BP-Mathematik in dimensionslosen Einheiten. Der Abstand zwischen den beiden Körpern ist gleich 1, ebenso die Gravitationskonstante. Sie umkreisen einen gemeinsamen Schwerpunkt auf Kreisbahnen mit einer Periode von$2 \pi$. Es ist einfacher zu visualisieren und zu berechnen, wenn Sie es in einem rotierenden Rahmen tun, sodass die beiden Massen fest sind. Der dritte Körper an Position$x, y, z$wird angenommen, dass es keine Gravitationswirkung auf die ersten beiden hat,

$$\mu = \frac{m_2}{m_1 + m_2}$$

$$x_1 = -\mu$$ $$x_2 = 1-\mu$$

$$r_1 = \sqrt{(x-x1)^2 + y^2 + z^2}$$ $$r_2 = \sqrt{(x-x2)^2 + y^2 + z^2}$$

Die Jacobi-Energie$C$ist eine Erhaltungsgröße in diesem Rotationssystem:

$$C = x^2 + y^2 + 2\frac{1-\mu}{r_1} + 2\frac{\mu}{r_2} - (\dot{x}^2 + \dot{y}^2 + \dot{z}^2)$$

bei dem die$x^2 + y^2$ist das Pseudopotential. Wenn Sie die geschwindigkeitsabhängigen Terme einstellen$(\dot{x}^2 + \dot{y}^2 + \dot{z}^2)$auf Null, erhalten Sie eine Oberfläche mit Nullgeschwindigkeit , diese Oberfläche, die in viele/die meisten Fragen zu drei Körperumlaufbahnen eingefügt wird. Diese Diagramme gelten nicht, wenn sich ein Objekt bewegt, und Sie können ihnen daher keine Umlaufbahnen überlagern!

Die Beschleunigung, die der dritte Körper in diesem rotierenden Rahmen spürt, hat beides, was erwartet wird$1/r^2$Kräfte und eine geschwindigkeitsabhängige Pseudokraft , die nicht real ist, aber die Tatsache berücksichtigt, dass sich der Rahmen dreht und nicht träge ist.

$$\ddot{x} = x + 2\dot{y} - \frac{(1-\mu)(x+\mu)}{r_1^3} - \frac{\mu(x-1+\mu)}{r_2^3}$$ $$\ddot{y} = y - 2\dot{x} - \frac{(1-\mu)y}{r_1^3} - \frac{\mu y}{r_2^3}$$ $$\ddot{z} = -\frac{(1-\mu) z}{r_1^3} - \frac{\mu z}{r_2^3}$$

---

Hier sind einige Berechnungen. Ich entschied mich$\mu = 0.001$was der Situation Jupiter und Sonne ziemlich nahe kommt. Ich wählte eine Reihe von Startpunkten am gegenüberliegenden Punkt aus$m_2$bei etwa$x=-1$ aber das habe ich wirklich nicht gemacht. Was ich wirklich getan habe, ist eine Reihe von Startgeschwindigkeiten zu wählen$-0.08 < \dot{y} < 0.08$und für jeden berechnete ich die Position auf der$x$Achse in der Nähe$x=-1$wo die Beschleunigung in der$x$Richtung war Null.

Das gibt den Lösungen ein kleines bisschen Startsymmetrie, aber Halo-Orbits sind holprig und wackelig und nicht immer so stabil, also ist dieser Aufwand nicht wirklich notwendig.

Ich propagierte jede Umlaufbahn, bis sie in denselben Bereich zurückkehrte, und stoppte sie, als sie die x-Achse kreuzte, wodurch eine Familie von Halbzyklen erzeugt wurde.

Um es kurz zu machen, die in der Antwort von @AtmosphericPrisonEscape gezeigte Methode zur Schätzung der Zykluszeit durch Berechnung der synodischen Periode im Trägheitsrahmen stimmt ziemlich gut mit diesen Halo-Umlaufbahnen überein, und das sollte nicht sehr überraschend sein!

oben: Halbzyklen einiger wackeliger Hufeisenbahnen

oben: Zeiten bis zum ersten Kreuzen der x-Achse der gleichen wackeligen Hufeisenbahnen, die zur Berechnung der Halbzykluszeiten verwendet werden.

oben: Zykluszeiten aus dieser Berechnung (schwarze Punkte) versus aus der Schätzmethode der synodischen Periode (rote Punkte). Gute qualitative Übereinstimmung. Auch die Startgeschwindigkeiten y an jedem Startpunkt in x.

unten: Python-Skript für diese Plots.

def x_acc(x, ydot):
    r1    = np.abs(x-x1)
    r2    = np.abs(x-x2)
    xddot = x + 2*ydot  -  ((1-mu)/r1**3)*(x+mu) - (mu/r2**3)*(x-(1-mu))
    return xddot

def C_calc(x, y, z, xdot, ydot, zdot):
    r1 = np.sqrt((x-x1)**2 + y**2 + z**2)
    r2 = np.sqrt((x-x2)**2 + y**2 + z**2)
    C = (x**2 + y**2 + 2.*(1-mu)/r1 + 2.*mu/r2 - (xdot**2 + ydot**2 + zdot**2))
    return C

def deriv(X, t): 
    x, y, z, xdot, ydot, zdot = X
    r1 = np.sqrt((x-x1)**2 + y**2 + z**2)
    r2 = np.sqrt((x-x2)**2 + y**2 + z**2)
    xddot = x + 2*ydot  -  ((1-mu)/r1**3)*(x+mu) - (mu/r2**3)*(x-(1-mu))
    yddot = y - 2*xdot  -  ((1-mu)/r1**3)*y      - (mu/r2**3)*y
    zddot =             -  ((1-mu)/r1**3)*z      - (mu/r2**3)*z
    return np.hstack((xdot, ydot, zdot, xddot, yddot, zddot))

# http://cosweb1.fau.edu/~jmirelesjames/hw4Notes.pdf

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint as ODEint
from scipy.optimize import brentq

halfpi, pi, twopi = [f*np.pi for f in (0.5, 1, 2)]

mu = 0.001

x1 = -mu
x2 = 1. - mu

x = np.linspace(-1.4, 1.4, 1201)
y = np.linspace(-1.4, 1.4, 1201)

Y, X = np.meshgrid(y, x, indexing='ij')
Z    = np.zeros_like(X)

xdot, ydot, zdot = [np.zeros_like(X) for i in range(3)]

C = C_calc(X, Y, Z, xdot, ydot, zdot)
C[C>8] = np.nan

if True:
    plt.figure()
    plt.imshow(C)
    plt.colorbar()
    levels = np.arange(2.9, 3.2, 0.04) 
    CS = plt.contour(C, levels,
                 origin='lower',
                 linewidths=2) 
    plt.show()

ydot0s   = np.linspace(-0.08, 0.08, 20)
x0ydot0s = []
for ydot0 in ydot0s:
    x0, infob =  brentq(x_acc, -1.5, -0.5, args=(ydot0), xtol=1E-11, rtol=1E-11,
                           maxiter=100, full_output=True, disp=True)
    x0ydot0s.append((x0, ydot0))

states = [np.array([x0, 0, 0, 0, ydot0, 0]) for (x0, ydot0) in x0ydot0s]

times  = np.arange(0, 150, 0.01)

results = []
for X0 in states:
    answer, info = ODEint(deriv, X0, times, atol = 1E-11, full_output=True)
    results.append(answer.T.copy())

resultz = []
for x0ydot0, thing in zip(x0ydot0s, results):
    y     = thing[1]
    check = y[2:]*y[1:-1] < 0
    zc    = np.argmax(y[2:]*y[1:-1] < 0) + 1
    if zc > 10:
        resultz.append((thing, zc, x0ydot0))

if True:
    plt.figure()
    hw = 1.6
    for j, (thing, zc, x0ydot0) in enumerate(resultz):
        x, y = thing[:2,:zc]
        plt.plot(x, y)
    plt.xlim(-hw, hw)
    plt.ylim(-hw, hw)
    plt.plot([x1], [0], 'ok')
    plt.plot([x2], [0], 'ok')
    plt.show()

if True:
    plt.figure()
    for j, (thing, zc, x0ydot0) in enumerate(resultz):
        x, y = thing[:2]
        plt.plot(times[:zc], y[:zc])
    plt.show()

if True:
    plt.figure()
    for j, (thing, zc, x0ydot0) in enumerate(resultz):
        x0, ydot0 = x0ydot0
        cycle_time = 2. * times[zc] / twopi
        ratio = abs(x0/x2)
        T_simple_model = twopi * abs(x0/x2)**1.5
        T_synodic_simple_model = 1. / (1. - twopi/T_simple_model) # https://astronomy.stackexchange.com/a/25002/7982
        plt.subplot(2, 1, 1)
        plt.plot(x0, cycle_time, 'ok')
        plt.plot(x0, abs(T_synodic_simple_model), 'or')
        plt.subplot(2, 1, 2)
        plt.plot(x0, ydot0, 'ok')
    plt.subplot(2, 1, 1)
    plt.xlabel('x0', fontsize=16)
    plt.ylabel('cycle times (periods)', fontsize=16)
    plt.subplot(2, 1, 2)
    plt.xlabel('x0', fontsize=16)
    plt.ylabel('ydot0', fontsize=16)
    plt.show()