{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Simulating a Yo-Yo"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*Modeling and Simulation in Python*\n",
    "\n",
    "Copyright 2021 Allen Downey\n",
    "\n",
    "License: [Creative Commons Attribution-NonCommercial-ShareAlike 4.0 International](https://creativecommons.org/licenses/by-nc-sa/4.0/)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# install Pint if necessary\n",
    "\n",
    "try:\n",
    "    import pint\n",
    "except ImportError:\n",
    "    !pip install pint"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# download modsim.py if necessary\n",
    "\n",
    "from os.path import basename, exists\n",
    "\n",
    "def download(url):\n",
    "    filename = basename(url)\n",
    "    if not exists(filename):\n",
    "        from urllib.request import urlretrieve\n",
    "        local, _ = urlretrieve(url, filename)\n",
    "        print('Downloaded ' + local)\n",
    "    \n",
    "download('https://github.com/AllenDowney/ModSimPy/raw/master/modsim.py')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# import functions from modsim\n",
    "\n",
    "from modsim import *"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Yo-yo\n",
    "\n",
    "Suppose you are holding a yo-yo with a length of string wound around its axle, and you drop it while holding the end of the string stationary. As gravity accelerates the yo-yo downward, tension in the string exerts a force upward. Since this force acts on a point offset from the center of mass, it exerts a torque that causes the yo-yo to spin.\n",
    "\n",
    "The following diagram shows the forces on the yo-yo and the resulting torque. The outer shaded area shows the body of the yo-yo. The inner shaded area shows the rolled up string, the radius of which changes as the yo-yo unrolls.\n",
    "\n",
    "![Diagram of a yo-yo showing forces due to gravity and tension in the\n",
    "string, the lever arm of tension, and the resulting\n",
    "torque.](https://github.com/AllenDowney/ModSim/raw/main/figs/yoyo.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this system, we can't figure out the linear and angular acceleration independently; we have to solve a system of equations:\n",
    "\n",
    "$\\sum F = m a $\n",
    "\n",
    "$\\sum \\tau = I \\alpha$\n",
    "\n",
    "where the summations indicate that we are adding up forces and torques.\n",
    "\n",
    "As in the previous examples, linear and angular velocity are related because of the way the string unrolls:\n",
    "\n",
    "$\\frac{dy}{dt} = -r \\frac{d \\theta}{dt} $\n",
    "\n",
    "In this example, the linear and angular accelerations have opposite sign.  As the yo-yo rotates counter-clockwise, $\\theta$ increases and $y$, which is the length of the rolled part of the string, decreases."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Taking the derivative of both sides yields a similar relationship between linear and angular acceleration:\n",
    "\n",
    "$\\frac{d^2 y}{dt^2} = -r \\frac{d^2 \\theta}{dt^2} $\n",
    "\n",
    "Which we can write more concisely:\n",
    "\n",
    "$ a = -r \\alpha $\n",
    "\n",
    "This relationship is not a general law of nature; it is specific to scenarios like this where there is rolling without stretching or slipping.\n",
    "\n",
    "Because of the way we've set up the problem, $y$ actually has two meanings: it represents the length of the rolled string and the height of the yo-yo, which decreases as the yo-yo falls.  Similarly, $a$ represents acceleration in the length of the rolled string and the height of the yo-yo.\n",
    "\n",
    "We can compute the acceleration of the yo-yo by adding up the linear forces:\n",
    "\n",
    "$\\sum F = T - mg = ma $\n",
    "\n",
    "Where $T$ is positive because the tension force points up, and $mg$ is negative because gravity points down."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Because gravity acts on the center of mass, it creates no torque, so the only torque is due to tension:\n",
    "\n",
    "$\\sum \\tau = T r = I \\alpha $\n",
    "\n",
    "Positive (upward) tension yields positive (counter-clockwise) angular acceleration.\n",
    "\n",
    "Now we have three equations in three unknowns, $T$, $a$, and $\\alpha$, with $I$, $m$, $g$, and $r$ as known parameters.  We could solve these equations by hand, but we can also get SymPy to do it for us."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sympy import symbols, Eq, solve\n",
    "\n",
    "T, a, alpha, I, m, g, r = symbols('T a alpha I m g r')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "eq1 = Eq(a, -r * alpha)\n",
    "eq1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "eq2 = Eq(T - m * g, m * a)\n",
    "eq2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "eq3 = Eq(T * r, I * alpha)\n",
    "eq3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "soln = solve([eq1, eq2, eq3], [T, a, alpha])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "soln[T]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "soln[a]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "soln[alpha]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The results are\n",
    "\n",
    "$T      = m g I / I^*   $\n",
    "\n",
    "$a      = -m g r^2 / I^* $\n",
    "\n",
    "$\\alpha = m g r / I^*    $\n",
    "\n",
    "where $I^*$ is the augmented moment of inertia, $I + m r^2$.\n",
    "\n",
    "You can also see [the derivation of these equations in this video](https://www.youtube.com/watch?v=chC7xVDKl4Q)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can use these equations for $a$ and $\\alpha$ to write a slope function and simulate this system."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Exercise:**  Simulate the descent of a yo-yo.  How long does it take to reach the end of the string?\n",
    "\n",
    "Here are the system parameters:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "Rmin = 8e-3     # m\n",
    "Rmax = 16e-3    # m\n",
    "Rout = 35e-3    # m\n",
    "mass = 50e-3    # kg\n",
    "L = 1           # m\n",
    "g = 9.8         # m / s**2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* `Rmin` is the radius of the axle.  `Rmax` is the radius of the axle plus rolled string.\n",
    "\n",
    "* `Rout` is the radius of the yo-yo body.  `mass` is the total mass of the yo-yo, ignoring the string.  \n",
    "\n",
    "* `L` is the length of the string.\n",
    "\n",
    "* `g` is the acceleration of gravity."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "1 / (Rmax)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Based on these parameters, we can compute the moment of inertia for the yo-yo, modeling it as a solid cylinder with uniform density ([see here](https://en.wikipedia.org/wiki/List_of_moments_of_inertia)).\n",
    "\n",
    "In reality, the distribution of weight in a yo-yo is often designed to achieve desired effects.  But we'll keep it simple."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "I = mass * Rout**2 / 2\n",
    "I"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And we can compute `k`, which is the constant that determines how the radius of the spooled string decreases as it unwinds."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "k = (Rmax**2 - Rmin**2) / 2 / L   \n",
    "k"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The state variables we'll use are angle, `theta`, angular velocity, `omega`, the length of the spooled string, `y`, and the linear velocity of the yo-yo, `v`.\n",
    "\n",
    "Here is a `State` object with the initial conditions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "init = State(theta=0, omega=0, y=L, v=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And here's a `System` object with `init` and `t_end` (chosen to be longer than I expect for the yo-yo to drop 1 m)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "system = System(init=init, t_end=2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Write a slope function for this system, using these results from the book:\n",
    "\n",
    "$ r = \\sqrt{2 k y + R_{min}^2} $ \n",
    "\n",
    "$ T      = m g I / I^*  $\n",
    "\n",
    "$ a      = -m g r^2 / I^* $\n",
    "\n",
    "$ \\alpha  = m g r / I^*  $\n",
    "\n",
    "where $I^*$ is the augmented moment of inertia, $I + m r^2$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution goes here"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Test your slope function with the initial conditions.\n",
    "The results should be approximately\n",
    "\n",
    "```\n",
    "0, 180.5, 0, -2.9\n",
    "```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution goes here"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Notice that the initial acceleration is substantially smaller than `g` because the yo-yo has to start spinning before it can fall.\n",
    "\n",
    "Write an event function that will stop the simulation when `y` is 0."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution goes here"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Test your event function:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "# Solution goes here"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then run the simulation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "# Solution goes here"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Check the final state.  If things have gone according to plan, the final value of `y` should be close to 0."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solution goes here"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "How long does it take for the yo-yo to fall 1 m?  Does the answer seem reasonable?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The following cells plot the results.\n",
    "\n",
    "`theta` should increase and accelerate."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [],
   "source": [
    "results.theta.plot(color='C0', label='theta')\n",
    "decorate(xlabel='Time (s)',\n",
    "         ylabel='Angle (rad)')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`y` should decrease and accelerate down."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [],
   "source": [
    "results.y.plot(color='C1', label='y')\n",
    "\n",
    "decorate(xlabel='Time (s)',\n",
    "         ylabel='Length (m)')\n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plot velocity as a function of time; is the acceleration constant?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [],
   "source": [
    "results.v.plot(label='velocity', color='C3')\n",
    "\n",
    "decorate(xlabel='Time (s)',\n",
    "         ylabel='Velocity (m/s)')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can use `gradient` to estimate the derivative of `v`.  How does the acceleration of the yo-yo compare to `g`?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = gradient(results.v)\n",
    "a.plot(label='acceleration', color='C4')\n",
    "decorate(xlabel='Time (s)',\n",
    "         ylabel='Acceleration (m/$s^2$)')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And we can use the formula for `r` to plot the radius of the spooled thread over time."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [],
   "source": [
    "r = np.sqrt(2*k*results.y + Rmin**2)\n",
    "r.plot(label='radius')\n",
    "\n",
    "decorate(xlabel='Time (s)',\n",
    "         ylabel='Radius of spooled thread (m)')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.16"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}