{
    "nbformat": 4,
    "nbformat_minor": 5,
    "metadata": {
        "kernelspec": {
            "display_name": "Python",
            "language": "python",
            "name": "python"
        },
        "language_info": {
            "name": "python"
        }
    },
    "cells": [
        {
            "id": "18957a32",
            "cell_type": "markdown",
            "source": "<div class='alert alert-warning'>\n\nSciPy's interactive examples with Jupyterlite are experimental and may not always work as expected. Execution of cells containing imports may result in large downloads (up to 60MB of content for the first import from SciPy). Load times when importing from SciPy may take roughly 10-20 seconds. If you notice any problems, feel free to open an [issue](https://github.com/scipy/scipy/issues/new/choose).\n\n</div>",
            "metadata": {}
        },
        {
            "id": "e1a4c6f8",
            "cell_type": "markdown",
            "source": "The second order differential equation for the angle `theta` of a\npendulum acted on by gravity with friction can be written\n```\n\ntheta''(t) + b*theta'(t) + c*sin(theta(t)) = 0\n\n```\nwhere `b` and `c` are positive constants, and a prime (') denotes a\nderivative. To solve this equation with `odeint`, we must first convert\nit to a system of first order equations. By defining the angular\nvelocity ``omega(t) = theta'(t)``, we obtain the system\n```\n\ntheta'(t) = omega(t)\nomega'(t) = -b*omega(t) - c*sin(theta(t))\n\n```\nLet `y` be the vector [`theta`, `omega`]. We implement this system\nin Python as:\n",
            "metadata": {}
        },
        {
            "id": "1eed2d4d",
            "cell_type": "code",
            "metadata": {},
            "execution_count": null,
            "source": "import numpy as np\ndef pend(y, t, b, c):\n    theta, omega = y\n    dydt = [omega, -b*omega - c*np.sin(theta)]\n    return dydt\n",
            "outputs": []
        },
        {
            "id": "144f0e07",
            "cell_type": "markdown",
            "source": "We assume the constants are `b` = 0.25 and `c` = 5.0:\n",
            "metadata": {}
        },
        {
            "id": "fc316cb9",
            "cell_type": "code",
            "metadata": {},
            "execution_count": null,
            "source": "b = 0.25\nc = 5.0",
            "outputs": []
        },
        {
            "id": "291a3c62",
            "cell_type": "markdown",
            "source": "For initial conditions, we assume the pendulum is nearly vertical\nwith `theta(0)` = `pi` - 0.1, and is initially at rest, so\n`omega(0)` = 0.  Then the vector of initial conditions is\n",
            "metadata": {}
        },
        {
            "id": "e290cc92",
            "cell_type": "code",
            "metadata": {},
            "execution_count": null,
            "source": "y0 = [np.pi - 0.1, 0.0]",
            "outputs": []
        },
        {
            "id": "1909b814",
            "cell_type": "markdown",
            "source": "We will generate a solution at 101 evenly spaced samples in the interval\n0 <= `t` <= 10.  So our array of times is:\n",
            "metadata": {}
        },
        {
            "id": "d67a4c89",
            "cell_type": "code",
            "metadata": {},
            "execution_count": null,
            "source": "t = np.linspace(0, 10, 101)",
            "outputs": []
        },
        {
            "id": "4eaf87f6",
            "cell_type": "markdown",
            "source": "Call `odeint` to generate the solution. To pass the parameters\n`b` and `c` to `pend`, we give them to `odeint` using the `args`\nargument.\n",
            "metadata": {}
        },
        {
            "id": "9d2f3416",
            "cell_type": "code",
            "metadata": {},
            "execution_count": null,
            "source": "from scipy.integrate import odeint\nsol = odeint(pend, y0, t, args=(b, c))",
            "outputs": []
        },
        {
            "id": "6b840677",
            "cell_type": "markdown",
            "source": "The solution is an array with shape (101, 2). The first column\nis `theta(t)`, and the second is `omega(t)`. The following code\nplots both components.\n",
            "metadata": {}
        },
        {
            "id": "f59a8f0c",
            "cell_type": "code",
            "metadata": {},
            "execution_count": null,
            "source": "import matplotlib.pyplot as plt\nplt.plot(t, sol[:, 0], 'b', label='theta(t)')\nplt.plot(t, sol[:, 1], 'g', label='omega(t)')\nplt.legend(loc='best')\nplt.xlabel('t')\nplt.grid()\nplt.show()",
            "outputs": []
        }
    ]
}