{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "408ce513",
   "metadata": {},
   "source": [
    "# Simple pendulum"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d3742b46-4699-42c8-9835-d398d09eae16",
   "metadata": {},
   "source": [
    "<img src=\"https://upload.wikimedia.org/wikipedia/commons/b/b2/Simple_gravity_pendulum.svg\" width=300>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8e129fa5",
   "metadata": {},
   "source": [
    "Let $L$ be the length of the rod, $g$ be the acceleration of gravity, and $m$ be the mass of the bob. The equation for the rotational motion of the bob is as follows:\n",
    "$$m L^2 \\, \\frac{\\mathrm{d}^2\\theta}{\\mathrm{d}t^2} = - m g \\, L \\sin\\theta,$$\n",
    "or\n",
    "$$\\frac{\\mathrm{d}^2\\theta}{\\mathrm{d}t^2} + \\omega^2 \\sin\\theta = 0,$$\n",
    "where we introduced the notation\n",
    "$$\\omega = \\sqrt{\\frac{g}{L}}.$$\n",
    "Let's introduce the dimensionless 'time' $$\\tau = \\omega t.$$\n",
    "Then the pendulum motion is described by the following universal second order nonlinear ordinary differential equation:\n",
    "\\begin{equation}\n",
    "\\frac{\\mathrm{d}^2\\theta}{\\mathrm{d}\\tau^2} + \\sin\\theta = 0.\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4afd1451-9fe7-4ea9-9cfe-a60c1ee75145",
   "metadata": {},
   "source": [
    "To solve Eq. (1) numerically, let's rewrite it as a system of two first order differential equations. Let\n",
    "$$\\nu = \\frac{\\mathrm{d}\\theta}{\\mathrm{d}\\tau}.$$\n",
    "Then, Eq. (1) can be written as follows:\n",
    "\\begin{eqnarray}\n",
    "    \\frac{\\mathrm{d}\\theta}{\\mathrm{d}\\tau} & = & \\nu ,\\\\\n",
    "    \\frac{\\mathrm{d}\\nu}{\\mathrm{d}\\tau} & = & - \\sin\\theta .\n",
    "\\end{eqnarray}"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3e8d920f-102d-456a-bf98-6116fb637815",
   "metadata": {},
   "source": [
    "If we introduce the vector of unknowns $\\vec{u} = (\\theta, \\nu)$ then the system of equations (2) can be rewritten as follows:\n",
    "\\begin{eqnarray}\n",
    "    \\frac{\\mathrm{d}u[1]}{\\mathrm{d}\\tau} & = & u[2] ,\\\\\n",
    "    \\frac{\\mathrm{d}u[2]}{\\mathrm{d}\\tau} & = & - \\sin(u[1]) .\n",
    "\\end{eqnarray}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7efeec1b",
   "metadata": {},
   "outputs": [],
   "source": [
    "using OrdinaryDiffEqTsit5"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3b80ff26",
   "metadata": {},
   "outputs": [],
   "source": [
    "using PyPlot"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "849e406b",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define the equation(s) of motion - RHS of Eq. (3)\n",
    "function pendulum!(dudt, u, p, tau)\n",
    "    dudt[1] = u[2]\n",
    "    dudt[2] = -sin(u[1])\n",
    "end"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "778ad18d-a59a-4366-be11-f2a616b78251",
   "metadata": {},
   "source": [
    "### Small oscillations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e808818c",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Initial conditions - small oscillations\n",
    "θ₀ = 0.1\n",
    "ν₀ = 0.0\n",
    "u₀ = [θ₀, ν₀];"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "01b418c2",
   "metadata": {},
   "outputs": [],
   "source": [
    "tspan = (0.0, 6π)\n",
    "p = 0.0 # dummy parameter\n",
    "prob = ODEProblem(pendulum!, u₀, tspan, p)\n",
    "sol = solve(prob, Tsit5(), reltol=1e-10, abstol=1e-10)\n",
    "sol.retcode"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f920f712",
   "metadata": {},
   "outputs": [],
   "source": [
    "np = 100\n",
    "t = range(tspan[1], tspan[2], np)\n",
    "sl = sol(t);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a63df54c",
   "metadata": {},
   "outputs": [],
   "source": [
    "ν = sl[2, :]\n",
    "θ = sl[1, :];"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5de32c6d",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot(t, θ, label=L\"$\\theta(t)$\")\n",
    "plot(t, ν, label=L\"$\\dot{\\theta}(t)$\")\n",
    "grid(true)\n",
    "xlabel(L\"$ω t$\")\n",
    "title(\"Simple pendulum: small-amplitude oscillations\")\n",
    "legend();"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "57329b0e-b066-48e7-a6be-06d916f41c40",
   "metadata": {},
   "source": [
    "### Large amplitude oscillations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f48aa068-4a34-457f-8f71-c1e614a2e34d",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Initial conditions - large oscillations\n",
    "θ₀ = π - 0.01\n",
    "ν₀ = 0.0\n",
    "u₀ = [θ₀, ν₀];"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8f5c56f7-03d0-4936-8ca4-406233929088",
   "metadata": {},
   "outputs": [],
   "source": [
    "tspan = (0.0, 20π)\n",
    "p = 0.0\n",
    "prob = ODEProblem(pendulum!, u₀, tspan, p)\n",
    "sol = solve(prob, Tsit5(), reltol=1e-10, abstol=1e-10)\n",
    "sol.retcode"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fd713f7a-f083-455f-a985-f7bdaecd5b21",
   "metadata": {},
   "outputs": [],
   "source": [
    "np = 1000\n",
    "t = range(tspan[1], tspan[2], np)\n",
    "sl = sol(t);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6c7ac7b3-a6b4-459f-9a34-544d99034090",
   "metadata": {},
   "outputs": [],
   "source": [
    "ν = sl[2, :]\n",
    "θ = sl[1, :];"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a85f2355-56f0-4794-be95-3670907c1bd9",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot(t, θ, label=L\"$\\theta(t)$\")\n",
    "plot(t, ν, label=L\"$\\dot{\\theta}(t)$\")\n",
    "grid(true)\n",
    "xlabel(L\"$ω t$\")\n",
    "title(\"Simple pendulum: large-amplitude oscillations\")\n",
    "legend();"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5ced8eed-5e65-4aec-b5f3-9dd2faa278f6",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Julia 1.10.5",
   "language": "julia",
   "name": "julia-1.10"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "1.10.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}