lectures/06-sympy/sympy-examples.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# SymPy examples"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "sources:\n",
    "http://docs.sympy.org/latest/tutorial/\n",
    "http://nbviewer.ipython.org/github/ipython/ipython/blob/master/examples/notebooks/SymPy%20Examples.ipynb"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from __future__ import division\n",
    "import sympy as sym\n",
    "import numpy as np\n",
    "# make things look pretty in the notebook\n",
    "from sympy.interactive import printing\n",
    "printing.init_printing()\n",
    "\n",
    "import math"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "SymPy provides support for symbolic math to python, similar to what you would do with Mathematica or Maple.  The major difference is that it acts just like any other python module, so you can use the symbolic math together in your own python projects with the rest of python functionality."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## SymPy types and basic symbolic manipulation\n",
    "\n",
    "Sympy defines its own types, you can convert them to python types, but you don't always want to (and will probably lose accuracy when you do).  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#timeit\n",
    "#np.sqrt(2)\n",
    "type(math.sqrt(2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "type(sym.sqrt(8))\n",
    "print(sym.sqrt(8))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "help(sym.sqrt(8))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can do symbolic math not just on numbers, but we can tell SymPy what to treat as a symbol, using `symbols()`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sympy import symbols\n",
    "x, y, z = sym.symbols(\"x y z\")\n",
    "from sympy.abc import pi, theta, phi, e"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "expr = x + 2*y + theta\n",
    "expr\n",
    "expr1 = e**x + y\n",
    "#f = x*exp(-x)+x*(1-x)\n",
    "#print(f)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "expr - 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "expr - y"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f = x*expr\n",
    "f"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "g = sym.expand(f)\n",
    "g"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sym.factor(g)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Using SymPy interactively"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "you can tell SymPy to setup some sane defaults that make an interactive session nicer"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "IPython console for SymPy 1.2 (Python 3.7.2-64-bit) (ground types: gmpy)\n",
      "\n",
      "These commands were executed:\n",
      ">>> from __future__ import division\n",
      ">>> from sympy import *\n",
      ">>> x, y, z, t = symbols('x y z t')\n",
      ">>> k, m, n = symbols('k m n', integer=True)\n",
      ">>> f, g, h = symbols('f g h', cls=Function)\n",
      ">>> init_printing()\n",
      "\n",
      "Documentation can be found at http://docs.sympy.org/1.2/\n",
      "\n"
     ]
    }
   ],
   "source": [
    "from sympy import init_session\n",
    "init_session()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## substitution\n",
    "\n",
    "SymPy provides methods to substitute values for symbols in symbolic expressions.  Note, the follow likely does not do what you expect:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAE0AAAAVBAMAAADx+n4ZAAAAMFBMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAMnZmzRC73UTviSKZVKsmUdgqAAAACXBIWXMAAA7EAAAOxAGVKw4bAAABlUlEQVQoFW2STSgEYRjHfzPs7Nqh1oWbBjciObhw2Au5mSJnFIXSOrkhTpRSTiuxSimnyUfKQVtyW7W5cVKScvH9FTWe16rZmfHUPD3///N/3+fjHQiaGSTACVPQFCbn/ihzxctVZDBzXQmPkEjPF2Ds26P74Jj4h0eoqPMP2h59AjeQ9ggVdfihQsuwYXPv+DJlRRcVErEk3MJ9wqcrbUGrr0vur1O9tttsSa488yuYQd9wXTeLNpXYg/gV/ZBiGtYt7UskclJMe2JnabO1UnDDUBZKkizYZJiAUXgWRTSvdNGsltedWokOiIuvSBF9a0TpRuBJmLJZceTkM9TcEDkVpz0Q2ZIuArqSrORqpB+xI+Viy3Jx72ORrlTVvcR0OONBKSzs37rncF2kMyypkkG3I6+oR5mHK5ljgAubcdVcV6E/2QCHudwtepJVOTI2tJ2QbVkstjbMVqWHqz4H295FE0tBj+u+YOS5k8Vx0Y3M74jz26QfFtBhmPznyUF+jqAZ4RKYA0GV7PUfzrDDOuTZg9YOP2LMZw8xaFSwAAAAAElFTkSuQmCC\n",
      "text/latex": [
       "$$\\sin{\\left (2 \\pi z \\right )}$$"
      ],
      "text/plain": [
       "sin(2⋅π⋅z)"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "expr = sin(z*2*pi)\n",
    "z = 0\n",
    "expr"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We've now redefined `z` to be a python type"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "int"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "type(z)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "to do substitution, we use the `subs()` method"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAE4AAAAVBAMAAAAazcUaAAAAMFBMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAMnZmzRC73UTviSKZVKsmUdgqAAAACXBIWXMAAA7EAAAOxAGVKw4bAAABpUlEQVQoFW2RMUgbYRSAv7skZ8zFNh20m566tVREJFOHgJC1ZytODnZQEEHiXNCsFsRC6aCUJIVOGcqhVTcJOgkphC5FF8VShC5RW21F4Xx/To2588H9/3vf/z3+d3fgD9MP4B4Ez4Je1PKY+b5+1lLALKcTdaCy697oZR0PwQax/3WgsrDt1debKrbgFyx6/GYNlW6y2/0d5G2qzi2oJR8aS4im4BCqicaDAbTurtTXHI8/LvdYchYv1IRZ9LzruiVe9gwWK/CFYcgwAzlLuxAl3Ks87YTi/Ke+R1DSD+IZ6GfOpsAUTMAfMZqkWdaSVtGdTgg5+sUDaf1B09lTlDcOJ2I0Z5VXlsdQ761hpBToIPJZpvB5tc/QLvNIPLTUOscqr47veGF17y6mwzeOlLHgRGx173f4ecczLLmsgG5HTpGfMpytEhM7yY7NpBou7c0X24f1cvkQPcUS7LVuMy/eJm/7nmTbFsfazl8n/4kTzcAL1/2LUeE3jHS1bFXEeyNPY0w3ll5lSrcvnvvqWhnvDVDDCSBYC7LIaJAxdg9bCbKQwxU5XmtS2k/HjgAAAABJRU5ErkJggg==\n",
      "text/latex": [
       "$$\\sin{\\left (2 \\pi x \\right )}$$"
      ],
      "text/plain": [
       "sin(2⋅π⋅x)"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "expr = sin(x*2*pi)\n",
    "expr"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAABwAAAAvBAMAAAACzbekAAAAMFBMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAMom7q0QizZkQVN1mdu8ilFydAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAA9klEQVQoFbWRPQ4BURSFz/gZP0NiByxAwQ60omAJU0gUChKJkmkk0xGtRlQ6WpUVSKzAEoiEiELGvc+7403PLe4535u8ecUH6CkGPEKoho2LVYlgohXBToQwge0PCnKW6aON2F0wvUQdOAsugMsMTVdzGWgATbqc5Se26nhP+9QH8kdG6wXY8yvgqI+OR5G5AQeFvtpTFzsuOU9hz8uWuKxguxTxZ2pGYZWQ50y+1rTR9X16mmY05D0OggcnNjUVshxX2j9TuQjXz18y/dDPTT+Eph9C0w+h9kNNhv18h/0Yo/x8+eNHWPsR1H40WtqPxtDPh0M/b1meRyRadn9nAAAAAElFTkSuQmCC\n",
      "text/latex": [
       "$$\\frac{\\sqrt{2}}{2}$$"
      ],
      "text/plain": [
       "√2\n",
       "──\n",
       "2 "
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a = expr.subs(x, 0.125)\n",
    "a"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that this is not a floating point number -- it is still a SymPy object.  To make it floating point, we can use evalf()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.707106781186548 <class 'sympy.core.numbers.Float'>\n"
     ]
    }
   ],
   "source": [
    "b = a.evalf()\n",
    "print(b, type(b))\n",
    "#a = expr.subs(x,0.3)\n",
    "#print(a.evalf())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This is still a SymPy object, because SymPy can do arbitrary precision "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABBwAAAAPBAMAAAC1uMEQAAAAMFBMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAiXaZIs1UEN1mu6tEMu+iw/3TAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAK+UlEQVRYCb1ZbYicVxk9uzu7s987m1gFoWTbCBYKZnCMH021oxFsq5gFv1M0U+lvd0orFITsauyHRujSVrQgdv0o1VrtakATY5sBEQXBbO2PKkXNj2LEVLOt3WyMTdbznHNn7u2mRfSHS/vuve95nnPPeZ5n3pmdAOmn8pWLV907/H1du9hcvJwAbvVd86yg2rw46v9+Z5OMVzbxykh29l+q/58TkVXn1ebDXw75Dy3cTPHynvu2vCkChzY2Ni4gNl7hdxwAY33Tl2KO8Inqlq0dhgaii/cBT6599E7exJMp+/SH7mqj0ti7wnsTdSNmE2Ui6iHv3PoIY6anV1LORJBNb2+mxJsvuzpWEd7YxcOFKGeg8RNg987vAwWlwgVLhnIq01d1bMLnKCbnyB6PeDJ5tqDsTOGZV2ziTRQhXbAFhkqXJDvL9pydpTsme5bqQr8QGVcJZNfOeogoU+OoJa0EhwmITQKVrQqWLQ39co8P47Yml303NBp3aKNV9dA+1l1Y5Us4jEsONd7aGgB+DCGGtRc8fPdnGQ4cTDx3nX0N8F70vch7w8tGzCZKJWakWscDLVzbGVtwzmd2rgODrernU2IdYy0T3TYz8YwR5VQ/hz82K8vY1kFBGeGmlAzl3IKRMzbhc0SZc2SPRxy0Z2cXzhSeecUmXlNIumJ0UxeVpHCW7Sk7S3dM4TlUF/qFiFIVs105y4go3cI4xyvBMmG2EOhs1Ta31CWR+5EF9C2xEPRVqWujFbCjDWODdTyMOl+EuB/4NEOJ+KK94DHe4U//QfNAE8AX7t9584rYEBGbKRNRF+k/g/lZfACjS3BOP8fhq8BrU+IiqrMmehz4lhHljJ3AaH2Mr8QTXW2kHIlwU0qGcna0qms2oXMUU+TUw15W6ezsTOEFr9jEaych3ZS6qYtKUjjL9pSdpTsG2bPqpZJIvxBRqpa2uwkRpRsX56QWBqVM6AgJdLZqK89Zv9xjdBlDzMIQMOaNVuq3sbk24RbYYrbnJeOgveA0Dr+/wTweh9NtHGti4sshnYjYTGmiHjL+LKZmpQLOCR/fA4KaiRzXGIcIXwO2GVHO1Cz6VifXMfZC0hYxChdsGcr5TYe5MqFzTJllCMkqnZ2dKbzgFZt4RVFQ6qYuKkl2VthTdpbuw9R0e456FfqFiFIlsN2IYQNZJyGiVON0TmphwDIhNul3ts5p9VpqbVFqTNUxfkbcuLS34btNPAOMfcJo3wqOP4S3GfHTQXvB8g7USAjyuA+HgWMr6BsM6UTEZkoT9RDic4uDJ/iLb0uREz7+ARxvKXH8XHtgkW9ng8vVczx4JSHMmV9C/3m+xYydT9qCMoUTlgzncN7/SXqa0DmOKWUQySq5ZHZ2pvCCN9jMK4pMqZtPSWUqSddZaS+0FNJ1WOHZre7qV9OTibBnu45hnV5ij6X3OWqCKFP3yCb9zk4V7LY05UT35jkOUSi+Vyx1N1yp38bWntje4f56Pj42tre44qD4or3gsXu2s5z9zRiHyD5y+du5Au4GXhfjEIjYTGmiHsIXwxcw9fEt+sAYOWEzJC8qEceefzO5IjxeCr80ohw+HfTxZDVpE6XDSWkZkdMBfrVIDprwOYopZASCnkoJQuFM4ZlXbOIVRUGpm0ZUkp6z0l5oiadDkh4xpWcXr6dfXbUJ2eMniNVU4HjZGLE9lj6doxYKdveCzfqZnWrba6lz1L09NUw8TzVsdgdpw5X6rW11bRF/5v5n/H/bWtNImgnuDY+2K/R2HWIcIvsh7OE1mCuLMQ6BiC2dEEQZwcSVj2H+Y5hsOScKUHkO+NOMKfs2PsgpC6LHgUfvMxI58QYcs8xAaTOlwgO2jMhZxG+Pcl7DhM9RTJEjez2VElQ6U3jmFZt5w0lBqZu6uCQ9Z6W90FJIl9TCs4qX9WscRGl7susCR78tw/ZY+nSOWkg4dU9s9sxipdr2WuocdS/6w/ryZw//84ar3jg8V90A3sG/Glf5hHzVozHQvadD7BMMfIfPWo2DsjG6wMjhOvoQ4xCIyH2CiDLCyLnW/CoGTjqnV5qaKW+8fqNjov2diWM3xTjUlFN9Fk/HO92BpM2UCg9KIpQROVwOMTBM+BzFFDICySq5nmuVzhRe8pJNvKIoKHUznRgl6Tkr7YWWQrqldl8C4UzF6+rXOIjS9mTXMTEO6bCwFy1M56gJMQ7uno5I+g/4pVbLLXWOujdfw7ifDpeQyxuu1HRvzwIPsBuzwL0YutAuxkF7w8A1zZGmxkHZGHiBHA3gvhgHIWIzpRIzwsjJ9akFjMSwMScKUGXTj88ocaSOuZMmqr7x1TueEqIc3HLZJ8/zr6J60mZKhQclYygjcla4PN2WCZ1jyixDSFZp8sJZUDon8ZJNvKIoKHWze+I1TREJ1oXbsCctWbope54Jq3hd/RoHUcqe7TomxqF7GO2BpU/nqAkBd01Mrif9LJZrm1tq/eoeP9tpsIAv6nNlbLhS043dznFYxPASwBfvnlYxDtoLfpDjsPIuaByYPTqLAbaJAiozMQ5CxGbKSLwpIxMrGD0zPIsRPviZo3GIj5JzLSVO8tl0xkTErmkLUQ63A7R8L39nSoULTjIi5w98x2jJhM5RTJEjJKs0+e1h/MFwpvDMC7MFr0pQUEqgLkrMzgp73ewkPTmJ3smZVJf6AwnKW1Ux2U3OukiyF41L56iFAcuE2OxZxdI5uaXKcfdGlzCi46r/Yhu1iZWa7i0feXw6zC+jusoPBot5HLwX/Bd+UMGvjx698ENE9uQCxhh8IyrvO3r02KmaELHposSvZWR+nePAB2Q8HZjT1OvhYR7bVuI81RwZVzhXe/ktCBHlcDu4hKEa3l1QKlxwksGc6kY7xoEm4kE88qJiihwhWaXJszNT8rQjQhJb8KokmZIhFKiLSqJwwdnee6ylJ92HZc9SrZtJf2r63lSCsLsJSYLUQpXRLYxEmRCbTKhYrmBuqbS5eyMnMDAbBsY5Mt7ESk339jg/O/DBXQeeYTGbeRy8F1wD/hpJp8zTx4lYwkQN/W2ulo2IzZSJqIvw9+T6OD87nEg54ePrwE4nTjaBb3PJsE+1K+eNKGfoJPas4Ang9V1tcTvCBVuGciiMT9MwMR7nmLKQEQh/TrkEykZ2Zsoeb9g83RavKApK3dSlFiURkWBdbE/ZWbpi/ESUZ6nWTevXoIjStQy7mxHbi8b5HLcwyigTYrNnFcu17bXUOeoe7sD+FX5jhyG+JWvjlT4wBobRWpV/WczR3i/auIdBO9rpor3g3ejnYxPgX9zBU6njDW1c22gc5r2ppYSITZdE1EVGFuOP7G9gP7+pVk74GG5V/+bE8Y9gZIZLhl/RHuX7ViDK6b8T30XlQGPXQlcbYxQu2DKUc7++dwgTOseUhQwh0h8CLSg7U3jm5deZ/OZAvKYI6abUTV1UksJZtqfsLN0xarqcSXWhX4goVQLZtbOM2J5aqHO0EiwTYpNAF8u1zS1Vjtzj5ssf02fVSvw1qY1WPz19cMZbbN21Ary/BUwcanQAIbp4H3ClcShmZOfGN6Hs3bveAuzb2OCYDew7N2NE5Lo4MSM/bzzCr4gaFKKcsR+c/RH/TepKHivKpxtXcxXhfdNcGVHO9NYOBvmva7NJmygVLlgylDPeuKpjEz5HMYWMsFeoVLaM25nCM6/YxGsKSResm7o4MTvL9qwlS1dM4Vmqs34holQJbFcxGTGlSq9ztBLs7oktBDrbFcwtVU6U+t/RluFrbcdLPgAAAABJRU5ErkJggg==\n",
      "text/latex": [
       "$$0.7071067811865475244008443621048490392848359376884740365883398689953662392310535194251937671638207864$$"
      ],
      "text/plain": [
       "0.7071067811865475244008443621048490392848359376884740365883398689953662392310\n",
       "535194251937671638207864"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a.evalf(100)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "want regular python types?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "c = float(b)\n",
    "print(c, type(c))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Python and SymPy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "x, y, z, t = symbols('x y z t')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "SymPy symbols are just objects and when you do operations on two sympy objects the result is a sympy object.  \n",
    "\n",
    "When you combine a sympy and python object, the result is also a sympy object.  \n",
    "\n",
    "But we need to be careful when doing fractions.  For instance doing `x + 1/3` will first compute `1/3` in python (giving `0` or `0.333...` depending on the division operator) and then add it to the sympy `x` symbol.  The `Rational()` function makes this all happen in sympy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAHcAAAAqBAMAAAB7B9mQAAAAMFBMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAMnZmzRC73UTviSKZVKsmUdgqAAAACXBIWXMAAA7EAAAOxAGVKw4bAAACZElEQVRIDe1VTWgTURD+dpNsN9nUxoPtTWO9KcZQSk+CC0KvXX/oSageKoogEQ85CDXgqWKpIIIRaSJ4KiKLv+BBFnsQYYXgRfSSqGhBhNRq1FLLOi/Zn+zLxsb02oGdnfne983OvN3HAl3agNalEIh8HeteDOzbFP/f1m9umLNfihN49wCozRe21xM5kZx0Ivf+uTBlg8o1F0RvEYo5mvAAFgUUdAnyHzfEUeAZYr89gEXhf52CprUF4BOQ94tDhp1f8uN8dhUoaKjqfvy6nb72w1wmq8AiUE348YN2yomFXYPqwzkM3LqfShIjXqzTLkIsWJZl4HDq0HwJuFdHAU48DmQwBcwlhVVihNOMJixjfub20FbAED/EM8BwsHhaQxFngdPAd2L00GPIG0JJ1HcCIV1c3UL13gSLe37uAROfBJaJEc0xmkmXxHZdgKQyYAdzL7LZfDZLXQI0kmUtIXKHJuPE9TezncYh60syP80cGTfzIxz51iQOs7bfQdHxCkuMPqtHtLZtU62PTWIpSb0WIWqRGugzG89VEaMSI6wOGffktxrOsIFHGzPHKsAT01yEqOImUN72EjMkek4XM058eWh3rj8/2b9yYuRXBZAzwJhl/YBUwhfg2GDvQolEF5iSjBM3QM+f80IvUqhk3Vzx07unvGU32u9GTUE8bSfOwYgYOMA64kzSOYClj3lMXMOEyoP0/zveimGSx+Qb6FN5kPIHrVgoqJuyM0orf11EPr8upS1BSL1vu9bBQjlgtzuQNSjRWsdUjigkIK1xWMfpRG0D4mgF3bcdSmMD7/mKSQfjL2a9nFyk2EiNAAAAAElFTkSuQmCC\n",
      "text/latex": [
       "$$\\sin{\\left (2 \\pi x \\right )} + \\frac{1}{3}$$"
      ],
      "text/plain": [
       "sin(2⋅π⋅x) + 1/3"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f = expr + Rational(1,3)\n",
    "f"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAARQAAAAVBAMAAAB1dQxbAAAAMFBMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAMnZmzRC73UTviSKZVKsmUdgqAAAACXBIWXMAAA7EAAAOxAGVKw4bAAADO0lEQVRIDb1VS2gTURQ9STpJ2jQ2CtadjXXnp5YqFUFoQOjW+KErobpQFMFWXGQhaMCFVCxVRLAibRRxEUSCX3AhwYoiRBgEEd20KlJwYbW2/mqI994ZvS8fYUDrgznv3XPe/cx98wEqR6SSAGpQ1ZvmgVldHTMcr+bmiYmc0cDRDCKF7pgSvKoqL7sqx3wk+zBtAO5e3QsDVDZI2aOKrByZI4Z/MDpjO3APDV9/mc5cl6ywbeswMzsQKBpg5dFlK3hSJIT4OCmMTGPAW2C4PHUg79rHnHk9cIVXa2xrxgB/Eb0JBU+KhBAfN4VOp4DRJKak/8qedZdPnfkAMM7lD6QxY0D4HJoSCp4UCSE+ms1ZhRPAJDAVKxc2uaZbCh3okC1U9BNNChhvJ1vBiyJ72Me3vDVxcwRLLlxvi1OUxgwBcBT+0VKplMfWti1ZSnpNWMApxeJSOCdwnCeF8CEyFTwp4s0+PUA/jgAjcd8chamTDL5pZAcvdiwE8v7Xjf3AOpJ4OKVEvgMtcbYXvY+Z4Gt7BSgYspKyMhQJIeRAEhnQ2e8DuNchmzOE8j7bn1sGBHL+uQVU3XNmaRilZISIFmlSwDj7K3hRZA/7hD6vBJeyB5imKPVpAhToCvKb5EMwwUQLw6NUajiVog5a1JWhOFPyjJtQP0ucAr8CrqykrAxF9hBpXaYnoqIUeXeX0sHRaIozDjDQcLoCelbkzukE+2wFXwzBogK8KLJHfHAL2z4apdRxg18iksMTfODkQzmLbqz8gLiFXcRapST6nv0Gu3eWSlFQ+c+KeDsy3+cbo5RgnE4lA3/SmgV9cnvSU2iggjq5KhpuVzYAD9gc4eYr1E/Q2SgYspKyMhTxdsgXSeznu+x2npWGCeBOoTAJfwLn6RgWP8Yg5bxPFw+3lJBtnaYPK9byJ0Uh0E6fBwVPiniLD050rEg3D+9u/rar88sEfRP6gc2l0gyCNt7RT6Y1OmZTCfLH0VKsS20xes0QLtDvUAEnC/Q7VPCiyB7Hh+KXjYNllmtEqEAZbldca56njbXiN7a7rPs7rLXn33PBXI2Yt2tw809ZO2vk2F2D+w/UjeocgVqdqt72l8xPI8fAb3J5t50AAAAASUVORK5CYII=\n",
      "text/latex": [
       "$$\\sin{\\left (2 \\pi x \\right )} + 0.333333333333333$$"
      ],
      "text/plain": [
       "sin(2⋅π⋅x) + 0.333333333333333"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "expr + 1/3"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## equality"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`=` is still the assignment operator of python (it does not mean symbolic equality), and `==` is still the logical test (exact structural equality).  There is a separate object, `Eq()` to specify symbolic equality.\n",
    "\n",
    "And testing for _algebraic_ equality is not always accomplished using `==`, since that tests for _structural equality_."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "False"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x + 1 == 4"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAFUAAAAQBAMAAACVTi62AAAAMFBMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEHarMolUmd1mIrvNRO9/G2jnAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAA0klEQVQoFWNgwAJYGjEFWQowxYAiwjUfMcU5DFDEeCdAuUxY1GqgqmVTwK2WZQ7xall5QGqFlEVCE8AG4jNXC6y2gMmIvQGr2ooOEGgGyXFNAKnlXMD0i3kCVrVgQQjBygBSy8LAfQDEZ+/o6HvR0QFyFQNmOMwCq2Vg4FcAyQMBbvdyKUDV3l/AKIBVLcK9vMbG9o8DGKQv2DOwgpXiMRcozwd0Xb5jHcM1ItTyb2BgEFdiq0xAVcu874cRRARBsp//rYDgIfkNWRAHG552cMgDAEuhMSwbiTTOAAAAAElFTkSuQmCC\n",
      "text/latex": [
       "$$x + 1 = 4$$"
      ],
      "text/plain": [
       "x + 1 = 4"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Eq(x + 1, 4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = (x + 1)**2\n",
    "b = x**2 + 2*x + 1    # these are algebraically equal"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "False"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a == b"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can use `simplify()` to test for algebraic equality"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "simplify(a - b)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = cos(x) + I*sin(x)\n",
    "a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "simplify(a)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## More substitution"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "note that substitution returns a new expression: SymPy expressions are immutable"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "expr = cos(x)\n",
    "expr.subs(x, 0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "expr"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "x"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "multiple substitutions, pass a list of tuples"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "expr = x**3 + 4*x*y - z\n",
    "expr"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "expr.subs([(x, 2), (y, 4), (z, 0)])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## simplifying"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "There is not unique definition of what the simplest form of an expression is.\n",
    "\n",
    "`simplify()` tries lots of methods for simplification"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAAgAAAAOBAMAAADgeEClAAAAIVBMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAADdcGRXAAAACnRSTlMAMt0Qq5nNdrvvxbMB0AAAAAlwSFlzAAAOxAAADsQBlSsOGwAAABtJREFUCB1jYGBUZmAwCVvMwMDAThbB0rnIHQBAqQoWQ8kakgAAAABJRU5ErkJggg==\n",
      "text/latex": [
       "$$1$$"
      ],
      "text/plain": [
       "1"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "simplify(sin(x)**2 + cos(x)**2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAC4AAAAOBAMAAACiOzMyAAAAMFBMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEHarMolUmd1mIrvNRO9/G2jnAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAAb0lEQVQYGWNgwAJYGrEIMjAI13zEKs7ARI64kLJIaAKqeWBzCpiM2BswxTkXMP1inoApzsLAfQAkyniiAwQcQGyIe/gVQGxkABG/v4BRAFkUol76gj0DK6ow2Jx8xzqGa1jExZXYKhNQxZn3/TACACosHU2+4BzhAAAAAElFTkSuQmCC\n",
      "text/latex": [
       "$$x - 1$$"
      ],
      "text/plain": [
       "x - 1"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "simplify( (x**3 + x**2 - x - 1)/(x**2 + 2*x + 1) )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "simplify(gamma(x)/gamma(x - 2))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "but sometimes it doesn't have your idea of what the simplest form is"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "simplify(x**2 + 2*x + 1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "instead factor may be what you want"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "factor(x**2 + 2*x + 1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### polynomial simplification"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "expand((x + 1)**2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "expand((x + 2)*(x - 3))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "expand( (x + 1)*(x - 2) - (x - 1)*x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "factor(x**2*z + 4*x*y*z + 4*y**2*z)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "factor_list(x**2*z + 4*x*y*z + 4*y**2*z)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "collect collects common powers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "expr = x*y + x - 3 + 2*x**2 - z*x**2 + x**3\n",
    "expr"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "collected_expr = collect(expr, x)\n",
    "collected_expr"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "cancel cancels"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = (x**2 + 2*x + 1)/(x**2 + x)\n",
    "a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cancel(a)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "trigsimp simplifies trigonometric identities"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "trigsimp(sin(x)**4 - 2*cos(x)**2*sin(x)**2 + cos(x)**4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "trigsimp(sin(x)*tan(x)/sec(x))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "the tutorial discusses some of the nuances of simplification of powers and special functions"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Calculus"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Calculus operations are simple in SymPy\n",
    "\n",
    "### derivatives"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "diff(cos(x), x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "diff(exp(x**2), x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "third derivative"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "diff(x**4, x, 3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "differentiate different variables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "expr = exp(x*y*z)\n",
    "diff(expr, x, y, z)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "unevaluated derivaties can be useful for building up ODEs and PDEs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "deriv = Derivative(expr, x, y, z)\n",
    "deriv"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "deriv.doit()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### integrals\n",
    "\n",
    "definite and indefinite integrals are supported"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "integrate(cos(x), x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "definite integral -- note the construction of the infinity"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "integrate(exp(-x), (x, 0, oo))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "double integral"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "integrate(exp(-x**2 - y**2), (x, -oo, oo), (y, -oo, oo))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "if it is unable to do the integral, it returns an Integral object"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "expr = integrate(x**x, x)\n",
    "print(expr)\n",
    "expr"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = x / sqrt(x**4 + 10*x**2 - 96*x - 71)   # example from Wikipedia Risch algorithm page)\n",
    "a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "integrate(a, x)     # this has a known solution, but SymPy fails to find it"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### limits"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "limit(sin(x)/x, x, 0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### series expansions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "expr = exp(sin(x))\n",
    "a = expr.series(x, 0, 10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(a)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "c = log(x).series(x, x0=1, n=6)\n",
    "c"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "simplify(c.removeO())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## solvers\n",
    "\n",
    "`solveset()` is the main interface to solvers in SymPy.  Note that it used to be `solve()`, but this has been replaced (see http://docs.sympy.org/latest/modules/solvers/solveset.html)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If no Eq() is done, then it is assumed to be equal to 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "solveset(x**2 - x, x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "you can restrict the domain of the solution (e.g. to reals).  Recall that Z is the set of integers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "solveset(sin(x) - 1, x, domain=S.Reals)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### linear systems\n",
    "\n",
    "`linsolve()` is the interface to linear systems"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "linsolve([x - y + 2, x + y - 3], [x, y])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "linsolve([x + y + z - 1, x + y + 2*z - 3 ], (x, y, z))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "roots will report if a solution is multiple by listing it multiple times"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "roots(x**3 - 6*x**2 + 9*x, x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "0 is 1 root, and 3 is 2 more roots"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Differential equations"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "you need an undefined function (f and g already are by our init_session() above, but we've probably reset these"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f, g = symbols('f g', cls=Function)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f(x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f(x).diff(x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "diffeq = Eq(f(x).diff(x, 2) - 2*f(x).diff(x) + f(x), sin(x))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "diffeq"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dsolve(diffeq, f(x))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Matrices"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "consider the Euler equations:\n",
    "\n",
    "$$q_t + A(q) q_x = 0$$\n",
    "\n",
    "where\n",
    "\n",
    "$$q = \\left ( \\begin{array}{c} \\rho \\\\ u \\\\ p \\end{array} \\right )\n",
    "\\qquad\n",
    "A(q) = \\left ( \\begin{array}{ccc} u  & \\rho     & 0 \\\\                          \n",
    "                                  0  &  u       & 1/\\rho \\\\                     \n",
    "                                  0  & c^2 \\rho & u \\end{array} \\right ) $$\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sympy.abc import rho\n",
    "rho, u, c = symbols('rho u c')\n",
    "A = Matrix([[u, rho, 0], [0, u, rho**-1], [0, c**2 * rho, u]])\n",
    "A"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A.row(0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The eigenvalues of the system are the speeds at which information propagates"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A.eigenvals()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can diagonalize it, such that\n",
    "$$ A = PDP^{-1}$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "P, D = A.diagonalize()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$D$ will be a matrix of the eigenvalues"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "D"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$P$ will be the matrix of right eigenvectors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "P"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Inverse"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A**-1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.7.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}