Bonjour
J'essaye de coder un algorithme qui résout les équations du problème à N corps. Cependant quand je trace la trajectoire avec la méthode d'Euler pour 2 corps j'obtiens des droites ce qui est manifestement faux. Si vous pouviez me dire où est l'erreur ça serait très sympathique car je n'ai pas cessé de chercher. Merci d'avance.(et ça ne fonctionne pas aussi pour les méthodes de Runge-Kutta et de Verlet)

Code : Sélectionner tout - Visualiser dans une fenêtre à part
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
 
 
def mullist(coeff,L):
    for k in range(0,len(L)):
        L[k]=coeff*L[k]
    return(L)
 
def somlist(L,M):
    N=[]
    for k in range(0,len(L)):
        N.append(L[k]+M[k])
    return(N)
 
def diflist(L,M):
    N=[]
    for k in range(0,len(L)):
        N.append(L[k]-M[k])
    return(N)
 
def normevecteur(V):
    normecarre=0
    for k in range(0,len(V)):
        normecarre+=V[k]**2
    return(sqrt(normecarre))
 
 
 
def Force_entre_deux_corps(m1,p1,m2,p2):#force de 2 sur 1 à un instant fixé
    G=39.4
    Force=[0,0,0]
    for k in range(0,3):
        Force[k]=(-G*m1*m2*(p2[k]-p1[k]))/(normevecteur(diflist(p2,p1))**3)
    return(Force)
 
def Force_entre_N_corps(j,m,position): #Résultante des forces exercées sur j à un instant fixé
    Resultante=[0,0,0]
    for k in range(0,len(position)):
        if k!=j:
            Force=Force_entre_deux_corps(m[j],position[j],m[k],position[k])
            Resultante[0]=Resultante[0]+Force[0]
            Resultante[1]=Resultante[1]+Force[1]
            Resultante[2]=Resultante[2]+Force[2]
    return((Resultante))
 
def Acceleration(m,position):
    L=[]
    for k in range(0,len(position)):
        L.append(mullist(1/m[k],Force_entre_N_corps(k,m,position))) #Calcul du vecteur acceleration pour chaque corps"
    return(L)
 
def fonction_de_passage(m,Y):
    f=[Y[1],Acceleration(m,Y[0])]
    return(np.array(f))
 
def suite_Euler_pas_fixe(m,t_min,t_max,n,position_initiale,vitesse_initiale):
    t=t_min
    temps=[t_min]
    h=(t_max-t_min)/n  #pas de discrétisation
    Y_corps=np.array([position_initiale,vitesse_initiale])
    Y_temps=[list(Y_corps)]
    while t<=t_max:
        Y_corps=Y_corps+h*fonction_de_passage(m,Y_corps)
        Y_temps.append(list(Y_corps))
        t+=h
        temps.append(t)
    return(temps,np.array(Y_temps))