@Julio “s uitstekend antwoord beschrijft een vliegroutehoek en legt uit dat deze is de hoek tussen de tangentiële richting (loodrecht op de radiale vector ten opzichte van het centrale lichaam) en de huidige snelheidsvector.
Ik heb eerst geprobeerd de hoek uit deze uitdrukking te halen, maar het is duidelijk verkeerd , aangezien $ \ arccos $ een even functie is en de hoek kan gaan van $ – \ pi / 2 $ naar $ \ pi / 2 $:
$$ \ arccos \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) – \ frac {\ pi} {2} \ \ \ \ text {(onjuist!)} $$
Ik heb banen geïntegreerd voor GM ($ \ mu $) en SMA ($ a $) van eenheid en startafstanden van 0.2 tot 1.8. Dat maakt de periode altijd $ 2 \ pi $. Als ik het resultaat van mijn functie plot, krijg ik te veel kronkels.
Welke uitdrukking kan ik gebruiken om het juiste gamma van de vliegbaanhoek te verkrijgen, uitgaande van toestandsvectoren?
Herziene python voor het foutieve deel zou worden gewaardeerd, maar zeker niet nodig voor een antwoord.
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() Opmerkingen
- mocht deze vraag worden TLDRed om aan te geven dat het een coderingsfout was, omdat het nog steeds lijkt te vragen wat ‘ er mis is met de formule
Antwoord
Dit is een probleem dat groepen mensen heeft geplaagd die zeer goed geïnformeerd zijn over orbitale dynamica, maar die leerden verschillende leerboeken te gebruiken: er zijn twee verschillende definities van “vliegroutehoek “!!
Naast $ \ gamma $, de hoek tussen de tangentiële richting en de snelheidsvector, is er $ \ beta $, de hoek tussen de radiale richting en de snelheidsvector. Mensen zeggen vaak “vliegbaanhoek” zonder te zeggen welke definitie ze “gebruiken . Verwarrend! (Ik heb net gemerkt dat het diagram in Julios antwoord ook $ \ beta $ laat zien)
Als u met $ \ beta $ werkt in plaats van $ \ gamma $, wordt $ \ beta $ gegeven door
$$ \ arccos \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) \ tag {1} $$
die van 0 (“recht omhoog”) naar $ \ pi $ ( “recht naar beneden”). Als u $ \ gamma $ gebruikt, is “recht omhoog” $ \ pi / 2 $ en “recht naar beneden” $ – \ pi / 2 $, dus als u $ \ beta $ naar $ \ gamma $ converteert, trekt u $ \ beta $ af $ \ pi / 2 $:
$$ \ gamma = \ pi / 2 – \ arccos \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) \ tag {2} $$
Dit komt overeen met
$$ \ gamma = \ arcsin \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) \ tag {3} $$
Ik ben niet bekend met de taal je hebt gebruikt voor je berekeningen en plots, dus ik heb “niet naar je algoritme gekeken om te zien waarom er” te veel wiebelen “zijn.
Opmerkingen
- Bedankt! Ik ‘ heb tags (getallen) toegevoegd aan vergelijkingen. Zou je zeggen dat er te veel wiebelen is, of is dat wiebelgedrag eigenlijk redelijk? Omdat je $ \ beta $ (eq 1) hetzelfde is als mijn foutieve $ \ gamma $ behalve een offset van een halve pi, moeten de bewegingen in mijn plot hetzelfde zijn als die in een goede plot van je $ \ beta $ (vergelijking 1).
- Het lijkt erop dat er te veel wiebelt. Ik ‘ controleer dat later.
- @uhoh, Eigenlijk is mijn vergelijking 1 slechts het negatief van je vergelijking. Er klopt nog iets anders. Natuurlijk weet je dat iets niet klopt, want alle geplotte $ \ gamma $ s zijn negatief of nul, wat ‘ niet kan zijn, behalve voor een binnenwaartse spiraal. Voor een Kepleriaanse excentrische baan moet $ \ gamma $ precies twee keer nul kruisen, bij periapsis en apoapsis, en monotoon zijn tussen de extremen, zowel voor de korte (extremum via periapsis naar het andere extremum) als lange (extremum via apoapsis naar het andere extremum). ) segmenten. Ik ‘ zal kijken of ik een voorbeeld kan maken van hoe de $ \ gamma $ -curve eruit zou moeten zien.
- Oeps, ik had hierboven moeten zeggen: ” Mijn eq. 2 is gewoon het negatief van jou. ” Ik moet uitloggen en naar bed gaan!
- @uhoh ” tangentieel ” aan een bol gecentreerd op de primaire, niet op de baan. Persoonlijk zeg ik ‘ d liever ” laterale snelheid ” maar mijn eerste professor in orbitale dynamiek op Stanford gebruikte ” tangential “.
Antwoord
Ik vond de fout in het script, het was te wijten aan mijn “homebrew” puntproduct.Ik had een extra kwadraat:
dotted = ((xx*vv)**2).sum(axis=0) # WRONG dotted = (xx*vv).sum(axis=0) # Correct Dus met behulp van deze plus @TomSpilker “s uitstekende verduidelijkingen Ik gebruik de volgende twee methoden om gamma te berekenen:
Methode 1:
$$ \ gamma_1 = \ arcsin \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) \ tag {3} $$
Methode 2:
Een brute-force alternatieve methode om dubbel te controleren:
$$ \ 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 $$
De modulo-bewerking is alleen echt nodig in het computerprogramma, aangezien elke theta afkomstig is van een afzonderlijke arctan2-bewerking:
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) Reacties
- Inderdaad de nieuwe $ \ gamma $ plot is wat ik had verwacht. Hoera! Goed speurwerk.