Untitled

 avatar
unknown
plain_text
2 years ago
7.2 kB
5
Indexable
import matplotlib.animation as animation
import math
import matplotlib.pyplot as plt
import numpy as np


def conditions_initiales_planetes(planete):
    periode = dictionnaire_planetes[planete]['periode']
    x0 = dictionnaire_planetes[planete]['périhélie']
    y0 = 0
    vx0 = 0
    vy0 = math.sqrt((G * masse_soleil * (1 + dictionnaire_planetes[planete]['excentricite'])) /
                    ((dictionnaire_planetes[planete]['longueur_demi_axe']) *
                     (1 - dictionnaire_planetes[planete]['excentricite'])))
    return x0, y0, vx0, vy0, periode


def methode_euler(planete):
    x, y, v_x, v_y, periode = conditions_initiales_planetes(planete)
    delta_t = periode / nb_valeur

    liste_x = [x]
    liste_y = [y]
    liste_v_x = [v_x]
    liste_v_y = [v_y]
    liste_a_x = []
    liste_a_y = []

    for i in range(8 * nb_valeur):
        r = math.sqrt(x ** 2 + y ** 2)
        a = -G * masse_soleil / r ** 3

        a_x = a * x
        a_y = a * y
        liste_a_x.append(a_x)
        liste_a_y.append(a_y)
        x = x + v_x * delta_t
        y = y + v_y * delta_t
        liste_x.append(x)
        liste_y.append(y)
        v_x = v_x + a_x * delta_t
        v_y = v_y + a_y * delta_t
        liste_v_x.append(v_x)
        liste_v_y.append(v_y)

    return liste_x, liste_y, liste_v_x, liste_v_y


def methode_euler_asymetrique(planete):
    x, y, v_x, v_y, periode = conditions_initiales_planetes(planete)
    delta_t = periode / nb_valeur

    liste_x = [x]
    liste_y = [y]
    liste_v_x = [v_x]
    liste_v_y = [v_y]
    liste_a_x = []
    liste_a_y = []
    liste_t = []

    for i in range(8*nb_valeur):
        r = math.sqrt(x ** 2 + y ** 2)
        a = -G * masse_soleil / r ** 3

        x = x + v_x * delta_t
        y = y + v_y * delta_t
        liste_x.append(x)
        liste_y.append(y)
        a_x = a * x
        a_y = a * y
        liste_a_x.append(a_x)
        liste_a_y.append(a_y)
        v_x = v_x + a_x * delta_t
        v_y = v_y + a_y * delta_t
        liste_v_x.append(v_x)
        liste_v_y.append(v_y)

        liste_t.append(i*delta_t)

    return liste_x, liste_y, liste_v_x, liste_v_y, liste_t, delta_t



def moyenne_d_une_liste(liste):
    somme_elmt = 0
    for elmt in liste:
        somme_elmt += elmt
    moyenne_liste = somme_elmt / len(liste)
    return moyenne_liste


def verifivation_conservation_energetique(liste_x, liste_y, liste_vx, liste_vy):
    liste_energie_potentielle = []
    liste_energie_cinetique = []
    liste_energie_totale = []
    for i in range(nb_valeur):
        energie_potentielle = (-G * masse_soleil * dictionnaire_planetes['Terre']['masse']) / \
                              abs(math.sqrt(liste_x[i] ** 2 + liste_y[i] ** 2))
        liste_energie_potentielle.append(energie_potentielle)
        energie_cinetique = 1 / 2 * dictionnaire_planetes['Terre']['masse'] * \
                            abs(math.sqrt(liste_vx[i] ** 2 + liste_vy[i] ** 2)) ** 2
        liste_energie_cinetique.append(energie_cinetique)
        energie_totale = energie_cinetique + energie_potentielle
        liste_energie_totale.append(energie_totale)
    moyenne_liste_energie_totale = moyenne_d_une_liste(liste_energie_totale)
    liste_energie_totale_veref = liste_energie_totale
    for i in range(len(liste_energie_totale_veref)):
        liste_energie_totale_veref[i] = abs(liste_energie_totale[i] / moyenne_liste_energie_totale - 1)
    conservation = True
    for elmt in liste_energie_totale_veref:
        if elmt > 0.001:
            conservation = False
    if conservation is False:
        print("L'energie n'est pas conservée")
    else:
        print("L'energie est conservée")
    return liste_energie_potentielle, liste_energie_cinetique, liste_energie_totale


def construction_soleil():
    rayon_soleil = 696342000
    x = []
    y = []
    for i in range(nb_valeur):
        x.append(rayon_soleil * math.cos(i * 2 * math.pi / nb_valeur))
        y.append(rayon_soleil * math.sin(i * 2 * math.pi / nb_valeur))
    plt.fill(x, y, color='orange')
    plt.show()
    return x, y


dictionnaire_planetes = {'Terre': {"masse": 5.972e24, 'périhélie': 147098291000, 'excentricite': 0.0167,
                                   'longueur_demi_axe': 149598023000, 'periode': 365.25 * 24 * 60 ** 2},
                         'Venus': {"masse": 4.867e24, 'périhélie': 107476170000, 'excentricite': 0.007,
                                   'longueur_demi_axe': 108209500000, 'periode': 224.7 * 24 * 60 ** 2},
                         'Mercure': {"masse": 3.285e23, 'périhélie': 46001009000, 'excentricite': 0.21,
                                     'longueur_demi_axe': 57909050000, 'periode': 88.0 * 24 * 60 ** 2},
                         'Mars': {"masse": 6.4185e23, 'périhélie': 206655000000, 'excentricite': 0.09339,
                                  'longueur_demi_axe': 227944000000, 'periode': 686.9 * 24 * 60 ** 2},
                         'Jupiter': {"masse": 1.898e27, 'périhélie': 740680000000, 'excentricite': 0.048,
                                     'longueur_demi_axe': 778340000000, 'periode': 4332 * 24 * 60 ** 2},
                         'Saturne': {"masse": 5.684e26, 'périhélie': 1349800000000, 'excentricite': 0.0539,
                                     'longueur_demi_axe': 1426700000000, 'periode': 10754 * 24 * 60 ** 2},
                         'Neptune': {"masse": 1.024e26, 'périhélie': 4459800000000, 'excentricite': 0.00859,
                                     'longueur_demi_axe': 4498400000000, 'periode': 60216 * 24 * 60 ** 2},
                         'Uranus': {"masse": 8.681e25, 'périhélie': 2735000000000, 'excentricite': 0.0047,
                                    'longueur_demi_axe': 2870700000000, 'periode': 30698 * 24 * 60 ** 2}
                         }

G = 6.67408e-11
masse_soleil = 1.989e30
nb_valeur = 400

# Conditions initiales :
liste_xT, liste_yT, liste_vxT, liste_vyT, liste_tT, delta_tT = methode_euler_asymetrique('Terre')

liste_xV, liste_yV, liste_vxV, liste_vyV, liste_tV, delta_tV = methode_euler_asymetrique('Venus')

liste_xM, liste_yM, liste_vxM, liste_vyM, liste_tM, delta_tM = methode_euler_asymetrique('Mercure')



plt.plot(liste_xM, liste_yM, color='grey')
plt.plot(liste_xT, liste_yT, color='blue')
plt.plot(liste_xV, liste_yV, color='brown')
#plt.plot(liste_xMa, liste_yMa, color='orange')
#plt.plot(liste_xJ, liste_yJ, color='red')
#plt.plot(liste_xU, liste_yU, color='green')
#plt.plot(liste_xS, liste_yS, color='purple')
#plt.plot(liste_xN, liste_yN, color='black')
construction_soleil()


Et = verifivation_conservation_energetique(liste_xT, liste_yT, liste_vxT, liste_vyT)
print(Et)


#animation

N=400


fig=plt.figure()
line,=plt.plot([liste_xT[0]],[0],marker="o")
plt.xlim(0,liste_xT[8*N])
plt.ylim(0,liste_yT[8*N])



def animate(i):

    x=liste_xT[i]
    y=liste_yT[i]
    line.set_data([x],[y])
    return line,


ani = animation.FuncAnimation(fig, animate, frames=N, interval=20, blit=True, repeat=False)
plt.show()
Editor is loading...