Skip to content

Latest commit

 

History

History
34 lines (28 loc) · 7.66 KB

File metadata and controls

34 lines (28 loc) · 7.66 KB

@def title="Double pendulum solver" @def hassim=true @def params = (g=(val=9.81, desc="Gravitational acceleration (in \(\mathrm{m}\cdot \mathrm{s}^{-2}\))."),l1=(val=1, desc="Length of pendulum 1 in metres."),l2=(val=1, desc="Length of pendulum 2 in metres."),m1b=(val=1, desc="Mass of pendulum bob 1 in kilograms."),m1r=(val=1, desc="Mass of pendulum rod 1 in kilograms."),m2b=(val=1, desc="Mass of pendulum bob 2 in kilograms."),m2r=(val=1, desc="Mass of pendulum rod 2 in kilograms."),b1b=(val=0.0, desc="Linear dissipation coefficient of pendulum bob 1."),b1r=(val=0.0, desc="Linear dissipation coefficient of pendulum rod 1."),c1b=(val=0.0, desc="Quadratic dissipation coefficient of pendulum bob 1."),c1r=(val=0.0, desc="Quadratic dissipation coefficient of pendulum rod 1."),b2b=(val=0.0, desc="Linear dissipation coefficient for pendulum bob 2."),b2r=(val=0.0, desc="Linear dissipation coefficient for pendulum rod 2."),c2b=(val=0.0, desc="Quadratic dissipation coefficient for pendulum bob 2."),c2r=(val=0.0, desc="Quadratic dissipation coefficient for pendulum rod 2."),tf=(val=60, desc="End time for the simulation in seconds."),theta10=(val=0, desc="Initial value of \(\theta_1\) in radians."),dtheta10=(val=0, desc="Initial value of \(\dot{\theta}_1\) in radians per second."),theta20=(val=0, desc="Initial value of \(\theta_2\) in radians."),dtheta20=(val=0, desc="Initial value of \(\dot{\theta}_2\) in radians per second."),epsilon=(val=1e-10,)) @def vars = ["theta1", "dtheta1", "theta2", "dtheta2"]; @def funcs = ["generateTable()","removeTable()", "generateTheta1Theta2PhasePlot(solveProblem(RKF45, readInputs()))","removeTheta1Theta2PhasePlot()", "generateDtheta1Theta1PhasePlot(solveProblem(RKF45, readInputs()))","removeDtheta1Theta1PhasePlot()", "generateTheta1Dtheta2PhasePlot(solveProblem(RKF45, readInputs()))","removeTheta1Dtheta2PhasePlot()", "generateTheta2Dtheta1PhasePlot(solveProblem(RKF45, readInputs()))","removeTheta2Dtheta1PhasePlot()","generateDtheta2Theta2PhasePlot(solveProblem(RKF45, readInputs()))","removeDtheta2Theta2PhasePlot()", "generateDtheta1Dtheta2PhasePlot(solveProblem(RKF45, readInputs()))","removeDtheta1Dtheta2PhasePlot()", "generatePendulumPlots(readInputs(), solveProblem(RKF45, readInputs()))","removePendulumPlots()", "generateTimePlot(solveProblem(RKF45, readInputs()))","removeTimePlot()", "generatePlots(readInputs())","removePlots()", "generateAnimation()","removeAnimation()", "generateTheta1PhaseAnimation()","removeTheta1PhaseAnimation()", "generateTheta2PhaseAnimation()","removeTheta2PhaseAnimation()", "generateAnimations()","removeAnimations()","generateAllOutputs()", "removeAllOutputs()"] @def labels = ["Tabulate the solution","Remove the solution table","Generate a \(\theta_2\) against \(\theta_1\) phase plot","Remove \(\theta_2\) against \(\theta_1\) phase plot","Generate a \(\dot{\theta}_1\) against \(\theta_1\) phase plot","Remove \(\dot{\theta}_1\) against \(\theta_1\) phase plot","Generate a \(\dot{\theta}_2\) against \(\theta_1\) phase plot","Remove \(\dot{\theta}_2\) against \(\theta_1\) phase plot","Generate a \(\dot{\theta}_1\) against \(\theta_2\) phase plot","Remove \(\dot{\theta}_1\) against \(\theta_2\) phase plot","Generate a \(\dot{\theta}_2\) against \(\theta_2\) phase plot","Remove \(\dot{\theta}_2\) against \(\theta_2\) phase plot","Generate a \(\dot{\theta}_2\) against \(\dot{\theta}_1\) phase plot","Remove \(\dot{\theta}_2\) against \(\dot{\theta}_1\) phase plot","Generate pendulum plots","Remove pendulum plots","Generate a \(\theta_1\), \(\dot{\theta}_1\), \(\theta_2\) and \(\dot{\theta}_2\) against time plot","Remove \(\theta_1\), \(\dot{\theta}_1\), \(\theta_2\) and \(\dot{\theta}_2\) against time plot","Generate all solution plots","Remove all plots","Generate an animation of the system","Remove system animation","Generate a \(\theta_1\) phase plot animation","Remove \(\theta_1\) phase plot animation","Generate a \(\theta_2\) phase plot animation","Remove \(\theta_2\) phase plot animation","Generate all animations","Remove all animations", "Generate all outputs", "Remove all outputs"] @def ids=["tableOutputs", "phasePlotTheta1Theta2","phasePlotDtheta1Theta1","phasePlotDtheta1Theta2","phasePlotDtheta2Theta1","phasePlotDtheta2Theta2","phasePlotDtheta2Dtheta1","pendulumPlot","pendulumTimePlot","timePlot","animation","animationTheta1Phase","animationTheta2Phase", "animation", "animation"]

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 pendulum. See this article for the derivation of the equations of motion.

<figure>
    <img src="/doublePendulum/Double pendulum.svg" width="500px"></img>
    <figcaption style="font-weight: bold; font-size: 18px;">Figure 1: Diagram of the double pendulum.</b></figcaption>
</figure>

The equations being solved in this webpage are

\begin{aligned} &\ddot{\theta}1 = \dfrac{-1}{\left(\dfrac{m{1r}}{12} + m_{1b}+m_{2b}\right)l_1} \left[m_{2b}l_2 ( \ddot{\theta}2\cos{(\theta_1-\theta_2)} +\dot{\theta}2^2\sin{(\theta_1-\theta_2)}) + g \cos{\theta_1}\left(\dfrac{m{1r}}{2} +m{2r} +m_{1b} + m_{2b}\right) -(b_{1b} + c_{1b} l_1 \dot{\theta}1)l_1 \dot{\theta}1 \right.\ &\left.+\left(b{2b}+c{2b}\sqrt{l_1^2 \dot{\theta}1^2 + l_2^2 \dot{\theta}2^2 +2l_1 l_2\dot{\theta}1 \dot{\theta}2 \cos{(\theta_1-\theta_2)}}\right)(l_1 \dot{\theta}1 + l_2 \dot{\theta}2 \cos{(\theta_1-\theta_2)}) +\left(b{1r} + \dfrac{c{1r}l_1 \dot{\theta}1}{2}\right) \dfrac{l_1 \dot{\theta}1}{4} \right.\ &\left.+\left(b{2r} + c{2r}\sqrt{l_1^2 \dot{\theta}1^2 + \dfrac{l_2^2 \dot{\theta}2^2}{4} + l_2 \dot{\theta}1 \dot{\theta}2 \cos{(\theta_1 -\theta_2)}}\right)\left(l_1 \dot{\theta}1 + \dfrac{l_2\dot{\theta}2 \cos{\left(\theta_1 - \theta_2\right)}}{2}\right)\right]\ &\ddot{\theta}2 = \dfrac{1}{\left(\dfrac{m{2r}}{12} + m{2b}\right)l_2 - \dfrac{m{2b}^2l_2\cos^2{(\theta_1-\theta_2)}}{\left(\dfrac{m{1r}}{12} + m{1b}+m{2b}\right)}}\left[\dfrac{m{2b}\cos{(\theta_1-\theta_2)}}{\left(\dfrac{m{1r}}{12} + m{1b}+m{2b}\right)}\left[m{2b}l_2\dot{\theta}2^2\sin{(\theta_1-\theta_2)} + g \cos{\theta_1}\left(\dfrac{m{1r}}{2} +m_{2r} +m_{1b} + m_{2b}\right) -(b_{1b} + c_{1b} l_1 \dot{\theta}1)l_1 \dot{\theta}1 \right.\right.\ &\quad\left.\left.+\left(b{2b}+c{2b}\sqrt{l_1^2 \dot{\theta}_1^2 + l_2^2 \dot{\theta}_2^2 +2l_1 l_2\dot{\theta}_1 \dot{\theta}_2 \cos{(\theta_1-\theta_2)}}\right)(l_1 \dot{\theta}1 + l_2 \dot{\theta}2 \cos{(\theta_1-\theta_2)})+\left(b{1r} + \dfrac{c{1r}l_1 \dot{\theta}1}{2}\right) \dfrac{l_1 \dot{\theta}1}{4} +\left(b{2r} + c{2r}\sqrt{l_1^2 \dot{\theta}_1^2 + \dfrac{l_2^2 \dot{\theta}_2^2}{4} + l_2 \dot{\theta}_1 \dot{\theta}_2 \cos{(\theta_1 -\theta_2)}}\right)\left(l_1 \dot{\theta}1 \right.\right.\right. \ &\quad\left.\left.\left.+ \dfrac{l_2\dot{\theta}2 \cos{\left(\theta_1 - \theta_2\right)}}{2}\right)\right] + m{2b}(l_1\dot{\theta}1^2\sin{(\theta_1-\theta_2)}-g\cos{\theta_2}) -\dfrac{1}{4}\left(b{2r} + c{2r}\sqrt{l_1^2 \dot{\theta}_1^2 + \dfrac{l_2^2 \dot{\theta}_2^2}{4} + l_1 l_2 \dot{\theta}_1 \dot{\theta}_2 \cos{(\theta_1 -\theta_2)}}\right)(2l_1 \dot{\theta}1 \cos{(\theta_1-\theta_2)}+ l_2 \dot{\theta}2)\right.\ &\quad\left.-\left(b{2b}+c{2b}\sqrt{l_1^2 \dot{\theta}_1^2 + l_2^2 \dot{\theta}_2^2 +2l_1l_2\dot{\theta}_1 \dot{\theta}_2 \cos{(\theta_1-\theta_2)}}\right)(l_1 \dot{\theta}_1 \cos{(\theta_1-\theta_2)} + l_2 \dot{\theta}_2)\right]. \end{aligned}

{{ insert template.html}}