@Julio « s excellente réponse décrit un angle de trajectoire de vol et lexplique est langle entre la direction tangentielle (perpendiculaire au vecteur radial au corps central) et le vecteur de vitesse actuel.
Jai dabord essayé dobtenir langle à partir de cette expression, mais cest évidemment faux , puisque $ \ arccos $ est une fonction paire et que langle peut aller de $ – \ pi / 2 $ à $ \ pi / 2 $:
$$ \ arccos \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) – \ frac {\ pi} {2} \ \ \ \ text {(incorrect!)} $$
Jai intégré des orbites pour GM ($ \ mu $) et SMA ($ a $) dunité et de distances de départ de 0,2 à 1,8. Cela rend la période toujours $ 2 \ pi $. Lorsque je trace le résultat de ma fonction, jobtiens trop de tremblements.
Quelle expression puis-je utiliser pour obtenir le gamma dangle de trajectoire de vol correct à partir de vecteurs détat?
Python révisé pour la partie erronée serait apprécié, mais certainement pas nécessaire pour une réponse.
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() Commentaires
- sur cette question être TLDRed pour indiquer quil sagissait dune erreur de codage, car il semble toujours demander ce que ‘ est incorrect avec la formule
Réponse
Cest un problème qui a tourmenté des groupes de personnes très bien informées sur la dynamique orbitale mais qui ont appris en utilisant différents manuels: il existe deux définitions différentes de « langle de trajectoire de vol « !!
En plus de $ \ gamma $, langle entre la direction tangentielle et le vecteur vitesse, il y a $ \ beta $, langle entre la direction radiale et le vecteur vitesse. Les gens disent souvent « angle de trajectoire de vol » sans dire quelle définition ils « utilisent . Confus! (Je viens de remarquer que le diagramme dans la réponse de Julio montre également $ \ beta $)
Si vous travaillez avec $ \ beta $ au lieu de $ \ gamma $, $ \ beta $ est donné par
$$ \ arccos \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) \ tag {1} $$
qui va de 0 (« tout droit ») à $ \ pi $ ( « vers le bas »). En utilisant $ \ gamma $, « straight up » est $ \ pi / 2 $ et « straight down » est $ – \ pi / 2 $, donc en convertissant $ \ beta $ en $ \ gamma $ il vous suffit de soustraire $ \ beta $ de $ \ pi / 2 $:
$$ \ gamma = \ pi / 2 – \ arccos \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) \ tag {2} $$
Ceci est équivalent à
$$ \ gamma = \ arcsin \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) \ tag {3} $$
Je ne connais pas la langue vous avez utilisé pour vos calculs et graphiques, donc je nai pas regardé votre algorithme pour voir pourquoi il y a « trop de tremblements ».
Commentaires
- Merci! Jai ‘ ajouté des balises (nombres) aux équations. Diriez-vous quil y a trop de tremblements ou est-ce que ce comportement est en fait raisonnable? Puisque votre $ \ beta $ (eq 1) est le même que mon $ \ gamma $ erroné sauf pour un décalage de moitié pi, alors les ondulations de mon tracé devraient être les mêmes que celles dun tracé approprié de votre $ \ beta $ (eq 1).
- Il me semble que trop de remue-méninges. Je ‘ vérifierai cela plus tard.
- @uhoh, En fait, mon eq 1 est juste le négatif de votre équation. Quelque chose ne va pas. Bien sûr, vous savez que quelque chose ne va pas parce que tous les $ \ gamma $ s tracés sont négatifs ou nuls, ce qui peut ‘ être sauf pour une spirale intérieure. Pour une orbite excentrique képlérienne, $ \ gamma $ doit croiser zéro exactement deux fois, au périapside et à lapoapside, et être monotone entre les extrema, à la fois pour le court (extremum passant par le périapsis jusquà lautre extremum) et le long (extremum par lapoapside vers lautre extremum ) segments. Je ‘ Je verrai si je peux dessiner un exemple de ce à quoi la courbe $ \ gamma $ devrait ressembler.
- Oups, jaurais dû dire ci-dessus, » Mon éq. 2 est juste le négatif du vôtre. » Je devrais me déconnecter et aller me coucher!
- @uhoh » tangentiel » à une sphère centrée sur la sphère primaire et non sur lorbite. Personnellement, je ‘ préfère dire » vitesse latérale » mais mon premier professeur de dynamique orbitale à Stanford utilisé » tangentiel « .
Réponse
Jai trouvé lerreur dans le script, elle était due à mon produit scalaire « homebrew ».Jai eu un quadrillage supplémentaire:
dotted = ((xx*vv)**2).sum(axis=0) # WRONG dotted = (xx*vv).sum(axis=0) # Correct Donc, en utilisant ce plus @TomSpilker « s excellentes clarifications Jai utilisé les deux méthodes suivantes pour calculer le gamma:
Méthode 1:
$$ \ gamma_1 = \ arcsin \ left (\ frac {\ mathbf {r \ centerdot v}} {| \ mathbf {r} | \ | \ mathbf {v} |} \ right) \ tag {3} $$
Méthode 2:
Une méthode alternative de force brute pour vérifier:
$$ \ 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 $$
Lopération modulo nest vraiment nécessaire que dans le programme informatique puisque chaque thêta provient dune opération arctan2 distincte:
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) Commentaires
- En effet le nouveau $ \ gamma $ plot est ce à quoi je mattendais. Hourra! Bonne recherche.