05_chemical_equations.ipynb

{
 "cells": [
  {
   "cell_type": "raw",
   "metadata": {},
   "source": [
    "---\n",
    "output-file: chemical_equations.html\n",
    "title: Exploring Chemical Equations\n",
    "\n",
    "---\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "19a01dbb-9e0f-492e-979a-eb43da88dbb5",
   "metadata": {},
   "outputs": [],
   "source": [
    "#| default_exp chem"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9a9d7a16-563c-4e82-ad91-c85ca16cd57b",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Populating the interactive namespace from numpy and matplotlib\n"
     ]
    }
   ],
   "source": [
    "%pylab inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fe7208a0-54f9-45ba-82ff-517834e39b6d",
   "metadata": {},
   "outputs": [],
   "source": [
    "#| export\n",
    "def ChemSimulation(eqnstr,verbose=True,**kwargs):\n",
    "    from pyndamics3 import Simulation\n",
    "    lines=eqnstr.strip().split('\\n')\n",
    "    \n",
    "    components=set([])\n",
    "    parameters=set([])\n",
    "    for line in lines:\n",
    "        lhs=[_.strip() for _ in line.split('--')[0].strip().split('+')]\n",
    "        rhs=[_.strip() for _ in line.split('->')[-1].strip().split('+')]\n",
    "        middle=line.split('--')[1].lstrip('\\t ->').split('->')[0].rstrip('\\t ->')\n",
    "\n",
    "        components=components.union(set(lhs+rhs))\n",
    "        parameters=parameters.union([middle])\n",
    "        \n",
    "        if verbose:\n",
    "            print(lhs,middle,rhs)\n",
    "        \n",
    "        \n",
    "    components=sorted(components)\n",
    "    parameters=sorted(parameters)\n",
    "    \n",
    "    if verbose:\n",
    "        print(\"Components\",components)\n",
    "        print(\"Parameters\",parameters)\n",
    "    \n",
    "    \n",
    "    diffeqs=[]\n",
    "    for c in components:\n",
    "        \n",
    "        if c=='ϕ':\n",
    "            continue\n",
    "        \n",
    "        eqn=\"%s' = \" % c\n",
    "\n",
    "        for line in lines:\n",
    "            lhs=[_.strip() for _ in line.split('--')[0].strip().split('+')]\n",
    "            rhs=[_.strip() for _ in line.split('->')[-1].strip().split('+')]\n",
    "            middle=line.split('--')[1].lstrip('\\t ->').split('->')[0].rstrip('\\t ->')\n",
    "\n",
    "\n",
    "            sign='0'\n",
    "            if c in lhs:\n",
    "                sign='-'\n",
    "                if len(rhs)>1:\n",
    "                    sign+='%d*' % len(rhs)\n",
    "\n",
    "\n",
    "            if c in rhs:\n",
    "                sign='+'\n",
    "                if len(lhs)>1:\n",
    "                    sign+='%d*' % len(lhs)\n",
    "\n",
    "            if sign=='0':\n",
    "                continue\n",
    "                \n",
    "            if c in rhs and c in lhs:  # A -> A +B doesn't change A\n",
    "                continue\n",
    "\n",
    "            plhs='*'.join(lhs)\n",
    "\n",
    "            eqn+=f\" {sign}{middle}*{plhs}\"\n",
    "\n",
    "        if eqn==\"%s' = \" % c:\n",
    "            eqn=\"%s' = 0\" % c\n",
    "            \n",
    "        diffeqs.append(eqn)  \n",
    "        \n",
    "    if verbose:\n",
    "        print(\"diffeqs\",diffeqs)\n",
    "        \n",
    "    sim=Simulation()\n",
    "        \n",
    "    for c,d in zip(components,diffeqs):\n",
    "        if not c in kwargs:\n",
    "            c0=0\n",
    "        else:\n",
    "            c0=kwargs[c]\n",
    "            \n",
    "            \n",
    "        sim.add(d,c0)\n",
    "        \n",
    "    for p in parameters:\n",
    "        if p in kwargs:\n",
    "            sim.params(**{p:kwargs[p]})\n",
    "        \n",
    "    return sim"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6d870697-34e6-44b5-9aaf-b2f7662ffa81",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "pyndamics3  version  0.0.27\n",
      "['D'] k1 ['D', 'M']\n",
      "['M'] k2 ['M', 'P']\n",
      "['M'] k3 ['ϕ']\n",
      "['P'] k4 ['ϕ']\n",
      "Components ['D', 'M', 'P', 'ϕ']\n",
      "Parameters ['k1', 'k2', 'k3', 'k4']\n",
      "diffeqs [\"D' = 0\", \"M' =  +k1*D -k3*M\", \"P' =  +k2*M -k4*P\"]\n"
     ]
    }
   ],
   "source": [
    "sim=ChemSimulation(\n",
    "\"\"\"\n",
    "D --k1--> D+M\n",
    "M --k2--> M+P\n",
    "M --k3--> ϕ\n",
    "P --k4--> ϕ\n",
    "\"\"\",D=1,M=0,ϕ=0,P=0,k1=.01,k3=0.00577,k4=0.0001925,k2=20*0.00577)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "560c8668-0c4c-4e7a-a465-17259a2f9bf3",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "['A'] k1 ['B']\n",
      "['B'] k_1 ['A']\n",
      "['A', 'B'] k2 ['C']\n",
      "['C'] k_2 ['A', 'B']\n",
      "Components ['A', 'B', 'C']\n",
      "Parameters ['k1', 'k2', 'k_1', 'k_2']\n",
      "diffeqs [\"A' =  -k1*A +k_1*B -k2*A*B +k_2*C\", \"B' =  +k1*A -k_1*B -k2*A*B +k_2*C\", \"C' =  +2*k2*A*B -2*k_2*C\"]\n",
      "\n",
      "A'=-k1*A +k_1*B -k2*A*B +k_2*C\n",
      "B'=+k1*A -k_1*B -k2*A*B +k_2*C\n",
      "C'=+2*k2*A*B -2*k_2*C\n",
      "k1=1\n",
      "k2=3\n",
      "k_1=2\n",
      "k_2=4\n",
      "\n"
     ]
    }
   ],
   "source": [
    "sim=ChemSimulation(\"\"\"\n",
    "A   -->k1->   B\n",
    "B   -->k_1->  A\n",
    "A+B -->k2->   C\n",
    "C   -->k_2->  A+B\n",
    "\"\"\",A=1,B=1,C=0,k1=1,k_1=2,k2=3,k_2=4)\n",
    "\n",
    "print()\n",
    "\n",
    "print(sim.equations())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "92c63fc7-84d8-48fb-a279-944497be7f31",
   "metadata": {},
   "outputs": [],
   "source": [
    "sim.params()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f5c81cb7-0270-4c97-b3ce-b7e81be1adc1",
   "metadata": {},
   "outputs": [],
   "source": [
    "assert sim.equations()==\"\"\"\n",
    "A'=-k1*A +k_1*B -k2*A*B +k_2*C\n",
    "B'=+k1*A -k_1*B -k2*A*B +k_2*C\n",
    "C'=+2*k2*A*B -2*k_2*C\n",
    "k1=1\n",
    "k2=3\n",
    "k_1=2\n",
    "k_2=4\n",
    "\"\"\".lstrip()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d596514d-20b5-4060-8dd3-230241d978e6",
   "metadata": {},
   "outputs": [],
   "source": [
    "for c in sim.components:\n",
    "    c.plot=1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a4fed08a-9433-4e1f-858a-8f203b1294fc",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAukAAAHyCAYAAABWPVDOAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8/fFQqAAAACXBIWXMAAAsTAAALEwEAmpwYAABKlElEQVR4nO3deZgddZn3//fda/YNSEQIhCCbghOwCQiiQRY33JdnFBV4BAY3cHtGx5U44uj4qLjiAyjgKKMi44aoEaFFQMAE4w9lC0vCnhASks7a2/f3R51OOp3upDvdfarOOe/XdfVVdarq1LnTRcKnv33XtyKlhCRJkqTiqMu7AEmSJEnbMqRLkiRJBWNIlyRJkgrGkC5JkiQVjCFdkiRJKhhDuiRJklQwDXkXUES77757mjVrVtk/d/369YwfP77sn6vy8jrXBq9zbfA6Vz+vcW3I8zovWrRoZUppj77bDen9mDVrFgsXLiz757a2tjJv3ryyf67Ky+tcG7zOtcHrXP28xrUhz+scEcv62267iyRJklQwhnRJkiSpYAzpkiRJUsEY0iVJkqSCMaRLkiRJBWNIlyRJkgrGkC5JkiQVjCFdkiRJKhhDuiRJklQwhnRJkiSpYAzpkiRJUsEY0iVJkqSCMaRLkiRJBWNIlyRJkgrGkC5JkiQVjCFdkiRJKpiGvAuQKsGiZav5zh8f4K7H10AEk5obWLupY8t6R1c3jfV122wbaP2p1etpvu36bd+3sR0imDimMVsnmNDcQNvmDkjBhOZ62jZ1AjC+uZ7Obqivg3WbOwlgfHMD6zb37G9gfc96UwMd3YmGumB9e7ZtXFMD69s7iVRa7+iElB3bc8z4pgY6u7tL7+sqva+eDb3XN2frY0vbo7TeWfq83tsGsz6Y9/V8du9j+62tfWttG3utd3Un6utim22DWR/K+3qOXbNuI01/+v3W/Snt0p9/NGobzT9/rdVW191B940LCllbpb6vaLWtXb+JxhsXFLK2In/fKq22zvbNvLvxIU4/dj+KIlJKeddQOC0tLWnhwoVl/9zW1lbmzZtX9s+tdjsL2Gs2tpNSMLapjrWbs8A6prGejs5u6gLWburcElQlSVL1+vzrD+NtR+1T1s+MiEUppZa+2x1JV0W78raH+fFfHqa9szsL3sC4xnqe2dRBdzeQEk9v6NjmPY8NsC5Jkmrbb/7+RNlD+kAM6Sq8npHwh55aR2N9HavXt9OZEps7umjb7Ai3JEkaGa84dM+8S9jCkK5C6R3II4LlazexttSLXQRTxjYwfkzjMHvS22geN3bI79vVzyvHZ5T7fZVQW+/rXLTaivx9q7Ta6jo30V3fXMjaKvV9Raut3P9mV8v3rdJq27xhPWee8LzCjKKDIV056wnl/3jsGda3d7Fm4+gF8oEC9mD+ko9tauB/H7vfiPzl9d6D2uB1rg1e5+rnNa4Nra2tzCtQQAdDusqs90j5uk2dPNm2eUTO+6xJzdTXRb9hu6mhjv915D6F+ulYkiRpRwzpKosrb3uYS//0AA+u3LDL59hjQhO7T2g2gEuSpKpnSNeoufK2h/nx7Q/z8Or1rN4w9DaWfaeNo7O7e0RbTSRJkiqBIV0jatGy1Vx9x6PceO8KHn1m06Df19MvvtfkMRwwYyJvOGJvXrDv1FGsVJIkqbgM6RoxX7j2br5z44ODPn7faeOYMq7RdhVJkqQ+DOkatitve5iv/P5eVq5r3+mxs3Ybx7HP2d2RckmSpB0wpGuXXXnbw1x43b2saNtxOH/WpGZmTBrjiLkkSdIgGdI1ZIuWrebfrv4b961Yv8PjnjN9gjd8SpIk7QJDuoZkMH3nc2dN5aOvOMR2FkmSpF1kSNegnf39hSy4a/mA+w3nkiRJIyP3kB4RbwJeAswB/gmYCPwwpfT2XTjX3sBngZcDuwFPAD8H5qeUVo9QyTXnytse5ku/u3vAuc6fu+dE/v11hxnOJUmSRkjuIR34JFk4Xwc8Chy8KyeJiP2BW4DpwC+Ae4C5wHnAyyPi2JTS0yNScQ15/38v4ld/e7LffXtPGcN7jj/AnnNJkqQRVoSQ/kGycH4/2Yj6Dbt4nm+TBfRzU0rf6NkYEV8pfcYFwDnDK7V2LFq2ms/+6h/87dE1/e5/3Zxnc+E/H17mqiRJkmpDXd4FpJRuSCktSSmlXT1HRMwGTgaWAt/qs/szwHrgHRExfpcLrSFX3vYwb7roln4DegDnvHi2AV2SJGkU5R7SR8hLS8sFKaXu3jtSSm3AzcA44OhyFzZY9/zlOjb89cfc85frcq1j0bLVfOJnd9LfT0zP3XMiP333MXzslYeUvS5JkqRaUoR2l5FwUGl53wD7l5CNtB8I/KEsFQ3BPX+5jv2veQsH0kXHNT/hHn7MwUeeWPY6Fi1bzbt/sLDfgH7Oi2cbziVJksqkWkL65NKy/wbqrdunDHSCiDgbOBtgxowZtLa2jlRtO7Xhr1dzIF3UBTSkLh686WqeXF/eS9P6cAeX39X/k0NPf24TR49bTmvrwNMvavDWrVtX1v++lA+vc23wOlc/r3FtKOJ1rpaQvjNRWg7Y955Suhi4GKClpSXNmzevDGVl7hnfSfc1P6KObrqpY/aL3sjBR5bv86+87WEuv+vO7bYHcMHrD3P2lhHW2tpKOf/7Uj68zrXB61z9vMa1oYjXuVp60ntGyicPsH9Sn+MK5eAjT+Qfe7wcgFsmvqysrS5X3vYwn/hZPwE9DOiSJEl5qZaQfm9peeAA+w8oLQfqWc/d2D2zfu/1dRPL9pkD3SRaF3DB6wzokiRJeamWdpeeudVPjoi63jO8RMRE4FhgI3BrHsUNxtjxWTjv3ry+bJ/5pd/ds11Af870CXzxjc/36aGSJEk5qqiR9IhojIiDS08X3SKl9ACwAJgFvLfP2+YD44Hvp5TKl4CHaPyErFOnY1Mbi5atHvXP++TP7uTWB1dts60uMKBLkiQVQO4j6RHxOuB1pZfPKi1fGBGXl9ZXppQ+UlrfC7gbWEYWyHt7D3AL8PWIOKF03FHA8WRtLp8Y+epHzsr2eqYBTd0bOfXSW/nhmUePWlj+91/9gx/c9vA22wL43OsOM6BLkiQVQO4hHZgDnNZn2+zSF2SB/CPsRErpgYhoAT4LvBx4JfAE8HVgfkpp1Y7en7f7ViUOBMazmY7Obm598OlRCcyLlq7iuzcv3Wabs7hIkiQVS+4hPaV0PnD+II9dytbpFPvb/whwxkjUVW4HzHwW/B3GxSbq6+s4evZuI/4Zi5at5pwfLNpu+7+8eLYBXZIkqUByD+nKHDRzBgDj2My/vuygER9FX7RsNW/5f7fQ1b3t9iNnTfVJopIkSQVTUTeOVrWmCQCMYxOTxjaO+Om/ef2S7QJ6fcDHXmFAlyRJKhpDelE0jQdgfGziV397fERneFm0bDU33PvUNtsi4N+9UVSSJKmQDOlF0TgOgMmsY/39t3DqpbeOWFA//5d/3+Z14MOKJEmSisyQXhQr7gKgmU5+2PR5Du26h1sffHrYp72o9X7ufGztNttOfO4MA7okSVKBGdKL4uE/k8jaUBrp5JiGe4Y9w8vCpav44m/v3WZbfcA5L9l/gHdIkiSpCJzdpShmHbdltZN6Xn7Km3nuMPrFFy1bzbv7TLdoH7okSVJlcCS9KGbOpb1xCgDndbyXZWOft8un6plu8al17dtsP/EQ21wkSZIqgSG9QDbUTQTgobQn5/1o8S7fOPr//vhAv9Mt2uYiSZJUGQzpBdKWmoDsgUYdXd27dOPoomWr+f1dy7fZZpuLJElSZTGkF0hD01gAxsUmCJg6rmnI5/jWDUtIvV473aIkSVLlMaQXSFNzFtLHs4mU4LPX/GNILS9X3vYw19+z7UOLnG5RkiSp8hjSC6SrfgwA49gEQEfn4FteFi1bzSd/fuc22+xDlyRJqkyG9ALZEtJjMwCN9XWDniv9otb76e7V52IfuiRJUuUypBdIT0h/1phOAN55zKxBhexFy1Zz3d0rttnmdIuSJEmVy5BeIN11WUjv3rwegO/d9NCgetK/vGDbp4rW2eYiSZJU0QzpBdIzkn5M/J0j4j46uxNX3/HoDt/zw1uXccsD2/atn3DIDNtcJEmSKpghvUAa29cAcGTdPfyw6fMcEffx04WPDjiavmjZaj75i79vs82bRSVJkiqfIb1Amjdn0yfWBTTSydF1d9PR1T3gaPqFv7+P5M2ikiRJVceQXiBtE58DQAI6aODW7kNIwE8XbT+afvtDT/On+1dus82bRSVJkqqDIb1A2iYdCEBM2Zev7Pl/uSNlrzv7zJe+aNlqzv6vRdu815tFJUmSqochvUC66rMnjjJmMrPmzNuyvRto29gBZE8VfdNFt/DMho5t3uvNopIkSdWjIe8CtFVnw/hsZfNantm4bQj/zo0P8tO/PsrKtvbt3ufNopIkSdXFkfQC6WwYl61sWsvRs3ejPrbd319A92ZRSZKk6mNIL5Cu+lJI39zGC/aZwgmHzNjh8QFc8LrDvFlUkiSpyhjSC6S7vgnqm6C7Azo38S8v2Z+GvsPpJXUBF7zegC5JklSN7EkvmuZJsGFlNpq+73R+fPYLufqOR7l/eRur1rczbXwTB8yYyBuO2NsWF0mSpCplSC+a+sZs+dBNcNgbeMG+Uw3jkiRJNcZ2lwKZtOYeaHsye/Hzc+CR2/MtSJIkSbkwpBfIlGf+Tva8UbK+9KV/yrUeSZIk5cOQXiDPTDkUonRJ6hpg1nH5FiRJkqRcGNILZO3kg2H/E7IXx7wfZs7NtyBJkiTlwpBeNNNmZ8vxe+RbhyRJknJjSC+a5onZctPafOuQJElSbgzpRTNmUrbcbEiXJEmqVYb0omk2pEuSJNU6Q3rR2O4iSZJU8wzpRTNmcrZ88k4fZiRJklSjDOlFs+bRbLnqAbjiNQZ1SZKkGmRIL5qn7t263tXuU0clSZJqkCG9aGa/ZOt6fZNPHZUkSapBhvSi2f+l2TLq4J2/8KmjkiRJNciQXjQNzdA4HlI3zHhu3tVIkiQpB4b0Iho7NVtuWJVvHZIkScqFIb2IxpVC+sbV+dYhSZKkXBjSi2isIV2SJKmWGdKLKKVs+ahzpEuSJNUiQ3rRPHI7LLslW//jl3yYkSRJUg0ypBfN0j9lM7sAdHf6MCNJkqQaZEgvmlnHQV19tl5X78OMJEmSapAhvWhmzoVjzs3WZx/vw4wkSZJqkCG9iPZuyZZ1DfnWIUmSpFwY0ovIKRglSZJqmiG9iLaEdJ84KkmSVIsM6UU0dlq2XPOYUzBKkiTVIEN6ET29JFt2rIcrXmNQlyRJqjGG9CJ65Lat613tzpUuSZJUYwzpRTTrOCCy9foG50qXJEmqMYb0Ipo5F/Y4OFt/+RecK12SJKnGGNKLatrsbDlu93zrkCRJUtkVJqRHxN4R8b2IeDwiNkfE0oi4MCKmDvE8r4qIBRHxaERsjIgHI+KqiHjhaNU+KsaXwvn6p/KtQ5IkSWVXiJAeEfsDi4AzgNuBrwIPAucBf46I3QZ5ni8C1wBHAL8FvgbcAbwWuDki3j7y1Y+SCdOzpSFdkiSp5hTlufPfBqYD56aUvtGzMSK+AnwQuAA4Z0cniIhnAR8BlgPPTymt6LXveOB64LPAD0a8+tEwfo9sueT3sP9L7UuXJEmqIbmPpEfEbOBkYCnwrT67PwOsB94REeN3cqp9yf48t/UO6AAppRuANmCPkai5LDY9ky0fW+hc6ZIkSTUm95AOvLS0XJBS6u69I6XUBtwMjAOO3sl5lgDtwNyI2OZuy4h4MTARuG5EKi6HNY9uXXeudEmSpJpShJB+UGl53wD7S4/f5MAdnSSltAr4KDADuCsiLo6I/4iInwALgN8D/zIC9ZbHrBdvXa9vcq50SZKkGlKEnvTJpeWaAfb3bJ+ysxOllC6MiKXA94Czeu26H7i8bxtMbxFxNnA2wIwZM2htbd3Zx424devWbfncho5mXgR01TXyt8POZ+0DG+CB8tekkdf7Oqt6eZ1rg9e5+nmNa0MRr3MRQvrOlB69SdrpgRH/Cnwe+DrwTeBJ4GDgP4AfRsSclNK/9vfelNLFwMUALS0tad68ecOvfIhaW1vZ8rnd3XBLPfXdHRzxqjOgobns9Wh0bHOdVbW8zrXB61z9vMa1oYjXuQjtLj0j5ZMH2D+pz3H9ioh5wBeBX6aUPpRSejCltCGldAfweuAx4MOlG1WLr64OxpS+JUsqp5VekiRJw1eEkH5vaTlQz/kBpeVAPes9Tiktb+i7I6W0gWz+9Trg8KEWmItHboeNq7P1n57h7C6SJEk1pAghvSdUnxwR29QTEROBY4GNwK07OU9PP8hA0yz2bG/flSLLbumf2NLh093h7C6SJEk1JPeQnlJ6gGz2lVnAe/vsng+MB76fUloPEBGNEXFw6SmlvfWk2LMjYq/eOyLiFWRhfxNwy8j+CUbJrOMg6rP1qHd2F0mSpBpSlBtH30MWnr8eEScAdwNHAceTtbl8otexe5X2LyML9j1+SjYP+onA3RHxM7IbRw8ha4UJ4GMppadH9U8yUmbOhX/6Z1j8Q3j+W3ziqCRJUg3JfSQdtoymtwCXk4XzDwP7k83S8sLBBOvSg5BeCXwQuIvsZtEPkz0E6VrgZSmlr41G/aNm7yOzZcSOj5MkSVJVKcpIOimlR4AzBnHcUrZOy9h3XwdwYemr8k16drZcdkt246ij6ZIkSTWhECPpGsCGVdly1YNwxWuc4UWSJKlGGNKL7On7t653tTvDiyRJUo0wpBfZASdvXa9vdIYXSZKkGmFIL7J9joKJz8rWX3uRPemSJEk1wpBedFNnZ8vxu+VbhyRJksrGkF50DaUHqS69Kd86JEmSVDaG9CJ75PatN4ve9BVnd5EkSaoRhvQiW/on6O7K1ru7nN1FkiSpRhjSi2zWcdmsLpA9ddTZXSRJkmqCIb3IZs6F1307W58ww9ldJEmSaoQhvegOPiVbrlsOy/6cby2SJEkqC0N60T15Z7ZM3fBfr/PmUUmSpBpgSC+63jeLdrV786gkSVINMKQX3azjIOqz9aj35lFJkqQaYEgvuplz4fn/K1uf8bx8a5EkSVJZGNIrwcRnZcsnFsMVr7EvXZIkqcoZ0ivBhpVb1+1LlyRJqnqG9Epw0Cu3rtc32ZcuSZJU5QzpleDAl0N9c7Z+/Cd8qJEkSVKVM6RXgkf/At3t2fr1/25PuiRJUpUzpFeCpX+ClLL1rg570iVJkqqcIb0SzDoO6hqy9bo6e9IlSZKqnCG9EsycC/P+LVuf8Kx8a5EkSdKoM6RXiin7ZMu1jzpXuiRJUpUzpFeK1Uu3rjtXuiRJUlUzpFeK2fOA2Pp67G55VSJJkqRRZkivFDPnwl4vyNZTN/z2Y7a8SJIkVSlDeiUZO6W0kmx5kSRJqmKG9Eqy/wlb1+ubnIpRkiSpShnSK8lhb8qWUQ8v/4+sBUaSJElVx5BeSXpmeEld8NuP2pMuSZJUpQzplaR3D3pnhz3pkiRJVcqQXklmHQd1Ddl64DSMkiRJVcqQXklmzoUjz8zWnYZRkiSpahnSK03zxK3rnZtteZEkSapChvRKM2nvXi+6bXmRJEmqQob0SrPx6V4vos9rSZIkVQNDeqWZdRzUNW597Ui6JElS1TGkV5qZc2HeR0svkjePSpIkVSFDekXqddk6N8HfrsyvFEmSJI04Q3ol2u84tl66BH+90tF0SZKkKmJIr0Qz58K+x2x93eXTRyVJkqqJIb1SPefEXi+cilGSJKmaGNIrVve2L59cnEsVkiRJGnmG9Eo16ziI+q2v7UuXJEmqGob0SjVzLhzwsq2v7UuXJEmqGob0SnaAfemSJEnVyJBeyTY9s+1r+9IlSZKqgiG9ks06Duoatr62L12SJKkqGNIr2cy5cPjbt77uavfpo5IkSVXAkF7p5pyKTx+VJEmqLob0SjdzLhz8yq2vneVFkiSp4hnSq0Hfp49uWptbKZIkSRo+Q3o12Lhq29e3fMOWF0mSpApmSK8GfZ8+mrq8gVSSJKmCGdKrwcy5cNArtt227ql8apEkSdKwGdKrxbHnQV2v0fR7fwMLL8+tHEmSJO06Q3q1mDkXjjht6+vUBb/+kL3pkiRJFciQXk3+6a1sc0lTF9x8YV7VSJIkaRcZ0qtJ3znTAe79raPpkiRJFaYwIT0i9o6I70XE4xGxOSKWRsSFETF1F851XERcHRFPlM71REQsiIhX7vzdFe7Y8yB6j6Z3+3AjSZKkClOIkB4R+wOLgDOA24GvAg8C5wF/jojdhnCuTwI3Ai8Gfgt8GfgVMBWYN6KFF9HMuXDMub02JB9uJEmSVGEa8i6g5NvAdODclNI3ejZGxFeADwIXAOfs7CQR8Wbg34HrgDeklNr67G8cyaILa8ykbV/f/DWYuh+0nJ5LOZIkSRqa3EfSI2I2cDKwFPhWn92fAdYD74iI8Ts5Tx3wRWAD8La+AR0gpdQxEjUX3qzjoK73z1/JmV4kSZIqSO4hHXhpabkgpdTde0cpaN8MjAOO3sl5jgH2A64FVkfEqyLioxFxXkS8cKSLLrSZc+GVXwZi67bUBb98n0FdkiSpAhQhpB9UWt43wP4lpeWBOznPkaXlcuAO4BrgC8CFwC0R8ceI2GMYdVaWltPh4Fdtu+2pe+GyVxrUJUmSCq4IPemTS8s1A+zv2T5lJ+eZXlqeAzwEnAjcBuxLdvPoy4CrGODm0Yg4GzgbYMaMGbS2tu608JG2bt26Ef3cSWOPYw7XEnRvGVNP3R2s/Pkn+cdhHx+xz9HQjPR1VjF5nWuD17n6eY1rQxGvcxFC+s5syZc7Oa6+1/FvSin9rfT6HxHxerKR+pdExAtTSn/u++aU0sXAxQAtLS1p3rx5wy58qFpbWxnZz50He42Baz5Az7cvgD2e/gvz9h+XtcWo7Eb+OquIvM61wetc/bzGtaGI17kI7S49I+WTB9g/qc9xA1ldWj7YK6ADkFLaCPyu9LK2kmnL6XDKhX02dtufLkmSVGBDCukRURcR74+IWyNiTUR09tp3eER8OyJ21jve172l5UDvO6C0HKhnve95nhlgf0+IHzu4sqpIy+lw8CnbbrM/XZIkqbAGHdIjogn4PdmNmPsDbWwzfQgPAf8bOHWINdxQWp5cmkax92dOBI4FNgK37uQ8NwKdwAGlWvs6tLRcOsT6qsOx50HUb7utuwOu+0w+9UiSJGlAQxlJ/z/A8cB8YAZwae+dKaVnyILyy4ZSQErpAWABMAt4b5/d84HxwPdTSusheyBRRBxcekpp7/OsBH5M1jbz6d77IuKkUl1ryJ5CWntmzoVXfYVtf64Clt0CvzeoS5IkFclQbhw9Fbg5pfRZgIjo70bOh4BX70Id7wFuAb4eEScAdwNHkf1QcB/wiV7H7lXav4ws2Pf2odL7PhERLwZuJ5vd5fVAF3BW6YeJ2tTzxNFrztt2+80XZsuT5pezGkmSJA1gKCPp+7HzlpNVwLShFlEaTW8BLicL2R8ma6n5OvDClNLTgzzPitL7vwrMBM4le1jSr4HjUkpXDbW2qtNyOhz7ge2333whXH1WmYuRJElSf4Yykr6Rnc9Vvg8D37i5QymlR4AzBnHcUrbr2dhm/yqyEfUP7UodNaFnxLxnBL3HnT+BtY/CifOdnlGSJClHQxlJX0x2c2d/N2USEZPJ+r6dLqQSnDS//xH1ZbfA914GCy8vd0WSJEkqGUpIv4SsheSHETGp946ImELWqjIV+M5IFadRdtJ8OOwt229P3VnfujeUSpIk5WLQ7S4ppf+OiBPJWlJeQ2ne8YhYCDwPaAa+lVK6djQK1Sh54yUw6dnbt75Atm3t49kxkiRJKpshPcwopfQusrnQ7wL2IOsNPwK4H3hXSun9I16hRt9J8+GUr9Hvfw53/gQuPMz2F0mSpDIaUkgHSCldnlI6HJgA7A1MTCkdllK6bMSrU/m0nA7v+h3se8z2+555OGt/uehFPqFUkiSpDIYc0nuklDamlB7veciQqsDMuXDGb/rvUwdYfid89yS47BWGdUmSpFG0yyFdVeyNl/Q/80uPZbdkYf1HpxrWJUmSRsGgbxyNiAcHeWhKKe2/i/WoKE6aDwe/Cn79IXjyzv6Pueea7GvfY5xbXZIkaQQNZSS9juxG0b5fU4BZpa+mIZ5TRTZzLpxzU3ZT6eR9Bj6uZ2T9Oy+Caz7o6LokSdIwDTpQp5RmpZT26+drGnAg8FvgAeCQ0SpWOWk5HT54545bYCAbcV/4PfvWJUmShmlERr1TSvcDbwD2AnwCTrU6aT686/dZG8zO9Iyuf+lAuOyVjrBLkiQNwaB70ncmpbQpIn4PvBX4t5E6rwpm5lz45yuzwH3zhfDIQli/fODj1y/PvpbdnI2yT5sNY6fC4e/MRuglSZK0nREL6SWdwLNG+Jwqop6wDtmDjm79Nqy8d+fvW1W6//ixRXD952DcVGhohvomg7skSVLJiIX0iNgdeD3wyEidUxWi5fTsa7Cj6z02PJV99XhsEVz/71lor2uAMZOgsx12PwCOPc/ZYyRJUs0YyhSMn97BOWYCrwUmY6tL7epvdH3TGlj35ODPsWHl9ttW3ptN9ThhBoyZCh3rIeqyEL9xDUTAsw4zyEuSpKoxlJH083eyfy3wuZTSf+56OaoaPaPrsHWEfeX9sLkN2h7ftXOuW5599eeZZVmQb54M46ZB50agHsZMhtQFDU1bA33vcN+z3jQOjnq37TaSJKkQIqU0uAMjXjLArm5gNXBPSqlzpArLU0tLS1q4cGHZP7e1tZV58+aV/XPLbuHl8NfvQ1d71s4ynOA+0hrGQsMYaBwHdENdE3RsyAJ98wRo35Ad1zwR2tcBAc2ToLsz66vf3JY9PWDslAF/KNiw+knGjR23tZ1nZz9A9Lc+lPeV4zPK/b4KqG2b61yw2or8fau02tZ2BJMaugtZW8W+r2C1lf3f7Cr5vlVabWvWb2LyS96Xy2BdRCxKKbX03T7okfSU0h9HtiTVrN6j7D16B/eevzx19VtvNC2Xzo3Z16bV2+9bv2Lreu8Wnv5+wFgz8K0ZYwE273KFqhBe59owMe8CNOr8u1wbJgFcc172oiC/VR/p2V2kXdNfcIesVeZvV8JT98H6lf3/ZNy+DjauKnfFuyzyLkBl4XWuDV7n6uc1rg1brvPdvyh+SI+IfXb1pCmlh3f1vdI2Zs4d3M2gPTeqdm4a+q+92lYMbjaaEZLwH/1a4HWuDV7n6uc1rg1brvMhr825kq12NJK+lKzmoUo7Oa808gYaiR+snptbn7hz1PvtNtqTXt1/ptKxG+1Jr4na2uxJr/o/U9n/za6S71ul1bY2x570gewoTH+fXQvpUuXpPX3kKLu9Vm4QrnFe59pwh9e56vl3uTb8tbWVeS3z8i5jGwOG9JTS6WWsQ5IkSVJJXd4FSJIkSdqWIV2SJEkqmCHf4BkRRwIvA/YCmvs5JKWU3jXcwiRJkqRaNeiQHhEBXA68HQi2n5Uo9dpuSJckSZJ20VDaXd4HvAP4L6CFLJBfCBwDfBxoA34EzB7ZEiVJkqTaMpR2l9OAe3tmfckG1nkmpXQrcGtE/A64Ffg9cNkI1ylJkiTVjKGMpB8EXN9n25aQn1L6K3AN8J4RqEuSJEmqWUMJ6QGs6fV6PTCtzzFLgIOHW5QkSZJUy4YS0h8jm9Glx4PAC/occwBZeJckSZK0i4YS0m9n21D+G2BuRHwqIp4XEe8FXkvWly5JkiRpF+0wpEfEayKi55irgfqI2K/0+j+BZcB84P8DvgE8A3xsdEqVJEmSasPOZnf5OfBYRFwGfDeldEjPjpTSqog4HDgL2B9YCnw/pfTEKNUqSZIk1YSdhfQ/AC8FPgl8PCJ+D1wM/DKl1JVSWgP831GuUZIkSaopO2x3SSmdRPZwoguAJ4GXAT8FHomICyLCBxdJkiRJI2ynN46mlJallD4F7AO8hmwu9N2BfwPui4gFEfGmiBjKg5EkSZIkDWDQs7uklLpTSteklF5LFtg/SdaHfiLwY7Le9S9GxAGjUqkkSZJUI4YyBeMWKaUnU0qfTyk9BzgJ+AkwCfgIcPcI1idJkiTVnJFoUfkj2ZNH9wPmjsD5JEmSpJq2yyE9Ig4CzgTeSdajHsBDwHdHpjRJkiSpNg0ppEfEGOAtZOH8WLJg3gH8D3BJSmnBiFcoSZIk1ZhBhfSImEP20KK3kfWeB/AAcClwWUppxWgVKEmSJNWaHYb0iPgXsnB+OFkwbweuAi5OKV0/+uVJkiRJtWdnI+kXlZb3AZcAV6SUVo5uSZIkSVJt21lI/2+yUfM/lqMYSZIkSTsJ6SmlU8tViCRJkqTMLj3MSJIkSdLoMaRLkiRJBWNIlyRJkgpml584quJbvGIxl/39MpauXUpjXSNt7W0A7Dl+T2ZPmc1r9n8Nc6bPybdISZIkbceQXoWuuvcqLvrbRTy18al+9z++/nEWrVjEVfddxZ7j9iQICMO7JElSURjSq8hV917FN//6TVZtXjXo9zyx4Ykt673D++SmyYxtGMvk5sl0dHcwa9Iszjj0DMO7JElSGRjSq8SX/vIlvn/X90fsfGva17CmfQ1PbngSgAfXPMj1j1zPpKZJjG0Yy5TmKbS1tzGmYQxvP+TtvPmgN4/YZ0uSJNU6Q3oV+NE9P9phQN99zO7sNnY32trbaO9qZ+WmXX9o7Nr2taxtX8vyDcu3bPvsrZ/lq4u+yrQx09jctZm6qGNi00R74CVJknaRIb3CXXXvVVxw2wX97nv2+Gdz5mFnbjfK3d8NpcMN720dbbR1tG3dsH7rau82mt3H7g5AfdQzqXESXXRtc1Or4V6SJMmQXtEWr1jM52793Hbbnz3+2XzxxV8cMNjOmT6Hr730a/2e71cP/IoHnnmAJ9ZnveoNdQ083PbwiNW8cuPWHwSWs3z7AwYI99Oap9FQ30BjXSObOjYREYxvHM+Gzg3URR0TGiewvmM9QTCxeSLr2tcBWejv6O4Y8AeBHa0P5X0DHWs7kCRJ2hWRUsq7hsJpaWlJCxcuLPvntra2Mm/evEEf/97r3suNj924zbaGaOCyl182oiPP/YX3iU0TeXrj08Mafa8lTXVNNNU30VzXzPrN62lsbKS5vpmu1EV91LO5azMAzfXNO16PbL07dVMXdQMe297VvmV9Z8du7tpMENtt7/u+MQ1j2NS5abv1sQ1j2di5kSAY0zCGztRJQzRss20w60N5Xzk+Y7i1rd24lqbGpkLWVuTvW6XVljoTNFDI2ir1fUWrrfff5aLVVuTvW6XV1r65nX854l9423PfNrqBoB8RsSil1LLddkP69iohpC9esZh3/uadJLZev9mTZzP/mPllbQ0ZaC72npHk4bbRSJIklcunj/502X/7PVBIL0y7S0TsDXwWeDmwG/AE8HNgfkpp9S6e8x1Azx2VZ6WULh2BUgvhlsdv2Sag11Nf9oAOA7fO9LajID9Qm4jhXpIkldt1D19XmBbVQoT0iNgfuAWYDvwCuAeYC5wHvDwijk0pPT3Ec84EvgGsAyaMbMX5+9uKv23z+p3Pe2dhb64cTJDvT99wX47e8pHsSbcdSJKkynLiPifmXcIWhQjpwLfJAvq5KaVv9GyMiK8AHwQuAM4Z7MkiIoDLgKeB/wE+MqLV5mzxisXc8sQtW14HwaTmSTlWNDp2NdwXSc8PGvesugfIwvtTa59izJgxud/UWsnvq4Tael/notVW5O9bpdUWm4PUnApZW6W+r2i1lfvf7Gr5vlVabZvWb+K0I04rzCg6FCCkR8Rs4GRgKfCtPrs/A5wNvCMiPpxSWs/gnAu8FJhXWlaV/777v7d5XR/1tMzYrpVJBdDfDxpDvUFYlcnrXBu8ztXPa1wbWltbmXfQvLzL2EZd3gWwNUQvSCl1996RUmoDbgbGAUcP5mQRcQjwBeBrKaUbd3Z8pVm8YjG/W/q7La/ro56PH/Xxwra6SJIkaeiKENIPKi3vG2D/ktLywJ2dKCIagP8CHgY+PvzSimfh8oV00bXl9RsPeGOhfjUjSZKk4cu93QWYXFquGWB/z/YpgzjXp4HDgRellDYOpYiIOJustYYZM2bQ2to6lLePiHXr1u30cx9es+2DhepW1uVSq3bdYK6zKp/XuTZ4nauf17g2FPE6FyGk70yUljuc0D0i5pKNnn85pfTnoX5ISuli4GLI5knPo/9sMH1v1954LTyTrddRx4xZM5h32I7fo2Kxv7E2eJ1rg9e5+nmNa0MRr3MR2l16RsonD7B/Up/jttOrzeU+4FMjV1rxPLnuSSCb0aWpvskbRiVJkqpQEUL6vaXlQD3nB5SWA/WsQzYP+oHAIcCmiEg9X2QzxABcUtp24XALzsviFYtZ/NRiAOqijn898l+9YVSSJKkKFaHd5YbS8uSIqOs9w0tETASOBTYCt+7gHJuB7w6w7wiyPvWbyH4gGHIrTFH86dE/bXnKaEqJNe0D/nJBkiRJFSz3kJ5SeiAiFpDNlf5esqeE9pgPjAf+X88c6RHRCOwPdKSUHiidYyNwZn/nj4jzyUL6FSmlS0frz1EO6zu2ThPfTTeTmwbqEJIkSVIlyz2kl7wHuAX4ekScANwNHAUcT9bm8olex+5V2r8MmFXeMvP14JoHt6zXUedIuiRJUpUqQk86pRHxFuBysnD+YbLR8q8DL0wpPZ1fdcWxuWsz4E2jkiRJ1a4oI+mklB4BzhjEcUvZOi3jYM57PnD+rtZVJA898xAAL5v1Mk495FRvGpUkSapShRhJ18798ZE/srp9NQA3PHLDTo6WJElSJTOkV4jewbyzu5OFyxfmWI0kSZJGkyG9Qmzq3LRlvbGu0X50SZKkKmZIrwCLVyzmt0t/C2SzuvgQI0mSpOpmSK8AC5cvpCt1bXnt1IuSJEnVzZBeAXq3tjTUNdjqIkmSVOUM6RXgoGkHAdn86P/nyP9jq4skSVKVM6RXgAVLFwCQSHx54ZdZvGJxvgVJkiRpVBnSK8DNj9+8Zb2ju8PpFyVJkqqcIb0CTG6aDGTtLk6/KEmSVP0a8i5AOxcRALxorxdx9vPPtiddkiSpyjmSXgHuefoeAI6YfoQBXZIkqQYY0gtu8YrFLH5qMQAX/e0ibxqVJEmqAYb0glu4fCGJBEBn6vSmUUmSpBpgSC+45+32vC3rTXVN3jQqSZJUAwzpBbfb2N0AmNo8lUtOvsSedEmSpBpgSC+4Wx6/BYC9J+5tQJckSaoRhvQCW7xiMV+742sA3PX0Xd40KkmSVCMM6QW2cPlCurq7AOhO3d40KkmSVCMM6QXWMqNly4OM6uvqvWlUkiSpRhjSC2zO9DlbZnc59/Bz7UmXJEmqEYb0guvs7gRg7rPm5lyJJEmSysWQXnBPb3wa2DoVoyRJkqqfIb3AulM3KzeuBODRtkdzrkaSJEnlYkgvsJsfu5luugF493XvdgpGSZKkGmFIL7CbH795y3pHd4dTMEqSJNUIQ3qBzZw4c8t6Y12jUzBKkiTVCEN6gU1ungzAAVMO4JKTL3EKRkmSpBphSC+wnpldjtrzKAO6JElSDTGkF9jTm5x+UZIkqRYZ0gtsyeolAKzvWJ9zJZIkSSonQ3pBLV6xmFseuwWAK/5xhdMvSpIk1RBDekEtXL5wyxzpXd1dTr8oSZJUQwzpBdUyo4UgAGiob3D6RUmSpBpiSC+oOdPnMLZhLABfP/7rzu4iSZJUQwzpBdXV3cWGzg0EwdF7Hp13OZIkSSojQ3pBtbW3ATCxaSL1dfU5VyNJkqRyMqQX1Jr2NcDWp45KkiSpdhjSC2rN5iykt3e1O/2iJElSjTGkF9RfV/wVgOUblnPWgrMM6pIkSTXEkF5Qd668c8t6R3eH86RLkiTVEEN6QU0fOx2AIGisa3SedEmSpBrSkHcB6t/EpokAHPmsI3n/4e93nnRJkqQa4kh6QfXM7jJv5jwDuiRJUo0xpBdUz+wuTsEoSZJUewzpBbUlpDcZ0iVJkmqNIb2gfJiRJElS7TKkF9SK9SsAeHz94zlXIkmSpHIzpBfQ4hWLeXLDkwB8+uZP+yAjSZKkGmNIL6C/PPmXLeud3Z0+yEiSJKnGGNIL6J/2+Kct6z7ISJIkqfYY0gto9pTZAIxtGMslJ1/iPOmSJEk1xpBeQBs6NgCw25jdDOiSJEk1yJBeQOs61gEwoWlCzpVIkiQpD4b0AlrfsR6A8Y3jc65EkiRJeTCkF5AhXZIkqbYZ0guop93FkC5JklSbDOkFtL49G0mf0GhPuiRJUi1qyLsAbe++1fcBW9teJEmSRtPmzZtZtWoVbW1tdHV15V1O2U2ePJm777572Oepr69n4sSJTJs2jebm5mGdy5BeMItXLObqJVcDsGDpAt568FudhlGSJI2azZs38/DDDzN16lRmzZpFY2MjEZF3WWXV1tbGxIkTh3WOlBIdHR2sXbuWhx9+mH322WdYQb0w7S4RsXdEfC8iHo+IzRGxNCIujIipg3z/bhFxZkT8LCLuj4iNEbEmIm6KiHdFRGH+rDuycPlCulL2E2xX6mLh8oU5VyRJkqrZqlWrmDp1KrvvvjtNTU01F9BHSkTQ1NTE7rvvztSpU1m1atWwzleI4BoR+wOLgDOA24GvAg8C5wF/jojdBnGaNwOXAEcBtwEXAlcDhwKXAj+JCvivrmVGC3Wlnyfq6+ppmdGSc0WSJKmatbW1MWnSpLzLqCqTJk2ira1tWOcoREgHvg1MB85NKb0upfSxlNJLycL6QcAFgzjHfcBrgL1TSqemlP4tpfS/gYOBR4A3Am8YnfJHzpzpczhi+hEAvPv577bVRZIkjaquri4aGxvzLqOqNDY2Dru3P/eQHhGzgZOBpcC3+uz+DLAeeEdE7HA+wpTS9SmlX6WUuvtsfxL4TunlvJGoebQ11GW3Chy6+6E5VyJJkmpBBTQbVJSR+H7mHtKBl5aWC/oJ2G3AzcA44OhhfEZHadk5jHOUTc+sLuMax+VciSRJkvJQhJB+UGl53wD7l5SWB+7KySOiAXhn6eVvd+Uc5dbzMCPnSZckSapNRZiCcXJpuWaA/T3bp+zi+b9AdvPotSml3w10UEScDZwNMGPGDFpbW3fx43bdunXraG1tZfW61QDcuehOHm14tOx1aHT1XGdVN69zbfA6V79auMaTJ08e9k2Ola6rq2vEvwebNm0a1n87RQjpO9PT1JOG/MaIc4EPA/cA79jRsSmli4GLAVpaWtK8efOG+nHD1trayrx58+i4sgO64IQXn8CkJu+2rjY911nVzetcG7zO1a8WrvHdd9897DnCK13fedIvuOACPvnJTwJwzz33cNBBBw301gGNGTOGww8/fJdrKkK7S89I+eQB9k/qc9ygRMR7ga8BdwHHp5SGN1llmXSn7i096eMbdnivrCRJkkZYSonvfve7W27+vOSSS3Kpowgh/d7ScqCe8wNKy4F61rcTER8Avgn8nSygP7nL1ZXZxs6NAIxtGEt9XX3O1UiSJNWWBQsW8NBDD3HaaacxY8YMrrjiCtrb28teRxFC+g2l5cl9nwoaEROBY4GNwK2DOVlEfJRsfvXFZAF9xciVOvrWtWc3jdZRx+IVi/MtRpIkqcb0jJyfddZZnHrqqaxcuZKf/exnZa8j95CeUnoAWADMAt7bZ/d8YDzw/ZTSeoCIaIyIg0tPKd1GRHyK7EbRRcAJKaWVo1n7aFi0fBEA6zvXc9aCswzqkiRJZbJ8+XJ++ctfcuCBB3LMMcdwxhlnAHDxxReXvZai3Dj6HuAW4OsRcQJwN3AUcDxZm8sneh27V2n/MrJgD0BEnAZ8FugC/gSc289E8ktTSpePyp9ghPx1xV+3rHd0d7Bw+UKfOipJkqrComWrufXBpzl69m68YN+peZezncsuu4yOjg5OP/10AA499FCOOOIIbrjhBu6//36e85znlK2WQoT0lNIDEdFCFrJfDrwSeAL4OjB/kDd97lda1gMfGOCYPwKXD6vYUTZ7yuwt6411jbTMaMmxGkmSVMtmfezXeZcwKEu/8KphnyOlxKWXXkpdXR3vfOc7t2w//fTTueOOO7j00kv5whe+MOzPGazc2116pJQeSSmdkVLaM6XUlFLaN6V0Xt+AnlJamlKKlNKsPtvPL23f0de8cv6ZdsXeE/YGYN+J+3LJyZc4ii5JklQG119/PQ888AAnnXQSe+2115btb3vb22hqauLyyy+no6NjB2cYWYUYSddWPbO7PGfqcwzokiQpVyMxQt1j0bLVnHrprXR0dtPYUMcPzzy6UC0vPX3nPa0uPXbbbTde/epXc/XVV/OLX/yCN73pTWWppzAj6cps6toEZFMwSpIkVYsX7DuVH555NB86+aDCBfSVK1fy85//HIC3vvWtRMQ2X1dffTVQ3htIHUkvmI0d2Uj6mIYxOVciSZI0sl6w79RChfMeV155Je3t7bzgBS9gzpw5/R7zy1/+kuuuu46HHnqI/fbbr99jRpIhvWAcSZckSSqvK664AoBvf/vbzJ07t99jPvWpT/G5z32OSy+9lAsuuGDUa7LdpWA2dG4AYEy9I+mSJEmjrbW1lSVLlnDYYYcNGNAB3vWudxERXHbZZXR2do56XYb0gtnUmY2kj2scl3MlkiRJ1a/nCaNnnnnmDo+bNWsWJ554Ik888QS/+tWvRr0u210Kpmd2F0fSJUmSRt8Pf/hDvvOd7zBx4sSdHrtgwYIyVJRxJL1gekbS7UmXJEmqXYb0gukJ6c7uIkmSVLsM6QXz5Pont1lKkiSp9hjSC+ShzQ/xt6f+BsC3Fn+LxSsW51uQJEmScmFIL5Alm5bQTTcAXd1dLFy+MOeKJEmSlAdDeoEcMOYAggCgoa6BlhktOVckSZKkPBjSC2S/5v3YY9weAFzwoguYM31OvgVJkiQpF4b0gkkpAXD49MNzrkSSJEl5MaQXzJaHGTkFoyRJUs0ypBdMzzzp4xrG5VyJJEmS8mJIL5Cu1EVn6qQ+6mmoa8i7HEmSJOXEkF4gm9NmAMY2jCUicq5GkiRJeTGkF0hHdwdgP7okSVKtM6QXSO+RdEmSJNUuQ3qBtKd2wJF0SZKkcouI7b6am5uZNWsWp512GnfffXdZ6/HuxAJp785CuiPpkiRJ+fjMZz6zZX3NmjXcfvvtfP/73+fqq6/mpptuYs6cOWWpw5BeID0j6WPrDemSJEl5OP/887fb9v73v59vfvObXHjhhVx++eVlqcN2lwKx3UWSJKl4Tj75ZACeeuqpsn2mIb1AekL64+seZ/GKxfkWI0mSJACuu+46AFpaWsr2mba7FMjjmx8HYMkzSzhrwVlccvIlzJk+J9+iJEmSRsojt8PSP8Gs42Dm3Lyr6Vfvdpe1a9fyl7/8hZtvvplTTjmFj3zkI2Wrw5BeIE90PLFlvaO7g4XLFxrSJUlSfs6fnHcFg3P+mhE71fz587fb9tznPpe3vvWtTJw4ccQ+Z2dsdymQaQ3TAAiCxrpGWmaU71cqkiRJgpTSlq9169Zx2223MWPGDE499VQ+8YlPlK0OR9ILZGJ99tNZy4wWzj3iXEfRJUlSvkZwhJpHbocrXgNd7VDfBKf9srAtLz3Gjx/P3Llz+Z//+R/23ntv/vM//5NzzjmHmTNnjvpnO5JeIJ2pE4C5e841oEuSpOoyc24WzF/6iYoI6L1NmTKFgw46iM7OTu64446yfKYj6QXSQQcAzfXNOVciSZI0CmbOrahw3tvq1asB6O7uLsvnOZJeIB0pC+lN9U05VyJJkqQeP//5z3nooYdobGzkmGOOKctnOpJeID3tLo6kS5Ik5aP3FIzr16/nrrvu4je/+Q0An//855kxY0ZZ6jCkF4ghXZIkKV+9p2Csr69njz324NWvfjXve9/7OOmkk8pWhyG9QGx3kSRJykdKKe8StmFPeoE4ki5JkiQwpBeKI+mSJEkCQ3qh9IR0R9IlSZJqmyG9QGx3kSRJEhjSC6VnJL2xrjHnSiRJkpQnQ3qBOJIuSZIkMKQXij3pkiRJAkN6oWzu3gzAktVLcq5EkiRJeTKkF8TiFYvZmDYC8KE/fojFKxbnW5AkSZJyY0gviIXLF25Z7+zu3Oa1JEmSaoshvSAOn374lvXGukZaZrTkWI0kSZLyZEgviEOmHQJAQzRwycmXMGf6nHwLkiRJUm4M6QXR3tUOwLjGcQZ0SZKkGmdIL4jNXdnMLk6/KEmSJEN6QbR3ZyPpTfVNOVciSZJUu+655x7e//73c+ihhzJ58mSampp49rOfzate9Sq++93vsmnTprLU0VCWT9FO9bS7OJIuSZKUj89+9rPMnz+f7u5ujj76aE477TQmTJjA8uXLaW1t5cwzz+Siiy5i4cLRn4XPkF4QtrtIkiTl5/Of/zyf+cxnmDlzJldddRVHHXXUdsdcc801fPnLXy5LPYb0gugZSbfdRZIkqbyWLVvG+eefT2NjI9deey2HHnpov8edcsopnHTSSWWpyZ70gnAkXZIkKR8/+MEP6Ojo4I1vfOOAAb1Hc3N5spohvSB6Qroj6ZIkqVotXrGYS++8lMUrFuddyjZuvfVWAE444YScK9nKdpeC8MZRSZJUNIddcVjeJQzKnafdOaz3P/nkkwDsvffeI1HOiHAkvSAcSZckScpHSgmAiMi5kq0cSS+ILTeO1hnSJUlSMQx3hLq3xSsWc9aCs+jo7qCxrpFLTr6kME9Z33PPPbnvvvt49NFH8y5lC0fSC8IbRyVJUjWbM30Ol5x8Ce87/H2FCugARx99NAB/+MMfcq5kq8KE9IjYOyK+FxGPR8TmiFgaERdGxNQ8zlNutrtIkqRqN2f6HM487MxCBXSAt7/97TQ2NnL11Vdz11137fDYzZs3l6WmQoT0iNgfWAScAdwOfBV4EDgP+HNE7FbO8+TBG0clSZLyse+++3L++efT3t7Oq171qgGfKPrb3/6WV7ziFWWpqSg96d8GpgPnppS+0bMxIr4CfBC4ADinjOcpO9tdJEmS8vPxj3+czs5O5s+fz5FHHskxxxxDS0sLEyZMYPny5dx4440sWbKElpaWstST+0h6RMwGTgaWAt/qs/szwHrgHRExvhznyYtPHJUkScrXpz/9af7+97/zvve9jzVr1nDZZZfxpS99iV//+tfsv//+XHrppdx0001lqaUII+kvLS0XpJS6e+9IKbVFxM1k4ftoYEfd/CN1nlw8ui67m3jlxpU5VyJJklS7DjnkEL7xjW/s/MBRlvtIOnBQaXnfAPuXlJYHluk8Zbd4xWKuf/h6AH50748K9xQuSZIklVcRRtInl5ZrBtjfs33KaJ4nIs4GzgaYMWMGra2tO/m4kbNgzQK6UhcAXd1d/OTPP+GZyc+U7fNVXuvWrSvrf1/Kh9e5Nnidq18tXOPJkyfT1taWdxm56urqGvHvwaZNm4b1304RQvrO9Dz6KY3meVJKFwMXA7S0tKR58+YN8+MGb8qKKSxYsID2rnaa6pt4ywvfUripiTRyWltbKed/X8qH17k2eJ2rXy1c47vvvpuJEyfmXUau2traRvx7MGbMGA4//PBdfn8R2l16RrgnD7B/Up/jRvs8ZTdn+hwuPflSTplyCpeefKkBXZIkqcYVYST93tJyoF7xA0rLgXrNR/o8uZgzfQ7PTH7GgC5JkqRCjKTfUFqeHBHb1BMRE4FjgY3ArWU6jyRJkpSr3EN6SukBYAEwC3hvn93zgfHA91NK6wEiojEiDi49XXSXzyNJkiQVVRHaXQDeA9wCfD0iTgDuBo4CjidrT/lEr2P3Ku1fRhbId/U8kiRJAlJKRMTOD9SgpDTc+U4KMJIOW0bBW4DLyUL1h4H9ga8DL0wpPV3O80iSJNWK+vp6Ojo68i6jqnR0dFBfXz+scxRlJJ2U0iPAGYM4bilbp1Pc5fNIkiQJJk6cyNq1a9l9993zLqVqrF27dthTOhZiJF2SJEn5mDZtGqtXr2blypW0t7ePSKtGLUop0d7ezsqVK1m9ejXTpk0b1vkKM5IuSZKk8mtubmafffZh1apVLF26lK6urrxLKrtNmzYxZsyYYZ+nvr6eiRMnss8++9Dc3DyscxnSJUmSalxzczN77rkne+65Z96l5KK1tXVYTwcdDba7SJIkSQVjSJckSZIKxpAuSZIkFYwhXZIkSSoYQ7okSZJUMIZ0SZIkqWAM6ZIkSVLBGNIlSZKkggkf/bq9iHgKWJbDR+8OrMzhc1VeXufa4HWuDV7n6uc1rg15Xud9U0p79N1oSC+QiFiYUmrJuw6NLq9zbfA61wavc/XzGteGIl5n210kSZKkgjGkS5IkSQVjSC+Wi/MuQGXhda4NXufa4HWufl7j2lC462xPuiRJklQwjqRLkiRJBWNIlyRJkgrGkJ6ziNg7Ir4XEY9HxOaIWBoRF0bE1Lxr0/BFxG4RcWZE/Cwi7o+IjRGxJiJuioh3RYR/B6tURLwjIlLp68y869HIiojjIuLqiHii9G/3ExGxICJemXdtGr6IeFXpej5a+nf7wYi4KiJemHdtGpqIeFNEfCMi/hQRa0v/Jv9gJ+85JiKujYhVEbEhIv6/iPhARNSXq26AhnJ+mLYVEfsDtwDTgV8A9wBzgfOAl0fEsSmlp3MsUcP3ZuAi4AngBuBhYAbwBuBS4BUR8ebkzSFVJSJmAt8A1gETci5HIywiPgn8O9mDT64h+/u9O3A4MA+4NrfiNGwR8UXgX4GngZ+TXefnAK8F3hgR70wp7TDkqVA+CfwT2b/HjwIH7+jgiHgtcDWwCfgxsAp4NfBV4Fiy/6+XhTeO5igifgecDJybUvpGr+1fAT4I/L+U0jl51afhi4iXAuOBX6eUunttfxZwOzATeFNK6eqcStQIi4gAfg/sB/wP8BHgrJTSpbkWphEREW8GfgJcB7whpdTWZ39jSqkjl+I0bKV/mx8DngKen1Ja0Wvf8cD1wEMppdk5laghKl23R4H7gZeQDZj9MKX09n6OnVQ6bjJwbEppYWn7GLJr/0LgrSmlH5Wjdn/VnpOImE0W0JcC3+qz+zPAeuAdETG+zKVpBKWUrk8p/ap3QC9tfxL4TunlvLIXptF0LvBS4Ayyv8eqEqX2tC8CG4C39Q3oAAb0ircvWTa6rXdAB0gp3QC0Ads9vl3FlVK6IaW0ZJC/sX4T2fX9UU9AL51jE9mIPMC7R6HMfhnS8/PS0nJBPwGuDbgZGAccXe7CVDY9/zPvzLUKjZiIOAT4AvC1lNKNedejEXcM2W9IrgVWl/qWPxoR59mrXDWWAO3A3IjYvfeOiHgxMJHstyiqTj3Z7Lf97LuR7Af0YyKiuRzF2JOen4NKy/sG2L+EbKT9QOAPZalIZRMRDcA7Sy/7+8dAFaZ0Tf+L7L6Dj+dcjkbHkaXlcuAO4LDeOyPiRrL2tafKXZhGRkppVUR8FPgKcFdE/JysN31/4DVkrWz/kl+FGmUDZrOUUmdEPAQ8D5gN3D3axRjS8zO5tFwzwP6e7VNGvxTl4AvAocC1KaXf5V2MRsSnyW4cfFFKaWPexWhUTC8tzwEeAk4EbiNrkfgy8DLgKmxhq2gppQsjYinwPeCsXrvuBy7v2wajqlKobGa7S3FFaemdvVUmIs4FPkw2m887ci5HIyAi5pKNnn85pfTnvOvRqOmZfi3IRsz/kFJal1L6B/B6spvTXmLrS2WLiH8FfgpcTjaCPh54AfAg8MOI+M/8qlPOyprNDOn56flpbPIA+yf1OU5VICLeC3wNuAs4PqW0KueSNEy92lzuAz6VczkaXatLywdTSn/rvaP025Oe34rNLWtVGjERMY/s5uBfppQ+lFJ6MKW0IaV0B9kPYo8BHy5N/qDqU6hsZkjPz72l5YED7D+gtByoZ10VJiI+AHwT+DtZQH8y34o0QiaQ/T0+BNjU6wFGiWymJoBLStsuzKtIjYief7efGWB/T4gfO/qlaJScUlre0HdHSmkD2dS5dWStbao+A2az0oDMfmSTPTxYjmLsSc9Pzz8AJ0dEXZ85tCeSTZi/Ebg1j+I0sko3In0BWAyclFJamW9FGkGbge8OsO8Isv+Z30T2j7+tMJXtRrL/QR8QEU0ppfY++w8tLZeWtSqNpJ5ZOwaaZrFne99rr+pwPXAq8HLgv/vsezHZrHs3ppQ2l6MYR9JzklJ6AFgAzALe22f3fLIeuO+nlJxnucJFxKfIAvoi4AQDenVJKW1MKZ3Z3xfwy9JhV5S2/TjPWjU8pb+7Pyb7Vfine++LiJPIbhxdgzM2VbI/lZZnR8RevXdExCvIBtA2kT0tXNXnp2RPmP3niGjp2Vh6mNHnSi8vKlcxPnE0RxGxP9lf9OnAL8im8zkKOJ6szeWYlNLT+VWo4YqI08huPuoie0x8f31sS1NKl5exLJVJRJxP1vLiE0erRERMJ3uOxXPIAt3tZLO7vJ7sZrK3pZSuyq9CDUfpgVW/I5u5pw34GfAkWTvbKWQ3Dn4gpfS13IrUkETE64DXlV4+i+yH6QfZ+gPZypTSR/oc/1OyH8Z+BKwim37zoNL2twzywUjDZkjPWUTMBD5L9quV3YAngJ8D872psPL1Cmk78seU0rzRr0blZkivThExjezpg68H9iILczcB/5FSskWxwkVEI9lvuP8ZeC5Zi8Mqsh/Ivp5SWpBjeRqiQfx/eFlKaVaf9xwLfAJ4ITCGbPrN75Fd/67RqXR7hnRJkiSpYOxJlyRJkgrGkC5JkiQVjCFdkiRJKhhDuiRJklQwhnRJkiSpYAzpkiRJUsEY0iVJkqSCMaRLkvoVEadHRIqI0/OuRZJqTUPeBUiSRl9EDPXJdWeMSiGSpEExpEtSbZjfz7YPAJOBrwHP9Nm3GHgIuBV4YhTrkiT1I1Ia6uCKJKkaRMRSYF9gv5TS0nyrkST1Zk+6JKlfA/WkR8TS0teEiPhqRDwSERsjYnFEvK50TENEfDwilkTEpoh4ICLet4PPellEXBsRKyNic+n4L0XElFH9Q0pSQdnuIknaFY3A74FpwC+AJuCtwNURcTLwHuAo4DfAZuDNwDci4qmU0o97nygiPk3WjrMKuAZYATwf+Ajwyoh4YUppbVn+VJJUEIZ0SdKueDZwBzAvpbQZICL+C7gRuAp4ADg0pfRMad9XgHuAjwFbQnpEHE8W0P8MvLLn+NK+04HLSvs/ONp/IEkqEttdJEm76gM9AR0gpfQnsptNpwIf7R24U0oPAjcDh0VEfa9znFtantX7+NJ7Lie7gfXUUahdkgrNkXRJ0q54JqX0QD/bHwf2Axb1s+8xoB54Vmkd4IVAB/DmiHhzP+9pAvaIiN1SSk8Pv2xJqgyGdEnSrlgzwPZOgJRSf/s7S8vGXtt2I/t/0Wd28nkTAEO6pJphSJck5WkNUJdSmpZ3IZJUJPakS5LydCswNSKel3chklQkhnRJUp6+WlpeEhHP7rszIsZHxNFlrkmScme7iyQpNymlP0TEx4D/AJZExLVkM8RMIHsa6kuAm4CX51elJJWfIV2SlKuU0hcj4may6RhfBLyWrFf9MeBi4Mocy5OkXERKKe8aJEmSJPViT7okSZJUMIZ0SZIkqWAM6ZIkSVLBGNIlSZKkgjGkS5IkSQVjSJckSZIKxpAuSZIkFYwhXZIkSSoYQ7okSZJUMIZ0SZIkqWD+f89JnLrvnBfkAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 864x576 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "<Figure size 864x576 with 0 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "sim.run(10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "78779261-b625-4841-b084-0028897d8d98",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "['D'] k1 ['D', 'M']\n",
      "['M'] k2 ['M', 'P']\n",
      "['M'] k3 ['ϕ']\n",
      "['P'] k4 ['ϕ']\n"
     ]
    }
   ],
   "source": [
    "sim=ChemSimulation(\n",
    "\"\"\"\n",
    "D --k1--> D+M\n",
    "M --k2--> M+P\n",
    "M --k3--> ϕ\n",
    "P --k4--> ϕ\n",
    "\"\"\",D=1,M=0,ϕ=0,P=0,k1=.01,k3=0.00577,k4=0.0001925,k2=20*0.00577)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4b71ce14-a821-4cb5-b1f4-e8b5c77110d8",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "D'=0\n",
      "M'=+k1*D -k3*M\n",
      "P'=+k2*M -k4*P\n",
      "k1=0.01\n",
      "k2=0.1154\n",
      "k3=0.00577\n",
      "k4=0.0001925\n",
      "\n"
     ]
    }
   ],
   "source": [
    "print(sim.equations())"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "34d5653a-8240-47dd-9216-0fd66618a722",
   "metadata": {},
   "source": [
    "## These are wrong -- need to deal with squares"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a005dcc8-3f8b-4ef5-9b00-21459a0e9877",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "['A'] k1 ['X']\n",
      "['B', 'X'] k2 ['Y', 'D']\n",
      "['2*X', 'Y'] k3 ['3*X']\n",
      "['X'] k4 ['E']\n",
      "Components ['2*X', '3*X', 'A', 'B', 'D', 'E', 'X', 'Y']\n",
      "Parameters ['k1', 'k2', 'k3', 'k4']\n",
      "diffeqs [\"2*X' =  -k3*2*X*Y\", \"3*X' =  +2*k3*2*X*Y\", \"A' =  -k1*A\", \"B' =  -2*k2*B*X\", \"D' =  +2*k2*B*X\", \"E' =  +k4*X\", \"X' =  +k1*A -2*k2*B*X -k4*X\", \"Y' =  +2*k2*B*X -k3*2*X*Y\"]\n",
      "\n",
      "2*X'=-k3*2*X*Y\n",
      "3*X'=+2*k3*2*X*Y\n",
      "A'=-k1*A\n",
      "B'=-2*k2*B*X\n",
      "D'=+2*k2*B*X\n",
      "E'=+k4*X\n",
      "X'=+k1*A -2*k2*B*X -k4*X\n",
      "Y'=+2*k2*B*X -k3*2*X*Y\n",
      "k1=1\n",
      "k2=1\n",
      "k3=1\n",
      "k4=1\n",
      "\n"
     ]
    }
   ],
   "source": [
    "sim=ChemSimulation(\"\"\"\n",
    "A   -->k1->   X\n",
    "B+X   -->k2->  Y+D\n",
    "2*X+Y -->k3->   3*X\n",
    "X   -->k4->  E\n",
    "\"\"\",A=1,B=2,D=3,X=0,Y=1,k1=1,k2=1,k3=1,k4=1)\n",
    "\n",
    "print()\n",
    "\n",
    "print(sim.equations())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7f395341-a415-48fb-abdb-add4d4053ebd",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "python3",
   "language": "python",
   "name": "python3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}