Double Pendulum Simulator

Coordinates:

$x\_1\; =\; L\_1\; \backslash sin\; \backslash theta\_1\; \backslash \backslash \; y\_1\; =\; L\_1\; \backslash cos\; \backslash theta\_1\; \backslash \backslash \; x\_2\; =\; x\_1\; +\; L\_2\; \backslash sin\; \backslash theta\_2\; =\; L\_1\; \backslash sin\; \backslash theta\_1\; +\; L\_2\; \backslash sin\; \backslash theta\_2\; \backslash \backslash \; y\_2\; =\; y\_1\; +\; L\_2\; \backslash cos\; \backslash theta\_2\; =\; L\_1\; \backslash cos\; \backslash theta\_1\; +\; L\_2\; \backslash cos\; \backslash theta\_2$

Lagrangian ( $\mathcal{L} = T - V$ ):

$T\; =\; \backslash frac\{1\}\{2\}\; m\_1\; (L\_1\; \backslash dot\{\backslash theta\}\_1)^2\; +\; \backslash frac\{1\}\{2\}\; m\_2\; \backslash left[\; (L\_1\; \backslash dot\{\backslash theta\}\_1)^2\; +\; (L\_2\; \backslash dot\{\backslash theta\}\_2)^2\; +\; 2\; L\_1\; L\_2\; \backslash dot\{\backslash theta\}\_1\; \backslash dot\{\backslash theta\}\_2\; \backslash cos(\backslash theta\_1\; -\; \backslash theta\_2)\; \backslash right]\; \backslash \backslash \; V\; =\; -m\_1\; g\; L\_1\; \backslash cos\; \backslash theta\_1\; -\; m\_2\; g\; (L\_1\; \backslash cos\; \backslash theta\_1\; +\; L\_2\; \backslash cos\; \backslash theta\_2)$

Equations of Motion (from Euler-Lagrange $\frac{d}{dt} \frac{\partial \mathcal{L}}{\partial \dot{\theta}_i} - \frac{\partial \mathcal{L}}{\partial \theta_i} = Q_i$ ):

Where $M$ is the mass matrix, $F$ contains gravitational and centrifugal/Coriolis terms, and $Q_{ext}$ are external generalized forces (damping, wind, Lorentz). $$\backslash ddot\{\backslash theta\}\_1\; =\; \backslash frac\{-g(2m\_1+m\_2)\backslash sin\backslash theta\_1\; -\; m\_2\; g\; \backslash sin(\backslash theta\_1-2\backslash theta\_2)\; -\; 2\backslash sin(\backslash theta\_1-\backslash theta\_2)m\_2(L\_2\backslash dot\{\backslash theta\}\_2^2\; +\; L\_1\backslash dot\{\backslash theta\}\_1^2\backslash cos(\backslash theta\_1-\backslash theta\_2))\}\{L\_1(2m\_1+m\_2-m\_2\backslash cos(2\backslash theta\_1-2\backslash theta\_2))\}\; +\; Q\text{'}\_\{1\}\; \backslash \backslash \; \backslash ddot\{\backslash theta\}\_2\; =\; \backslash frac\{2\backslash sin(\backslash theta\_1-\backslash theta\_2)\backslash left(L\_1\backslash dot\{\backslash theta\}\_1^2(m\_1+m\_2)\; +\; g(m\_1+m\_2)\backslash cos\backslash theta\_1\; +\; L\_2\backslash dot\{\backslash theta\}\_2^2\; m\_2\; \backslash cos(\backslash theta\_1-\backslash theta\_2)\backslash right)\}\{L\_2(2m\_1+m\_2-m\_2\backslash cos(2\backslash theta\_1-2\backslash theta\_2))\}\; +\; Q\text{'}\_\{2\}$$ $Q'_{i}$ represents the effect of external torques after inverting the mass matrix.

External Forces (Cartesian Components & Generalized Torques $Q_i$):

$\backslash vec\{F\}\_\{g,j\}\; =\; (0,\; -m\_j\; g)\; \backslash \backslash \; \backslash vec\{F\}\_\{d,j\}\; =\; -d\_j\; \backslash vec\{v\}\_j\; \backslash \backslash \; \backslash vec\{F\}\_\{wind\}\; =\; (F\_W\; \backslash cos\; \backslash phi\_W,\; F\_W\; \backslash sin\; \backslash phi\_W)\; \backslash \backslash \; \backslash vec\{F\}\_\{B,j\}\; =\; q\_j\; (\backslash vec\{v\}\_j\; \backslash times\; \backslash vec\{B\})\; =\; q\_j\; (v\_\{j,y\}\; B\_z,\; -v\_\{j,x\}\; B\_z)\; \backslash \backslash \; Q\_\{i,ext\}\; =\; \backslash sum\_\{j=1,2\}\; (\backslash vec\{r\}\_j\; \backslash times\; (\backslash vec\{F\}\_\{d,j\}\; +\; \backslash vec\{F\}\_\{wind\}\; +\; \backslash vec\{F\}\_\{B,j\}))\_z$

Total Energy ( $E = T + V$ ):

$E\; =\; \backslash text\{Kinetic\; Energy\; (T)\}\; +\; \backslash text\{Potential\; Energy\; (V)\}$

Numerical Integration (4th Order Runge-Kutta):

$\backslash vec\{y\}\_\{n+1\}\; =\; \backslash vec\{y\}\_n\; +\; \backslash frac\{h\}\{6\}(\backslash vec\{k\}\_1\; +\; 2\backslash vec\{k\}\_2\; +\; 2\backslash vec\{k\}\_3\; +\; \backslash vec\{k\}\_4)\; \backslash quad\; \backslash text\{where\}\; \backslash quad\; \backslash vec\{y\}\; =\; (\backslash theta\_1,\; \backslash theta\_2,\; \backslash dot\{\backslash theta\}\_1,\; \backslash dot\{\backslash theta\}\_2)$

Energy Conservation: In an ideal system with no damping ($d_1=d_2=0$) and no external forces (wind/Lorentz forces = 0), the total mechanical energy $E = T + V$ should remain constant over time. This is because gravity is a conservative force. However, numerical integration methods like RK4 introduce small computational errors at each step, which can cause the calculated energy to drift slightly over long simulations. The "Energy Variation" metric reflects these numerical fluctuations or actual energy changes if non-conservative forces are active. For the default parameters ($L_1=L_2=1$, $m_1=m_2=1$, $g=9.80665$, $\theta_1=120^\circ$, $\theta_2=-30^\circ$, $\dot{\theta}_1=\dot{\theta}_2=0$), the initial total energy is approximately $1.314$ Joules.