# Importation de la classe capsule
from Orbiter import Capsule,Stage
%matplotlib inline
plt.rc('xtick',labelsize=14)
plt.rc('ytick',labelsize=14)
# Constante célestes:
G=6.67408e-11 # Constante de gravitation -
Rt=6380000 # Rayon terrestre - m
Rl=1736000 # Rayon luniare - m
Ms=1.9884e30 # Masse du soleil - kg
Mt=5.9722e24 # Masse de la terre - kg
Ml=7.3477e22 # Masse Lune - kg
al=384e6 # Distance Terre-Lune - m
at=149.598262e9 # Distance Terre-Soleil - m
Le module Orbiter permet la résolution simplifiée des équations du mouvement d'un objet spatial soumis à une force de gravité et à une force de poussée. Il tient également compte des effets atmosphériques au décollage et des flux thermique lors d'une rentrée atmosphérique. Le module contient 3 classes d'utilisation différentes:
Vous pouvez trouver en ligne la documentation complète du module Orbiter. Ce notebook présente l'utilisation des deux première class. L'utilisation de la class Lanceur est présenté dans le notebook consacré au lanceur SaturnV.
L'utilisation de la class capsule permet de se familiariser avec les trajectoires spatiales et les manoeuvres. Ici on se propose de placer une capsule sur une orbite initiale circulaire à 400km d'altitude. La class prend comme argument initiaux, l'altitude et la vitesse $r\dot{\theta}$ initiale:
zs=5000e3
us=np.sqrt(G*Mt/(Rt+zs))
cnap=Capsule(z0=zs,u=us) # cnap est le nom de notre capsule
print(cnap)
Si on veut regarder l'évolution de la capsule pendant une periode entière de révolution, on intégre la trajectoire sur une periode avec la méthode run():
# Periode de l'orbite à zs:
T1=2*np.pi*(Rt+zs)/us
cnap.run(dt=1,ntmax=T1)
La capsule a donc parcouru un tour complet ($\theta=360$), on peut tracer la trajectoire à l'aide de la méthode plotraj():
cnap.plotraj('polar')
plt.ylim(0,20000)
On va maintenant changer d'orbite: Calculons le ΔV nécessaire pour bifurquer vers une orbite elliptique de demi-grand axe a2=20000km ayant la terre comme foyer et calculons la nouvelle période T2:
a1=Rt+cnap.z
a2=20000e3
dv=cnap.u*(np.sqrt(1+(G*Mt/cnap.u**2)*((1/a1)-(1/a2)))-1)
T2=2*np.pi*np.sqrt(a2**3/(G*Mt))
print(dv)
On notera ici que l'on peut accéder à la vitesse de notre capsule en faisant cnap.u.
On va maintenant réaliser la manoeuvre à l'aide de la méthode mvr() qui prend comme argument le $\Delta V$ et l'angle de la manoeuvre (ici 0):
cnap.mvr(dv,0) # Manoeuvre
Notez que la vitesse à changée. On intègre maintenant pendant une periode $T_2$ complète et on trace la nouvelle trajectoire:
cnap.run(dt=1,ntmax=T2) # Intégration pendant la periode orbitale T2
# On trace la trajectoire:
cnap.plotraj('polar')
plt.ylim(0,30000)
On peut réaliser les simulations précédentes en utilisant un niveau supérieure de modélisation. La class Stage permet de spécifier les caractéristiques propulsives comme l'ISP, le débit massique de carburant $q_m$ et la paramètre constructif $k_s$ représentant le ratio entre la masse des structures et la masse d'ergols.
On va donc initialiser notre véhicule spatial sur la même orbite que précedemment. Pour les paramètres propulsifs, nous allons d'abord estimer la masse $\delta m$ de carburant nécessaire pour réaliser la manoeuvre précédente de $\Delta V=1161$m/s. Pour cela on utilise la formule de Tsiolkovsky vue dans le notebook précédent:
$$ \Delta V =g_0 ISP \ln{\left(\dfrac{m_i}{m_f}\right)}=g_0 ISP \ln{\left(\dfrac{m}{m-\delta m}\right)} $$
soit
$$ \dfrac{\delta m}{m_e} =(k_s+1)\left(1-e^{\frac{-\Delta V}{g_0 ISP}}\right) $$
Ainsi pour $k_s=0.24$ et $ISP=445$ (typiquement 3ème étage d'Ariane5) on a:
ks=0.24
ISP=445
dv=1161
(ks+1)*(1-np.exp(-dv/(9.81*ISP)))
Ainsi, en prenant 1 tonne d'ergols au départ (avant manoeuvre) on consommera donc:
1000*0.28956982786072571
soit 289.5 kg. En prenant comme référence le débit du moteur vinci: $qm=280$kg/s, la manoeuvre durera donc:
tm=289.5698278607257/280
print(tm)
On peut à présent initialiser notre objet Stage et réaliser la simulation:
cnap=Stage(ISP=445,IC=0.24,qm=280,mc=1000,z0=zs,u0=us,gamma0=0)
print(cnap)
On remarque que par défaut, le moteur est éteint. On peut donc réaliser une orbite sans rien faire comme dans la première partie:
cnap.run(dt=1,ntmax=T1-1)
cnap.plotraj('polar')
plt.ylim(0,20000)
Maintenant nous allons réaliser la manoeuvre en allumant le moteur pendant le temps $t_m$ calculé précedemment. Attention, la précision d'intégration est très importante ici !! Il faudra donc adapter le pas de temps pour atteindre la précision voulue:
cnap.engine_start()
cnap.run(dt=0.001,ntmax=tm/0.001)
cnap.engine_stop()
Puis on réalise une periode complète sur la nouvelle orbite:
cnap.run(dt=1,ntmax=T2)
cnap.plotraj('polar')
plt.ylim(0,30000)
On pourra refaire tous les calculs en changeant les caractéristiques propulsives (changer le moteur Vinci).
from IPython.core.display import HTML
style=open('notebooks.css', "r").read()
HTML(style)