import matplotlib.pyplot as plt
import math as m
import scipy.optimize as op 
#déclaration des listes
t,M,u,theta,R,X,Y = [],[],[],[],[],[],[]

 # Données d’astronomie de Mercure
T_rev = 0.240 #période de révolution (an) 
e = 0.206 # excentricité
a = 0.387# demi grand axe en (UA)
N = 40 # Nombre de positions

# résolution des équations de Kepler 
# détermination de la position de l'astre autour de son orbite
for i in range(N):
   t.append(i*T_rev /N)
   M.append(2*m.pi/T_rev *t[i]) 
   u.append(float(op.fsolve (lambda x:x-e*m.sin(x)-M[i],0) ))
# Calcul des coordonnées polaires
   theta.append(2*m.atan((m.sqrt((1+e)/ (1-e))*m.tan(u[i]/2))))
   R.append(a*(1-e**2)/(1+e*m.cos(theta[i]))) 
# calcul des corrdonnées cartésiennes
   X.append(R[i]*m.cos(theta[i])) 
   Y.append(R[i]*m.sin(theta[i]))
#affichage de l’orbite 
plt.grid(True) 
plt.xlabel("distance (U.A)") 
plt.ylabel("distance (U.A)") 
plt.axis('equal')
plt.plot(X,Y,"bo") 
plt.plot(0,0,"go")
plt.show()