Double Pendulum Simulator

Is a double pendulum chaotic?

Yes, once it swings high enough. A double pendulum is a pendulum hung from the bob of another, and released from large angles its motion is chaotic: two releases that differ by far less than anything you could measure end up swinging completely differently within seconds. This simulator solves its exact equations of motion twice at once, for the pendulum and for a ghost whose upper arm starts a millionth of a degree further round, and at the default release, both arms at 120°, the two lower bobs are visibly apart after 12.1 seconds.

What grows is the gap between them, and it grows roughly exponentially, δ(t) ≈ δ₀ e^(λt), where δ₀ is the starting gap and λ the growth rate. That is why a tiny head start is no protection. At the default release the gap doubles about every half second, so 17 nanometres becomes 20 centimetres in about 23 doublings.

The equations of motion

Measure both angles from straight down and call them θ₁ and θ₂, with ω₁ and ω₂ their rates of change, and write d = θ₁ − θ₂. For two point masses on light rigid rods, Lagrange’s equations give the two angular accelerations exactly:

• α₁ = [−m₂L₂ω₂² sin d − (m₁ + m₂)g sin θ₁ − m₂ cos d (L₁ω₁² sin d − g sin θ₂)] / [L₁(m₁ + m₂ sin² d)]
• α₂ = [(m₁ + m₂)(L₁ω₁² sin d − g sin θ₂) + cos d (m₂L₂ω₂² sin d + (m₁ + m₂)g sin θ₁)] / [L₂(m₁ + m₂ sin² d)]

They come from the Lagrangian Landau and Lifshitz give for the coplanar double pendulum, L = ½(m₁ + m₂)L₁²ω₁² + ½m₂L₂²ω₂² + m₂L₁L₂ω₁ω₂ cos d + (m₁ + m₂)gL₁ cos θ₁ + m₂gL₂ cos θ₂, and nothing in them is approximated: there is no small-angle assumption anywhere. Many references print the same equations with a denominator of 2m₁ + m₂ − m₂ cos(2θ₁ − 2θ₂), which is twice m₁ + m₂ sin² d. The simulator integrates them with fourth-order Runge-Kutta at a step of 1 ms or less and checks itself against the one quantity the true motion must keep fixed, the total energy E = ½(m₁ + m₂)L₁²ω₁² + ½m₂L₂²ω₂² + m₂L₁L₂ω₁ω₂ cos d + (m₁ + m₂)gL₁(1 − cos θ₁) + m₂gL₂(1 − cos θ₂), measured from the hanging rest position. The Energy drift readout is how far the computed value has strayed from its start.

Worked example: how long until the ghost parts company?

Take the defaults: 1 m arms, 1 kg bobs, g = 9.81 m/s², both arms released from rest at 120°, and the ghost’s upper arm a millionth of a degree further round.

1. The energy. Each bob is lifted from its lowest point, so E = (m₁ + m₂)gL₁(1 − cos θ₁) + m₂gL₂(1 − cos θ₂) = 2 × 9.81 × 1.5 + 9.81 × 1.5 = 44.1 J. That is more than both flip thresholds, 2m₂gL₂ = 19.62 J for the lower arm and 2(m₁ + m₂)gL₁ = 39.24 J for the upper, so either arm may go over the top and the Enough energy to flip readout says Both arms.
2. The starting gap. Turning the upper arm by δ carries the whole lower arm with it, so the lower bob moves by 2L₁ sin(δ/2): 2 × 1 m × sin(0.0000005°) = 1.75 × 10⁻⁸ m, about 17 nanometres.
3. The line to cross. Visibly apart means more than a tenth of the full reach, 0.1 × (1 m + 1 m) = 0.2 m.
4. The moment. Stepping both pendulums, the gap first passes 0.2 m at 12.077 s, which the Visibly apart after readout shows as 12.1 s.
5. The growth rate. λ = ln(0.2 m / 1.75 × 10⁻⁸ m) / 12.077 s = 16.25 / 12.077 s = 1.35 per second, the Growth rate λ readout. The gap therefore doubled on average every ln 2 / λ = 0.515 s.

Precision buys time only slowly

Change the ghost’s head start and watch the Visibly apart after readout. At the default release:

• a head start of 1°: visibly apart after 3.37 s
• 0.001°: after 8.55 s
• 0.000001°, a millionth of a degree: after 12.1 s
• 0.000000001°, a billionth of a degree: after 17.1 s

Each thousandfold improvement in how precisely the start is known buys only 3.5 to 5.2 more seconds of agreement. That is exponential growth seen from the other side: if the gap grows as e^(λt), shrinking the starting gap a thousandfold delays any given gap by only ln(1000)/λ, which is 4.9 s at the λ of about 1.4 per second that a long run of this release settles on. Knowing the angles to a billionth of a degree, far finer than any protractor or camera could manage, would still predict this pendulum for under 20 seconds.

That is what separates chaos from ordinary uncertainty. In an ordinary calculation an uncertainty in an input carries through to the answer in proportion, which is what the error propagation calculator works out: halve the uncertainty and the answer is twice as precise. Here halving it buys about half a second. Gravity can behave the same way: the single planet in the orbit simulator follows a regular ellipse, and a second body can make the problem chaotic.

When the double pendulum is not chaotic

Small swings are perfectly regular. With both angles small, sin θ ≈ θ, the equations become linear, and the motion is a mix of two normal modes whose angular frequencies are

• ω² = g/(2m₁L₁L₂) × [(m₁ + m₂)(L₁ + L₂) ∓ √((m₁ + m₂)((m₁ + m₂)(L₁ + L₂)² − 4m₁L₁L₂))]

For 1 m arms and 1 kg bobs that gives periods of 2.62 s, with the arms swinging in step and the lower one √2, about 1.41, times as far, and 1.09 s, with the arms swinging opposite ways. Each mode on its own is simple harmonic motion, the kind the mass spring simulator shows, and the Small-swing periods readout gives both for whatever you set. Release the pendulum from 5° and 7°, close to that √2 ratio, and it swings almost purely in the slow mode: the fast mode’s amplitude is only 0.5 percent of the slow mode’s. Release it from 30° and 30° and the ghost keeps company with it: after 20 seconds the gap is 21 nanometres, hardly more than the 17 it started at. A single arm is never chaotic, whatever the amplitude, which is why the simple pendulum simulator can give an exact period for any release.

Energy draws a sharper line. An arm can only go over the top if the energy covers lifting it there, which takes more than 2m₂gL₂ for the lower arm and more than 2(m₁ + m₂)gL₁ for the upper. Each is a plain mgh, the sum the gravitational potential energy calculator does, with m the mass lifted and h twice the arm’s length. Released from 90° and 90°, the pendulum has 29.43 J against thresholds of 19.62 J and 39.24 J, so the lower arm may flip but the upper arm never will, however long it swings. With equal masses and equal arms, neither arm can flip whenever 2 cos θ₁ + cos θ₂ is 1 or more.

Why the swing on screen is not the true motion for ever

Numerical error grows the same way the ghost’s offset does. At the default release the 1 ms Runge-Kutta step is accurate to about 7 × 10⁻⁹ rad over the first five seconds, but runs at that step and at half of it put the lower bob a tenth of the reach apart after 15.6 seconds. Beyond that the swing on screen still keeps its energy to a few parts in a billion, but it is no longer the swing that began from exactly the angles you set. The billionth-of-a-degree ghost above parts after that horizon, which is why its moment depends on the step: 17.1 s at 1 ms, 17.0 s at a quarter of it. No finite step escapes this, and a smaller one only postpones it by the same slow, logarithmic amount the ghost shows. The model notes below give the measured figures.

The same sensitivity is why the simulator carries its own sine and cosine. JavaScript lets each browser engine round Math.sin in its own way, and measured on a Mac the engines behind Chrome and Safari disagree in the last binary digit of about 4 percent of sines. With the built-in functions, the default swing put the lower bob more than 10 cm apart in the two engines after 25 seconds, so one shared link would have shown two readers two different pendulums. The simulator’s functions use only arithmetic that every engine rounds identically, so a link reproduces the same run everywhere.

Reading the two plots

The upper plot is the gap between the lower bobs on a logarithmic scale, so exponential growth is a straight line, rising λ/ln 10 decades a second. The gap climbs in a ragged ramp from its starting size to the dashed line at a tenth of the reach and then levels off, because two lower bobs on pendulums of the same size can never be more than twice the reach apart. The lower plot is the lower arm’s phase portrait, its angular velocity against its angle. Regular motion draws a tidy band that it keeps retracing. Chaotic motion wanders over a whole region, and every break in the line is the arm going over the top, where its angle jumps between 180° and −180°.

What this model leaves out

• Friction and air. The pivots are frictionless and there is no drag, so the energy is constant and the chaos never dies away. A real double pendulum loses energy to its bearings and the air, so its swings shrink until the motion is small and regular.
• Arms with mass. The rods are massless and the bobs are points. Arms with mass change the equations through their moments of inertia and move the flip thresholds, but not the chaos: Shinbrot, Grebogi, Wisdom and Yorke, reporting on a real double pendulum in 1992, measured small uncertainties growing exponentially at about 7.5 per second on their apparatus.
• The third dimension. Both arms stay in one vertical plane. A real joint can let the lower arm wobble out of that plane, which a planar model cannot represent.
• A perfect ghost. Nobody can release a real pendulum twice from angles a billionth of a degree apart, or even from the same angle twice. That impossibility is the lesson: the ghost is the only way to watch two nearly identical releases side by side.

Common mistakes

• Thinking chaos means randomness. Nothing here is random. The same start always gives the same swing, which is why a shared link reproduces a run exactly. The difficulty is that the same start has to be exactly the same.
• Expecting a better computer to predict it for ever. Precision buys time logarithmically: a start known a thousand times more precisely buys about 5 more seconds.
• Assuming every double pendulum is chaotic. Small swings are regular: released from 30° and 30°, the ghost is still within 21 nanometres of the pendulum after 20 seconds.
• Treating it as one long pendulum. Swinging in step at small angles, the default pendulum has a period of 2.62 s, not the 2.84 s of a single 2 m pendulum, because half its mass is the upper bob, which hangs only 1 m down, where a pendulum swings in 2.01 s.
• Reading the energy drift as the chaos. The drift is the integrator’s error, a few parts in a billion over a minute at the default settings. The chaos is in the gap, which grows more than ten millionfold.
• Measuring the lower angle from the upper arm. Here, as in most textbooks, both angles are measured from straight down, so θ₂ is not the bend at the middle joint.