Méthode d'Euler ordre 2 trajectoire satellite

Bonjour,

Je souhaite programmer la représentation de la trajectoire d'un satellite autour de la Terre en Python, je me concentre pour le moment sur la trajectoire du satellite et non sur l'aspect graphique que je programmerai avec Tkinter par la suite. J'ai obtenu deux équations différentielles d'ordre 2, une selon x et une autre selon y en coordonnées cartésiennes à l'aide d'un PFD, qui sont les suivantes :

dx²/dt²+GM_T*x/((x²+y²)^3/2)=0
dy²/dt²+GM_T*y/((x²+y²)^3/2)=0

(j'ai simplifié par la masse m du satellite et j'ai négligé les frottements pour le moment)

avec G la constante de gravitation et M_T la masse de la Terre.

J'ai appliqué la méthode d'Euler et j'ai le programme suivant :

from scipy import*
from pylab import*
 
def euler(x0,y0,n,tf):
    T=[0]
    X=[0]
    Y=[0]
    h=tf/float(n)
    x,y,G,M_T=x0,y0,6.67*10**(-11),5.972*10**24
    t=0
    while t<tf:
        x+=h*(-G*M_T*x)/((x**2+y**2)**3/2)
        y+=h*(-G*M_T*y)/((x**2+y**2)**3/2)
        t+=h
        X.append(x)
        Y.append(y)
        T.append(t)
    plot(Y,X)
    show()

Mon programme ne fonctionne pas (il renvoi une droite) et je comprends bien sûr que mon problème est dû au fait que je ne prends pas en compte le temps t dont dépendent x et y mais je ne vois pas comment faire. Pourriez-vous m'éclairer ?