.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "api/gallery/chaos/continuous_systems/plot_double_pendulum.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_api_gallery_chaos_continuous_systems_plot_double_pendulum.py: Double Pendulum =============== The double pendulum is a simple mechanical system -- two point masses :math:`m_1`, :math:`m_2` on massless rigid rods of length :math:`l_1`, :math:`l_2`, hanging from a fixed pivot under gravity :math:`g`, with no friction -- whose motion is nonetheless chaotic for large enough swings. In terms of the two rod angles :math:`\theta_1, \theta_2` from vertical and their angular velocities :math:`\omega_1 = \dot\theta_1`, :math:`\omega_2 = \dot\theta_2`, the equations of motion are .. math:: \dot\omega_1 &= \frac{-g(2m_1+m_2)\sin\theta_1 - m_2 g \sin(\theta_1 - 2\theta_2) - 2\sin(\theta_1-\theta_2)\, m_2 \left(\omega_2^2 l_2 + \omega_1^2 l_1 \cos(\theta_1-\theta_2)\right)} {l_1 \left(2m_1 + m_2 - m_2 \cos(2\theta_1 - 2\theta_2)\right)} \\ \dot\omega_2 &= \frac{2 \sin(\theta_1-\theta_2)\left(\omega_1^2 l_1 (m_1+m_2) + g(m_1+m_2)\cos\theta_1 + \omega_2^2 l_2 m_2 \cos(\theta_1-\theta_2)\right)} {l_2 \left(2m_1 + m_2 - m_2 \cos(2\theta_1 - 2\theta_2)\right)} with :math:`\dot\theta_1=\omega_1`, :math:`\dot\theta_2=\omega_2`. This example integrates a trajectory with :meth:`physicskit.chaos.systems.continuous.DoublePendulum.trajectory`, traces the path of the outer bob, and checks how well the integrator conserves total mechanical energy using :func:`physicskit.chaos.utils.metrics.energy_drift`. .. GENERATED FROM PYTHON SOURCE LINES 30-41 .. code-block:: Python import matplotlib.pyplot as plt import numpy as np from physicskit.chaos.systems.continuous import DoublePendulum from physicskit.chaos.utils.metrics import energy_drift from physicskit.chaos.visualizers.dynamic_plots import animate_double_pendulum system = DoublePendulum(m1=1.0, m2=1.0, l1=1.0, l2=1.0, g=9.81) state0 = np.array([np.pi / 2, np.pi / 2, 0.0, 0.0]) .. GENERATED FROM PYTHON SOURCE LINES 42-46 Animation: the pendulum swinging, arms and all ------------------------------------------------------ Rendered kinematically -- the two rigid rods and bobs swinging in real space -- rather than only as an abstract (theta, omega) trajectory. .. GENERATED FROM PYTHON SOURCE LINES 46-50 .. code-block:: Python anim = animate_double_pendulum(system, state0=state0, dt=0.005, n_steps=10000, skip=20) plt.show() .. container:: sphx-glr-animation .. raw:: html .. GENERATED FROM PYTHON SOURCE LINES 51-55 To save the animation to a file instead of (or in addition to) displaying it interactively, use e.g.:: anim.save("double_pendulum_animation.gif", writer="pillow", fps=30) .. GENERATED FROM PYTHON SOURCE LINES 57-61 Integrate --------- Start from a near-horizontal, at-rest configuration -- energetic enough to be chaotic. .. GENERATED FROM PYTHON SOURCE LINES 61-63 .. code-block:: Python t, states = system.trajectory(state0=state0, n_steps=8000, dt=0.005) .. GENERATED FROM PYTHON SOURCE LINES 64-68 Outer bob trajectory --------------------- Convert the generalized coordinates (theta1, theta2) to the Cartesian position of the second (outer) bob. .. GENERATED FROM PYTHON SOURCE LINES 68-77 .. code-block:: Python th1, th2 = states[:, 0], states[:, 1] x2 = system.l1 * np.sin(th1) + system.l2 * np.sin(th2) y2 = -system.l1 * np.cos(th1) - system.l2 * np.cos(th2) fig, ax = plt.subplots(figsize=(6, 6)) ax.plot(x2, y2, lw=0.4, color="indigo") ax.set_aspect("equal") ax.set_title("Double pendulum: outer bob trajectory") .. image-sg:: /api/gallery/chaos/continuous_systems/images/sphx_glr_plot_double_pendulum_002.png :alt: Double pendulum: outer bob trajectory :srcset: /api/gallery/chaos/continuous_systems/images/sphx_glr_plot_double_pendulum_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none Text(0.5, 1.0, 'Double pendulum: outer bob trajectory') .. GENERATED FROM PYTHON SOURCE LINES 78-83 Energy conservation check -------------------------- The double pendulum is a conservative system, so total mechanical energy should stay (nearly) constant; the small residual drift below reflects the RK4 integrator's local truncation error rather than any physical damping. .. GENERATED FROM PYTHON SOURCE LINES 83-92 .. code-block:: Python drift = energy_drift(t, states, system.energy) fig2, ax2 = plt.subplots(figsize=(7, 4)) ax2.plot(t, drift) ax2.set_xlabel("t") ax2.set_ylabel("relative energy drift") ax2.set_title(f"Double pendulum energy drift (max |drift| = {np.max(np.abs(drift)):.2e})") plt.show() .. image-sg:: /api/gallery/chaos/continuous_systems/images/sphx_glr_plot_double_pendulum_003.png :alt: Double pendulum energy drift (max |drift| = 5.51e-05) :srcset: /api/gallery/chaos/continuous_systems/images/sphx_glr_plot_double_pendulum_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 20.089 seconds) .. _sphx_glr_download_api_gallery_chaos_continuous_systems_plot_double_pendulum.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_double_pendulum.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_double_pendulum.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_double_pendulum.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_