1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
|
import numpy as np
import matplotlib.pyplot as plt # importations des modules numpy et matplotlib.pyplot
#données initiales
xM=0 # variable position x de Mars
yM=0 # variable position y de Mars
xT=0
yT=0
Ms= 2*10**30
G= 6.6743* 10**(-11) # constante de gravitation
Rm=3389.5 * 10**(3) # rayon de Mars
Rt=6371* 10**3 # Rayon de la Terre
Mt=5.972* 10**(24) # Masse de la Terre
Mm=6.42* 10**(23) # Masse de Mars
dSM=207* 10**9 # distance Soleil-Mars
dST=150* 10**9 # distance Soleil-Terre
phi=0 # variable angle de déphasage de Mars
def rk4(x, dx, y, deriv):
# /* d1, d2, d3, d4 = estimations des derivees
# yp = estimations intermediaires des fonctions */
ddx = dx/2. # /* demi-pas */
d1 = deriv(x,y) # /* 1ere estimation */
yp = y + d1*ddx
d2 = deriv(x+ddx,yp) #/* 2eme estimat. (1/2 pas) */
yp = y + d2*ddx
d3 = deriv(x+ddx,yp) #/* 3eme estimat. (1/2 pas) */
yp = y + d3*dx
d4 = deriv(x+dx,yp) # /* 4eme estimat. (1 pas) */
#/* estimation de y pour le pas suivant en utilisant
# une moyenne ponderee des derivees en remarquant
# que : 1/6 + 1/3 + 1/3 + 1/6 = 1 */
return y + dx*( d1 + 2*d2 + 2*d3 + d4 )/6
def integrationEDO(start, end, step, v_ini, derivee):
'''
Application de la méthode de RK4 pour chaque temps
'''
n = v_ini.size # nombre d'equa-diff du 1er ordre
# Création du tableau temps
interval = end - start # Intervalle
num_points = int(interval / step) + 1 # Nombre d'éléments
t = np.linspace(start, end, num_points) # Tableau temps t
print('Pas de temps = ',t[1]-t[0],' s') # affiche intervalle de temps
# Initialisation du tableau v
v = np.empty((n, num_points))
# tableau y qui contient les valeurs actuelles de v au temps t[i]
y = np.empty(n)
# Condition initiale
v[:, 0] = v_ini
y[:] = v_ini
distVm=np.empty(num_points)
# Boucle for
for i in range(num_points - 1):
if (y[0]-xT)**2 + (y[1]-yT)**2 < Rt**2:
print("Le vaisseau atterit sur la Terre") #vérification que le vaisseau n'a pas crashé sur Terre
break
if (y[0]-xM)**2 + (y[1]-yM)**2 < Rm**2:
print("Le vaisseau atterit sur Mars") #vérification que le vaisseau n'a pas crashé sur Mars
break
y = rk4(t[i], step, y, derivee) # modifie le tableau y
v[:, i + 1] = y
distVm[i]= ((y[0]-xM)**2 + (y[1]-yM)**2)**(1/2)
# Argument de sortie
return t[:i], v[:,:i], i, distVm[:i]
def derivee_u (t, u):
global xM
global yM
global xT
global yT
global phi
#position de la Terre
xT= dST * np.cos(2*np.pi*t/(365.25*3600*24))
yT= dST * np.sin(2*np.pi*t/(365.25*3600*24))
#positions de Mars
xM= dSM * np.cos(2*np.pi*t/(686.7*3600*24) + np.radians(phi))
yM= dSM * np.sin(2*np.pi*t/(686.7*3600*24) + np.radians(phi))
# Initialisation de la dérivée
du = np.empty(u.size)
# Déclaration des équations
dxt= u[0]-xT
dyt= u[1]-yT
dxm= u[0]-xM
dym= u[1]-yM
du[0] = u[2]
du[1]= u[3]
du[2] = -(G *Mt * dxt)/((dxt**2 + dyt**2)**(3/2)) - (G *Mm * dxm)/((dxm**2 + dym**2)**(3/2))
du[3]= -(G *Mt * dyt)/((dxt**2 + dyt**2)**(3/2)) - (G *Mm * dym)/((dxm**2 + dym**2)**(3/2))
return du
def decoupe(secondes):
""" convertit les secondes en jours/heures/minutes/secondes"""
minutes = secondes // 60
secondes %= 60
heures = minutes // 60
minutes %= 60
jours = heures // 24
heures %= 24
return(f"{jours}jour(s), {heures}heure(s), {minutes}minute(s), {secondes}seconde(s)")
def fonction(tmax, dt, v0, angle, phi_local):
""" affiche la trajectoire du vaisseau (les positions x et y) en fonction du temps ,
trouve les valeurs du temps et des coordonnées à l'aide des condition initales imposés et de la fonction integrationEDO"""
global phi
phi=phi_local # pour que la variable phi soit utilisable dans les autres fonctions
# Conditions initiales
u_ini = np.array([dST, Rt, v0*np.cos(np.radians(angle)) - dST * (2*np.pi/(365.25*3600*24))*np.sin(2*np.pi/(365.25*3600*24)) , v0*np.sin(np.radians(angle))+ dST * (2*np.pi/(365.25*3600*24))*np.cos(2*np.pi/(365.25*3600*24))])
# Méthode d'Euler
t, u, i, distVm = integrationEDO(0, tmax, dt, u_ini, derivee_u)
d = min(list(distVm)) # valeur minimum de de la distance vaisseau-Mars
D = list(distVm) # Liste des distances vaisseau-Mars
#position de la Terre
xT= dST * np.cos(2*np.pi*t/(365.25*3600*24))
yT= dST * np.sin(2*np.pi*t/(365.25*3600*24))
#positions de Mars
xm= dSM * np.cos(2*np.pi*t/(686.7*3600*24) + np.radians(phi))
ym= dSM * np.sin(2*np.pi*t/(686.7*3600*24) + np.radians(phi))
#dimensions du cadre
plt.figure(figsize=(7,7))
#traçage trajectoire du vaisseau
plt.plot(u[0,:],u[1,:],'r',label='Trajectoire du vaisseau')
#traçage Terre
plt.plot(xT, yT,'k',label='Tracé de la Terre')
#traçage trajectoire de la Mars
plt.plot(xm,ym,'b',label='Trajectoire de Mars')
#légendes
plt.title('Trajectoire de retour libre')
plt.xlabel('Position en x')
plt.ylabel('Position en y')
plt.legend()
plt.grid()
plt.figure()
#traçage distance vaisseau Mars dans le temps
plt.plot(t,distVm)
plt.title('Distance Vaisseau Mars en fonction du temps')
plt.xlabel('Temps(s)')
plt.ylabel('Distance Vaisseau Mars (m)')
#plt.ylim(1500* 10**3,8000* 10**3)
#traçage force attraction mars sur vaisseau
#plt.plot(t, )
plt.grid()
plt.figure()
#traçage vitesse vaisseau dans le temps
plt.plot(t, (u[2,:]**2+u[3,:]**2)**(1/2) ,'b')
plt.title('Vitesse du vaisseau en fonction du temps')
plt.ylabel('Vitesse Vaisseau (m.s-1)')
plt.xlabel('Temps(s)')
plt.grid()
plt.figure()
print("distance minimale = ",(distVm.min()-Rm)*10**(-3), 'km', 'à', decoupe(D.index(d)*dt))
plt.show()
print('durée voyage=',decoupe(i*dt)) |
