{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# A Simple Model of the Glycolysis of Human Erythrocytes\n",
    "\n",
    "This is a model for the glycolysis of human erythrocytes which takes into account ATP-synthesis and -consumption.\n",
    "This model is based on the model introduced in the following publication.\n",
    "\n",
    "* Rapoport, T.A. and Heinrich, R. (1975) \"Mathematical analysis of multienzyme systems. I. Modelling of the glycolysis of human erythrocytes.\", Biosystems., 7, 1, 120-129. \n",
    "* Heinrich, R. and Rapoport, T.A. (1975) \"Mathematical analysis of multienzyme systems. II. Steady state and transient control.\", Biosystems., 7, 1, 130-136."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "from ecell4 import *\n",
    "util.decorator.ENABLE_RATELAW = True"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false,
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "with reaction_rules():\n",
    "    2 * ATP > 2 * A13P2G + 2 * ADP | (3.2 * ATP / (1.0 + (ATP / 1.0) ** 4.0))\n",
    "    A13P2G > A23P2G | 1500\n",
    "    A23P2G > PEP | 0.15\n",
    "    A13P2G + ADP > PEP + ATP | 1.57e+4\n",
    "    PEP + ADP > ATP | 559\n",
    "    AMP + ATP > 2 * ADP | (1.0 * (AMP * ATP - 2.0 * ADP * ADP))\n",
    "    ATP > ADP | 1.46"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "y0 = {\"A13P2G\": 0.0005082, \"A23P2G\": 5.0834, \"PEP\": 0.020502,\n",
    "      \"AMP\": 0.080139, \"ADP\": 0.2190, \"ATP\": 1.196867}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEPCAYAAABFpK+YAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XmUFPW5//H3M4woSyMDOiwugdyAK0ZICCSINiiKipKo\nwXjiBlmOniBorkpwiYMxv98FvYlrfiYgyygoXg2KG66MSLioBGQTSURBRBGiDFtQUZ7fH93T9MB0\nT81S3TNTn9c5dWaqqqvq6e85PNR861vfx9wdERGJjoJ8ByAiIrmlxC8iEjFK/CIiEaPELyISMUr8\nIiIRo8QvIhIxoSZ+M2trZo+Z2Soze9vM+oZ5PRERqV5hyOe/C3jW3S8ws0KgVcjXExGRalhYL3CZ\n2cHAEnf/ZigXEBGRWgmzq6crsNnMppjZYjObaGYtQ7yeiIgEEGbiLwR6AX9y917ATuA3IV5PREQC\nCLOP/0PgQ3d/M7n+GPskfjPTREEiIrXg7lbbY0O743f3jcB6M+ue3HQasLKKz2lx55Zbbsl7DA1l\nUVuoLdQW2Ze6CntUz1XAdDNrDqwBhod8PRERqUaoid/dlwK9w7yGiIjUjN7cbSDi8Xi+Q2gw1BZ7\nqS32UlvUn9DG8Qe6uJnn8/oiIo2RmeEN8eGuiIg0TEr8IiIRo8QvIhIxSvwiIhGjxC8iEjFK/CIi\nEaPELyISMUr8IiIRo8QvIhIxSvwiIhGjxC8iEjFK/CIiEaPELyISMUr8IiIRo8QvIhIxSvwiIhGj\nxC8iEjFK/CIiEaPELyISMUr8IiIRo8QvIhIxSvwiIhGjxC8iEjFK/CIiEaPELyISMUr8IiIRUxj2\nBcxsLbAN+BrY7e7fC/uaIiKSWeiJH3Ag7u6fVbWz1597YWYAGFZpX8X2ffelb6/quEyfC3JMdccF\nle38jVV9tIuI5F8uEj+QOQtOPGciAI5X2u6+dz19X/r2qo7L9Lkgx1R3XFDZzt9Y1Ue7iEj9OJmT\n63S8hf0P2szeA7aS6Or5s7tPTNvnSigiIjVjZrh7rf8Ez8Udfz93/9jMDgVeNLN33P21HFxXRESq\nEHrid/ePkz83m9ks4HtAKvGXlJSkPhuPx4nH42GHJCLSqJSVlVFWVlZv5wu1q8fMWgLN3H27mbUC\nXgDGufsLyf3q6hERqaGG3tXTAZiVHA1SCEyvSPoiIpIfoT/czXpx3fGLiNRYXe/49eauiEjEKPGL\niESMEr+ISMQo8YuIRIwSv4hIxCjxi4hEjBK/iEjEKPGLiERMtYnfzIaZWZvk7zeb2Swz6xV+aCIi\nEoYgd/w3u/s2MzsJOBV4APh/4YYlIiJhCZL4v07+HAJMdPengebhhSQiImEKMknbBjP7CzAI+C8z\nOwg9GxCJPJXizI0w5jOrdpK25HTKZwDL3f2fZtYJ6FEfs2xqkjaRxis5UVi+w2jSMrVx6JO0uftO\nYDNwUnLTV8C7tb2giIgEd+eddzJx4sTqP1gDQUb1lADXA2OTm5oDD9ZrFCIiUqXi4mJ27txZr+cM\n0lf/I2AosBPA3TcAsXqNQkREciZI4v/C3fdUrCT7/EVEpJEKkvj/x8z+DLQ1s18CLwOTwg1LRETC\nEuTh7u3A48mlO4kXuu4OOzARkbqKx+O0a9eOL7/8MrVt7ty5DBgwgLZt29K1a9f9jhkwYADFxcW0\nadOGY445ptKD1alTp9KsWTNisRgHH3wwPXv25JlnngFg4cKFDBo0iPbt21NcXMywYcPYuHFjpXMv\nWrSIIUOG0K5dO4qKijjuuOO46aabKC8vD6kFqhZoPL67v+Du1yaXF8MOSkSkrtauXcsbb7xBcXEx\ns2fPTm1v3bo1P//5z7n99turPO7uu+9mw4YNbNu2jWnTpnHVVVexevXq1P5+/fqxfft2ysvL+dnP\nfsawYcMoLy+nvLycK664gnXr1rFu3TpisRjDhw9PHbdgwQIGDBhA//79Wb16NVu2bGHOnDkUFhay\ndOnS8BqiChlf4DKzHUCmQbru7m3CCUlEpO5KS0s57bTT6NOnD9OmTeOCCy4AoHfv3vTu3ZuXXnqp\nyuN69OhRab1169a0abM33VWMqzczhg8fzqhRo3jvvfcYPHhwpeN+9atfEY/HU+vXX389I0aMYMyY\nMaltRxxxBCUlJXX5mrWS8Y7f3Vu7eyzDoqQvIg1aaWkpF154IcOGDeP5559n06ZNgY8dMmQILVq0\nIB6PM3nyZDp16rTfZ7766ismTZpELBajW7du++2fN28exx9/PAA7d+5k4cKFnH/++bX/QvUoyDj+\nI6tachGciDReZvWz1Mb8+fPZsGED5557Lt26dePYY49lxowZgY9/+umn2bFjB6WlpVx++eV88MEH\nqX0LFy6kqKiITp06MXPmTGbNmkUsVnmE+7Jly/jd736X6k7asmULe/bsoWPHjqnPXH/99RQVFdG6\ndWt+//vf1+6L1lKQPv5ngWeSy8vAe8BzYQYlIo2fe/0stTFt2jROP/30VEL+8Y9/zLRp02p0jmbN\nmnHBBRfQp08fZs2aldret29ftmzZwubNm1mwYAEDBw6sdNy7777LWWedxd13302/fv0AKCoqoqCg\ngI8//jj1uQkTJrBlyxZ+9KMf8fXXX5NL1U7S5u7Hp68n5+L/VWgRiYjUwa5du3j00UfZs2dPqovm\niy++oLy8nGXLlnHCCSfU6Hy7d++mVatgry+tW7eOQYMG8dvf/paf/vSnqe2tWrWiT58+PP7445xy\nyimVjnH3nM95VONZNt19MdAnhFhEROrsiSeeoLCwkFWrVrF06VKWLl3KqlWr6N+/P6Wlpbg7n3/+\nObt378bd+eKLL1LDPVevXs1zzz3Hrl272L17Nw899BCLFi3i9NNPr/a6GzZsYODAgYwcOZJf/vKX\n++2fMGECkydPZvz48annDR9++CFr167N+UynQfr4/zNtuc7MHgY25CA2EZEaKy0tZcSIERx++OEU\nFxdTXFxMhw4dGDlyJDNmzOCVV16hZcuWnH322axfv54WLVqkRuS4O+PGjaNDhw507NiRSZMm8cwz\nz3DkkYnHmmaWMUlPmjSJ999/n5KSEmKxGLFYrNJooH79+vHKK68wb948jjrqKIqKijjzzDMZMGAA\nV111VfgNkybItMwl7B3W+RWwFnjc3T+v88U1LbNIo6VpmcNnZkyfPp1NmzZx9dVXV9pel2mZg/Tx\nl9T25CIi0vAE6ep50czapq23M7Png17AzJqZ2RIze6q2QYqISP0J8nD3UHdPTSTh7p8BHWpwjdHA\n22R+C1hERHIoULF1M/tGxYqZdQH2ZPx0GjM7HDiLxGyeKtApItIABCm2fiPwmpnNS66fDOw/Vqlq\nfwSuAzTFg4hIAxHk4e4cM/sOe8fuX+3u/6ruODMbAmxy9yVmFs/0ufQJiuLxeKVJjUREBDZu3Fiv\nk7kFGc5ZAPwU6Orutybn6eno7m9Uc9z/AS4hMQT0IBJ3/Y+7+6Vpn9FwTpFGSsM5wxfWcM4gffx/\nAr4PXJRc35HclpW73+DuR7h7V+AnwCvpSV9ERPIjSB9/H3fvaWZLIDGqx8wOqMW1dGsgItIABLnj\n/9LMmlWsmNmhBBzVU8HdX3X3c2sanIhIXVRVevH222+nR48etGnThm9+85vccccdlY5R6cWEe4BZ\nQHGy3/5vwP8NNSoRkTrKVHoR4MEHH6S8vJw5c+Zw7733MnPmzNS+KJRerPbhLoCZHQOcmlx92d1X\n1cvF9XBXpNFq6A93b731VhYtWkSfPn1YuHAhTz1V9eQBo0ePxt25++6799v3xhtvcNZZZ7F8+XI6\nderE1KlTeeCBB3jttdeARGWtWCzGokWL6NWrV6VjFy9eTDweZ9u2bQCcdNJJfOc73+Guu+4K/B1y\n/nA3OTVDOzNrB3wCPJxcPkluExFpsIKUXnT3SiUSKzT10ovZHu4uJkuxdeCb9R+OiDQVNq5+Xtb3\nW2r+V0V66cVYLJYqvZh+1wx73yNK75KBROnFr7/+mlmzZnH55Zfz1ltvpaZmrii9WFhYSLdu3bKW\nXqzoYspUenHixIns3r2bsWPHcuONN9b4e9ZWxsTv7l1yFoWINDm1Sdj1JVPpxfTEf++99/LQQw/x\n2muvccAB+w9UrCi9+MADDzBr1ixGjx4NJEovVnT1VKW60ovdu3cHEoVZJkyYwCWXXNLwSi8CmNlQ\nElM1OPCqu2umTRFpkIKUXpw8eTITJkxg3rx5dO7cOev5mmLpxWoTv5n9F9AbmE5iorVRZvYDdx8b\ndnAiIjVVUXpx6dKlNG/eHEgk12HDhlFaWkrPnj258cYbmTt3Ll26dKl07OrVq3nvvfeIx+MUFhYy\nc+ZMFi1axOTJk6u9bpDSi2eccQaHHXYYw4cPp7i4OFV6seKvgJyp+N8m0wIsB5qlrTcDlld3XJAl\ncXkRaYwa6r/fwYMH+7XXXrvf9kcffdQ7duzoXbt29ebNm3vr1q1Ty5VXXunu7qtWrfI+ffp4LBbz\ndu3a+SmnnOLz589PnWPq1Knev3//Kq9bUlLiZlbpvLFYrNJnXn/9dT/rrLO8bdu23rZtWz/++OP9\npptu8s8++6zKcwI+ffp0/+Mf/7jfdq9D7g0yV88yYIC7f5pcbw/Mdfealaqv+txe3fVFpGFq6MM5\nm4K8lV4k8bLWYjMrS66fAvymthcUEZH8CjIt88Nm9iqJfn4Hxrj7xmoOExGRBipj4jezXvts+jD5\ns7OZdXb3xeGFJSIiYcl2x78IWAF8mmH/gPoPR0REwpYt8f8a+DHwb2AmMMvdt+ckKhERCU3GuXrc\n/U537weMAg4HXjaz/zGzE3MWnYiI1Ltqp2V29zXAk8ALJB7wHhV2UCIiEp5sD3f/g0TJxKHAByS6\ne37v7rtyFJuIiIQgWx//P0m8tfsEsA04ErjSzIzEW2N/yEF8IiJSz7J19dwK/JVEmcXWySWW9lNE\npEGrqvTi5ZdfzoEHHkibNm1o06YNPXr04IYbbkgVTIHsJRabgmwPd0vcfVymJZdBiojUVKbSi2bG\nmDFj2LZtG//617+YMmUKCxcupF+/fvz73/9Ofa6qEotbt27Nx1epd0Fq7oqINDqlpaWcdtppXHLJ\nJUybNq3Svoo5hpo3b853v/tdZs+ezaeffsqUKVP2+4yZMXz4cHbt2sWaNWty9wVCpMQvIk3SvqUX\nN2/enPGzrVu3ZtCgQVUWWKmuxGJjlK3m7ujkz5NyF46INBlm9bPUQnrpxW7dunHssccyffr0rMd0\n6tSJzz77LLVeUWKxU6dOzJw5s8oSi41Vtjv+Ecmf9+QiEBFpYtzrZ6mFTKUXE2FVfc4NGzbQvn37\n1Hrfvn3ZsmULmzdvZsGCBQwcOLBWsTRE2YZzvm1m/wQOM7Pl++zz+piPX0SkvmUqvbh161aWLVuG\nmWH7/CWxY8cOXnrpJW6++eZ8hJxz2YqtX2RmHUm8sXsOibKLIiINWnWlFyvWIfEfwooVKxgzZgzt\n27dn+PDheYs7l7I+3HX3jck7+4/ZO5Z/g7uvzUFsIiI1VlpayogRIzj88MMpLi6muLiYDh06MHLk\nSKZPn85XX33FhAkTaNOmDYcccgiXXXYZvXv3ZsGCBbRo0QKgyr8KmpIgpRfjwDRgXXLTkcBl7v5q\nnS+u0osijZZKL4Yvn6UX/wCc7u6rkxfsDjwC7FuoZT9mdhDwKnAg0Bx40t3H1jZYERGpuyCJv7Ai\n6QO4+z/MLMhxuPvnZjbA3f+dPGa+mZ3k7vNrG7CIiNRNkAT+dzObBDxE4gHvT0lU5wrE3SvegW4O\nNAM+y/JxEREJWZA3d68EVpEoyHIVsDK5LRAzKzCzt4BPgLnu/nZtAhURkfpR7R2/u38O/HdyqTF3\n3wOcaGYHA8+bWdzdyyr2l5SUpD4bj8eJx+O1uYyISJO1cePGSrmyrqod1VOfzOxmYJe735Fc16ge\nkUZKo3rCF9aonlAnaTOzQ8ysbfL3FsAgYEmY1xQRkeyyJn4za2Zmd9Th/J2AV5J9/K8DT7n7y3U4\nn4iI1FHWPn53/9rMTrJa9sm4+3ICjPcXEZHcCdLV8xbwpJldYmbnJ5fzwg5MRKSuMpVeLCgoqFSV\nC+Caa66hoKAgNYtnUy6/GCTxH0Ri7P1AYEhyOSfMoERE6ipb6cXu3bunJmyDRLGVRx99lG9961uV\n5uhpquUXq0387n55chmevuQiOBGR2spWevGcc85h/vz5lJeXAzBnzhy+/e1v06FDh0qfa6rlF6tN\n/GZ2lJm9bGYrk+snmNlN4YcmIlJ72UovHnTQQQwdOpRHHnkk9dlLL70047maWvnFIFM2TASuA+5P\nri8HHgZuCysoEWn8rKysXs7jtXipM730YiwWS5VeTB8Lf+mll3Lddddx0UUXMW/ePEpLS7nvvvsq\nnaei/GJhYSHdunVrMuUXgyT+lu7+ekW/l7u7me0ONywRaexqk7DrS6bSi1dffTXujpnRr18/Nm/e\nzG233cY555zDQQcdtN95+vbtW2UB9sYuSOLfbGbfqlgxswtIFGYREWlwgpRerHDxxRdz6623UlZP\nf500FkES/0jgL8DRZvYR8D6JGTpFRBqcmpReHDVqFCeffDL9+/fPW7z5EGSStjXAqWbWCihw9+3h\nhyUiUjvppRfTjRw5klGjRjFo0KDUXX9RUREDBgyo8jxNufxikNKLhwC3ACcBDrwG3Orun9b54pqk\nTaTR0iRt4cvnJG2PAJuA84ALgM3AzNpeUERE8itIH39Hd/9d2vptZnZhWAGJiEi4gtzxv2BmFyUr\naRUkk/4LYQcmIiLhyHjHb2Y7SPTpA1wNPJj8vQDYCfxnuKGJiEgYMiZ+d2+dy0BERCQ3gvTxY2Yn\nAF3SP+/ufw0pJhERCVG1id/MpgA9gJXAnrRdSvwiIo1QkDv+PsBxGnAvItI0BBnV8yZwbNiBiIhI\nbgRJ/FOA/zWzf5jZ8uSyLOzARETqKr304hVXXEEsFiMWi3HggQfSvHnz1PrZZ5/NunXrKCgoSG3r\n2rUr48ePz/dXCEWQrp4HgIuBFVTu4xcRabAqSi8eeeSRzJ49m/vvv5/770+UFRk3bhxr1qypVH5x\n7dq1AGzdupWCggIWLlzIqaeeyoknnsgZZ5yRj68QmiB3/Jvcfba7v+fuayuWsAMTEamLbKUX3b3a\neYb69u3Lcccdx8qVK8MMMy+CJP4lZjYj+fbu+cnlvNAjExGpg31LL27atCnQcRX/Kfztb39j5cqV\n9OzZM+RIcy9QBS7gS+D0fbZrOKeIZFRmZfVynrjHa3xMVaUXZ8yYUWmGy0wOOeQQzIyOHTsyfvz4\njNM2N2ZB5uO/PAdxiEgTU5uEXV+ylV6szqeffkpBQZDOkMYr6Atc6RzA3UeEEpGISB1kKr1YXl7O\nsmXLOOGEE/IcYf4F6ep5hr2TtbUAfgR8FFpEIiJ1UF3pxTvuuCPPEeZfkK6ex9LXzWwG8LfQIhIR\nqYNspRdHjx7N+PHjM5ZVbKqlFvdVbenF/Q4wOxp42t2/FeCzRwClQDGJvxr+4u53p+3XTBAijZRK\nL4YvrNKLQfr40+fld+ATYEzA8+8GrnH3t8ysNfB3M3vR3VfVKloREamzIF09tZ6X3903AhuTv+8w\ns1VAZ0CJX0QkT7JV4Doy24Hu/kFNLmRmXYCewOs1OU5EROpXtjv+Z9nbxZPu0OTSLOhFkt08jwGj\n3X1H+r6SkpLU7/F4nHg8HvS0IiKRsHHjxkq5sq4CP9xN3rH/BjgNuMvd7wl43AHA08Bz7n7nPvv0\ncFekkdLD3fCF9XC32tfTzKy7mU0F5gB/B46pQdI3ErN7vr1v0hcRkfzImPjNrIeZPQI8DrxMogrX\nRHffXYPz9yMxpfMAM1uSXAbXLWQREamLbH38bwEfkuim+R7wvbSXG9zdR1V3cnefT7AZQEVEJEey\nJf6fJX/u24lnVWwTEZFGIuPduLtPTS7T9lmmuvu0TMeJiORTly5daNmyJbFYjI4dOzJ8+HB27txJ\nPB6nRYsWqdKKsViMoUOHAlBWVpYqu9imTRuOPvpopk6dmt8vEiJ1w4hIk2JmPP3002zfvp3Fixez\naNEibrvtNsyM++67j+3bt6eWJ598MnXcYYcdxvbt29m2bRvjx4/nF7/4BatWNc13TZX4RaTJ6ty5\nM2eeeSYrVqyo0XFDhw6lqKhIiV9EpLGoeL9g/fr1PPvss/Tq1avS9mz27NnDrFmzKC8vp0ePHqHG\nmS9BJmk7CvgT0NHdjzOzE4Bz3f220KMTkUarrKx+pjiOx2s2lsTd+eEPf0hhYSEHH3wwQ4YMYezY\nsbz66quMGjWKa6+9NvXZUaNGMW7cOAA++ugjioqKKCgo4Bvf+AYPPfQQ3bp1q5fv0NAEKcQyEbgO\nuD+5vhx4GFDiF5GMapqw64uZ8eSTTzJw4MD9tt9zzz2MGFF18cDOnTuzfv36XISYd0G6elq6e2pi\nteQcCzV5iUtERBqQIIl/s5mliq6Y2QXAx+GFJCISDs0tlBAk8Y8E/gwcbWYfAdcAV4YalYhICEaO\nHFlpHH/v3r1T+6JSdhGCFWJZA5xqZq2AAnffHn5YIiK18/7771e5fe7cuRmPicfjfPBBjUqMNGpB\nRvUcBJwPdAGaJWfcdHe/NeTYREQkBEFG9TwJlJOYkvnzcMMREZGwBUn8h7n7GaFHIiIiORHk4e6C\n5EtbIiLSBGQrtr4S2EOitu5wM3sf+CK5291d/xmIiDRC2bp6OgMnkph/X0REmohsiX+tu6/LWSQi\n0uhEaex7U5It8R9qZr+m6jt+d/c/hBSTiDQCFW/B3nnnnRQXF+c5GqmJbIm/GRDLVSAi0ji1atWK\nTZs25TuMJq1Vq1b1ej7LNHeFmS1x9571erX9r+GaO0NEpGbMDHevdT+bCrGIiERMtjv+9u7+aagX\n1x2/iEiN1fWOP2PizwUlfhGRmqtr4g8yZUOotv7v1nyHICISKXlP/Gt+vSbfIYiIRIq6ekREGhmN\n6hERkRoJNfGb2WQz+8TMlod5HRERCS7sO/4pwOCQryEiIjUQauJ399eALWFeQ0REakZ9/CIiEaPE\nLyISMXkfx19SUpL6PR6PE4/H8xaLiEhDVFZWRllZWb2dL/Rx/GbWBXjK3XtUsU/j+EVEaqhBj+M3\ns4eBBUB3M1tvZsPDvJ6IiFRPb+6KiDQyDfqOX0REGh4lfhGRiFHiFxGJGCV+EZGIUeIXEYkYJX4R\nkYhR4hcRiRglfhGRiFHiFxGJGCV+EZGIUeIXEYkYJX4RkYhR4hcRiRglfhGRiFHiFxGJGCV+EZGI\nUeIXEYkYJX4RkYhR4hcRiRglfhGRiFHiFxGJGCV+EZGIUeIXEYkYJX4RkYgpzHcAHHVUviMQEYkU\nc/f8XdzM/Z138nZ9EZHGyI4+Gne3Wh+f98Sfx+uLiDRGZlanxK8+fhGRiAk18ZvZYDN7x8z+aWZj\nwryWiIgEE9rDXTNrBtwLnAZsAN40s9nuvir9c+etWBFWCI3K5jff5NDevfMdRoOgtthLbbGX2qL+\nhDmq53vAu+6+FsDMHgGGApUS/8UdOoQYQuMx8+23uXDIkHyH0SCoLfZSW+yltthrVh2PDzPxHwas\nT1v/EOiz74fOO/TQEENoPJa1aqW2SFJb7KW22EttUX/CTPyBhuuUldX6wXSTsnYtlJWNy3cYDYLa\nYi+1xV5qi/oT2nBOM+sLlLj74OT6WGCPu49P+4zGcoqI1EKDHMdvZoXAauBU4CPgDeCifR/uiohI\nboXW1ePuX5nZSOB5oBnwgJK+iEj+5fXNXRERyb28vbkb5Ze7zOwIM5trZivNbIWZjUpub2dmL5rZ\nP8zsBTNrm+9Yc8XMmpnZEjN7KrkeybYws7Zm9piZrTKzt82sT4TbYmzy38hyM5thZgdGpS3MbLKZ\nfWJmy9O2Zfzuybb6ZzKnnl7d+fOS+NNe7hoMHAtcZGbH5COWPNkNXOPuxwF9gV8lv/9vgBfdvTvw\ncnI9KkYDb7N3NFhU2+Iu4Fl3PwY4AXiHCLaFmXUBfgH0cvceJLqLf0J02mIKifyYrsrvbmbHAheS\nyKWDgT+ZWdbcnq87/tTLXe6+G6h4uSsS3H2ju7+V/H0HiZfaDgPOBaYlPzYN+GF+IswtMzscOAuY\nBFSMVIhcW5jZwUB/d58Miedk7r6VCLYFsI3EDVLL5ECRliQGiUSiLdz9NWDLPpszffehwMPuvjv5\nwuy7JHJsRvlK/FW93HVYnmLJq+SdTU/gdaCDu3+S3PUJEJXXmv8IXAfsSdsWxbboCmw2sylmttjM\nJppZKyLYFu7+GfDfwAckEn65u79IBNsiTabv3plEDq1QbT7NV+LXE2XAzFoDjwOj3X17+r7kfNVN\nvp3MbAiwyd2XsPduv5KotAWJUXa9gD+5ey9gJ/t0ZUSlLczsP4CrgS4kEltrM7s4/TNRaYuqBPju\nWdslX4l/A3BE2voRVP4fq8kzswNIJP0H3f2J5OZPzKxjcn8nYFO+4suhHwDnmtn7wMPAQDN7kGi2\nxYfAh+7+ZnL9MRL/EWyMYFt8F1jg7p+6+1fAX4HvE822qJDp38S++fTw5LaM8pX4FwHdzKyLmTUn\n8WBidp5iyTkzM+AB4G13vzNt12zgsuTvlwFP7HtsU+PuN7j7Ee7elcTDu1fc/RKi2RYbgfVm1j25\n6TRgJfAUEWsLEg+1+5pZi+S/l9NIPPyPYltUyPRvYjbwEzNrbmZdgW4kXpjNzN3zsgBnkniz911g\nbL7iyNN3P4lEf/ZbwJLkMhhoB7wE/AN4AWib71hz3C6nALOTv0eyLYBvA28CS0nc5R4c4ba4nsR/\nfMtJPMw8ICptQeKv34+AL0k8Dx2e7bsDNyRz6TvAGdWdXy9wiYhEjEoviohEjBK/iEjEKPGLiESM\nEr+ISMQo8YuIRIwSv4hIxCjxS2SYWfvk1M9LzOxjM/sw+ft2M7s33/GJ5IrG8UskmdktwHZ3/0O+\nYxHJNd3xS5QZgJnF0wrAlJjZNDObZ2Zrzew8M7vDzJaZ2XPJKYIxs++YWZmZLTKzORVzqIg0Bkr8\nIvvrCgwLbOsFAAAA1klEQVQgMf/5QySKX5wA7ALOTk6wdw9wvrt/l0TRjN/nK1iRmgqt2LpII+XA\nc+7+tZmtAArc/fnkvuUkpgnuDhwHvJSYP4xmJOZVEWkUlPhF9vclgLvvMbPdadv3kPg3Y8BKd/9B\nPoITqSt19YhUVmUxmH2sBg41s76QqK2QrHsq0igo8UuUedrPqn6H/SsZuSfqRF8AjDeziqm1vx9m\noCL1ScM5RUQiRnf8IiIRo8QvIhIxSvwiIhGjxC8iEjFK/CIiEaPELyISMUr8IiIRo8QvIhIx/x91\nHtoUu0bmuAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f7fe1b99350>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "run_simulation(100, y0=y0)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "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.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}