{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Enzyme Kinetics (Part 2)\n",
    "\n",
    "We have already implemented a model in which an enzyme E aids a reaction in which a substrate S is consumed to produce a product P.\n",
    "\n",
    "\\begin{eqnarray*}\n",
    "\\frac{dS}{dt} &=& -k_1\\cdot S\\cdot E + k_{-1}\\cdot SE\\\\\n",
    "\\frac{dE}{dt} &=& -k_1S\\cdot E + k_{-1}SE + k_2\\cdot SE\\\\\n",
    "\\frac{dP}{dt} &=& k_2\\cdot SE\\\\\n",
    "\\frac{dSE}{dt} &=& k_1\\cdot S\\cdot E - k_{-1}\\cdot SE - k_2\\cdot SE\\\\\n",
    "\\end{eqnarray*}\n",
    "\n",
    "where SE is the compound of S and E.\n",
    "\n",
    "Our goal today is to rewrite the model with ODEs for S and P, but not E or ES.\n",
    "\n",
    "First, we rerun the model from Tuesday."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import scipy.integrate\n",
    "\n",
    "def enzyme_model( t, y, params ):\n",
    "    S = y[0]\n",
    "    E = y[1]\n",
    "    P = y[2]\n",
    "    SE = y[3]\n",
    "\n",
    "    k1f = params[0]\n",
    "    k1b = params[1]\n",
    "    k2 = params[2]\n",
    "\n",
    "    dydt = np.zeros( (4,) )\n",
    "    dydt[0] = -k1f*S*E + k1b*SE # S\n",
    "    dydt[1] = -k1f*S*E + k1b*SE + k2*SE #E\n",
    "    dydt[2] = k2*SE #P\n",
    "    dydt[3] = k1f*S*E - k1b*SE - k2*SE #SE\n",
    "    return dydt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Simulating when Enzyme concentration is low\n",
    "A more realistic scenario is one in which the enzyme concentration is low.\n",
    "\n",
    "Use k1f = 2, k1b=2, k2 = 3, S(0) = 100, E(0) = 1, SE(0) = 0, P(0) = 0.\n",
    "\n",
    "Simulate with from t=0 to t=100 with a timestep of 0.1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "k1f = 2\n",
    "k1b = 1\n",
    "k2 = 3\n",
    "params = (k1f, k1b, k2)\n",
    "\n",
    "S0 = 100\n",
    "E0 = 1\n",
    "SE0 = 0\n",
    "P0 = 0\n",
    "yinit = (S0,E0,P0,SE0)\n",
    "\n",
    "em = lambda t,y: enzyme_model(t,y,params)\n",
    "sol = scipy.integrate.solve_ivp( fun=em, t_span=(0,100), y0=yinit, \n",
    "                                 t_eval=np.arange(0,100,1) )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plot the simulation over time."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEWCAYAAACJ0YulAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy86wFpkAAAACXBIWXMAAAsTAAALEwEAmpwYAAA5yElEQVR4nO3dd3hUZdrH8e+dSSdACAmd0EJoEgNEFBSNUhS7uyhrWWBXwVVQ7G33te1aV5SmWFDBsljWhgUsCGtBUUqQ3luooSSUkP68f5yTOMYkTJKZOcnM/bmuuTJz6u/MwNxz2vOIMQallFIKIMTpAEoppeoOLQpKKaXKaFFQSilVRouCUkqpMloUlFJKldGioJRSqowWBaW8SERmiMi//Li+dBHJ9Nf6PCUiiSJyVERcTmdR1aNFQVVIRLaKyHH7P3bpY6rTufzNfh8GOZ2jriv/PhljthtjYowxxU7mUtUX6nQAVaddZIz5yukQyjtEJNQYU+R0DlW36Z6CqjYRGSUi34nIUyJySES2iMhQe1y/cnsXefavyBYikisiTd2W01tEskQkzF7m9yLyjIhki8hmEelvD98hIvtEZKTbvBH2+reLyF4ReV5EoqrI3FZE3rfXd6B0r0dEOonI1/aw/SLypojE2uNeBxKBj+1tucse/q6I7BGRHBH5RkR6VLLOlSJykdvrMHsdvURkarn3qUhEHrSnayUi79lZt4jIzW7LiLIPUR0SkdXAKSf4rIyIjBWRDcAGe9iFIpJhv88LRSTFbfp7RGSTiBwRkdUiclm55Y0WkTVu43tX9D6JSHt73aFu2zRbRA6KyEYRGe22zAdF5B0Rec1e7ioRSatqu5QPGWP0oY/fPYCtwKBKxo0CCoHRgAu4AdgFSLnpwoD/AY/Zrz8DbnAb/wwwxW2ZRcBf7GX+C9gOPAtEAEOAI0CM27yzgTigIfBx6XoqyOsCltvzNAAigTPscUnAYHsdCcA3wMSq3gfgr/Y6I4CJQIbbuBnAv+zndwFvu427BFhRQb5UIAvohfVDbQlwPxAOdAQ2A+fa0z4OfGtvd1tgJZBZxedogC/t6aPsdewDTrXfl5H2NkbY018OtLJzDAeOAS3dxu3EKkRiv3ftKnqfgPb2ukPt198Az9nvfen2nmOPexDIA863Mz0G/Oj0/4FgfTgeQB9182H/Jz8KZLs9RtvjRgEb3aaNtr8AWpRbxjTgEyDEfj0c+N5+7gL2AH3dlrnBbd6e9jKbuw07YH+hiP1l1cltXD9gSyXb0s/+Egr1YLsvBZaVex8qLI72+Fg7Z2P79Qx+LQqtsApZI/v1f4G7ys2fYK/jT/brU4Ht5aa5F3jVfr4ZOM9t3BhOXBTOKfeZ/LPcNOuAsyqZPwO4xH7+OTC+in8vFRYFrOJVDDR0G/8YMMN+/iDwldu47sBxp/8PBOtDzymoqlxqKj+nsKf0iTEmV0QAYkqHicj1QDpwqjGmxB78EfC8iHQAugA5xpif3Ja51+35cXvZ5YfFYH2RRgNL7PWCVShc9rrnAAPs4ddj7dVsMxUcTxeR5sAke/qGWL+QD1WyzdhX0zyC9as5ASjdtnggx31aY8wuEfke+KOIfAAMBca7LSsMq1D8xxjzlj24HdBKRLLdFuXC2jsAq9DscBu3rbKsbtynbweMFJGb3IaF28tFREYAt2F9qYP1fsfbz9sCmzxYX3mtgIPGmCNuw7YB7oeI9rg9zwUi9RyIM7QoKK8TkQHAP7EO0RwuHW6MyRORd4BrgK7A6zVcxX6sAtHDGLOz/EhjzNByefoBiZV8yTyK9Yu2pzHmoIhcCrhfZVW+GeGrsA4DDcL6ddwYq4gIFZsJXIf1f+2HcnmnAIeBf7gN24G1x9O5kuXtxvpyXmW/TqxkOnfu27ADeMQY80j5iUSkHfASMNDOWiwiGfy6bTuATh6so7xdQJyINHQrDIlYh6JUHaMnmpVXiUhb4B1ghDFmfQWTvIZ1qOhialgU7D2Pl4BnRKSZvd7WInJuJbP8hPVl+riINBCRSBE53R7XEOswWY6ItAbuLDfvXqzj+rhNn491KCsaq6hU5UOgN9YewmulA+09qbOAq932pEqzHhGRu+2Tyi4ROUlESk8ovwPcKyJNRKQN4P6L3xMvAX8TkVPF0kBELhCRhljnWwzWoTZE5C/ASW7zTgfuEJE+9rxJdiGB379PZYwxO4CFwGP2e58CXAu8Uc3syg+0KKiqlF5NUvr4wIN5BgLNgf+6zVf6qxZjzPdYh1yWGmM8OfRRmbuBjcCPInIY+ArrkNTvGOta+YuwToxuBzKxzm8APIT1pZ0DfAq8X272x4B/2Ffq3IH1xb4N61fuauDHqkIaY44D7wEdyi37Sqwv0V1u79N9dtYLsc6dbMHaK5qOtUdSmnebPe4LqllYjTGLsS4QmIq1h7MRq0hjjFkNTAB+wPqS7wl87zbvu1iHzv6Dda7kQ6wT2PD796m8K7EOSe0CPgAeqOLQpHKQGKOd7Cj/EpGvsY6jT3c6iz+IyP1AsjHmGqezKHUiek5B+ZV9GKQ31nH5gCcicViHSv7sdBalPKGHj5TfiMhMrMM8t5S7EiUg2Tdo7QDmGGO+cTqPUp7Qw0dKKaXK6J6CUkqpMvX6nEJ8fLxp37690zGUUqpeWbJkyX5jTEJF4+p1UWjfvj2LFy92OoZSStUrIlLp5eB6+EgppVQZLQpKKaXKaFFQSilVRouCUkqpMloUlFJKlfFZURCRV8TqQnGl27A4EflSRDbYf5vYw0VEJtvd9P0iIr19lUsppVTlfLmnMAM4r9ywe4B5dlvx8+zXYHU+0tl+jMHqHUoppZSf+ew+BWPMNyLSvtzgS7B64wKr85EFWE0gXwK8Zqw2N34UkVgRaWmM2e2LbMt3ZDN31R7uPq+rLxavHJBbmMu2w9s4lH+I7LxsjhQcobCkkIKSAopKiigxJRjsLgcr6A9Gm3tR9U1623ROij/pxBNWk79vXmvu9kW/B6vdfYDW/LbLwEx72O+KgoiMwdqbIDHRk06nfu+XzGymLdjE+Se1pGebxieeQdVJP+/5mS+3fUnGvgzWHVpHyW/6qqk+qbTzNKXqnmbRzQKiKJQxxhgRqfbPM2PMi8CLAGlpaTX6eXdJr9Y88tka/vPTdh5r07Mmi1AOytiXwZRlU/hpz09EhUaREp/C6J6j6RLXhbjIOJpENqFReCPCQsIId4UTKqGICCESgiC49euslCrH30Vhb+lhIRFpCeyzh+/E6ne2VBt82H9ro8gwLkppxeyMnfz9gm7ERNTr1j6CxvGi4zyw8AHmbJlDXGQc9/S9h2HJw4hwRTgdTamA4e9LUmcDI+3nI4GP3IaPsK9COg3I8dX5hFJXnprIsYJiPl6+y5erUV5y4PgBrvv8OuZumcuNJ9/InD/M4epuV2tBUMrLfPYTWURmYZ1UjheRTOAB4HHgHRG5Fquf2SvsyT8DzsfqLzYX+IuvcpXq1TaWri0aMuun7VzZt2bnJpR/bM3Zyg1f3UDW8SyeSX+Gge0GOh1JqYDly6uPrqxk1O/+R9tXHY31VZaKiAhX9k3kgdmrWLkzh5Na6wnnumj30d385fO/UFxSzMvnvszJCSc7HUmpgBbUdzRf2qs1EaEhzPppu9NRVAWOFhzlxnk3kl+Uz6vnvaoFQSk/COqi0DgqjAtTWvFRxi6O5BU6HUe5KSop4o5v7mBrzlYmpE+gU2wnpyMpFRSCuigAjOrfnqP5RUz/dovTUZSbJ39+ku93fs/fT/s7/Vr1czqOUkEj6ItCzzaNGXpSC6Z/u5kDR/OdjqOAuVvnMmvtLEZ0H8Gw5GFOx1EqqAR9UQC4fUgXjhcW8+z8TU5HCXo7j+7k4YUPk5KQwi19bnE6jlJBR4sCkNQshmF92vDGj9vIPJTrdJygVVhSyF3f3IXB8OSZTxIWEuZ0JKWCjhYF2/hBySAw6asNTkcJWtMypvFL1i880O8BWse0djqOUkFJi4KtdWwUfz6tHe8tzWTD3iNOxwk6GfsymL5iOpclXcZ5Hcq3uK6U8hctCm7Gnp1EdHgoE75Y73SUoJJXlMf/ff9/tGzQkrv73u10HKWCmhYFN3ENwhk9oCNzV+0hY0e203GCxrMZz7L18FYe7P8gDcIaOB1HqaCmRaGcawd0oGmDcJ6cu9bpKEEhY18Gr61+jWHJw/R+BKXqAC0K5cREhDL27CQWbjrAdxv2Ox0noOUX53P/wvtpFt2M2/vc7nQcpRRaFCp09WmJtI6N4om5a7WbRh+avmI6W3K28GC/B4kJj3E6jlIKLQoVigh1cevgZFbszOGzFXucjhOQNmdvZvqK6Zzf4XxOb32603GUUjYtCpW4rFdrOjeLYcKX6ygqrl3fv+q3SkwJD/3wENGh0dx1yl1Ox1FKudGiUAlXiHDnuV3YnHWM95ZmOh0noHy48UOW7lvK7Wm30zSqqdNxlFJutChUYXD35vRKjGXiVxvIKyx2Ok5AOHD8ABMWT6B3s95cmnSp03GUUuVoUaiCiHDXuV3ZnZPH6z9sczpOQHhmyTPkFubyQL8HCBH956dUXaP/K0+gX6emnJmcwLMLNnJYO+KplWX7lvHRpo8Y0WMEHWM7Oh1HKVUBLQoeuOvcLmTnFvLSN5udjlJvFZUU8a8f/0WLBi24PuV6p+MopSqhRcEDJ7VuzAUpLXn5uy1kHdGOeGpi1tpZrD+0nrtPuZvosGin4yilKqFFwUO3D04mv6iEZ+dvdDpKvZOVm8WzGc9yeuvTGZg40Ok4SqkqaFHwUMeEGK5Ia8ubi7ax46B2xFMdTy95moLiAu7rex8i4nQcpVQVtChUw/iBnQkR4ZmvtGltTy3bt4xPNn/CqB6jSGyU6HQcpdQJaFGohhaNIxnVvz0fLNvJ2j2HnY5T5xWXFPPookdpHt2c63pe53QcpZQHtChU0w3pnYiJCOWpz3Vv4UTe2/Aeaw+u5Y60O/TkslL1hBaFaoqNDudvZ3XiqzV7WbLtoNNx6qyc/BwmL5vMKS1O4dz25zodRynlIS0KNfCX09sTHxPBE3PWadPalZi6bCpHCo5w9yl368llpeoRLQo1EB0eyviBSfy09SAL1mc5HafOWX9oPe+sf4fLky+nS1wXp+MopapBi0INDT8lkcS4aJ6cu46SEt1bKGWM4cmfniQmLIZxqeOcjqOUqiYtCjUUHhrC7UOSWbP7MB//ssvpOHXG1zu+ZtGeRYxNHUtsZKzTcZRS1eRIURCRW0VklYisFJFZIhIpIh1EZJGIbBSRt0Uk3Ils1XFRSiu6tWzEhC/WU1CkHfHkF+fz1M9PkRSbxBVdrnA6jlKqBvxeFESkNXAzkGaMOQlwAX8CngCeMcYkAYeAa/2drbpCQoS7zu3C9oO5vL14h9NxHPfG6jfIPJrJnafcSWhIqNNxlFI14NTho1AgSkRCgWhgN3AO8F97/EzgUmeiVU96lwT6to9j8rwN5BYUOR3HMfuP7+elFS+R3jad/q36Ox1HKVVDfi8KxpidwFPAdqxikAMsAbKNMaXfqplA64rmF5ExIrJYRBZnZTl/5Y+IcPfQLmQdyefV77c6HccxU5dNJb84nzvS7nA6ilKqFpw4fNQEuAToALQCGgDneTq/MeZFY0yaMSYtISHBRymrp0+7OAZ1a87z/9tEdm6B03H8bs2BNby/4X2u6noV7Rq1czqOUqoWnDh8NAjYYozJMsYUAu8DpwOx9uEkgDbATgey1did53bhaH4R0xZscjqKXxljePLnJ2kc0ZjrT9bOc5Sq75woCtuB00QkWqxbXQcCq4H5wDB7mpHARw5kq7EuLRpyWa/WzFi4ld05x52O4zfzts9j8d7FjEsdR6PwRk7HUUrVkhPnFBZhnVBeCqywM7wI3A3cJiIbgabAy/7OVlu3DkqmxBgmz9vgdBS/KCguYMLiCSTFJvHH5D86HUcp5QWOXDdojHkAeKDc4M1AXwfieE3buGiuPrUdr/+4jesGdKRTQozTkXxq1tpZZB7N5IVBL+glqEoFCL2j2cvGnZNERGgIT38R2E1rH8w7yAvLX+CM1mfQv7VegqpUoNCi4GXxMRFcN6Ajn67YzS+Z2U7H8ZlpGdPILcrVS1CVCjBaFHxg9IAOxDUI54m5a52O4hObszfz7vp3GZY8jE6xnZyOo5TyIi0KPtAwMoyxZyfx/cYDfLdhv9NxvG7CkglEhUZxY+qNTkdRSnmZFgUfufrURFrHRvHk52sDqiOeH3b9wDeZ3zAmZQxxkXFOx1FKeZkWBR+JDHNxy6DO/JKZw5yVe5yO4xXFJcU8tfgpWse05qpuVzkdRynlA1oUfOgPvdvQuVkMT32+jqLi+t+09uxNs1l/aD239L6FCFeE03GUUj6gRcGHXCHCned2YfP+Y7y7JNPpOLWSW5jLlGVTSElI4dz25zodRynlI1oUfGxw9+b0Toxl4lfrySssdjpOjc1YNYOs41ncmXYnVuskSqlApEXBx0SEu8/ryt7D+cxYuNXpODWyL3cfM1bNYEi7IaQ2S3U6jlLKh7Qo+MGpHZuS3iWB5+ZvJCe30Ok41TZ12VQKSwq5pc8tTkdRSvmYFgU/uevcrhzOK+L5b+pX09rrDq7jw40fclXXq2jbsK3TcZRSPqZFwU+6t2rEJamtePX7Lew9nOd0HI89veRpGoY3ZEzKGKejKKX8QIuCH90+uAtFxYZJ9aRp7e93fs/CXQu5PuV6Gkc0djqOUsoPtCj4UWLTaK46NZG3f97Blv3HnI5TpeKSYiYsmUCbmDZc2fVKp+MopfxEi4Kf3XROZyJCQ3jqi3VOR6nS7E2z2XBoA7f0uYUwV5jTcZRSfqJFwc8SGkZw7Rkd+PSX3azIzHE6ToVyC3OZumwqKQkpDGk3xOk4Sik/0qLggDFndqRJdBhPfl43m9Z+bfVr7Du+jzvS7tAb1ZQKMloUHFDatPa3G/bz/ca61bT2/uP7eXXlqwxKHESvZr2cjqOU8jMtCg655rR2tGocyRNz61bT2tMyplFQXKA3qikVpLQoOCQyzMWtg5P5JTOHz1bUjaa1N+ds5r0N73F5l8tp16id03GUUg4I9WQiEUkARgPt3ecxxvzVN7GCwx96t+Glbzfz1BfrGNKjOWEuZ2v0xCUTiQqN4m8n/83RHEop53hUFICPgG+Br4D629RnHWM1rd2V0a8t5t3FmVx1aqJjWZbsXcL8HfMZ33u89qimlK2wsJDMzEzy8upPKwTuIiMjadOmDWFhnl9W7mlRiDbG3F2zWKoqg7o1o0+7Jkz8aj2X9WpNVLjL7xmMMUxYPIHm0c25pts1fl+/UnVVZmYmDRs2pH379vXuSjxjDAcOHCAzM5MOHTp4PJ+nxys+EZHzaxZNVaW0ae19R/J5deEWRzJ8vu1zVuxfwbhe44gMjXQkg1J1UV5eHk2bNq13BQGs75amTZtWey/H06IwHqsw5InIEftxuNopVYX6dojjnK7NmLZgE9m5BX5dd2FxIZOWTCK5STIXdbzIr+tWqj6ojwWhVE2ye1QUjDENjTEhxphI+3lDY0yjaq9NVequ87pwNL+IaQv827T2O+vfIfNoJrf1uQ1XiP8PXSmlTuyRRx6hR48epKSkkJqayqJFi3y2Lk/PKSAiFwNn2i8XGGM+8U2k4NS1RSMuS23NjIVbGXV6e1o2jvL5Oo8UHOH55c9zWsvT6N+qv8/Xp5Sqvh9++IFPPvmEpUuXEhERwf79+yko8N0RBY/2FETkcaxDSKvtx3gRecxnqYLUrYOTMQYmfumfprVfXvEy2fnZ3Nbntnq9i6xUINu9ezfx8fFEREQAEB8fT6tWrXy2Pk/3FM4HUo0xJQAiMhNYBtzrq2DBqG1cNFeflsjMhVsZfWYHkpo19Nm69hzbwxtr3uDCjhfSrWk3n61HqUDx0MerWL3Lu6dSu7dqxAMX9ahymiFDhvDwww+TnJzMoEGDGD58OGeddZZXc7irzt1SsW7Pa9XjiojEish/RWStiKwRkX4iEiciX4rIBvtvk9qso74ad3YS0eGh/Ptz3zatPXXZVEpMCTf1usmn61FK1U5MTAxLlizhxRdfJCEhgeHDhzNjxgyfrc/TPYXHgGUiMh8QrHML99RivZOAucaYYSISDkQD9wHzjDGPi8g99vKD7t6IpjERjB7QkWe+Ws/S7Yfonej92rju4Dpmb5rNiO4jaBXju91QpQLJiX7R+5LL5SI9PZ309HR69uzJzJkzGTVqlE/W5enVR7OA04D3gfeAfsaYt2uyQhFpjFVUXraXXWCMyQYuAWbak80ELq3J8gPBdQM6EB8TzhNzfNNY3sSlE4kJj2F0ymivL1sp5V3r1q1jw4ZfzzNmZGTQrp3v2iarsiiISFf7b2+gJZBpP1rZw2qiA5AFvCoiy0Rkuog0AJobY3bb0+wBmleSaYyILBaRxVlZWTWMULc1iAjl5oGdWbTlIAvWe3cbF+1exHc7v2N0z9Ha77JS9cDRo0cZOXIk3bt3JyUlhdWrV/Pggw/6bH0nOnx0GzAGmFDBOAOcU8N19gZuMsYsEpFJlDsUZYwxIlLhT2RjzIvAiwBpaWl1p81pL/vTKYlM/3YLT8xZy1mdEwgJqf3VQSWmhKeXPE3LBi25qttVXkiplPK1Pn36sHDhQr+tr8o9BWPMGPvpUGPM2e4PrCuSaiITyDTGlN598V+sIrFXRFoC2H/31XD5ASE8NITbhySzds8RPlq+0yvLnLtlLqsPrGZcr3FEuCK8skylVGDx9OqjispUjUqXMWYPsENEutiDBmLd+zAbGGkPG4nVMmtQuyilFT1aNWLCF+vJL6pd47QFxQVMXjaZ5CbJXNDhAi8lVEoFmhOdU2ghIn2AKBHpJSK97Uc61hVDNXUT8KaI/AKkAo8CjwODRWQDMMh+HdRCQqzG8jIPHec/i7bXalnvrHuHnUd3cmufW7U5C6VUpU50TuFcYBTQBnjabfgRrEtIa8QYkwGkVTBqYE2XGagGdI6nf6emTPl6I8P6tKFhpOftopc6UnCEF355gVNbnMrprU73QUqlVKA40TmFmfb5g1HlzilcbIx5308Zg1pp09oHjxXw0jeba7SMV1e+SnZ+Nrem3arNWSilquTRzWvGmPdE5AKgBxDpNvxhXwVTvzq5bSwXpLRk+ndbuKZfO5o19LzPg73H9vL66tcZ2n4oPZo6d/ONUqp+8LRBvOeB4VjnAgS4HNCe3f3ojiFdKCgqYcq8jdWab9ryaRSZIm7qrc1ZKFUfuVwuUlNTyx6PP+7b062eNnPR3xiTIiK/GGMeEpEJwBxfBlO/1SG+AX/q25ZZP23n2jM60D6+wQnn2Zy9mQ82fsCVXa+kbcO2fkiplPK2qKgoMjIy/LY+Ty9JLe3PLVdEWgGFWHc4Kz+6eWBnwlwhPPWFZ43lTVw6kajQKMakjDnxxEophed7Ch+LSCzwb2Ap1t3ML/kqlKpYs4aRXDegA1O+3sj1Z+bQs03lzVQs27eM+Tvmc1Ovm4iLjPNjSqUC1Jx7YM8K7y6zRU8YWvXhoOPHj5Oamlr2+t5772X48OHezeHmhEVBREKwWi/NBt4TkU+ASGNMjs9SqUqNObMjb/y4jSfmruWN606tcBpjDE8vfpqEqASu6XaNnxMqpbzJ34ePTlgUjDElIvIs0Mt+nQ/k+zqYqljDyDBuOqczD3+ymm83ZDGgc8Lvpvl6x9dkZGVwf7/7iQ6rzT2GSqkyJ/hFHyg8PacwT0T+KHqRe51w9WmJtGkSxeNz1lJS8ts2AYtKipi0dBLtG7XnsqTLHEqolKqvPC0K1wPvAvkiclhEjoiId/ulUx6LCHVx+5BkVu06zMe/7PrNuI82fsSWnC3c0vsWQkM8PWWklKqrSs8plD7uuac2/ZudmKc3r/mus2BVI5ec3JoXv9nChC/WM/SkloSHhnC86DjPZTzHyQknc05iTVo1V0rVNcXFtWsMs7o8vXltnifDlP9YjeV1YfvBXP6zaBsAb655k33H93FrH23OQilVM1XuKYhIJFZrqPEi0gTrbmaARkBrH2dTJ3BWcgL9OjZl8tcbGXRSDC+veJn0Nun0ad7H6WhKqXrqRHsK1wNLgK7239LHR8BU30ZTJyIi3DPUaizvji+fIbcol/G9xzsdSylVj1W5p2CMmQRMEpGbjDFT/JRJVcPJbWMZeFIYi45+xtD2F5LUJMnpSEqpeszTE81TRKQ/0N59HmPMaz7KpaqhQYt5sEMoOTjE6ShKqXrO0xPNrwNPAWcAp9iPijrJUX627uA6FuycQ+fIoXy0+Bhb9h9zOpJSqh7z9EL2NKC7McaccErlV5OWTiImPIanhtzCRWuW8tTn63j26t5Ox1JKeYnL5aJnz54UFRXRrVs3Zs6cSXS071oq8PTmtZVAC5+lUDXy856f+Xbnt4zuOZpOTZtx3YCOfLpiNxk7sp2OppTyktK2j1auXEl4eDjPP/+8T9fnaVGIB1aLyOciMrv04ctgqmrGGJ5Z8gzNo5tzZdcrAauxvKYNwnl8zhp0p06pwDNgwAA2bqxeR1vV5enhowd9GUJV31fbv2LF/hU83P9hIkOt7jljIkK5eWBnHpi9igXrsji7azOHUyoVOJ746QnWHlzr1WV2jevK3X3v9mjaoqIi5syZw3nnnefVDOV5tKdgjPkfsBUIs5//jNWvgnJAUUkRk5dOplPjTlzc6eLfjLuybyLtmkbzxNy1FJfo3oJS9V1p20dpaWkkJiZy7bXX+nR9Hu0piMhoYAwQB3TCupv5eWCg76Kpynyw8QO2Ht7K5LMn4wpx/WZceGgIdwzpwk2zlvHBsp0M69PGoZRKBRZPf9F7W13tjnMscDpwGMAYswHQYxMOyC3MZVrGNHo160V62/QKp7mgZ0t6tm7M01+sI6/Qv41pKaXqN0+LQr4xpqD0hYiEYnXJqfzszTVvknU8q8pG70JChHuHdmVXTh6v/bDVvwGVUvWap0XhfyJyHxAlIoOx+lb42HexVEWy87J5ZeUrpLdNp1ezXlVO2z8pnrOSE3h2/iZycgv9lFAp5W1Hjx716/o8LQr3AFnACqxG8j4D/uGrUKpiL614yWr0rpdnjd7dfV5XDucV8twC317CppQKHJ4WhSjgFWPM5caYYcAr9jDlJ7uO7mLW2llc3Olijxu9696qEZeltubVhVvZmX3cxwmVUoHA4z6a+W0RiAK+8n4cVZlnM55FEMamjq3WfLcNSQYDT3+x3kfJlFKBxNOiEGmMKTuwZT/3XeMb6jfWH1rPx5s+5qpuV9GiQfVaG2nTJJpRp7fn/WWZrNmt3WorparmaVE4JiJlrayJSB9Aj0f4yeSlk4kJi+G6ntfVaP4b0zvRMCKUJ+Z6925MpVTg8bQo3AK8KyLfish3wNvAuNqsWERcIrJMRD6xX3cQkUUislFE3haR8NosP1As2buE/2X+j7/2/CuNIxrXaBmx0eGMPTuJBeuyWLhxv5cTKqUCiafNXPyM1SXnDcDfgG7GmCW1XPd4YI3b6yeAZ4wxScAhwLf3ctcDpY3eJUQlcHW3q2u1rJH929OqcSSPzVlLiTZ/oVS98sgjj9CjRw9SUlJITU1l0aJFpKen06VLF1JTU0lNTWXYsGFeWZenDeKB1bFOe3ue3iJS457XRKQNcAHwCHCbWHdhnQNcZU8yE6sRvmk1WX6gmL9jPsuzlnN/v/uJCq3dxV6RYS5uH9KF299dzicrdnPxya28lFIp5Us//PADn3zyCUuXLiUiIoL9+/dTUGDdS/zmm2+Slubd/s6c6nltInAXUGK/bgpkG2OK7NeZWO0rVZRljIgsFpHFWVlZtYhQtxWXFDN56WTaN2rPZUmXeWWZl/ZqTdcWDfn352vJL9LmL5SqD3bv3k18fDwREREAxMfH06qV737U+b3nNRG5ENhnjFkiIunVnd8Y8yLwIkBaWlrAHgeZvWk2m3I28XT604SGVGeHrnKuEOHe87sx8pWfePPH7fz1jA5eWa5SwWDPo4+Sv8a7F2tEdOtKi/vuq3KaIUOG8PDDD5OcnMygQYMYPnw4Z511FgBXX301UVHWUYTBgwfz73//u9aZPP22Ke15bXet12g1rHexiJwPRAKNgElArIiE2nsLbYCdXlhXvZRXlMezGc/SM74ngxIHeXXZZ3aO54ykeKZ8vYFhaW1oFBnm1eUrpbwrJiaGJUuW8O233zJ//nyGDx/O448/Dvjm8JGnRaG057WfgPzSgcaYiyufpWLGmHuBewHsPYU7jDFXi8i7wDDgLWAk8FF1lx0o3lr7Fntz9/LoGY9W2uhdTYkI9wztyoVTvmPagk3cfV5Xry5fqUB1ol/0vuRyuUhPTyc9PZ2ePXsyc+ZMn62rLvW8djfwloj8C1gGvOyHddY5hwsO89KKlzi91en0bdnXJ+s4qXVjLk1txSvfbWFEv3a0bKwtlihVV61bt46QkBA6d+4MQEZGBu3atWPlypU+WV91el5bCzS0H2vsYbVijFlgjLnQfr7ZGNPXGJNkt7GUf6L5A9GrK1/lcMFhxvf2rNG7mrp9SBeMNn+hVJ139OhRRo4cSffu3UlJSWH16tU8+OCDgHVOofSS1EGDvHOo2dOe164A/g0sAASYIiJ3GmP+65UUCoB9uft4Y/UbDO0wlG5Nu/l0XW3jrOYvXvp2M9cO6EDXFo18uj6lVM306dOHhQsX/m74ggULfLI+T+9o/jtwijFmpDFmBNAX+D+fJApizy9/nqKSIm5Kvckv6xubnkSjyDAe+0ybv1BKWTwtCiHGmH1urw9UY17lga05W3l/w/sMSx5G20Zt/bLOxtFhjDs7if+tz+K7Ddr8hVLK8y/2uSLyuYiMEpFRwKdYHe0oL5maMZVwVzjXn3y9X9c7on872jSJ4tHP1mjzF0qpqouCiCSJyOnGmDuBF4AU+/ED9g1kqvZWHVjF51s/Z0T3EcRHxft13RGhLu48twurdx/mo+VBe2uIUpXywj27jqlJ9hPtKUwEDtsLf98Yc5sx5jbgA3uc8oKJSyYSGxHLqB6jHFn/RSmt6Nm6MU99vp68Qm3+QqlSkZGRHDhwoF4WBmMMBw4cIDIyslrznejqo+bGmBUVrGyFiLSv1ppUhX7Y9QM/7v6RO9PuJCY8xpEMISHCved35aqXFjFj4Vb+dlYnR3IoVde0adOGzMxM6ms7a5GRkbRp06Za85yoKMRWMU7veKolYwwTl06kZYOWDO863NEs/TvFc07XZjw7fyPD09rSpIF2Z6FUWFgYHToEVxthJzp8tFhERpcfKCLXAbXtTyHofbHtC1YfWM3Y1LFEuCKcjsM9Q7tyLL+IyV9vcDqKUsohJ9pTuAX4QESu5tcikAaEA95pzzlIFZYUMmXZFJJik7iw44VOxwEguXlDhp/Sljd+3Mao/u1p17SB05GUUn5W5Z6CMWavMaY/8BCw1X48ZIzpZ4zZ4/t4gevDjR+y7fA2bu51M64Ql9Nxytw6KJkwVwhPzl3ndBSllAM8aubCGDMfmO/jLEHjeNFxpmVMIzUhlfS26U7H+Y1mjSIZc2ZHJn61gb9uO0Sfdk2cjqSU8iO9K9kBb655k6zjWdzS5xavN43tDaMHdCShYQSPfLq6Xl6Kp5SqOS0KfpaTn8MrK1/hzDZn0qd5H6fjVKhBRCi3D05m6fZs5q7Uo4RKBRMtCn728sqXOVpwlJt73ex0lCpdntaWLs0b8vjctRQUlZx4BqVUQNCi4Ed7j+3lP2v+wwUdL6BLXBen41TJZd/Qtu1ALq//uM3pOEopP9Gi4EfTlk+j2BQzNnWs01E8clZyAgM6xzN53gZycgudjqOU8gMtCn6yJWcLH278kCuSr6BNw+rddu4UEeHeod04nFfI1Pl6Q5tSwUCLgp9MXWY1jT0mZYzTUaqle6tGXN6nDTMXbmP7gVyn4yilfEyLgh+s2r+KL7Z9wcgeI2ka1dTpONV22+AuuEKEJ+ZqD21KBTotCn4wcelEmkQ0YWT3kU5HqZEWjSMZfWZHPl2xmyXbDjodRynlQ1oUfKy0aezRKaMdaxrbG64/syPNGkbwz0/W6A1tSgUwLQo+ZIxh0tJJVtPYXZxtGru2GkSEcse5XcjYkc3Hv+x2Oo5Syke0KPjQl9u+ZNWBVYxNHUu4q/73T/DH3m3o3rIRT8xZqz20KRWgtCj4SFFJEVOWTaFT4051pmns2nKFCP+4oBs7s4/z6vdbnY6jlPIBLQo+8uHGD9l6eCs3965bTWPXVv+keAZ1s3po23803+k4Sikv06LgA3lFeUxbPo2UhBTObnu203G87t7zu5FXWMzTX653OopSysu0KPjArLWz2Je7j1t731onm8aurU4JMVxzWjve+mk7a/ccdjqOUsqLtCh42eGCw0xfMZ0zWp9BWos0p+P4zC2DOtMwMoxHPtVLVJUKJFoUvGzGyhkcLjjM+N7jnY7iU7HR4dw8sDPfbtjP/HX7nI6jlPISvxcFEWkrIvNFZLWIrBKR8fbwOBH5UkQ22H/rXT+QWblZvLHmDYZ2GErXuK5Ox/G5P5/Wjg7xDfjXp2soLNY+F5QKBE7sKRQBtxtjugOnAWNFpDtwDzDPGNMZmGe/rlde+OUFCosLGZc6zukofhEeGsLfz+/G5qxjvKF9LigVEPxeFIwxu40xS+3nR4A1QGvgEmCmPdlM4FJ/Z6uN7Ye389769/hj8h9JbJTodBy/GditGWckxTPxqw0cOlbgdBylVC05ek5BRNoDvYBFQHNjTGn7CXuA5pXMM0ZEFovI4qysLP8E9cDUjKmEucK4PuV6p6P4lYjwfxd250heIRO/0ktUlarvHCsKIhIDvAfcYoz5zXWNxrqcpcJLWowxLxpj0owxaQkJCX5IemJrD65lzpY5XNPtGhKi60Ymf+rSoiFXn9qONxZtZ/3eI07HUUrVgiNFQUTCsArCm8aY9+3Be0WkpT2+JVBvLmmZtHQSjcIbMeqkUU5Hccytg5NpEO7in5+s1ktUlarHnLj6SICXgTXGmKfdRs0GSjscGAl85O9sNfHznp/5bud3XNfzOhqFN3I6jmPiGoQzflAy327Yz7w19aaeK6XKcWJP4XTgz8A5IpJhP84HHgcGi8gGYJD9uk4rbRq7WVQzrux6pdNxHDeiXzuSmsXwz09Xk1+kragqVR85cfXRd8YYMcakGGNS7cdnxpgDxpiBxpjOxphBxpg638XXgh0LWJ61nL+l/o3I0Ein4zguzBXC/Rd2Z9uBXF75bqvTcZRSNaB3NNdQcUkxk5dNpn2j9lyWdJnTceqMM5MTGNStOVO/3sC+w3lOx1FKVZMWhRr6bMtnbMzeyLhe4wgNCXU6Tp3yfxd2o7DY8PjctU5HUUpVkxaFGigoLuDZjGfp3rQ7g9sNdjpOndOuaQOuHdCB95fuZOn2Q07HUUpVgxaFGnh3/bvsPLqT8b3GEyL6FlZk3NlJNG8UwQMfraK4RC9RVaq+0G+0asotzOXFX16kb4u+9GvVz+k4dVaDiFD+fkF3VuzM4e2fdzgdRynlIS0K1fTa6tc4mHeQ8b3HB2QHOt50UUpLTu0Qx5Ofr9V2kZSqJ7QoVMOhvEPMWDWDgYkDSUlIcTpOnSciPHRJD47kFfHUF+ucjqOU8oAWhWqYvmI6x4uOc1Ovm5yOUm90bdGIEf3a8Z+ftrMiM8fpOEqpE9Ci4KE9x/bw1tq3uLjTxXSK7eR0nHrl1sHJNG0QwT8+XKEnnZWq47QoeGja8mkYDDeefKPTUeqdRpFh/OOCbizPzGHWT9udjqOUqoIWBQ9sztnMhxs/ZHiX4bSMael0nHrpktRW9OvYlCfnrmX/0Xyn4yilKqFFwQNTl00l0hXJ6JTRTkept0SEf156EscLi3n0szVOx1FKVUKLwgms3L+SL7d9yageo4iLjHM6Tr2W1CyGMWd25P2lO/lx8wGn4yilKqBF4QQmLZ1Ek4gmjOgxwukoAWHc2Z1pGxfFfe+vIK9Qm9dWqq7RolCFH3b9wI+7f2R0ymgahDVwOk5AiAp38ehlPdm8/xhTv97odBylVDlaFCphjGHy0sm0bNCSK7pc4XScgDKgcwJ/6N2a5/+3ibV7Dp94BqWU32hRqMS87fNYeWAlN5x8AxGuCKfjBJx/XNCdRlFh3POe3rugVF2iRaECRSVFTF42mY6NO3JRp4ucjhOQ4hqEc/+F3cnYkc2MhVudjqOUsmlRqMDHmz5mS84Wbu51s3ag40OXpLZiYNdmPDl3LZuyjjodRymFFoXfyS/O57nlz9EzvifnJJ7jdJyAJiI89oeeRIW7uP2d5RQVlzgdSamgp0WhnLfXvs2eY3u0aWw/adYokn9echIZO7J54ZvNTsdRKuhpUXBztOAoL614iX4t+3Fqy1OdjhM0Ljq5FRektGTiV+tZvUuvRlLKSVoU3MxcPZPs/GzG9x7vdJSg889LTiI2Opxxs5ZyLL/I6ThKBS0tCrYDxw/w2qrXGNxuMD3iezgdJ+jENQhn0p9S2br/GPd9sAJj9DJVpZygRcE2fcV08orzGNdrnNNRglb/TvHcNjiZjzJ28R9tYlspR2hRAHYd3cXb697m0qRL6di4o9NxgtqN6UmcmZzAQ7NXa09tSjlAiwLwXMZzCMINJ9/gdJSgFxIiTByeSnxMOH+d+TPbD+Q6HUmpoBL0RWFT9iY+3vwxf+r6J1o0aOF0HIV1fmHmX/tSWFzCn19ZxL4jeU5HUipoBH1RmLJsClGhUVzX8zqnoyg3nZs35JVRp7DvcD4jX/mZw3mFTkdSKigEdVH4JesX5m2fx6geo2gS2cTpOKqc3olNeP7Pfdi47wjDX/iRLfuPOR1JqYAXnEXBGExxEROXTiQuMo4R3bUDnbrqrOQEXhqRxu6c41w05Ts+W7Hb6UhKBTSpS9eDi8h5wCTABUw3xjxe1fRpaWlm8eLF1V/R6tl8/92j/C38KE8ev4Du3+/CFFfSC5gnTV1UNk0Vs0plI92XVelyK5nmN4Pl14EVTS/y6/Qiv+b5zXgpm85antuwELGGSQiEhCAhvz63xtnDXS5wuZCQEAh1IaFhiMuFhIchYfYjPAKJjCAkMpKQ6GhCYmIIaRCDKzaW0CaxSHg4O7OPM/bNpazfsYdruoVweWchKTwbObYPjh+C3IOQfwSK8qC4AIrywZSAKbb/GqDu/FtXqtYG3A49LqvRrCKyxBiTVtG4OtMEqIi4gGeBwUAm8LOIzDbGrPb2ukpCQnk95zD/+tLQfvtsijp2IKzF708ye1QwK5ukqnkrG+c+3JgKpzPuK3QfXW5eU8EyyuY1btO7T2PM76dxf2Cs98QAJSVlw01JifXFW2Ks4moMmBLreUkJpqgYSooxxcXW88oKcCVCIkIIjYbJUQWERuYTtrGYsOhicqOLCY0qwURHEtI4ltCYxkhYJBIaCeHRIC6rUEmIXeyk8kKrVH0T3tAni60zRQHoC2w0xmwGEJG3gEsArxeF7xf8zA0zDS5XCC36HCS2x0okakcFU3rw5V2pcl/wnkxX6ToqW9YJhpf/dVw2jXEbbH47vPx89hf8r9O5/61Zq6a/zi6YYigpEUyRUFJs/S0uDKGkUCgujqSoJIbiwkiK8sMoPAp5Bwop3lzRZaol5LmOUhiSR6ErjBJXKEaEkpAQrGIAxi4KlX4aWjBUPSIj9nJ2Z+8vty4VhdaA+zdzJvC7VulEZAwwBiAxMbFGK3K1a8uOns0Z/MwsIgq2w7LXoaS6h488Oazk4fSerOM301RnePlfx+6HiSp7Xm4ZUvrFKuX+uv0KF5fb4Sa3X+ghoRDi+vW5KwwRF4SGI65wcIXjCo2A0CgIi4TwBhDRCMJjrNcVMAUFFO7bR+GuXRzbvY99O/dxcHcWx7NzKM4voCS/AFNQYO3BlO7FVLTn5fbSvVToQSZVHyTExfpkuXXmnIKIDAPOM8ZcZ7/+M3CqMabSdidqfE5BKaWCWFXnFOrS1Uc7gbZur9vYw5RSSvlJXSoKPwOdRaSDiIQDfwJmO5xJKaWCSp05p2CMKRKRccDnWJekvmKMWeVwLKWUCip1pigAGGM+Az5zOodSSgWrunT4SCmllMO0KCillCqjRUEppVQZLQpKKaXK1Jmb12pCRLKAbTWcPR7Y78U49UUwbncwbjME53YH4zZD9be7nTEmoaIR9boo1IaILK7sjr5AFozbHYzbDMG53cG4zeDd7dbDR0oppcpoUVBKKVUmmIvCi04HcEgwbncwbjME53YH4zaDF7c7aM8pKKWU+r1g3lNQSilVjhYFpZRSZYKyKIjIeSKyTkQ2isg9TufxBRFpKyLzRWS1iKwSkfH28DgR+VJENth/mzid1dtExCUiy0TkE/t1BxFZZH/eb9tNswcUEYkVkf+KyFoRWSMi/YLks77V/ve9UkRmiUhkoH3eIvKKiOwTkZVuwyr8bMUy2d72X0Skd3XXF3RFQURcwLPAUKA7cKWIdHc2lU8UAbcbY7oDpwFj7e28B5hnjOkMzLNfB5rxwBq3108AzxhjkoBDwLWOpPKtScBcY0xX4GSs7Q/oz1pEWgM3A2nGmJOwmtz/E4H3ec8Azis3rLLPdijQ2X6MAaZVd2VBVxSAvsBGY8xmY0wB8BZwicOZvM4Ys9sYs9R+fgTrS6I11rbOtCebCVzqSEAfEZE2wAXAdPu1AOcA/7UnCcRtbgycCbwMYIwpMMZkE+CftS0UiBKRUCAa2E2Afd7GmG+Ag+UGV/bZXgK8Ziw/ArEi0rI66wvGotAa2OH2OtMeFrBEpD3QC1gENDfG7LZH7QGaO5XLRyYCdwEl9uumQLYxpsh+HYifdwcgC3jVPmw2XUQaEOCftTFmJ/AUsB2rGOQASwj8zxsq/2xr/f0WjEUhqIhIDPAecIsx5rD7OGNdjxww1ySLyIXAPmPMEqez+Fko0BuYZozpBRyj3KGiQPusAezj6JdgFcVWQAN+f5gl4Hn7sw3GorATaOv2uo09LOCISBhWQXjTGPO+PXhv6e6k/XefU/l84HTgYhHZinVY8BysY+2x9uEFCMzPOxPINMYssl//F6tIBPJnDTAI2GKMyTLGFALvY/0bCPTPGyr/bGv9/RaMReFnoLN9hUI41omp2Q5n8jr7WPrLwBpjzNNuo2YDI+3nI4GP/J3NV4wx9xpj2hhj2mN9rl8bY64G5gPD7MkCapsBjDF7gB0i0sUeNBBYTQB/1rbtwGkiEm3/ey/d7oD+vG2VfbazgRH2VUinATluh5k8EpR3NIvI+VjHnl3AK8aYR5xN5H0icgbwLbCCX4+v34d1XuEdIBGr2fErjDHlT2LVeyKSDtxhjLlQRDpi7TnEAcuAa4wx+Q7G8zoRScU6uR4ObAb+gvWjL6A/axF5CBiOdbXdMuA6rGPoAfN5i8gsIB2reey9wAPAh1Tw2drFcSrWYbRc4C/GmMXVWl8wFgWllFIVC8bDR0oppSqhRUEppVQZLQpKKaXKaFFQSilVRouCUkqpMloUlPKAiDQVkQz7sUdEdtrPj4rIc07nU8pb9JJUpapJRB4EjhpjnnI6i1LepnsKStWCiKS79dvwoIjMFJFvRWSbiPxBRJ4UkRUiMtdudgQR6SMi/xORJSLyeXVbsVTKl7QoKOVdnbDaXLoYeAOYb4zpCRwHLrALwxRgmDGmD/AKEHB31Kv6K/TEkyilqmGOMaZQRFZgNaMy1x6+AmgPdAFOAr60WiTAhdXss1J1ghYFpbwrH8AYUyIihebXk3YlWP/fBFhljOnnVEClqqKHj5Tyr3VAgoj0A6t5cxHp4XAmpcpoUVDKj+wuYIcBT4jIciAD6O9oKKXc6CWpSimlyuieglJKqTJaFJRSSpXRoqCUUqqMFgWllFJltCgopZQqo0VBKaVUGS0KSimlyvw/jw42qSFpPe4AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt \n",
    "\n",
    "plt.plot( sol.t, sol.y.T )\n",
    "plt.xlabel( \"Time\")\n",
    "plt.ylabel( \"Concentration\" )\n",
    "plt.legend( [\"S\",\"E\",\"P\",\"SE\"])\n",
    "plt.title( \"Enzyme-catalyzed reaction\");"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We are interested in learning how to express dP/dt in terms of S and have learned from the tutorial that the Michaelis-Menton equation can do that for us.\n",
    "\n",
    "The curve dP/dt vs S can be expressed as\n",
    "\n",
    "$v_{max}\\frac{S}{K_m+S}$\n",
    "\n",
    "where\n",
    "\n",
    "$v_{max} = E(0)\\cdot k_2$\n",
    "\n",
    "and \n",
    "\n",
    "$K_m = \\frac{k_{-1}+k_2}{k_1}$\n",
    "\n",
    "Plot the analytic curve (from the Michaelis-Mention model for enzyme kinetics) on the same plot as the approximated dP/dt (you can just replot the curve from above and the plot the analytic curve as well)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Estimate the curve dP/dt from the simulation output\n",
    "P = sol.y[2,:].T\n",
    "S = sol.y[0,:-1].T\n",
    "velocity_of_P = np.divide(np.diff(P), (sol.t[1]-sol.t[0]))\n",
    "\n",
    "# Now compute it according to the Michaelis-Mention equation\n",
    "# Write code here!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Model Enzyme Kinetics with the Michaelis-Menton Equation\n",
    "\n",
    "Implement a model that uses the M-M equation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "def MM_model( t, y, params ):\n",
    "    # Write code here\n",
    "    dydt = np.zeros( (2,) )\n",
    "    dydt[0] = # Write code here\n",
    "    dydt[1] = # Write code here\n",
    "    return dydt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using the same initial conditions as above, and computing Km and vmax according to\n",
    "\n",
    "$v_{max} = E(0)\\cdot k_2$\n",
    "\n",
    "and \n",
    "\n",
    "$K_m = \\frac{k_{-1}+k_2}{k_1}$\n",
    "\n",
    "simulate the M-M model for the same time steps as the full model. Save the solution in a variable named sol2 so that we can plot data from both sol and sol2."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Write code here"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plot the new model's simulation results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Write code here"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now compare P from one model to the other by plotting them both"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Write code here"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Awesome, the full model and the M-M model are nearly indistinguishable. When we are modeling a process using M-M kinetics, we need two parameters, $v_{max}$ and $K_m$. Above, we saw that the two parameters depend on a larger set of parameters and the initial concentration of E. Typically, we will be assigning values directly to $v_{max}$ and $K_m$, skipping the step of determining separate parameters for forward and backward reaction rates. We should develop and intuition for how $v_{max}$ and $K_m$ affect dP/dt and P(t).\n",
    "\n",
    "We saw in the tutorial that $v_{max}$ is simply the maximum value of dP/dt.\n",
    "\n",
    "$K_m$ is often referred to as the activation threshold because it identifies the concentration of S that causes dP/dt to be at half it maximum. \n",
    "\n",
    "Make an array of S values and compute the M-M equationn for a range of values of Km. The trend is that curves with larger Km are less steep."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Copy-paste from notebook with code"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Run the simulation with a small Km and compare it to a simulation with a large Km"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "# To do on your own after class."
   ]
  }
 ],
 "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.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}