@Julios „ ausgezeichnete Antwort beschreibt einen Flugbahnwinkel und erklärt dies ist der Winkel zwischen der Tangentialrichtung (senkrecht zum Radialvektor zum Zentralkörper) und dem aktuellen Geschwindigkeitsvektor.
Ich habe zuerst versucht, den Winkel aus diesem Ausdruck zu ermitteln, aber er ist offensichtlich falsch , da $ \ arccos $ eine gerade Funktion ist und der Winkel von $ – \ pi / 2 $ bis $ \ pi / 2 $ gehen kann:
$$ \ arccos \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) – \ frac {\ pi} {2} \ \ \ \ text {(falsch!)} $$
Ich habe Umlaufbahnen für GM ($ \ mu $) und SMA ($ a $) mit einer Einheit und Startentfernungen von 0,2 bis 1,8 integriert. Das macht den Zeitraum immer $ 2 \ pi $. Wenn ich das Ergebnis meiner Funktion zeichne, bekomme ich zu viele Wackelbewegungen.
Welchen Ausdruck kann ich verwenden, um den richtigen Flugbahnwinkel Gamma ausgehend von Zustandsvektoren zu erhalten?
Überarbeitete Python für den fehlerhaften Teil wäre willkommen, aber für eine Antwort sicherlich nicht erforderlich.
def deriv(X, t): x, v = X.reshape(2, -1) acc = -x * ((x**2).sum())**-1.5 return np.hstack((v, acc)) import numpy as np import matplotlib.pyplot as plt from scipy.integrate import odeint as ODEint halfpi, pi, twopi = [f*np.pi for f in (0.5, 1, 2)] T = twopi time = np.linspace(0, twopi, 201) a = 1.0 rstarts = 0.2 * np.arange(1, 10) vstarts = np.sqrt(2./rstarts - 1./a) # from vis-viva equation answers = [] for r, v in zip(rstarts, vstarts): X0 = np.array([r, 0, 0, v]) answer, info = ODEint(deriv, X0, time, full_output= True) answers.append(answer.T) gammas = [] for a in answers: xx, vv = a.reshape(2, 2, -1) dotted = ((xx*vv)**2).sum(axis=0) rabs, vabs = [np.sqrt((thing**2).sum(axis=0)) for thing in (xx, vv)] gamma = np.arccos(dotted/(rabs*vabs)) - halfpi gammas.append(gamma) if True: plt.figure() plt.subplot(4, 1, 1) for x, y, vx, vy in answers: plt.plot(x, y) plt.plot(x[:1], y[:1], ".k") plt.plot([0], [0], "ok") plt.title("y vs x") plt.subplot(4, 1, 2) for x, y, vx, vy in answers: plt.plot(time, x, "-b") plt.plot(time, y, "--r") plt.title("x (blue) y (red, dashed)") plt.xlim(0, twopi) plt.subplot(4, 1, 3) for x, y, vx, vy in answers: plt.plot(time, vx, "-b") plt.plot(time, vy, "--r") plt.title("vx (blue) vy (red), dashed") plt.xlim(0, twopi) plt.subplot(4, 1, 4) for gamma in gammas: plt.plot(time, gamma) plt.title("gamma?") plt.xlim(0, twopi) plt.show() Kommentare
- sollte diese Frage stellen TLDRed sein, um anzuzeigen, dass es sich um einen Codierungsfehler handelt, da immer noch gefragt wird, was ‚ mit der Formel
falsch ist
Dies ist ein Problem, das Gruppen von Menschen geplagt hat, die sich mit Orbitaldynamik sehr gut auskennen, aber mit verschiedenen Lehrbüchern gelernt haben: Es gibt zwei verschiedene Definitionen für „Flugbahnwinkel“ „!!
Zusätzlich zu $ \ gamma $, dem Winkel zwischen der Tangentialrichtung und dem Geschwindigkeitsvektor, gibt es $ \ beta $, den Winkel zwischen die radiale Richtung und der Geschwindigkeitsvektor. Die Leute sagen oft „Flugbahnwinkel“, ohne zu sagen, welche Definition sie verwenden . Verwirrend! (Ich habe gerade bemerkt, dass das Diagramm in Julios Antwort auch $ \ beta $ zeigt)
Wenn Sie mit $ \ beta $ anstelle von $ \ gamma $ arbeiten, wird $ \ beta $ durch
$$ \ arccos \ left (\ frac {\ mathbf {r \ centerdot v}} angegeben. {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) \ tag {1} $$
von 0 („gerade nach oben“) bis $ \ pi $ ( „gerade nach unten“). Wenn Sie $ \ gamma $ verwenden, ist „gerade nach oben“ $ \ pi / 2 $ und „gerade nach unten“ ist $ – \ pi / 2 $. Wenn Sie also $ \ beta $ in $ \ gamma $ konvertieren, subtrahieren Sie einfach $ \ beta $ von $ \ pi / 2 $:
$$ \ gamma = \ pi / 2 – \ arccos \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) \ tag {2} $$
Dies entspricht
$$ \ gamma = \ arcsin \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) \ tag {3} $$
Ich bin mit der Sprache nicht vertraut Sie haben für Ihre Berechnungen und Diagramme verwendet, daher habe ich Ihren Algorithmus nicht untersucht, um festzustellen, warum „zu viele Wackelbewegungen“ auftreten.
Kommentare
- Danke! Ich ‚ habe Tags (Zahlen) zu Gleichungen hinzugefügt. Würden Sie sagen, dass es zu viele Wackelbewegungen gibt, oder ist dieses Wackelverhalten tatsächlich vernünftig? Da Ihr $ \ beta $ (Gleichung 1) bis auf einen Versatz von einem halben pi mit meinem fehlerhaften $ \ gamma $ identisch ist, sollten die Wackelbewegungen in meinem Plot dieselben sein wie in einem richtigen Plot Ihres $ \ beta $ (Gl. 1).
- Sieht für mich nach zu vielen Wackelbewegungen aus. Ich ‚ werde das später überprüfen.
- @uhoh, eigentlich ist meine Gleichung 1 nur das Negative Ihrer Gleichung. Etwas anderes stimmt nicht. Natürlich wissen Sie, dass etwas nicht stimmt, weil alle eingezeichneten $ \ gamma $ s negativ oder null sind, was ‚ nur eine Einwärtsspirale sein kann. Für eine Keplersche exzentrische Umlaufbahn sollte $ \ gamma $ bei Periapsis und Apoapsis genau zweimal Null kreuzen und zwischen den Extrema sowohl für die kurze (Extremum durch Periapsis zum anderen Extremum) als auch für die lange (Extremum durch Apoapsis zum anderen Extremum) monoton sein ) Segmente. Ich ‚ werde sehen, ob ich ein Beispiel dafür zeichnen kann, wie die $ \ gamma $ -Kurve aussehen soll.
- Hoppla, ich hätte oben sagen sollen, “ Meine Gl. 2 ist nur das Negative von Ihnen. “ Ich sollte mich abmelden und ins Bett gehen!
- @uhoh “ tangential “ zu einer Kugel, die auf der Primärseite zentriert ist, nicht zur Umlaufbahn. Persönlich würde ich ‚ lieber “ Quergeschwindigkeit “ sagen, aber mein erster Professor für Orbitaldynamik in Stanford verwendet “ tangential „.
Antwort
Ich habe den Fehler im Skript gefunden, der auf mein „Homebrew“ -Punktprodukt zurückzuführen ist.Ich hatte eine zusätzliche Quadratur:
dotted = ((xx*vv)**2).sum(axis=0) # WRONG dotted = (xx*vv).sum(axis=0) # Correct Verwenden Sie also dieses Plus @TomSpilker „s ausgezeichnete Erläuterungen Ich habe die folgenden zwei Methoden zur Berechnung des Gammas verwendet:
Methode 1:
$$ \ gamma_1 = \ arcsin \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) \ tag {3} $$
Methode 2:
Eine alternative Brute-Force-Methode zur Überprüfung:
$$ \ theta_r = \ arctan2 (y, x) $$
$$ \ theta_v = \ arctan2 (vy, x) $$
$$ \ theta_ {tanj} = \ theta_r + \ frac {\ pi} {2} $$
$$ \ gamma_2 = \ theta_ {tanj} – \ theta_v $$
$$ \ gamma_ {2mod} = \ mod (\ gamma_2 + \ pi, 2 \ pi) – \ pi $$
Die Modulo-Operation wird nur im Computerprogramm wirklich benötigt, da jedes Theta aus einer separaten arctan2-Operation stammt:
gammas_1, gammas_2 = [], [] for a in answers: xx, vv = a.reshape(2, 2, -1) dotted = (xx*vv).sum(axis=0) rabs, vabs = [np.sqrt((thing**2).sum(axis=0)) for thing in (xx, vv)] gamma_1 = np.arcsin(dotted/(rabs*vabs)) # Per Tom Spilker"s answer Eq. 3 theta_r = np.arctan2(xx[1], xx[0]) theta_v = np.arctan2(vv[1], vv[0]) theta_tanj = theta_r + halfpi gamma_2 = theta_tanj - theta_v gamma_2 = np.mod(gamma_2 + pi, twopi) - pi gammas_1.append(gamma_1) gammas_2.append(gamma_2) plt.figure() plt.subplot(2, 1, 1) for gamma_1 in gammas_1: plt.plot(time, gamma_1) plt.title("gammas_1", fontsize=16) plt.subplot(2, 1, 2) for gamma_2 in gammas_2: plt.plot(time, gamma_2) plt.title("gammas_2", fontsize=16) Kommentare
- In der Tat das neue $ \ Gamma $ Plot ist das, was ich erwartet habe. Hurra! Gute Ermittlungen.