# Capacité numérique n°8 : non isochronisme d'un pendule

import numpy as np
import matplotlib.pyplot as plt
from math import *
from scipy.integrate import odeint

## Script 1 : évolutions temporelles de l'angle en fonction de t pour différentes amplitudes initiales

# Donnée du problème
p_propre = 1   #  pulsation propre de 1 rad/s

# Calculs
def pendule(X,t):
    theta = X[0] # première ligne de X(t)
    omega = X[1] # seconde ligne de X(t)
    X_point = [omega,-p_propre**2*np.sin(theta)] # dû à l'équa diff
    return X_point

t = np.linspace(0,20,200)

for theta0 in np.array([0.25*i for i in range (1,11)]):
    CI = np.array([theta0,0])
    X = odeint(pendule,CI,t)
    theta = X[:,0]
    plt.plot(t,theta)


# Représentations graphiques
plt.figure(1)
plt.title('Elongation angulaire pour 10 amplitudes initiales')
plt.ylabel('theta en rad')
plt.xlabel('temps en s')
plt.grid()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import quad

## Script 2 Tracé de T en fonction de theta0
from scipy.integrate import quad

xmax= theta0
xmin = 0
def periode(theta):
    return 8**(1/2)/p_propre*((np.cos(theta)-np.cos(theta0))**(1/2))
res= quad(periode, xmin, xmax)


liste_theta0=[]
liste_periode=[]
for i in range(0,90,1):
    xmin=0
    theta0=i/180*np.pi
    xmax=theta0
    liste_theta0.append(theta0*180/(np.pi))
    res= quad(periode, xmin, xmax)
    liste_periode.append(res[0])
    plt.figure(2)
plt.plot(liste_theta0,liste_periode)
plt.title('Période du pendule simple en fonction de theta(0)')
plt.xlabel("angle theta à t=0")
plt.ylabel("Période du pendule simple (s)")
plt.grid()
plt.show()