Double pendulum motion
A simple pendulum exhibits regular oscillations. For minor deviations from its equilibrium position, its movements are harmonic, typified by sine or cosine functions. During nonlinear oscillations, the pendulum’s period varies with its amplitude, yet the motion remains systematic. Thus, the simple pendulum’s behavior under small oscillation conditions accurately represents its fundamental dynamics.
Conversely, a double pendulum displays distinct behavior. At low energy levels, it exhibits beat phenomena during small oscillations. As the energy input increases, the nature of the pendulum’s oscillations shifts dramatically to chaotic motion. Even though a set of ordinary differential equations, which are fully deterministic, can model a double pendulum, it still exhibits unpredictable chaos. This mirrors the behavior seen in the Lorenz system, where a deterministic model comprising just three equations also leads to chaos.
Considering a pendulum of mass m_1 which has another pendulum of mass m_2 attached; the length of the first pendulum rod is L_1 and the second pendulum is L_2. As generalized coordinates, \theta is the angle of the pendulum against the vertical, and \alpha is the angle relative to the first rod (so relative to the angle of the second pendulum relative to the vertical is \theta + \alpha).
A few assumptions:
- the rod has no mass with length L_1 and L_2;
- the mass (represented by the ball of finite radius) is m_1 and m_2;
- all pivots are assumed to be without friction.
The choice of this angle simplify slightly the equations of motions and that could be explained: if there is no gravity, there is a rotational symmetry and therefore an associated conservation law; choosing \alpha relative to \theta and not absolute in that case would be independent when apply a rotation, in a similar case of the coordinates chosen in the previous example of a mass on an inclined plane.
The next step is to find the Lagrangian of the system, written as:
\mathcal L = T - V
The kinetic energy T is the sum of the kinetic energy of the two pendulums; rather than compute them directly in polar form, it is easier to refer temporarily to Cartesian coordinates and then apply a transformation.
T = \frac{m}{2} \left( \dot x^2 + \dot y^2 \right)
The relation for the first pendulum is:
\begin{aligned} & x_1 = L_1\sin(\theta) \\ & y_1 = L_1 \cos(\theta) \end{aligned}
The relation for the second pendulum is:
\begin{aligned} & x_2 = L_1\sin(\theta) + L_2 \sin(\alpha + \theta) \\ & y_2 = L_1 \cos(\theta) + L_2 \cos(\alpha + \theta) \end{aligned}
The velocity in Cartesian Coordinates is the time derivative of these quantity, and the kinetic energy can be calculated. For the first pendulum:
\begin{aligned} & \dot x_1 = L_1\cos(\theta) \dot \theta \\ & \dot y_1 = - L_1 \sin(\theta) \dot \theta \end{aligned}
which gives for the kinetic energy:
T_1 = \frac{m_1}{2} \left( \dot x_1^2 + \dot y_1^2 \right) = \frac{m_1}{2} \left( \cos(\theta)^2 + \sin(\theta)^2 \right)\dot L_1^2 \theta^2 = \frac{m_1}{2} L_1^2 \dot \theta^2
For the second pendulum:
\begin{aligned} & \dot x_2 = L_1\cos(\theta) \dot \theta + L_2 \cos(\alpha + \theta)(\dot \alpha + \dot \theta)\\ & \dot y_2 = - L_1 \sin(\theta) \dot \theta - L_2 \sin(\alpha + \theta)(\dot \alpha + \dot \theta) \end{aligned}
which gives for the kinetic energy:
\begin{aligned} T_2 & = \frac{m_2}{2} \left( \dot x_2^2 + \dot y_2^2 \right) = \frac{m_2}{2} \left[ L_1^2 \dot \theta^2 + L_2^2 (\dot \alpha + \dot \theta)^2 + 2\,L_1L_2 (\cos(\theta)\cos(\alpha + \theta) + \sin(\theta)\sin(\alpha + \theta))(\dot \alpha + \dot \theta)\dot \theta \right] \\ & = \frac{m_2}{2} \left[L_1^2 \dot \theta^2 + L_2^2(\dot \alpha + \dot \theta)^2 + 2\,L_1L_2\cos(\alpha)(\dot \alpha + \dot \theta)\dot \theta \right] \end{aligned}
The above was obtained using the identity:
\cos(\alpha + \theta) = \cos(\alpha)\cos(\theta) - \sin(\alpha)\sin(\theta)
which gives:
\begin{aligned} & \cos(\theta)\cos(\alpha + \theta) + \sin(\theta)\sin(\alpha + \theta) \\ &\quad = \cos(\theta)(\cos(\alpha)\cos(\theta) - \sin(\alpha)\sin(\theta)) + \sin(\theta)(\sin(\alpha)\cos(\theta) + \cos(\alpha)\sin(\theta)) \\ &\quad = \cos(\alpha)\cos(\theta)^2 - \sin(\alpha)\sin(\theta)\cos(\theta) + \sin(\alpha)\sin(\theta)\cos(\theta) + \cos(\alpha)\sin(\theta)^2 \\ &\quad = \cos(\alpha)\cos(\theta)^2 + \cos(\alpha)\sin(\theta)^2 \\ &\quad = \cos(\alpha) \end{aligned}
The total kinetic energy is:
T = T_1 + T_2 = \frac{m_1}{2} L_1^2 \dot \theta^2 + \frac{m_2}{2} \left[L_1^2 \dot \theta^2 + L_2^2(\dot \alpha + \dot \theta)^2 + 2\,L_1L_2\cos(\alpha)(\dot \alpha + \dot \theta)\dot \theta \right]
If there is no potential energy, and therefore L = T there is a cyclic coordinate (\theta) which does not appear into the Lagrangian and therefore there is associated a conservation law (\dot P_\theta = 0):
\frac{\mathrm d}{\mathrm dt} \left( \frac{\partial \mathcal L}{\partial \dot \theta}\right) = \frac{\mathrm d}{\mathrm dt} \left(m_1L_1^2\dot \theta + m_2L_1^2\dot \theta + m_2 L_2^2\left[(\dot \alpha + \dot \theta) + \frac{L_1}{L_2}\cos(\alpha)(\dot \alpha + 2\,\dot \theta)\right]\right) = 0
This law express the conservation of angular momentum.
If there is a gravitational field, then there is a potential gravitational energy and the angular momentum is not conserved; the potential energy is just a function of the vertical distance:
V = V_1 + V_2 = -m_1g\,L_1\cos(\theta) -m_2g(L_1\cos(\theta) + L_2 \cos(\alpha + \theta))
The Lagrangian is therefore:
\begin{aligned} \mathcal L = T - V = & \frac{m_1}{2} L_1^2 \dot \theta^2 + \frac{m_2}{2} \left[L_1^2 \dot \theta^2 + L_2^2(\dot \alpha + \dot \theta)^2 + 2\,L_1L_2\cos(\alpha)(\dot \alpha + \dot \theta)\dot \theta \right] \\ & + m_1g\,L_1\cos(\theta) + m_2g(L_1 \cos(\theta) + L_2 \cos(\alpha + \theta)) \end{aligned}
It is possible to compute the coefficient for Euler-Lagrange equations:
\begin{array}{lcl} \dfrac{\mathrm d}{\mathrm dt} \left( \dfrac{\partial \mathcal L}{\partial \dot \theta}\right) & = & \dfrac{\mathrm d}{\mathrm dt} \left(\,m_1L_1^2\dot \theta + \,m_2L_1^2\dot \theta + m_2 L_2^2\left[(\dot \alpha + \dot \theta) + \dfrac{L_1}{L_2}\cos(\alpha)(\dot \alpha + 2\,\dot \theta)\right] \right) = m_1L_1^2\ddot \theta + \,m_2L_1^2\ddot \theta \\ & & + m_2 L_2^2\left[ (\ddot \alpha + \ddot \theta) \dfrac{L_1}{L_2}\left( \cos(\alpha)(\ddot \alpha + 2\,\ddot \theta) -\sin(\alpha)(\dot \alpha + 2\,\dot \theta) \dot \alpha \right) \right] \\[10pt] \dfrac{\partial \mathcal L}{\partial \theta} & = & -\,L_1m_1g\sin(\theta) -\,L_1m_2g\sin(\theta) - m_2g\,L_2\sin(\alpha + \theta)\\[10pt] \dfrac{\mathrm d}{\mathrm dt} \left( \dfrac{\partial \mathcal L}{\partial \dot \alpha}\right) & = & \dfrac{\mathrm d}{\mathrm dt} \left( m_2L_2^2(\dot \alpha + \dot \theta) + m_2L_1L_2\cos(\alpha)\dot \theta \right) = m_2L_2^2(\ddot \alpha + \ddot \theta) + m_2L_1L_2\left( \cos(\alpha)\ddot \theta - \sin(\alpha)\dot \theta \dot \alpha \right) \\[10pt] \dfrac{\partial \mathcal L}{\partial \alpha} & = & -m_2g\,L_2\sin(\alpha + \theta) -m_2L_1L_2\sin(\alpha)(\dot \alpha + \dot \theta)\dot \theta \end{array}
It is finally possible to write the equations of motion:
\begin{aligned} &\dfrac{\mathrm d}{\mathrm dt} \left( \dfrac{\partial \mathcal L}{\partial \dot \theta}\right) = \dfrac{\partial \mathcal L}{\partial \theta} \\ &m_1L_1^2\ddot \theta + \,m_2L_1^2\ddot \theta + m_2 L_2^2\left[ (\ddot \alpha + \ddot \theta) + \dfrac{L_1}{L_2}\left( \cos(\alpha)(\ddot \alpha + 2\,\ddot \theta) -\sin(\alpha)(\dot \alpha + 2\,\dot \theta) \dot \alpha \right) \right]\\ &\quad \quad = -\,L_1m_1g\sin(\theta) -\,L_1m_2g\sin(\theta) - m_2g\,L_2\sin(\alpha + \theta) \\ &\dfrac{\mathrm d}{\mathrm dt} \left( \dfrac{\partial \mathcal L}{\partial \dot \alpha}\right) = \dfrac{\partial \mathcal L}{\partial \alpha} \\ & m_2L_2^2(\ddot \alpha + \ddot \theta) + m_2L_1L_2\left( \cos(\alpha)\ddot \theta - \sin(\alpha)\dot \theta \dot \alpha \right) \\ & \quad \quad = -m_2g\,L_2\sin(\alpha + \theta) -m_2L_1L_2\sin(\alpha)(\dot \alpha + \dot \theta)\dot \theta \\ & m_2L_2^2(\ddot \alpha + \ddot \theta) + m_2L_1L_2\left( \cos(\alpha)\ddot \theta \right) \\ & \quad \quad = -m_2g\,L_2\sin(\alpha + \theta) -m_2L_1L_2\sin(\alpha)\dot \theta^2 \end{aligned}
In this case there is no longer another conserved quantity except the energy.
Numerical simulation
I employed a fourth-order Runge-Kutta (RK4) method to solve the non-linear differential equations derived from Lagrangian mechanics and simulate the double pendulum dynamics.
This method is chosen for its balance between computational efficiency and accuracy in dealing with systems exhibiting chaotic behavior. The initial conditions, including angles and angular velocities, are set to explore a range of dynamics from regular periodic motion to chaotic behavior. The physical parameters —masses and lengths of the pendulums— are selected to illustrate the sensitivity of the system’s evolution to initial conditions.
The RK4 integration provides a detailed time evolution of the system, allowing for an in-depth analysis of the dynamics of the double pendulum. This approach not only reinforces the theoretical understanding of chaotic systems but also demonstrates the practical application of numerical methods in complex physical systems.