This is a simulation of a double pendulum. For large motions it is a chaotic system, but for small motions it is a simple linear system.
You can change parameters in the simulation such as mass, gravity, and length of rods (enable the "show controls" checkbox). You can drag the pendulum with your mouse to change the starting position.
(If you don't see the simulation try
instructions for enabling Java.)
For small angles, a pendulum behaves like a linear system (see Simple Pendulum). When the angles are small in the Double Pendulum, the system behaves like the linear Double Spring. In the graph, you can see similar Lissajous curves being generated. This is because the motion is determined by simple sine and cosine functions.
For large angles, the pendulum is nonlinear and the phase graph becomes much more complex. You can see this by dragging one of the masses to a larger angle and letting go.
We regard the pendulum rods as being massless and rigid. We regard the pendulum masses as being point masses. The derivation of the equations of motion is shown below, using the direct Newtonian method.
Kinematics of the Double Pendulum
Kinematics means the relations of the parts of the device, without regard to forces. In kinematics we are only trying to find expressions for the position, velocity, & acceleration in terms of the variables that specify the state of the device.
 x = horizontal position of pendulum mass
 y = vertical position of pendulum mass
 θ = angle of pendulum (0 = vertical downwards, counterclockwise is positive)
 L = length of rod (constant)
We place the origin at the pivot point of the upper pendulum. We regard
y as increasing upwards. We indicate the upper pendulum by subscript 1, and the lower by subscript 2. Begin by using simple trigonometry to write expressions for the positions
x_{1}, y_{1}, x_{2}, y_{2} in terms of the angles
θ_{1}, θ_{2}.
x_{1} = L_{1} sin θ_{1}
y_{1} = −L_{1} cos θ_{1}
x_{2} = x_{1} + L_{2} sin θ_{2}
y_{2} = y_{1} − L_{2} cos θ_{2}
The velocity is the derivative with respect to time of the position.
x_{1}' = θ_{1}' L_{1} cos θ_{1}
y_{1}' = θ_{1}' L_{1} sin θ_{1}
x_{2}' = x_{1}' + θ_{2}' L_{2} cos θ_{2}
y_{2}' = y_{1}' + θ_{2}' L_{2} sin θ_{2}
The acceleration is the second derivative.

x_{1}'' = −θ_{1}'^{2} L_{1} sin θ_{1} + θ_{1}'' L_{1} cos θ_{1}
 (1) 

y_{1}'' = θ_{1}'^{2} L_{1} cos θ_{1} + θ_{1}'' L_{1} sin θ_{1}
 (2) 

x_{2}'' = x_{1}'' − θ_{2}'^{2} L_{2} sin θ_{2} + θ_{2}'' L_{2} cos θ_{2}
 (3) 

y_{2}'' = y_{1}'' + θ_{2}'^{2} L_{2} cos θ_{2} + θ_{2}'' L_{2} sin θ_{2}
 (4) 
Forces in the Double Pendulum
upper mass
lower mass
We treat the two pendulum masses as point particles. Begin by drawing the free body diagram for the upper mass and writing an expression for the net force acting on it. Define these variables:
 T = tension in the rod
 m = mass of pendulum
 g = gravitational constant
The forces on the upper pendulum mass are the tension in the upper rod
T_{1}, the tension in the lower rod
T_{2}, and gravity
−m_{1} g. We write separate equations for the horizontal and vertical forces, since they can be treated independently. The net force on the mass is the sum of these. Here we show the net force and use Newton's law
F = m a.

m_{1} x_{1}'' = −T_{1} sin θ_{1} + T_{2} sin θ_{2}
 (5) 

m_{1} y_{1}'' = T_{1} cos θ_{1} − T_{2} cos θ_{2} −
m_{1} g
 (6) 
For the lower pendulum, the forces are the tension in the lower rod
T_{2}, and gravity
−m_{2} g.

m_{2} x_{2}'' = −T_{2} sin θ_{2}
 (7) 

m_{2} y_{2}'' = T_{2} cos θ_{2} − m_{2} g
 (8) 
In relating these equations to the diagrams, keep in mind that in the example diagram
θ_{1} is positive and
θ_{2} is negative, because of the convention that a counterclockwise angle is positive.
Direct Method for Finding Equations of Motion
Now we do some algebraic manipulations with the goal of finding expressions for θ_{1}'', θ_{2}'' in terms of θ_{1}, θ_{1}', θ_{2}, θ_{2}'. Begin by solving equations (7), (8) for T_{2} sin θ_{2} and T_{2} cos θ_{2} and then substituting into equations (5) and (6).

m_{1} x_{1}'' = −T_{1} sin θ_{1} − m_{2} x_{2}''
 (9) 

m_{1} y_{1}'' = T_{1} cos θ_{1} − m_{2} y_{2}'' − m_{2} g −
m_{1} g
 (10) 
Multiply equation (9) by cos θ_{1} and equation (10) by sin θ_{1} and rearrange to get

T_{1} sin θ_{1} cos θ_{1} = −cos θ_{1} (m_{1} x_{1}'' +
m_{2} x_{2}'')
 (11) 

T_{1} sin θ_{1} cos θ_{1} = sin θ_{1} (m_{1} y_{1}'' + m_{2} y_{2}'' +
m_{2} g +
m_{1} g)
 (12) 
This leads to the equation

sin θ_{1} (m_{1} y_{1}'' + m_{2} y_{2}'' + m_{2} g +
m_{1} g) =
−cos θ_{1} (m_{1} x_{1}'' + m_{2} x_{2}'')
 (13) 
Next, multiply equation (7) by cos θ_{2} and equation (8) by sin θ_{2} and rearrange to get

T_{2} sin θ_{2} cos θ_{2} = −cos θ_{2} (m_{2} x_{2}'')
 (14) 

T_{2} sin θ_{2} cos θ_{2} = sin θ_{2} (m_{2} y_{2}'' + m_{2} g)
 (15) 
which leads to

sin θ_{2} (m_{2} y_{2}'' + m_{2} g) = −cos θ_{2} (m_{2} x_{2}'')
 (16) 
Next we need to use a program such as Mathematica to solve equations (13) and (16) for θ_{1}'', θ_{2}'' in terms of θ_{1}, θ_{1}', θ_{2}, θ_{2}'. Note that we also include the definitions given by equations (14), so that we have 2 equations (13, 16) and 2 unknowns (θ_{1}'', θ_{2}''). The result is somewhat complicated, but is easy enough to program into the computer.
θ_{1}'' =

−g (2 m_{1} + m_{2}) sin θ_{1} − m_{2} g sin(θ_{1} − 2 θ_{2}) − 2 sin(θ_{1} − θ_{2}) m_{2} (θ_{2}'^{2} L_{2} + θ_{1}'^{2} L_{1} cos(θ_{1} − θ_{2}))

L_{1} (2 m_{1} + m_{2} − m_{2} cos(2 θ_{1} − 2 θ_{2}))

θ_{2}'' =

2 sin(θ_{1} − θ_{2}) (θ_{1}'^{2} L_{1} (m_{1} + m_{2}) + g(m_{1} + m_{2}) cos θ_{1} + θ_{2}'^{2} L_{2} m_{2} cos(θ_{1} − θ_{2}))

L_{2} (2 m_{1} + m_{2} − m_{2} cos(2 θ_{1} − 2 θ_{2}))

These are the equations of motion for the double pendulum.
Numerical Solution
The above equations are now close to the form needed for the RungeKutta method. The final step is convert these two 2nd order equations into four 1st order equations. Define the first derivatives as separate variables:
 ω_{1} = angular velocity of top rod
 ω_{2} = angular velocity of bottom rod
Then we can write the four 1st order equations:
θ_{1}' = ω_{1}
θ_{2}' = ω_{2}
ω_{1}' =

−g (2 m_{1} + m_{2}) sin θ_{1} − m_{2} g sin(θ_{1} − 2 θ_{2}) − 2 sin(θ_{1} − θ_{2}) m_{2} (ω_{2}^{2} L_{2} + ω_{1}^{2} L_{1} cos(θ_{1} − θ_{2}))

L_{1} (2 m_{1} + m_{2} − m_{2} cos(2 θ_{1} − 2 θ_{2}))

ω_{2}' =

2 sin(θ_{1}−θ_{2}) (ω_{1}^{2} L_{1} (m_{1} + m_{2}) + g(m_{1} + m_{2}) cos θ_{1} + ω_{2}^{2} L_{2} m_{2} cos(θ_{1} − θ_{2}))

L_{2} (2 m_{1} + m_{2} − m_{2} cos(2 θ_{1} − 2 θ_{2}))

This is now exactly the form needed to plug in to the RungeKutta method for numerical solution of the system.
This web page was first published February 2002.