Double elastic pendulum problem solver

This webpage uses the Runge-Kutta-Fehlberg fourth-order method with fifth-order error checking (RKF45) to approximate the solution to the problem of the double elastic pendulum.

The equations that will be integrated here are derived in this article.

Figure 1: Diagram of the double elastic pendulum.

The ordinary differential equation system being solved is

[r¨1r¨2θ¨1θ¨2]=[1m2cosΔm1+m20m2r2sinΔm1+m2cosΔ1r1sinΔ00m2sinΔ(m1+m2)r11m2r2cosΔ(m1+m2)r1sinΔr20r1cosΔr21]1[r1θ˙12gsinθ1+m2m1+m2(r2θ˙22cosΔ+2r˙2θ˙2sinΔ)+Qr1k1(r1l1)m1+m2r2θ˙22gsinθ2+r1θ˙12cosΔ2r˙1θ˙1sinΔ+Qr2k2(r2l2)m22r˙1θ˙1r1gcosθ1r1m2(m1+m2)r1[2r˙2θ˙2cosΔr2θ˙22sinΔ]+Qθ1(m1+m2)r122r˙2θ˙2r2gcosθ2r22r˙1θ˙1cosΔr2r1θ˙12sinΔr2+Qθ2m2r22].\begin{aligned} \begin{bmatrix} \ddot{r}_1 \\ \ddot{r}_2 \\ \ddot{\theta}_1 \\ \ddot{\theta}_2 \end{bmatrix} &= \begin{bmatrix} 1 & \dfrac{m_2\cos{\Delta}}{m_1+m_2} & 0 & -\dfrac{m_2r_2\sin{\Delta}}{m_1+m_2} \\ \cos{\Delta} & 1 & r_1\sin{\Delta} & 0 \\ 0 & \dfrac{m_2\sin{\Delta}}{(m_1+m_2)r_1} & 1 & \dfrac{m_2r_2\cos{\Delta}}{(m_1+m_2)r_1} \\ -\dfrac{\sin{\Delta}}{r_2} & 0 & \dfrac{r_1\cos{\Delta}}{r_2} & 1 \end{bmatrix}^{-1} \begin{bmatrix} r_1\dot{\theta}_1^2-g\sin{\theta_1} + \dfrac{m_2}{m_1+m_2}\left(r_2\dot{\theta}_2^2\cos{\Delta} + 2\dot{r}_2\dot{\theta}_2\sin{\Delta}\right) + \dfrac{Q_{r_1}-k_1(r_1-l_1)}{m_1+m_2}\\ r_2\dot{\theta}_2^2-g\sin{\theta_2} + r_1\dot{\theta}_1^2\cos{\Delta} - 2\dot{r}_1\dot{\theta}_1\sin{\Delta} + \dfrac{Q_{r_2}-k_2(r_2-l_2)}{m_2} \\ -\dfrac{2\dot{r}_1\dot{\theta}_1}{r_1} - \dfrac{g\cos{\theta_1}}{r_1} - \dfrac{m_2}{(m_1+m_2)r_1}\left[2\dot{r}_2\dot{\theta}_2\cos{\Delta} -r_2\dot{\theta}_2^2\sin{\Delta}\right] + \dfrac{Q_{\theta_1}}{(m_1+m_2)r_1^2} \\ -\dfrac{2\dot{r}_2\dot{\theta}_2}{r_2}- \dfrac{g\cos{\theta_2}}{r_2} - \dfrac{2\dot{r}_1\dot{\theta}_1\cos{\Delta}}{r_2} - \dfrac{r_1\dot{\theta}_1^2\sin{\Delta}}{r_2} + \dfrac{Q_{\theta_2}}{m_2r_2^2} \end{bmatrix}. \end{aligned}
Simulation parameter form.
Parameter Value Explanation
Gravitational acceleration in metres (m) per second (s) squared (ms2\mathrm{m}\cdot \mathrm{s}^{-2}).
Rest length (m) of pendulum 1.
Rest length (m) of pendulum 2.
Mass (kilograms or kg) of pendulum bob 1.
Mass (kg) of pendulum bob 2.
Spring coefficient for the first pendulum.
Spring coefficient for the second pendulum.
Linear dissipation coefficient for pendulum 1.
Linear dissipation coefficient for pendulum 2.
Quadratic dissipation coefficient for pendulum 1.
Quadratic dissipation coefficient for pendulum 2.
End time (s) for the simulation.
Initial value of r1r_1 (m).
Initial value of r˙1\dot{r}_1 in ms1\mathrm{m}\cdot \mathrm{s}^{-1}.
Initial value of r2r_2 (m).
Initial value of r˙2\dot{r}_2 in ms1\mathrm{m}\cdot \mathrm{s}^{-1}.
Initial value of θ1\theta_1 in radians (r).
Initial value of θ˙1\dot{\theta}_1 in rs1\mathrm{r}\cdot \mathrm{s}^{-1}.
Initial value of θ2\theta_2 in r.
Initial value of θ˙2\dot{\theta}_2 in rs1\mathrm{r}\cdot \mathrm{s}^{-1}.
Error tolerance.
Tolerance type, can be either absolute (0) or relative (1).
Initial step size.
Minimum allowed step size.
Time increment for skipping ahead in animation.
Time you want to skip ahead to in animation when you press the skip button.
Width (in px) of Plotly windows used for plotting and animation below.
Height (in px) of Plotly windows used for plotting and animation below.
Proportion of animation time passed per real time. tScale=1.0t_{\mathrm{Scale}}=1.0 means animation and real time match. tScale<1.0t_{\mathrm{Scale}}<1.0 means the animation is going more slowly than real time. tScale>1.0t_{\mathrm{Scale}}>1.0 means it is going more rapidly.