Double Pendulum Chaos Lab

Machining the rig...

Presets

Live readouts

Graphs

Angles th1, th2 vs time

Phase portrait: th1 vs w1

Kinetic, potential, total energy

Twin separation (log) with fitted exponential

Poincare section (th2, w2 when th1 crosses 0 upwards)

Benettin Lyapunov estimate vs time

Numerical checks

running...

How to use

Pick a preset, or set masses, lengths, gravity and the starting state. Orange is the real pendulum, cyan is its twin started 10^-n rad away. Watch the twin separation graph rise exponentially, then saturate. Drag to orbit, wheel to zoom, camera buttons for presets. Space pauses, R resets.

Models and baked PBR textures made in Blender; Three.js r160 animates every part from the shared physics module.

The science

Two point masses m1, m2 hang from rigid, massless rods of length l1, l2. With angles measured from the downward vertical the Lagrangian is L = T - V with

T = ½ (m1+m2) l1² w1² + ½ m2 l2² w2² + m2 l1 l2 w1 w2 cos(th1 - th2)
V = -(m1+m2) g l1 cos th1 - m2 g l2 cos th2

The Euler-Lagrange equations give the coupled second-order system integrated here (with an optional viscous joint term -b w):

th1'' = [-g(2m1+m2) sin th1 - m2 g sin(th1-2th2) - 2 sin(th1-th2) m2 (w2² l2 + w1² l1 cos(th1-th2))] / [l1 (2m1+m2 - m2 cos(2th1-2th2))]
th2'' = [2 sin(th1-th2) (w1² l1 (m1+m2) + g(m1+m2) cos th1 + w2² l2 m2 cos(th1-th2))] / [l2 (2m1+m2 - m2 cos(2th1-2th2))]

Small angles. Linearising around the rest state gives two normal modes. For m1 = m2 and l1 = l2 = l: w² = (g/l)(2 -+ sqrt2), with th2/th1 = +-sqrt2. The checks below measure those periods from the full nonlinear integrator and compare with the closed form.

Chaos. Away from the linear regime the system is non-integrable: two states a distance d0 apart separate roughly as d(t) = d0 e^(lambda t). The largest Lyapunov exponent lambda is estimated with the Benettin method (a shadow trajectory is renormalised to d0 every 0.05 s and the log growth is averaged). The twin pendulum drawn in cyan starts 10^-n rad away in th1; its separation saturates once it reaches the size of the accessible phase space, so the exponential law only holds on the rising part of the plot. Typical lambda for strongly chaotic starts is of order 1 to 3 per second (a typical value, it depends on the energy). Initial angles with too little energy to flip (roughly E below the level needed to lift the lower mass above the pivot) give lambda near 0 and a quasi-periodic Poincare section drawn as closed curves.

Flip-time map. The 2D page shows a coloured square that integrates 64x64 starting angle pairs from rest and records when either rod first passes the upright position. The fractal boundary between fast-flipping and non-flipping regions is the visible footprint of sensitive dependence.

Numerics. With a symplectic-like small step RK4 conserves total energy to about 1e-7 relative over 20 s even for a chaotic orbit. Explicit Euler gains energy steadily and RK2 drifts visibly; switch the integrator to see it. Simplifications: point masses, rigid massless rods, frictionless pivot unless damping is set, planar motion, no air drag. The time-reversal check integrates forward, flips the velocities and integrates back; recovering the start proves the integrator is accurate over that horizon (it is not a symmetry of the chaos itself, which would eventually amplify rounding errors).