notebooks/sir/octave.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## SIR model using Octave and LSODE"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Octave requires a function that returns the gradient of the system, given the current state of the system, `x` and the time, `t`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "function xdot = sir_eqn(x,t)\n",
    "    % Parameter values\n",
    "    beta=0.1;\n",
    "    mu=0.05;\n",
    "\n",
    "    % Define variables\n",
    "    s = x(1);\n",
    "    y = x(2);\n",
    "    r = x(3);\n",
    "\n",
    "    % Define ODEs\n",
    "    ds=-beta*s*y;\n",
    "    dy=beta*s*y-mu*y;\n",
    "    dr=mu*y;\n",
    "\n",
    "    % Return gradients\n",
    "    xdot = [ds,dy,dr];\n",
    "endfunction"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We then set up the time axis over which to integrate and the initial conditions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "t = linspace(0, 200, 2001)+.1;\n",
    "x0=[0.99,0.01,0];"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The function lsode can be used to solve ODEs of this form using Hindmarsh’s ODE solver [LSODE](https://computation.llnl.gov/casc/nsde/pubs/u113855.pdf)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "x = lsode(\"sir_eqn\",x0, t);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The following saves the times and the states."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [],
   "source": [
    "out=[transpose(t),x];\n",
    "save -ascii sir_octave.out out;"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " 1.00000000e-01 9.90000000e-01 1.00000000e-02 0.00000000e+00\n",
      " 2.00000000e-01 9.89900759e-01 1.00491165e-02 5.01241223e-05\n",
      " 3.00000000e-01 9.89801042e-01 1.00984640e-02 1.00494189e-04\n",
      " 4.00000000e-01 9.89700844e-01 1.01480441e-02 1.51112049e-04\n",
      " 5.00000000e-01 9.89600163e-01 1.01978578e-02 2.01979039e-04\n",
      " 6.00000000e-01 9.89499005e-01 1.02479029e-02 2.53092305e-04\n",
      " 7.00000000e-01 9.89397349e-01 1.02981882e-02 3.04462323e-04\n",
      " 8.00000000e-01 9.89295211e-01 1.03487075e-02 3.56081707e-04\n",
      " 9.00000000e-01 9.89192584e-01 1.03994631e-02 4.07953352e-04\n",
      " 1.00000000e+00 9.89089464e-01 1.04504565e-02 4.60079036e-04\n"
     ]
    }
   ],
   "source": [
    "!head sir_octave.out"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualisation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjAAAAGkCAIAAACgjIjwAAAJNmlDQ1BkZWZhdWx0X3JnYi5pY2MA\nAHiclZFnUJSHFobP933bCwvssnRYepMqZQHpvUmvogJL7yxLEbEhYgQiiog0RZCggAGjUiRWRLEQ\nFBSxoFkkCCgxGEVUUPLDOxPn3vHHfX49884755yZA0ARBQBARQFSUgV8Pxd7TkhoGAe+IZKXmW7n\n4+MJ3+X9KCAAAPdWfb/zXSjRMZk8AFgGgHxeOl8AgOQCgGaOIF0AgBwFAFZUUroAADkLACx+SGgY\nAHIDAFhxX30cAFhRX30eAFj8AD8HABQHQKLFfeNR3/h/9gIAKNvxBQmxMbkc/7RYQU4kP4aT6edi\nz3FzcOD48NNiE5Jjvjn4/yp/B0FMrgAAwCEtfRM/IS5ewPmfoUYGhobw7y/e+gICAAh78L//AwDf\n9NIaAbgLANi+f7OoaoDuXQBSj//NVI8CMAoBuu7wsvjZXzMcAAAeKMAAFkiDAqiAJuiCEZiBJdiC\nE7iDNwRAKGwAHsRDCvAhB/JhBxRBCeyDg1AD9dAELdAOp6EbzsMVuA634S6MwhMQwhS8gnl4D0sI\nghAROsJEpBFFRA3RQYwQLmKNOCGeiB8SikQgcUgqkoXkIzuREqQcqUEakBbkF+QccgW5iQwjj5AJ\nZBb5G/mEYigNZaHyqDqqj3JRO9QDDUDXo3FoBpqHFqJ70Sq0ET2JdqFX0NvoKCpEX6ELGGBUjI0p\nYboYF3PAvLEwLBbjY1uxYqwSa8TasV5sALuHCbE57COOgGPiODhdnCXOFReI4+EycFtxpbga3Alc\nF64fdw83gZvHfcHT8XJ4HbwF3g0fgo/D5+CL8JX4Znwn/hp+FD+Ff08gENgEDYIZwZUQSkgkbCaU\nEg4TOgiXCcOEScICkUiUJuoQrYjexEiigFhErCaeJF4ijhCniB9IVJIiyYjkTAojpZIKSJWkVtJF\n0ghpmrREFiWrkS3I3uRo8iZyGbmJ3Eu+Q54iL1HEKBoUK0oAJZGyg1JFaadco4xT3lKpVGWqOdWX\nmkDdTq2inqLeoE5QP9LEado0B1o4LYu2l3acdpn2iPaWTqer023pYXQBfS+9hX6V/oz+QYQpoifi\nJhItsk2kVqRLZETkNYPMUGPYMTYw8hiVjDOMO4w5UbKouqiDaKToVtFa0XOiY6ILYkwxQzFvsRSx\nUrFWsZtiM+JEcXVxJ/Fo8ULxY+JXxSeZGFOF6cDkMXcym5jXmFMsAkuD5cZKZJWwfmYNseYlxCWM\nJYIkciVqJS5ICNkYW53txk5ml7FPsx+wP0nKS9pJxkjukWyXHJFclJKVspWKkSqW6pAalfokzZF2\nkk6S3i/dLf1UBiejLeMrkyNzROaazJwsS9ZSlidbLHta9rEcKqct5ye3We6Y3KDcgryCvIt8uny1\n/FX5OQW2gq1CokKFwkWFWUWmorVigmKF4iXFlxwJjh0nmVPF6efMK8kpuSplKTUoDSktKWsoByoX\nKHcoP1WhqHBVYlUqVPpU5lUVVb1U81XbVB+rkdW4avFqh9QG1BbVNdSD1Xerd6vPaEhpuGnkabRp\njGvSNW00MzQbNe9rEbS4Wklah7XuaqPaJtrx2rXad3RQHVOdBJ3DOsOr8KvMV6Wualw1pkvTtdPN\n1m3TndBj63nqFeh1673WV9UP09+vP6D/xcDEINmgyeCJobihu2GBYa/h30baRjyjWqP7q+mrnVdv\nW92z+o2xjnGM8RHjhyZMEy+T3SZ9Jp9NzUz5pu2ms2aqZhFmdWZjXBbXh1vKvWGON7c332Z+3vyj\nhamFwOK0xV+WupZJlq2WM2s01sSsaVozaaVsFWnVYCW05lhHWB+1Ftoo2UTaNNo8t1WxjbZttp22\n07JLtDtp99rewJ5v32m/6GDhsMXhsiPm6OJY7DjkJO4U6FTj9MxZ2TnOuc153sXEZbPLZVe8q4fr\nftcxN3k3nluL27y7mfsW934Pmoe/R43Hc09tT75nrxfq5e51wGt8rdra1LXd3uDt5n3A+6mPhk+G\nz6++BF8f31rfF36Gfvl+A/5M/43+rf7vA+wDygKeBGoGZgX2BTGCwoNaghaDHYPLg4Uh+iFbQm6H\nyoQmhPaEEcOCwprDFtY5rTu4bircJLwo/MF6jfW5629ukNmQvOHCRsbGyI1nIvARwRGtEcuR3pGN\nkQtRblF1UfM8B94h3qto2+iK6NkYq5jymOlYq9jy2Jk4q7gDcbPxNvGV8XMJDgk1CW8SXRPrExeT\nvJOOJ60kByd3pJBSIlLOpYqnJqX2pymk5aYNp+ukF6ULMywyDmbM8z34zZlI5vrMHgFLkC4YzNLM\n2pU1kW2dXZv9ISco50yuWG5q7uAm7U17Nk3nOef9tBm3mbe5L18pf0f+xBa7LQ1bka1RW/u2qWwr\n3Da13WX7iR2UHUk7fiswKCgveLczeGdvoXzh9sLJXS672opEivhFY7std9f/gPsh4YehPav3VO/5\nUhxdfKvEoKSyZLmUV3rrR8Mfq35c2Ru7d6jMtOzIPsK+1H0P9tvsP1EuVp5XPnnA60BXBaeiuOLd\nwY0Hb1YaV9YfohzKOiSs8qzqqVat3le9XBNfM1prX9tRJ1e3p27xcPThkSO2R9rr5etL6j8dTTj6\nsMGloatRvbHyGOFY9rEXTUFNAz9xf2pplmkuaf58PPW48ITfif4Ws5aWVrnWsja0Latt9mT4ybs/\nO/7c067b3tDB7ig5BaeyTr38JeKXB6c9Tved4Z5pP6t2tq6T2VnchXRt6prvju8W9oT2DJ9zP9fX\na9nb+aver8fPK52vvSBxoewi5WLhxZVLeZcWLqdfnrsSd2Wyb2Pfk6shV+/3+/YPXfO4duO68/Wr\nA3YDl25Y3Th/0+LmuVvcW923TW93DZoMdv5m8lvnkOlQ1x2zOz13ze/2Dq8ZvjhiM3LlnuO96/fd\n7t8eXTs6/CDwwcOx8DHhw+iHM4+SH715nP146cn2cfx48VPRp5XP5J41/q71e4fQVHhhwnFi8Ln/\n8yeTvMlXf2T+sTxV+IL+onJacbplxmjm/Kzz7N2X615OvUp/tTRX9KfYn3WvNV+f/cv2r8H5kPmp\nN/w3K3+XvpV+e/yd8bu+BZ+FZ+9T3i8tFn+Q/nDiI/fjwKfgT9NLOcvE5arPWp97v3h8GV9JWVn5\nBy6ikLxSF1/9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAAHXRFWHRTb2Z0d2FyZQBHUEwgR2hvc3Rz\nY3JpcHQgOS4xOJQFEHMAABsRSURBVHic7d29bxv3HcDhH4vOFNLZ7JAC5kCvKSCudgBrbQBzDQep\nIw1IQxd7oJcEkIFqFQfNdFCvNpB0pYCk3cRBQzuUmRuI/wA7XMowFEnx5Y739jwQWlpvPDHUffT9\n3ZGsTCaTAABp+03aGwAAIQgSABkhSABkgiABkAmCBEAm7ClIX3zxxX6uCICcSjxI7969e/bs2c3N\nTdJXBECuVfbzOKR2u311dbWHKwIgpxxDAiATfpv2BoR6vZ72JgAQp9vb2y2+Kv0ghW03na1VKnta\nqmWWmz0Vbvb923rMsGQHQCbsaUJyRgNAvCqVStqb8LO4ZtBMLNmF+zerERvgIVlYjYyxi9kI0v3b\nVKKSlIU7cQm52VPhZs+RbATpvgcT5U4GUCxZDdJ9cwXSJ4BiyU+Q5qzokzgB5FBugzRnNkLiBJBD\nRQnSrIVxUiaAbCtikGZNO2RsAsi2ogdpytgEEJNvvvnm66+//uSTT+L9tqUJ0qz7Y5MyAaztz3/+\n8z/+8Y9PP/30p59++t3vfhfXty33c9lNJj+/VSo/vwGwhm+++SaEEO+QVO4gTSkTUADTPViMb4t8\n++23P/300x/+8IcXL17EuPmC9GtzZQLIkekeLMa3e/79739/+umnX3311b/+9a+Tk5MYN1+QljAw\nASzy3Xffff3119HleJfsSnlSw/rmTn9w7gNQeicnJ3/5y18+//zz2L+zIK0nSpEsAYTw1VdfTS8X\n7uUn8mI2S0GZAOIkSJuzjgeQAEHagXU8gPgI0s5kCSAOghQTWQLYjSDFSpYAtiVICZAlgM15pobE\nTJ/rAYA1CFLCPC0ewHos2SXPCh7AGgRpX2QJYCVLdvvlwBLAEoKUBgeWAO4RpJQYlQB+TZBSZVQC\n+D9BSptRCSCEIEhZYVQCSk+QMsOoBOREjK8SO0uQMkaTgLISpOyxfAeUkmdqyKTp0zp4TgdgbZUQ\n/x+yk7C/vZAgZVg0KmkSsJ59xiMJluyyzfIdUBqClHnOvgPKQZByQpOAzJgkcyhBkPJDk4BCE6Rc\ncUgJKC5ByhuHlICCEqR80iSgcAQptzQJKBZByjNNAgpEkHJOk4CiEKT80ySgEASpEJwODuSfIBWF\n08GBnBOkYtEkILcEqXA0CcgnQSoiTQJySJAKSpOAvBGk4tIkIFcEqdA0CcgPQSo6TQJy4rdpbwDJ\ni5qUzCs8AiVUqVSePXs2/efl5eXJyUkM3zahV6JdX71ev729TXcbSkGToFgqldR24LNXXalUPvnk\nk//+97/Tj269V7dkVxrW7oBkfPrpp7F8H0t2ZWLtDgotib85l+0wPv/88+nlb7/9NpbrEqSS0SQo\nrn3+Zk8jVKlU/vnPf84eUtqaJbvysXYHZJIglZImAfH57rvvYvk+glRWmgRsa+7svq+++iqWbytI\nJaZJQJYIEgCZEM9ZduPxuNPpRBe63W6j0Zh+aDgcnp+fhxCOjo5evHgRy9URGyfdAZkRz4TU7/eb\nzebV1VW32+31erMfOj8/v7i4uLq66vf7sVwXMbNwB2RDPEEaDAatViuE0Gg0RqPR7IeePHkyGo3G\n43G1Wo3luoifJgEZENsDY6e9mQtPrVbrdDq1Wi2EMBwOZ1fzpioze8PUn1uvpKzdATuY3Y0/fvx4\nu28S20kN4/F47kKk3++/f//+6urq7OxsbjVvajIjru1hY+YkYFux7MbjCVKz2YwOEQ2Hw2gYIpc0\nCUhPPEt2rVar0+kMBoPoLLvondEzkJ+dnUUn4B0cHJydncVydSTI2h3kR6VYf0F6PSTuESRgB14P\nifhYuAPSIEgsoknA3gkSS2gSsF+CxHKaBOyRILGSJgGb2GWH4SXMeYgTwYF7kvhL1YQEwCqVyoK3\nEMJksvhtayYk1mBIghJYNvTs7VdfkFiPJkGxLMxPur/igsTaNAly635+MvirLEhsQpMgD3KRn/sE\nCSD35gqUi/zcJ0hsyJAEacvpAPQgQWJzmgR7Nxuhov7yCRJb0SRIWDFW4TYiSGxLkyBuZRiDVhAk\ngDSVPEKzBIkdGJJgKyK0kCCxG02C9YjQgwSJnWkSLCFCGxEk4qBJMGPaIb8TGxEkgBgYhnYnSMTE\nkEQpGYZiJEjER5MoDR1KgiARK02i0HQoUYIE8AAd2g9BIm6GJIpCh/ZMkEiAJpFnOpQWQSIZmkQO\nRSlyt02LIAFlZyTKCEEiMYYkMs9IlCmCRJI0iUwyEmWTIJEwTSJLjERZJkhA8RmJckGQSJ4hifQY\niXJEkNgLTWLvpCh3BIl90ST2RYpySpCAgnCgKO8EiT0yJJEMI1ExCBL7pUnESoqKRJDYO00iDlJU\nPIIE5IwUFZUgkQZDEluRomITJFKiSWxCisrgN2lvACUWNQkeEv3pokaFZ0ICsstgVCqCRKos3LGE\nFJWQIJE2TeLXpKi0BAnICikqOSc1kAHObsCZCwgSWaFJJVapWLUlBEt2ZIiDSeVjjY5ZggSkw58f\nzLFkR5ZYuCsHa3QsZEIiYyzcFZo1OlYwIZE95qSCch4dq5mQgMQZjFiHCYlMMiQViMGINQkSWaVJ\n+efkBTZiyQ5IhBSxKRMSGWZIyieDEdsRJLJNk/LGESO2JkhkniblhMGIHTmGBMRAitidCYk8MCRl\nmxoRC0EiJzQpkyzTESNLdsCWpIh4mZDID0NSZhiMSEI8E9J4PO50OtGFbrfbaDTufyiEcHFxUa1W\nY7lGSspzgWeA/wIkJJ4Jqd/vN5vNq6urbrfb6/VmP3R+fn5ycnJ1dfXkyZN+vx/L1VFq5qRUqRHJ\niSdIg8Gg1WqFEBqNxmg0mv3Qx48fB4NBu90OIRwfH8dydUAq1IhExXYMaboWN7cod3d3F0K4urqq\nVqtv375d+LWVGXFtD0VmSNo7B41YLZbdeGxBGo/HcxemTk9PQwjHx8c3NzcLv3YyI67toeA0aY88\nGxAPimU3Hk+Qms1mdHxoOBzWarXZDx0dHV1fX4cQrq+v5z4EO9GkvTAYsTfxnGXXarU6nc5gMIjO\nsoveWa/Xb29vj4+PO51OlKKLi4tYrg7YAy/zyp5VUl8li7qV7jaQY/6AT4bbla1tvVf3wFhyzsJd\nAtSIVAgS+adJsVIj0uK57ICfOWhEukxIFIIhaWfO7SZ1gkRRaNIOLNORBYJEgWjSVtSIjBAkKDU1\nIjsEiWIxJG1CjcgUQaJwNGkNniyVDBIkikiTVnJ6N9nkcUhQLgYjMsuEREEZkhZRI7JMkCguTfo1\nNSLjFgepXq8ve3VXyBNN+j81IvsWB+nk5OTm5ub+a78CeaRG5MLikxouLy9DCJ999tn0PV6yiLyK\nhqQS74/L/dOTJ4uDdHJysuftgASVuEll/bnJpcVLdsfHxyGEm5ub09PTEEL0v5BjpTyYpEbky+IJ\n6fz8/OPHj3d3dyGEfr8fNIkCKNOc5KGv5NHiCanf779//z66/P79+6hJQC6oETm19HFI1Wo1unB3\ndxeNSpB7JVi48zp75NfSkxqePXsWQmi324PBwDkOFEehF+6K+5NRCouDdHp6OplMPnz48J///Of4\n+NgBJAqloE0q4s9EuSwO0rt373q9XnS51+v9/ve/f/HixR63CtiMGlEAi48hXV5eHh0d/fDDDz/8\n8MPR0VH0OFkojmIdTFIjimFxkEajUavVqlar1Wq11WqNRqM9bxYkrihNUiMKY3GQms1mv98fj8fj\n8bjf7zupgWLKf5PUiCKZP4ZUr9enlz98+DC97LwGiinPJzjkdsNhsfkgGYYonXw2KYebDA+YD5JJ\nCLJPjSikxad9v337du7MOi8/QZHlakjKz5bCZpa+HtLR0VGtVtvz1kBqctKkPGwjbGlxkA4ODlqt\n1uHh4Z63BtKU+SZle+tgV4uD9Pz58/Pz82azOX2PY0uUQoablNXtgthUJovu47Mnf0eSO4ZUr9cd\noCJDMrnjz+RGwWJb79WXLtl1u93nz5/vtlWQQ9kbkjK2OZCUxUGq1Wq9Xm84HE7fY8mOEslSkzKz\nIZC4xUG6ubmZ/m9EkCiXbDQp509sBJtZHCQHdSAjTUr7+mF/lr6EOZDus69moIawV4snpH2eZQeZ\nltKcpEaU0OIgzT7F6mAw8JQNsE9qRDktDtLsKQyNRuP8/Hxf2wPZs98hSY0orYePIY1GI68YS9nt\n62CSGlFmax1Dmn0OISip5OckNaLkHj6GVK1WW63WvrYHMizJJqkRPHwMCfhFMk1SIwj3g3T/hO+I\n077hZ3E3SY0gMh+k2cW6u7u7fr8fHEOCOfE1SY1gaj5I0WLdeDzu9Xr9fr/ZbJ6cnHilPpgXR5PU\nCGYtOIbU6/V6vV61Wr24uPAKFLDUbk1SI5gzH6Rnz56NRqNarXZ0dDQcDqevQOE0B1ggG0/ACsUw\n/4qx+z+pwSvGkm9bBUnFKLDYXjFWG2Azmw9JagQLefkJ2NkmTyykRrCMIEEc1muSGsEKggQxeahJ\nagSrCRLEZ3mT1AgeJEgQq0VNUiNYhyBB3Pb14klQMIIECZhpkvEI1iRIkIzJJFQqagTrEyRISiVM\nJsHaHaxLkCARP89GjifB2gQJ4verlTpNgvUIEsRswXGjqEmyBCsJEsRp6VkM0Xs1CZaLJ0jj8bjd\nbrfb7S+++GL6Ekqz3r171263Y7kuyKvJxPIdrBBPkKIXO7+6uup2u71eb+6jw+Hww4cPsVwRZNla\nJ3lrEiwRT5AGg0Gr1QohNBqN0Wg0+6HxePz69euzs7NYrggya4OHHGkSLBLbMaRqtTp3IfL69etu\nt9toNFZ8bWVGXNsD+7TxA2A1iWKJZTc+/4qxWxuPx1GKxuPx9J3X19eDweDu7i6EMBwO3759e3p6\nev9rJx7LTp5t+XQMm7/ULGTW7G68Xq9v903iCVKz2ez3+8fHx8PhsFarTd9/eHj4/fffR5fb7fbC\nGkGu7dQUTYIZ8QSp1Wp1Op3BYDAej7vdbvTOer1+e3sby/eHwtIk+L9K6stlukV+xZaSaNk97V9G\niMXWe3UPjIUtxTnYeIgSCBJsJ5FlNk2i3AQJNpbgQR9NosQECTaT+CkImkRZCRJsYE8nxHl2cEpJ\nkCCTnOZA+QgSrCuFxwtpEmUiSLCW1B69avmO0hAkeFjKz6Vg+Y5yECR4QFae2UeTKDpBglWyUqOI\nJlFoggS54pASxSVIsFS2xqMph5QoKEGCxTJaoymjEoUjSLBA1msUMSpRLIIE8/JRoymjEkUhSPAr\nOatRxKhEIQgSFIVRiZwTJPhFLsejWUYl8kyQ4Ge5r9GUUYl8EiQIoUg1ihiVyCFBgsLVaMqoRK78\nNu0NgJQVtkaR6GeLmlTkn5MiMCFBCVjBIw8EiVIr+Hg0xwoe2WbJjvIqV40iVvDIMEGipMpYoylZ\nIpMEiTIqdY2mZImMcQyJ0lGjX3G+A5khSIDzHcgES3aUi/FoKSt4pE2QKBE1epgskR5BoizUaAOy\nRBoEiVJQo23MZikoE4kTJIpPjXYyve0MTCRMkCg4NYqNdTwSJkjAJmSJxAgSRWY8SorDSyRAkCgs\nNUqcw0vESpAoJjXaKwMTcRAkCkiN0mFgYjeCRNGoUfoMTGxFkIBkzA1MQZl4gCBRKMajLFIm1iNI\nFIcaZZ0ysZIgURBqlCfKxCKCRBGoUV4pEzMEidxToyJQJgSJvFOjorlfpiBOZSFI5NjsLouimY2Q\nsakcBIl8s4MqBWNTOQgSeWWxrowWjk1BnApCkMglNWJpnII+5ZUgkT9qxLy5O4ThKZ8EiZxRIx5m\neMonQSJP1IiNrRie7n+UVAkSuaFGxECfMkyQyAc1IhGr+3T/E0iSIJEDHgDLntzPj0TtkSCRD3YC\npGOdRC38NDYnSGSdxTqyZeHdcdkU7767CUEi09SIfFhWKaHahCCRXWpEvi27+64I1YqvKgFBIqPU\niMJacc9e3arVX5t/gkQWqREl9eD9fp1TTnP7yyNIZI4awVLrFGv9x0lk7DdNkMgWNYKdrP/7s1G6\nNv3mW/lNot8dNqJGsD+TycZvUcMefNtWPBPSeDzudDrRhW6322g0ph+6vr4+Pz+vVqsHBwdnZ2eP\nHj2K5RopHjWCrFvzV7Re3+7bxzMh9fv9ZrN5dXXV7XZ7vd7sh169etXtdq+uro6Pjy8vL2O5OopH\njYB4gjQYDFqtVgih0WiMRqPZD7158yYamA4ODsbjcSxXR8GoERBiPIZUrVbnLkQODw9DCNfX151O\n5/j4eOHXVmbEtT3khRpBAcSyG48tSNPpZ24MGo/H7Xb7w4cPFxcXs8eWZk1mxLU95IIaQTHEshuP\nJ0jNZrPf74cQhsNhrVab/VC73T45Oel2u05nYI4aAbMqsQwlC8+yq9frf//7358+fdpsNqNPe/Lk\nyenp6dzX1uv129vb3beBfFEjKKqt9+rxBGkXglRCagQFtvVe3QNj2Tc1AhYSJPZKjYBlBIn9USNg\nBUFiT9QIWE2Q2Ac1Ah4kSCROjYB1CBLJUiNgTYJEgtQIWJ8gkRQ1AjYiSCRCjYBNCRLxUyNgC/G8\nhDlMZbBGlZDg62xNQsZ+WsgtQSJOqddoYXuSa0YlVPZ8jVBggkRsUqnRXA/2XIJlV5fuVkFOCRLx\n2FuNcrGvn9uqXGwzpE6Q2FWlEkJItkazO/Q87s31CdYhSOwk0cFouuMu2C57RZ8K9pPCRgSJ7SVR\noxLunWd/zBL++DAlSGwp3hoVdRjalDhRZoLENuKqkQ6tsDBObigKTJDYTFynMER7WLvXNU1vKGMT\nBSZIbGD3wchf+jsyNlFggsS6dqmRXWcSjE0UjCCxlq1rZGluD4xNFIMg8YDtDhrZLabl/tjkPwF5\nIUisssVgZCTKCGUidwSJpTaqkb1eZjnURF4IEgtstExnJMoLh5rIOEFi3pqDkT1arlnQI4MEiV9Z\np0ZGoiJRJrJDkPjZOst0UlRgykTqBIkQHhqM7KFKxUkQpEWQym71YGQkKjMnQbBnglRqKwYjKWKW\nBT32QJBKatlgZHfDaspEcgSpjBYORkYiNuJQE7ETpHJZOBhJEbtwqIm4CFKJzA1G9h3EzoIeuxCk\nUpgbjIxEJM2CHlsQpIKbTZE/Wtk/C3qsT5CKbLpGZyQiC4xNrCZIxRQNRmHy///3O0+WGJtYSJCK\nZpqi4DecPDA2MSVIxSFF5NrCsSm4M5eJIBXBNEWTMAl+e8k/cSonQco3KaLwxKk8BCmXKqESKpPw\ny0OL/GZSCuJUbIKUJz//BlYmIUzWeZVxKDBxKh5ByoeZFD38EuNQNsviFPQpVwQp03751ZIiWM9c\ngQxPOSJImTP3+7P6FV2B1QxPOSJIWTH3ePWoQxUpgvisGJ7uf5T9E6SU3X/eFCMR7Mfc8KRPqROk\ndCzrUJAiSMP9/OjT/gnS/iw8uKpDkE2r1/fufwK7E6RkrTjDx9Ic5MiDI9TCz2EjghS/1aeZGomg\nGNZJ1MJPYxlBiseDj3XQISi8hYkySK1PkLa0zgHPysyn6BCU0OI9w6JBatknl4ogbWDNh3wbhoAV\nllXKip8gLbXRoG0YAnaxbPeybJxa8SX5JUi/2HSpV4SApK0IVfFaVdIgbT0aixCQBasWbFa2avXX\npqv4QdpxWbby668WISDjHty/rc7Vmt8kCcUJUlwnrigQUGzrFGudaK3zrTaSpyA9eANtd9MoEMCs\nNfely7r1ODze7np/s92XJSH62Va8hRAmYbLiba1rufd9575LGVQqa/3tQ7zc7Klwsydnl13xQolP\nSMPh8PXr19VqNYRwcXERXZhz/6mvY7HwfliS5ADkTuITUq/X63a7V1dXzWaz3+8v/Jwdu7pspFrY\nbgCyKfEgjUajRqMRQmi1WoPBYLtvsnotb9kSHgA5kviS3XSNbuFiXWSdNV6BASi2xIM0Ho/nLtz3\n+HH9we9Tf/hTWNfjx4/rbtC9c7Onws2eI4kHqVarDYfDRqPR7/ebzeb9T7i9vU16GwDIvsok4bWw\ndc6yA4DEgwQA68jQA2MBKDNBAiAT0nwuu+jwUgih2+1Gj1UiOe12O7rw5MmT09PT8Xjc6XSi9zi2\nl4ToFr66uppejt4f3dru/Alpt9vRbR5CePv27WAwiO7bZ2dntVot+q8wHo/d7HG5vr4+Pz+vVqsH\nBwdnZ2ePHj2au29vtquZpKfT6dzd3d3d3XU6nRQ3oyS+/PLL2X9eXl72+/3JZNLv9y8vL1PaqMI6\nPz//7LPPprf5/VvbnT92Nzc3f/rTnx4/fjx9z/n5+WAwmP7z8vIyuvFvbm7c7HF5+vTpzc3NZDK5\nubl59erV5N59e6NdTZpLdnd3d9VqtVqt3t3dpbgZJTEajV6+fPny5csff/wxhDAYDJ4/fx5CePHi\nxdbPoMEyp6en33///fSf929td/7YNRqNv/3tb7OPLbm5uRkMBu12++3btyGEwWDQarWizxyNRqlt\naLG8efMmmjUPDg6iB5vO3bc32tU4hlQKP/7445MnT7rd7vHx8atXr6J3WqbbJ7f2/tVqtWazGa3g\nRU1a54lj2Mjh4WEI4fr6utPpHB8fL/yc9W/tPL0eElt79OjRX//61xDC7Lr5eDyuVqsrnkGDGLm1\n96/b7UYXTk9Po2Oo0X+FsPKJY9hIdIioVqtdXFw8evRo2eeseedPc0KKRrzxeHxwcJDiZpTBcDjs\n9Xph5vew2Wx+/PgxhPDx48eFz6BBjO7f2u78e/Dy5cvoDn99ff3kyZPpCw4Mh8NarZb21hVEu90+\nOTnpdrvTGs3dtzfa1aT5wFgnGu3T69evh8NhtVo9Ozvb+NQXtjI948tZdnsze5bdu3fv+v1+dAJY\nNC05yy5eP/7449OnT6eZic7g3eUsO8/UAEAmOKkBgEwQJAAyQZAAyARBAiATBAmATBAkADJBkADI\nBEGCXV1fX9fvCSHU6/XoKdSAdXhgLOxqPB4Ph8MQwpdfftlqtY6OjkIIh4eH19fXtVpt2RN8AXNM\nSLCrarV6eHgYPe3xwcHB9PKXX34ZPXlavV5/+fLls2fPopmp3W7X6/Uvvvgieqa1d+/e/fGPf6zX\n669fv/akn5SZIME+jEajN2/eHB0dXV5eNpvNi4uLm5ubXq83HA5fvXp1fHz8ww8/fPz48fz8PO0t\nhdQIEuxDs9k8PDyMnmT6+Pg4esmyEEL0RMhHR0fVavX58+fRP6GcvB4SpO/p06dpbwKkz4QEaYpG\npffv39/e3p6cnESvsQ3lJEiQpkaj8ebNm+g0h5ubG0GizJz2DUAmmJAAyARBAiAT/gdf7POFnmx6\npwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plot(t,x(:,1),\"-r\",t,x(:,2),\"-g\",t,x(:,3),\"-b\")\n",
    "xlim([0 200])\n",
    "xlabel(\"Time\",\"fontweight\",\"bold\")\n",
    "ylabel(\"Number\",\"fontweight\",\"bold\")\n",
    "h = legend(\"S\",\"I\",\"R\");\n",
    "legend(h,\"show\")"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Octave",
   "language": "octave",
   "name": "octave"
  },
  "language_info": {
   "file_extension": ".m",
   "help_links": [
    {
     "text": "GNU Octave",
     "url": "https://www.gnu.org/software/octave/support.html"
    },
    {
     "text": "Octave Kernel",
     "url": "https://github.com/Calysto/octave_kernel"
    },
    {
     "text": "MetaKernel Magics",
     "url": "https://github.com/calysto/metakernel/blob/master/metakernel/magics/README.md"
    }
   ],
   "mimetype": "text/x-octave",
   "name": "octave",
   "version": "4.2.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}