Everything this simulator computes, term by
term — with the equation behind each and where it shows up on
screen. Angles \(\theta_1\), \(\theta_2\) are measured from
vertical; rods have length \(L_1\), \(L_2\); bobs have mass
\(m_1\), \(m_2\); \(g\) is gravity.
Two pendulums, one chaotic system
A double pendulum is just a pendulum hanging off another
pendulum: a bob on a rod, with a second bob on a second rod hinged
to the first. Its state is fully described by two angles and
their rates of change, \((\theta_1,\dot\theta_1,\theta_2,\dot\theta_2)\).
The bob positions follow directly:
$$x_1=L_1\sin\theta_1,\quad y_1=-L_1\cos\theta_1$$
$$x_2=x_1+L_2\sin\theta_2,\quad y_2=y_1-L_2\cos\theta_2$$
Simple to describe, yet this is one of the simplest mechanical
systems known to behave chaotically — which is exactly why it's
worth simulating a swarm of them at once.
The Lagrangian
Rather than juggle tension forces in two rods, it's far cleaner
to write the system's kinetic and potential energy and let the
Euler–Lagrange equations produce the dynamics automatically:
$$T=\tfrac12(m_1+m_2)L_1^2\dot\theta_1^2+\tfrac12 m_2L_2^2\dot\theta_2^2+m_2L_1L_2\dot\theta_1\dot\theta_2\cos(\theta_1-\theta_2)$$
$$V=-(m_1+m_2)gL_1\cos\theta_1-m_2gL_2\cos\theta_2$$
With \(\mathcal{L}=T-V\), the equations of motion fall out of:
$$\frac{d}{dt}\frac{\partial \mathcal{L}}{\partial \dot\theta_i} - \frac{\partial \mathcal{L}}{\partial \theta_i} = 0$$
The equations of motion
Working through that recipe gives the closed-form angular
accelerations this simulator evaluates every substep
(\(\Delta=\theta_1-\theta_2\)):
$$\ddot\theta_1=\frac{-g(2m_1+m_2)\sin\theta_1-m_2g\sin(\theta_1-2\theta_2)-2\sin\Delta\,m_2(\dot\theta_2^2L_2+\dot\theta_1^2L_1\cos\Delta)}{L_1(2m_1+m_2-m_2\cos2\Delta)}$$
$$\ddot\theta_2=\frac{2\sin\Delta\,(\dot\theta_1^2L_1(m_1+m_2)+g(m_1+m_2)\cos\theta_1+\dot\theta_2^2L_2m_2\cos\Delta)}{L_2(2m_1+m_2-m_2\cos2\Delta)}$$
The \(\cos\Delta\) coupling terms are what make the two arms'
motion inseparable — neither pendulum moves independently of the
other, which is the root of the chaos.
Why tiny differences explode: Lyapunov exponents
Two trajectories that start a tiny distance \(\delta_0\) apart in
state space typically diverge exponentially fast, at a rate set by
the system's Lyapunov exponent \(\lambda\):
$$|\delta(t)|\approx|\delta_0|\,e^{\lambda t}$$
Every pendulum in the swarm starts within the spread
control's angle of its neighbor — sometimes thousandths of a
degree apart. The rainbow fan you watch spreading apart over a few
seconds is this equation, made visible. It's not
randomness: every pendulum runs the exact same deterministic
equations above, with only the tiniest difference in starting
condition.
What's simplified
The simulator integrates with 4th-order Runge–Kutta at a
small fixed timestep — extremely accurate, but not exact;
floating-point error itself becomes a (tiny) extra source of the
same sensitivity being demonstrated, which is fitting here but
worth noting.
No air resistance is physically modeled. The damping
control is a simple linear drag on each joint's angular velocity,
not a first-principles aerodynamic model — useful for taming the
chaos, not for matching a real pendulum's friction.
Rods are treated as massless and perfectly rigid, bobs as point
masses. A real double pendulum (two coupled rigid bodies with
distributed mass) has different moments of inertia and needs a
fuller Lagrangian. And once damping > 0, the system is no
longer energy-conserving — it settles toward hanging straight
down.