{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Three faces of classical mechanics\n",
    "\n",
    "- 1687 - edition of Principia Mathematics by Isaac Newton \n",
    "- 1788 - edition of the Mécanique analytique by Joseph Louis Lagrange  (Giuseppe Lodovico Lagrangia)\n",
    "- 1833 - formulation of mechanics by William Rowand Hamilton \n",
    "\n",
    "\n",
    "There are three approaches to describe classical mechanical systems. Historically, the first theory was formulated by Newton and is called Netwon Mechanics. The second is Lagrange Mechanics and is convenient for description of systems with constraints. The third is  Hamilton Mechanics. Lagrange mechanics is a basis for the  modern theory of elementary particles. In turn, the Hamilton approach is used  in formulation of quantum physics. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Harmonic Oscillator \n",
    "\n",
    "We demostrate three approaches to classical mechanical systems by considering  one of the simplest example:  a one-dimensional harmonic oscillator (a particle of mass $m$) moving along the x-axis and characterised by the position $x(t)$ at time $t$. We present it in a trivialized way. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### The Newton mechanics \n",
    "\n",
    "\n",
    "In the Newton description we have to know all forces $F$ which act on the particle of mass $m$ and its  dynamics is determined by the Newton second law (the equation of motion): \n",
    "\n",
    "$$m a =F,$$\n",
    "where $a=a(t)$ is an acceleration of the particle which is a time-derivative of the particle velocity $v=v(t)$, which in turn in a time-derivative of the particle position (coordinate) $x=x(t)$: \n",
    "\n",
    "$$a(t) = \\dot{v}(t) = \\frac{dv(t)}{dt}, \\qquad v(t) = \\dot{x}(t) = \\frac{dx(t)}{dt}, \\qquad F = F(x, v, t).$$\n",
    "Taking into account the above relations the Newton equation $ma=F$ is in fact a second-order differential equation  \n",
    "in the form \n",
    "\n",
    "$$m \\frac{d^2 x}{dt^2} = F(x, \\dot{x}, t)$$\n",
    "Therefore two initial conditions have to be imposed. Usually, it is  the initial position $x_0=x(0)$ and the initial velocity $v_0=v(0)$ of the particle. \n",
    "\n",
    "In a general case, the force $F$ depends on the particle position $x$, its velocity $\\dot x =v$  and explicitly on time $t$ (in the case of an external time-dependent driving as e.g. $A \\cos(B t)$). When $F=F(x)$ then it is conservative system. In such a case we can define  a potential energy $U(x)$ of the system by the relations:  \n",
    "\n",
    "$$ U(x) = - \\int F(x) dx \\qquad \\mbox{or} \\qquad  F(x) = - \\frac{dU(x)}{dx}$$\n",
    "\n",
    "For the harmonic oscillator, the force $F$ is proportional to the partice displacement $x$, namely, \n",
    "\n",
    "$$F = F(x) =-kx$$\n",
    "where $k$ is a positive parameter (e.g. $k$ is a measure of the stiffness of the spring for the mass-spring oscillator). The corresponding potential energy $E_p = U(x)$ is \n",
    "\n",
    "$$ U(x) = \\frac{1}{2} k x^2 $$\n",
    "i.e. it is a parabola. The Newton equation for the harmonic oscillator is rewritten in the standadrd form as \n",
    "\n",
    "$$ \\ddot x + \\omega_0^2 x =0, \\qquad \\mbox{where} \\qquad \\omega_0^2 = \\frac{k}{m}$$\n",
    "Its analysis in presented in the next notebook. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### The Lagrange mechanics \n",
    "\n",
    "In the Lagrange approach, the mechanical system is described by the scalar function $L$ called the Lagrange function. In the case of  conservative systems it is a difference between the kinetic energy $E_k$ and the potential energy $E_p$ of the system, i.e.,  \n",
    "\n",
    "$$ L = L(x, v) = E_k - E_p = \\frac{m v^2}{2} - U(x) $$\n",
    "For the harmonic oscillator it reads \n",
    "\n",
    "$$ L = \\frac{m v^2}{2} -  \\frac{1}{2} k x^2 $$ \n",
    "Dynamics of the system is determined by the Euler-Lagrange equation: \n",
    "\n",
    "$$ \\frac{d}{dt} \\frac{\\partial L}{\\partial v} - \\frac{\\partial L}{\\partial x} = 0$$\n",
    "It is a counterpart of the Newton equation.  \n",
    "\n",
    "$$ \\frac{\\partial L}{\\partial v}  = mv = p $$ \n",
    "is the (Newton) momentum $p$ of the particle and \n",
    "\n",
    "$$\\frac{\\partial L}{\\partial x} = -U'(x) =F$$\n",
    "is the force $F$ acting on the particle.  As a result, the Euler-Lagrange equation leads to the equation of motion\n",
    "\n",
    "$$ m \\frac{dv}{dt} = -U'(x) = F$$\n",
    "and is the same as the Newton equation because $v = \\dot x$. For the harmonic oscillator the potential is \n",
    "$U(x)= kx^2/2$, the force is $f(x)=-U'(x)=-kx$ the the Newton equation is \n",
    "\n",
    "$$ m \\frac{dv}{dt} = -kx .$$\n",
    "\n",
    "It has the same form as in the Newton mechanics. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### The Hamilton mechanics\n",
    "\n",
    "In the Hamilton approach, the mechanical system is described by the scalar function $H$ called the Hamilton  function. In the case of  conservative systems,  it is a total energy of the system, i.e.,  \n",
    "\n",
    "$$ H = H(x, p) = E_k + E_p = \\frac{p^2}{2m} + U(x) $$\n",
    "where the kinetic energy $E_k$ is expressed by the CANONICAL MOMENTUM $p$ !\n",
    "\n",
    "For the harmonic oscillator it reads \n",
    "\n",
    "$$ H = \\frac{p^2}{2m} +  \\frac{1}{2} k x^2 $$ \n",
    "Dynamics of the system is determined by the Hamilton equations:  \n",
    "\n",
    "$$ \\frac{dx}{dt} = \\frac{\\partial H}{\\partial p}, \\qquad  \\frac{dp}{dt} = -\\frac{\\partial H}{\\partial x} $$\n",
    "It is a counterpart of the Newton equation or the Euler-Lagrange equation. For the harmonic oscillator we obtain \n",
    "\n",
    "$$ \\frac{\\partial H}{\\partial p}  =  \\frac{p}{m} $$ \n",
    "is the velocity $v$ of the particle and \n",
    "\n",
    "$$\\frac{\\partial H}{\\partial x} = kx = - F$$\n",
    "is proportional to the force $F$ acting on the oscillator. As a result, the Hamilton equations lead to the equation of motion in the form \n",
    "\n",
    "$$ \\frac{dx}{dt} = \\dot x =\\frac{p}{m},  \\qquad \\frac{dp}{dt} = \\dot p = -kx$$\n",
    "and are equivalent to the Newton equation or the Euler-Lagrange equation. Indeed, if we differentiate the first equation we get $\\ddot x = \\dot p /m$. Next we insert to it the second equation for $\\dot p$ and finally we obtain $\\ddot x = -k x/m = -\\omega_0^2 x$.  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### REMARK: \n",
    "\n",
    "We want to stress that in a general case, the construction of the Lagrange function or the Hamilton function can be a more complicated task. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercise\n",
    "\n",
    "Solve the equation of motion for the harmonic oscillator with the mass $m=2$ and the eigen-frequency $\\omega_0 =1$. Find  $x(t)$ and $v(t)$ for a given initial conditions: $x(0) =1$ and $v(0) =0.5$. Next,  depict the time-evolution of the kinetic energy $E_k(t)$, potential energy $E_p(t)$  and the total energy $E(t) =E_k(t) + E_p(t)$: \n",
    "\n",
    "$$E_k(t) = \\frac{1}{2} m v^2(t), \\qquad E_p(t) = \\frac{1}{2} m \\omega_0^2 x^2(t)$$\n",
    "What is the main conclusion regarding the total $E(t)$. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Solution\n",
    "\n",
    "In SageMath we can write easile solve any system of ODE numerically. First we write Newton equation in a following form:\n",
    "\n",
    "\\begin{equation}\n",
    "\\label{intro_eq:odes4num}\n",
    "\\left\\{\n",
    "\\begin{array}{l}\n",
    "\\displaystyle\\frac{dx}{dt} = v \\\\\n",
    "\\displaystyle\\frac{dv}{dt} =  - \\omega_0^2 x\n",
    "\\end{array}\n",
    "\\right.\n",
    "\\end{equation}\n",
    "Then we use ``desolve_odeint`` to obtain a numerical solution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m = 2 \n",
    "omega0 = 1\n",
    "var('x,v')\n",
    "times = srange(0,4,0.01)\n",
    "xv = desolve_odeint([ v, -omega0^2*x ], [1,0.5] , times, [x,v])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can compute $E_k$ and $E_p$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Ek = 1/2*m*xv[:,1]^2\n",
    "Ep = 1/2*m*omega0^2*xv[:,0]^2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And plot the results:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "p_Ek = line( zip( times, Ek),\\\n",
    "            legend_label=r'$E_k(t)$', figsize=(6,2))\n",
    "p_Ep = line( zip( times,Ep ),\\\n",
    "            color='red',legend_label=r'$E_p(t)$')\n",
    "p_Etot = line( zip(times,Ek + Ep),\\\n",
    "             color='green',legend_label=r'$E_{tot}(t) = E_p(t)+E_k(t)$')\n",
    "p_Ek + p_Ep + p_Etot"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Damped harmonic oscillator \n",
    "\n",
    "\n",
    "### The Newton mechanics \n",
    "\n",
    "\n",
    "A system which interacts with its environment is  dissipative (it losses its energy) due to friction. It can be  conveniently described by the Newton mechanics. For small velocity of the particle, the friction force is proportional to its velocity (the Stokes force) $F=-\\gamma v = -\\gamma \\dot x$, where $\\gamma$ is a friction  or damping coefficient. The  Newton equation for the damped harmonic oscillator has the form \n",
    "\n",
    "$$ m\\ddot x = -\\gamma \\dot x -k x.$$ \n",
    "Its analysis in presented in the next notebook. \n",
    "\n",
    "\n",
    "### The Lagrange mechanics \n",
    "\n",
    "In the Lagrange approach, in this case the Lagrange function is not a difference between the kinetic energy and the potential energy but is constructed in such an artificial  way in order to obtain the correct equation of motion presented above. Let us propose the following  function: \n",
    "\n",
    "$$ L = L(x, v, t)  = e^{\\gamma t/m} \\left[\\frac{m v^2}{2} -  \\frac{1}{2} k x^2\\right] $$ \n",
    "Then in the Euler-Lagrange equation: \n",
    "\n",
    "\n",
    "$$ \\frac{\\partial L}{\\partial v}  = mv e^{\\gamma t/m} $$ \n",
    "Its time derivative is \n",
    "\n",
    "$$\\frac{d}{dt} \\frac{\\partial L}{\\partial v} = m \\dot v e^{\\gamma t/m} + \\gamma  v e^{\\gamma t/m}$$ \n",
    "The second part of the Euler-Lagrange equation is \n",
    "\n",
    "$$\\frac{\\partial L}{\\partial x} = -kx e^{\\gamma t/m} $$\n",
    "As a result, the final form of the Euler-Lagrange equation is \n",
    "\n",
    "$$ m\\ddot x  + \\gamma \\dot x + k x = 0$$ \n",
    "and is the same as in the  Newton approach.  \n",
    "\n",
    "\n",
    "### The Hamilton mechanics\n",
    "\n",
    "In the Hamilton approach, the Hamilton function is in the form  \n",
    "$$ H = H(x, p, t) =  \\frac{p^2}{2m} e^{-\\gamma t/m} +  \\frac{1}{2} k x^2 e^{\\gamma t/m}$$ \n",
    "The partial derivatives are \n",
    "$$ \\frac{\\partial H}{\\partial p}  =  \\frac{p}{m} e^{-\\gamma t/m} $$ \n",
    "and \n",
    "$$\\frac{\\partial H}{\\partial x} = kx e^{\\gamma t/m} $$\n",
    "As a result, the Hamilton equations lead to the equation of motion in the form \n",
    "\n",
    "$$ \\frac{dx}{dt} = \\frac{p}{m} e^{-\\gamma t/m},  \\qquad \\frac{dp}{dt} =  -kx e^{\\gamma t/m} $$\n",
    "One can show that they are equivalent to the Newton equation or the Euler-Lagrange equation. \n",
    "However, there is one important remark: From the first Hamilton equation it follows that \n",
    "\n",
    "$$ m v = p e^{-\\gamma t/m}$$\n",
    "The left side is the Newton momentum of the particle. In the right side, $p$ is the canonical momentum. \n",
    "The above equations of motion are investigated in detail in the next notebook. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercise\n",
    "\n",
    "Solve the equation of motion for the damped harmonic oscillatora with the mass $m=2$, the friction coefficient $\\gamma = 1$ and the eigen-frequency $\\omega_0 =1$. Find  $x(t)$ and $v(t)$ for a given initial conditions: $x(0) =1$ and $v(0) =0.5$. Next depict the time-evolution of the kinetic energy $E_k(t)$, potential energy $E_p(t)$  and the total energy $E(t) =E_k(t) + E_p(t)$: \n",
    "\n",
    "$$E_k(t) = \\frac{1}{2} m v^2(t), \\qquad E_p(t) = \\frac{1}{2} m \\omega^2 x^2(t)$$\n",
    "\n",
    "Compare the results with the frictionless case, $\\gamma =0$. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Solution\n",
    "\n",
    "In this case the Newton equation can be rewritten in the form:\n",
    "\n",
    "\\begin{equation}\n",
    "\\label{intro_eq:ode4num}\n",
    "\\left\\{\n",
    "\\begin{array}{l}\n",
    "\\displaystyle\\frac{dx}{dt} = v \\\\\n",
    "\\displaystyle\\frac{dv}{dt} =  - \\omega_0^2 x - \\gamma_0 v, \\quad \\gamma_0=\\gamma/m\n",
    "\\end{array}\n",
    "\\right.\n",
    "\\end{equation}\n",
    "\n",
    "Then we use ``desolve_odeint`` to obtain a numerical solution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m = 2 \n",
    "omega0 = 1\n",
    "gamma_ = 1\n",
    "var('x,v')\n",
    "times = srange(0,4,0.01)\n",
    "xv = desolve_odeint([ v, -x-gamma_*v ], [1,0.5] , times, [x,v])\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can compute $E_k$ and $E_p$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Ek = 1/2*m*xv[:,1]^2\n",
    "Ep = 1/2*m*omega0^2*xv[:,0]^2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And plot the results:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "p_Ek = line( zip( times, Ek),\\\n",
    "            legend_label=r'$E_k(t)$', figsize=(6,2))\n",
    "p_Ep = line( zip( times,Ep ),\\\n",
    "            color='red',legend_label=r'$E_p(t)$')\n",
    "p_Etot = line( zip(times,Ek + Ep),\\\n",
    "             color='green',legend_label=r'$E_{tot}(t) = E_p(t)+E_k(t)$')\n",
    "p_Ek + p_Ep + p_Etot"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\newpage"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "SageMath 8.6",
   "language": "",
   "name": "sagemath"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.15"
  },
  "nbTranslate": {
   "displayLangs": [
    "en",
    "pl"
   ],
   "hotkey": "alt-t",
   "langInMainMenu": true,
   "sourceLang": "pl",
   "targetLang": "en",
   "useGoogleTranslate": true
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}