In classical mechanics the double pendulum demonstrates the advantage of solving mechanics problems using the Euler-Lagrange equations derived from the calculus of variations, \[{d \over dt} {\partial L \over \partial \dot{f}} - {\partial L \over \partial f} = 0 ,\] over Newton's Second Law, \[ -{ \partial U \over \partial f } = m {\ddot{f}} .\] The problem consists of a pendulum suspended from a fixed point with another pendulum attached to it. The relevant parameters consist of \(l \) and \(\ell \), the length of the top and bottom pendulum respectively, \(m \) and \( M \), the mass of the top and bottom pendulum respectively, as well as \( \theta \) and \( \phi \) where \( \theta \) is the angle made by the vertical from the fixed point the top pendulum is attached to and the rod of the pendulum and \( \phi \) is the angle made by the vertical from the position of the first mass and the rod of the second pendulum as shown in the figure below:
Specifying the \(x\) and \(y\) coordinates as follows,
\[ x_{m} = l\sin{\theta}, \quad x_{M} = \ell \sin{\phi} + l\sin{\theta}\]
\[ y_{m} = l\cos{\theta}, \quad y_{M} = \ell \cos{\phi} + l\cos{\theta} ,\]
and the respective time derivatives as
\[ \dot{x}_{m} = \dot{\theta} l\cos{\theta}, \quad \dot{x}_{M} = \dot{\phi} \ell \cos{\phi} + \dot{\theta} l\cos{\theta}\]
\[ \dot{y}_{m} = -\dot{\theta} l\sin{\theta}, \quad \dot{y}_{M} = -\dot{\phi} \ell \sin{\phi} - \dot{\theta} l\sin{\theta} .\]
The kinetic energy's can be written as,
\[ T_{m} = \frac{1}{2}m(\dot{x}_{m}^2 + \dot{y}_{m}^2) \]
\[ T_{M} = \frac{1}{2}M(\dot{x}_{M}^2 + \dot{y}_{M}^2) ,\]
and the potential energy's as,
\[ U_{m} = mgy_{m} \]
\[ U_{M} = Mgy_{M} .\]
Thus the lagrangian (\( L = T - U \)) of the system can be found
\[ L = \frac{1}{2}m(l \dot{\theta})^2 + \frac{1}{2}M[ (l \dot{\theta})^2 + 2(l \dot{\theta})(\ell \dot{\phi})
\cos{(\phi - \theta) + (\ell \dot{\phi})^2}] - [-mgl \cos{\theta} - Mg(l \cos{\theta} + \ell \cos{\phi} )].
\]
Using the Euler-Lagrange equation the equations of motion for the two masses can be found as a system of two coupled
non-linear 2nd order differential equations
\[ \ddot{\theta} = \frac{-g(2m+M)\sin{\theta} - Mg\sin{(\theta - 2\phi}) - 2\sin{(\theta - \phi)}M(\dot{\phi}^2 \ell + \dot{\theta}^2l\cos{(\theta - \phi))} }{l(2m + M - M\cos{(2\theta - 2\phi)}) } \]
\[ \ddot{\phi} = \frac{2\sin{(\theta - \phi)}(\dot{\theta}^2l(m + M) + g(m+ M)\cos{\theta} + \dot{\phi}^2 \ell M\cos{(\theta - \phi)}) }{\ell(2m + M - M\cos{(2\theta - 2\phi)})} .\]
To numerically evaluate this set of ODE's we'll need to reduce it from a system of two 2nd order ODE's to a system of four 1st order
ODE's by letting \( u = \theta, v = \dot{\theta}, w = \phi, s = \dot{\phi} \), which gives us the system
\[ \dot{u} = v \]
\[ \dot{v} = \frac{-g(2m+M)\sin{u} - Mg\sin{(u - 2w}) - 2\sin{(u - w)}M(\dot{w}^2 \ell + \dot{u}^2l\cos{(u - w))} }{l(2m + M - M\cos{(2u - 2w)}) } \]
\[ \dot{w} = s \]
\[ \dot{s} = \frac{2\sin{(u - w)}(\dot{u}^2l(m + M) + g(m+ M)\cos{u} + \dot{w}^2 \ell M\cos{(u - w)}) }{\ell(2m + M - M\cos{(2u - 2w)})} .\]