diff --git a/double_pendule.py b/double_pendule.py new file mode 100644 index 0000000..5ddf3c9 --- /dev/null +++ b/double_pendule.py @@ -0,0 +1,70 @@ +import legrandchien as lgc +import matplotlib.pyplot as plt +import matplotlib.animation as animation +import numpy as np + +g = lgc.Const(-10) +m1 = lgc.Const(1) +m2 = lgc.Const(1) +r1 = lgc.Const(1) +r2 = lgc.Const(1) + +theta1 = lgc.Var("theta1") +theta2 = lgc.Var("theta2") + +x1 = r1 * lgc.cos(theta1) +y1 = r1 * lgc.sin(theta1) + +x2 = x1 + r2 * lgc.cos(theta2) +y2 = y1 + r2 * lgc.sin(theta2) + +T = 0.5 * m1 * (x1.diff()*x1.diff() + y1.diff()*y1.diff()) +T += 0.5 * m2 * (x2.diff()*x2.diff() + y2.diff()*y2.diff()) + +V = m1 * g * x1 +V += m2 * g * x2 + +L = T - V + +init = { + "theta1" : -3 * np.pi/4, + "d_theta1" : 0, + "theta2" : -np.pi/2, + "d_theta2" : 0 +} + +data = L.solve(init, 10, 0.01) + +# Extraire les positions +times = sorted(data.keys()) +positions = [data[t] for t in times] + +x1_positions = [r1.value * np.sin(p["theta1"]) for p in positions] +y1_positions = [r1.value * np.cos(p["theta1"]) for p in positions] + +x2_positions = [x1_positions[i] + r2.value * np.sin(p["theta2"]) for i, p in enumerate(positions)] +y2_positions = [y1_positions[i] + r2.value * np.cos(p["theta2"]) for i, p in enumerate(positions)] + +# Créer la figure et l'axe +fig, ax = plt.subplots(figsize=(8, 8)) +ax.set_xlim(-3, 3) +ax.set_ylim(-3, 3) +ax.set_xlabel("X (m)") +ax.set_ylabel("Y (m)") +ax.grid(True) +ax.invert_yaxis() + +line1, = ax.plot([], [], 'o-', lw=2, markersize=8, label='Pendule 1') +line2, = ax.plot([], [], 'o-', lw=2, markersize=8, label='Pendule 2') +trace, = ax.plot([], [], 'r-', alpha=0.3, lw=1, label='Trace') + +def animate(frame): + line1.set_data([0, x1_positions[frame]], [0, y1_positions[frame]]) + line2.set_data([x1_positions[frame], x2_positions[frame]], [y1_positions[frame], y2_positions[frame]]) + trace.set_data(x2_positions[:frame+1], y2_positions[:frame+1]) + return line1, line2, trace + +ax.legend() +anim = animation.FuncAnimation(fig, animate, frames=len(x1_positions), interval=50, blit=True, repeat=True) + +plt.show() \ No newline at end of file