Pendule simple

$$ \ddot{\theta} + \omega_0^2\sin{\theta} = 0 $$

ce qui correspond à

$$ \frac{\mathrm{d}\theta}{\mathrm{d}t} = \dot{\theta} $$ $$ \frac{\mathrm{d}\dot{\theta}}{\mathrm{d}t} = -\omega_0^2\sin{\theta} $$

ou sous forme vectorielle:

$$ \frac{\mathrm{d}}{\mathrm{d}t}\left( \begin{array}{l} \theta \\ \dot{\theta} \end{array} \right) = \left( \begin{array}{l} \dot{\theta} \\ -\omega_0^2\sin{\theta} \end{array} \right) $$

import matplotlib
matplotlib.use('nbagg')

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

def equation_vectorielle(theta_thetadot, t):
    theta, theta_dot = theta_thetadot
    return [theta_dot, -np.sin(theta)]
t = np.linspace(0, 5*np.pi, 1000)
theta_thetadot_initial = (np.pi/3, 0)
solution = odeint(equation_vectorielle, theta_thetadot_initial, t)
theta, theta_dot = solution.T

Ec = 1./2. * theta_dot**2
Ep = 1 - np.cos(theta)
Etot = Ec + Ep

plt.figure(1)
plt.plot(t, Ec, 'b-', label=r'$E_c$')
plt.plot(t, Ep, 'g-', label=r'$E_p$')
plt.plot(t, Etot, 'r-', label=r'$E_m$')
plt.legend()
plt.xlabel(r'$t$')
plt.ylabel(r'$E$')
plt.show()

plt.figure()
for theta_0 in np.linspace(np.pi/3, np.pi, 10):
    theta_thetadot_initial = (theta_0, 0)
    solution = odeint(equation_vectorielle, theta_thetadot_initial, t)
    theta, theta_dot = solution.T
    plt.plot(theta, theta_dot, '-')
plt.xlabel(r'$\theta$')
plt.ylabel(r'$\dot{\theta}$')
plt.show()