@Julio “s eccellente risposta descrive un angolo di traiettoria di volo e spiega che è langolo tra la direzione tangenziale (perpendicolare al vettore radiale al corpo centrale) e il vettore velocità corrente.
Ho prima cercato di ottenere langolo da questa espressione, ma è ovviamente sbagliato , poiché $ \ arccos $ è una funzione pari e langolo può andare da $ – \ pi / 2 $ a $ \ pi / 2 $:
$$ \ arccos \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) – \ frac {\ pi} {2} \ \ \ \ text {(errato!)} $$
Ho integrato orbite per GM ($ \ mu $) e SMA ($ a $) di unità e distanze iniziali da 0,2 a 1,8. Ciò rende il periodo sempre $ 2 \ pi $. Quando traccio il risultato della mia funzione, ottengo troppe oscillazioni.
Quale espressione posso usare per ottenere il corretto angolo gamma del percorso di volo a partire dai vettori di stato?
Python rivisto per la parte errata sarebbe apprezzato, ma certamente non necessario per una risposta.
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() Commenti
- dovrebbe porre questa domanda essere TLDRed per indicare che si trattava di un errore di codifica, poiché sembra ancora chiedere cosa ‘ non va con la formula
Risposta
Questo è un problema che ha afflitto gruppi di persone molto informate sulle dinamiche orbitali ma che hanno imparato utilizzando libri di testo diversi: ci sono due diverse definizioni di “angolo di traiettoria di volo “!!
Oltre a $ \ gamma $, langolo tra la direzione tangenziale e il vettore di velocità, cè $ \ beta $, langolo tra la direzione radiale e il vettore di velocità. Le persone spesso dicono “angolo della traiettoria di volo” senza dire quale definizione “stanno usando . Confuso! (Ho appena notato che il diagramma nella risposta di Julio mostra anche $ \ beta $)
Se lavori con $ \ beta $ invece di $ \ gamma $, $ \ beta $ è dato da
$$ \ arccos \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) \ tag {1} $$
che va da 0 (“verso lalto”) a $ \ pi $ ( “dritto giù”). Usando $ \ gamma $, “straight up” è $ \ pi / 2 $ e “straight down” è $ – \ pi / 2 $, quindi convertendo $ \ beta $ in $ \ gamma $ devi semplicemente sottrarre $ \ beta $ da $ \ pi / 2 $:
$$ \ gamma = \ pi / 2 – \ arccos \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) \ tag {2} $$
Questo è equivalente a
$$ \ gamma = \ arcsin \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) \ tag {3} $$
Non ho familiarità con la lingua hai usato per i tuoi calcoli e grafici, quindi non ho esaminato il tuo algoritmo per vedere perché ci sono “troppe oscillazioni”.
Commenti
- Grazie! Ho ‘ ho aggiunto tag (numeri) alle equazioni. Diresti che ci sono troppe oscillazioni o questo comportamento dimenante è effettivamente ragionevole? Dato che il tuo $ \ beta $ (eq 1) è lo stesso del mio errato $ \ gamma $ tranne che per un offset di metà pi greco, allora le oscillazioni nel mio grafico dovrebbero essere le stesse di quelle in un grafico corretto del tuo $ \ beta $ (eq 1).
- A me sembrano troppe oscillazioni. ‘ lo controllerò più tardi.
- @uhoh, In realtà, la mia eq 1 è solo il negativo della tua equazione. Qualcosaltro non va. Ovviamente sai che qualcosa non va perché tutti i $ \ gamma $ tracciati sono negativi o zero, il che ‘ non può essere tranne che per una spirale interna. Per unorbita eccentrica kepleriana $ \ gamma $ dovrebbe incrociare lo zero esattamente due volte, in corrispondenza del periasse e dellapoapsi, ed essere monotono tra gli estremi, sia per il corto (estremo attraverso il periapsi allaltro estremo) che per quello lungo (estremo attraverso lapopasso allaltro estremo ) segmenti. ‘ vedrò se riesco a disegnare un esempio di come dovrebbe apparire la curva $ \ gamma $.
- Oops, avrei dovuto dire sopra, ” La mia eq. 2 è solo il tuo negativo. ” Dovrei disconnettermi e andare a letto!
- @uhoh ” tangenziale ” a una sfera centrata sul primario, non allorbita. Personalmente ‘ preferirei dire ” velocità laterale ” ma il mio primo professore di Orbital Dynamics a Stanford ha utilizzato la ” tangenziale “.
Risposta
Ho trovato lerrore nello script, era dovuto al mio prodotto dot “homebrew”.Avevo una quadratura in più:
dotted = ((xx*vv)**2).sum(axis=0) # WRONG dotted = (xx*vv).sum(axis=0) # Correct Quindi utilizzando questo plus @TomSpilker “s eccellenti chiarimenti Ho” utilizzato i seguenti due metodi per calcolare la gamma:
Metodo 1:
$$ \ gamma_1 = \ arcsin \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) \ tag {3} $$
Metodo 2:
Un metodo alternativo a forza bruta per ricontrollare:
$$ \ 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 $$
Loperazione modulo è realmente necessaria solo nel programma del computer poiché ogni theta proviene da unoperazione arctan2 separata:
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) Commenti
- In effetti il nuovo $ \ Il grafico gamma $ è quello che mi aspettavo. Evviva! Buona investigazione.