ECR.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "***Note: this is the ECR.ipynb notebook. The\n",
    "PDF version \"Enzyme-catalysed Reactions\"\n",
    "is available [here](ECR.pdf).***"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Introduction\n",
    "The standard model of an enzyme-catalysed reaction is given by \n",
    "<cite data-cite=\"KeeSne09\">(Keener and Sneyd, 2009)</cite>, Sec. 1.4 is given by the two coupled reactions:\n",
    "\\begin{align}\n",
    "A + E &\\Leftrightarrow C \\\\\n",
    "C &\\rightarrow B + E \n",
    "\\end{align}\n",
    "where A is the substrate, B the product, E the enzyme and C a complex. The first reaction is reversible and the second is not. This is not thermodynamically feasible so <cite data-cite=\"KeeSne09\">(Keener and Sneyd, 2009)</cite>, Sec. 1.4.5 and \n",
    "<cite data-cite=\"GawCra14\">(Gawthrop and Crampin, 2014)</cite>\n",
    "consider the reversible reaction pair:\n",
    "\\begin{align}\n",
    "A + E &\\Leftrightarrow C \\\\\n",
    "C &\\Leftrightarrow B + E \n",
    "\\end{align}\n",
    "Such a reaction pair is thermodynamically feasible and so has a bond graph representation given by <cite data-cite=\"GawCra14\">(Gawthrop and Crampin, 2014)</cite> Sec. 5. The corresponding flow $v$ is given by:\n",
    "\\begin{align}\n",
    " v &= \\bar{\\kappa} \\frac{K_c e_0}{K_m +  \\sigma_v} \\delta_v\\\\\n",
    " \\text{where }\n",
    " \\bar{\\kappa} &= \\frac{\\kappa_1 \\kappa_2}{{\\kappa_1+\\kappa_2}}\\\\\n",
    " \\delta_v &= K_A x_A - K_B x_B\\\\\n",
    " K_m &= \\frac{K_c}{K_e}\\\\\n",
    " \\text{and }\n",
    " \\sigma_v &= \\frac{\\kappa_1 K_A x_A + \\kappa_2 K_B x_B}{\\kappa_1+\\kappa_2}\n",
    "\\end{align}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Implement the theory\n",
    "The flow equations can be expressed in Python as:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "def ECR_flow(x_A,x_B,e0=1,\n",
    "             K_A = 1,K_B=1,K_C=1,K_E=1,\n",
    "             kappa_r1 = 1,kappa_r2=1):\n",
    "    \"\"\"Theoretical flows in Enzyme-catalysed Reactions\n",
    "    \"\"\"\n",
    "    \n",
    "    kappa_bar = (kappa_r1*kappa_r2)/(kappa_r1+kappa_r2)\n",
    "    delta = K_A*x_A - K_B*x_B\n",
    "    sigma = (kappa_r1*K_A*x_A + kappa_r2*K_B*x_B)/(kappa_r1 + kappa_r2)\n",
    "    K_m = K_C/K_E\n",
    "    \n",
    "    v = kappa_bar*e0*K_C*delta/(K_m+sigma)\n",
    "    \n",
    "    return v"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Import some python code\n",
    "The bond graph analysis uses a number of Python modules:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "## Some useful imports\n",
    "import BondGraphTools as bgt\n",
    "import numpy as np\n",
    "import sympy as sp\n",
    "import matplotlib.pyplot as plt\n",
    "import IPython.display as disp\n",
    "\n",
    "## Stoichiometric analysis\n",
    "import stoich as st\n",
    "\n",
    "## SVG bg representation conversion\n",
    "import svgBondGraph as sbg\n",
    "\n",
    "## Set quiet=False for verbose output\n",
    "quiet = True"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# The bond graph representation\n",
    "The bond graph representation of the (reversible) enzyme-catalysed reaction is\n",
    "<cite data-cite=\"GawCra14\">(Gawthrop and Crampin, 2014)</cite>:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/svg+xml": [
       "<svg height=\"234pt\" viewBox=\"406 -1275 7287 3714\" width=\"460pt\" xmlns=\"http://www.w3.org/2000/svg\" xmlns:xlink=\"http://www.w3.org/1999/xlink\">\n",
       "<g fill=\"none\">\n",
       "<!-- Text -->\n",
       "<g transform=\"translate(2835,2025) rotate(-270)\">\n",
       "<text fill=\"#000000\" font-family=\"Helvetica\" font-size=\"216\" font-style=\"normal\" font-weight=\"bold\" text-anchor=\"middle\" x=\"0\" xml:space=\"preserve\" y=\"0\">Re:r1</text>\n",
       "</g><!-- Text -->\n",
       "<g transform=\"translate(5085,2025) rotate(-270)\">\n",
       "<text fill=\"#000000\" font-family=\"Helvetica\" font-size=\"216\" font-style=\"normal\" font-weight=\"bold\" text-anchor=\"middle\" x=\"0\" xml:space=\"preserve\" y=\"0\">Re:r2</text>\n",
       "</g><!-- Text -->\n",
       "<text fill=\"#000000\" font-family=\"Helvetica\" font-size=\"216\" font-style=\"normal\" font-weight=\"bold\" text-anchor=\"middle\" x=\"675\" xml:space=\"preserve\" y=\"990\">Ce:A</text>\n",
       "<!-- Text -->\n",
       "<text fill=\"#000000\" font-family=\"Helvetica\" font-size=\"216\" font-style=\"normal\" font-weight=\"bold\" text-anchor=\"middle\" x=\"1800\" xml:space=\"preserve\" y=\"2115\">1</text>\n",
       "<!-- Text -->\n",
       "<text fill=\"#000000\" font-family=\"Helvetica\" font-size=\"216\" font-style=\"normal\" font-weight=\"bold\" text-anchor=\"middle\" x=\"675\" xml:space=\"preserve\" y=\"2115\">0</text>\n",
       "<!-- Text -->\n",
       "<text fill=\"#000000\" font-family=\"Helvetica\" font-size=\"216\" font-style=\"normal\" font-weight=\"bold\" text-anchor=\"middle\" x=\"4050\" xml:space=\"preserve\" y=\"-1035\">Ce:E</text>\n",
       "<!-- Text -->\n",
       "<text fill=\"#000000\" font-family=\"Helvetica\" font-size=\"216\" font-style=\"normal\" font-weight=\"bold\" text-anchor=\"middle\" x=\"4050\" xml:space=\"preserve\" y=\"90\">0</text>\n",
       "<!-- Text -->\n",
       "<text fill=\"#000000\" font-family=\"Helvetica\" font-size=\"216\" font-style=\"normal\" font-weight=\"bold\" text-anchor=\"middle\" x=\"4050\" xml:space=\"preserve\" y=\"2115\">0</text>\n",
       "<!-- Text -->\n",
       "<text fill=\"#000000\" font-family=\"Helvetica\" font-size=\"216\" font-style=\"normal\" font-weight=\"bold\" text-anchor=\"middle\" x=\"4050\" xml:space=\"preserve\" y=\"990\">Ce:C</text>\n",
       "<!-- Text -->\n",
       "<text fill=\"#000000\" font-family=\"Helvetica\" font-size=\"216\" font-style=\"normal\" font-weight=\"bold\" text-anchor=\"middle\" x=\"7425\" xml:space=\"preserve\" y=\"990\">Ce:B</text>\n",
       "<!-- Text -->\n",
       "<text fill=\"#000000\" font-family=\"Helvetica\" font-size=\"216\" font-style=\"normal\" font-weight=\"bold\" text-anchor=\"middle\" x=\"7425\" xml:space=\"preserve\" y=\"2115\">0</text>\n",
       "<!-- Text -->\n",
       "<text fill=\"#000000\" font-family=\"Helvetica\" font-size=\"216\" font-style=\"normal\" font-weight=\"bold\" text-anchor=\"middle\" x=\"6300\" xml:space=\"preserve\" y=\"2115\">1</text>\n",
       "<!-- Line -->\n",
       "<polyline points=\" 6075,1800 4275,0 4275,225\" stroke=\"#000000\" stroke-width=\"15px\"/>\n",
       "<!-- Line -->\n",
       "<polyline points=\" 3825,0 2025,1800 2250,1800\" stroke=\"#000000\" stroke-width=\"15px\"/>\n",
       "<!-- Line -->\n",
       "<polyline points=\" 900,2025 1575,2025 1350,2250\" stroke=\"#000000\" stroke-width=\"15px\"/>\n",
       "<!-- Line -->\n",
       "<polyline points=\" 6525,2025 7200,2025 6975,2250\" stroke=\"#000000\" stroke-width=\"15px\"/>\n",
       "<!-- Line -->\n",
       "<polyline points=\" 675,1800 675,1125 900,1350\" stroke=\"#000000\" stroke-width=\"15px\"/>\n",
       "<!-- Line -->\n",
       "<polyline points=\" 4050,-225 4050,-900 4275,-675\" stroke=\"#000000\" stroke-width=\"15px\"/>\n",
       "<!-- Line -->\n",
       "<polyline points=\" 4050,1800 4050,1125 4275,1350\" stroke=\"#000000\" stroke-width=\"15px\"/>\n",
       "<!-- Line -->\n",
       "<polyline points=\" 7425,1800 7425,1125 7650,1350\" stroke=\"#000000\" stroke-width=\"15px\"/>\n",
       "<!-- Line -->\n",
       "<polyline points=\" 3150,2025 3825,2025 3600,2250\" stroke=\"#000000\" stroke-width=\"15px\"/>\n",
       "<!-- Line -->\n",
       "<polyline points=\" 4275,2025 4950,2025 4725,2250\" stroke=\"#000000\" stroke-width=\"15px\"/>\n",
       "<!-- Line -->\n",
       "<polyline points=\" 5400,2025 6075,2025 5850,2250\" stroke=\"#000000\" stroke-width=\"15px\"/>\n",
       "<!-- Line -->\n",
       "<polyline points=\" 2025,2025 2700,2025 2475,2250\" stroke=\"#000000\" stroke-width=\"15px\"/>\n",
       "</g>\n",
       "</svg>"
      ],
      "text/plain": [
       "<IPython.core.display.SVG object>"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "## Enzyme-catalysed reaction\n",
    "\n",
    "disp.SVG('RE_abg.svg')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This graphical representation may be converted to [bond-graph tools](https://pypi.org/project/BondGraphTools/) format using "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "sbg.model('RE_abg.svg')\n",
    "import RE_abg"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Reactions\n",
    "The reactions coreponding to this bond graph are generated as follows:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "\\begin{align}\n",
       "\\ch{A + E &<>[ r1 ] C }\\\\\n",
       "\\ch{C &<>[ r2 ] B + E }\n",
       "\\end{align}\n"
      ],
      "text/plain": [
       "<IPython.core.display.Latex object>"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "s = st.stoich(RE_abg.model(),quiet=quiet)\n",
    "disp.Latex(st.sprintrl(s,chemformula=True))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Stoichiometric analysis"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Stoichiometric matrix"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "\\begin{align}\n",
       "N &=\n",
       "\\left(\\begin{matrix}-1 & 0\\\\0 & 1\\\\1 & -1\\\\-1 & 1\\end{matrix}\\right)\n",
       "\\end{align}\n"
      ],
      "text/plain": [
       "<IPython.core.display.Latex object>"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "disp.Latex(st.sprintl(s,'N'))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Conserved moieties"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "\\begin{align}\n",
       "G &=\n",
       "\\left(\\begin{matrix}1 & 1 & 1 & 0\\\\-1 & -1 & 0 & 1\\end{matrix}\\right)\n",
       "\\end{align}\n"
      ],
      "text/plain": [
       "<IPython.core.display.Latex object>"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "disp.Latex(st.sprintl(s,'G'))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Chemostats and Pathways\n",
    "Chemostats are discussed by <cite data-cite=\"GawCra16\">(Gawthrop and Crampin, 2016)</cite> and pathways by <cite data-cite=\"GawCra17\">(Gawthrop and Crampin, 2017)</cite>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "\\begin{align}\n",
       "K &=\n",
       "\\left(\\begin{matrix}1\\\\1\\end{matrix}\\right)\n",
       "\\end{align}\n"
      ],
      "text/plain": [
       "<IPython.core.display.Latex object>"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "chemostats = ['A','B']\n",
    "sc = st.statify(s,chemostats=chemostats)\n",
    "disp.Latex(st.sprintl(sc,'K'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "\\begin{align}\n",
       "\\ch{A &<>[ pr1 ] B }\n",
       "\\end{align}\n"
      ],
      "text/plain": [
       "<IPython.core.display.Latex object>"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sp = st.path(s,sc)\n",
    "disp.Latex(st.sprintrl(sp,chemformula=True))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Theoretical flows\n",
    "The theoretical flows are generated as a function of x_A (with constant x_B) for two values of $K_B$ and two values of $K_E$. $K_B=0$ corresponds to the irreversible case."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZQAAAEQCAYAAACX5IJuAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzs3Xd8FNXawPHf2U2vpJAeIBBAqpRAsKBURbFdxStWEBD1gl5RUKyggIIvNhRQ4aKiCHYFFBDUIIIgnUCAkIRAeu+bTbac948JGDAJSUiyITlfPvvZKWdmnh2SfTJz5pwjpJQoiqIoysXS2ToARVEUpWVQCUVRFEVpECqhKIqiKA1CJRRFURSlQaiEoiiKojQIlVAURVGUBmHThCKEWCGEyBRCHK5mvRBCLBJCxAkhDgkh+lVaN04IcaLiNa7polYURVGqYusrlI+BUTWsvwHoXPGaDCwFEEJ4A7OASGAgMEsI4dWokSqKoig1smlCkVL+DuTWUORWYKXU7ATaCCECgeuBzVLKXCllHrCZmhOToiiK0shsfYVyIcFAUqX55Ipl1S1XFEVRbMTO1gE0NiHEZLTbZTg7O/cPDQ2t136sVis6XXPPv82HOl91o85X3VQ+XxKJVVo5809Kqb1X/LPKqqcl8h9lAazSena/Z5Y1FoE45/2caVHN8vOma5w/s1iCj50PeqGvV5yxsbHZUsq2FyrX3BNKClA5A4RULEsBhpy3PKqqHUgpPwQ+BIiIiJB79uypVyBRUVEMGTLkguUUjTpfddPaz5eUkmJTMbnGXPKMedqrLI/8snzyy/IpLCukoKyAwvJCCssLySzIxKQ3UWIqOZsAaiIQONk54WznfPblpHfC0c5RW653xkHvgJOdEw46Bxz0DjjqHXHUO+Kg1+btdfbatM4Be709dsLu7HI7nR32Ovuzy+10duh1+rPT58/rhA690COEuGDs1Z2vknILOcVl5JSUk1dYSkFGNoaMLMqyczDn5mHNz4PCAvRFhegK8xm5eAHtw0PqdTwhxKnalGvuCWUtMFUIsQatAr5ASpkmhNgEvFqpIv464FlbBakoStWklBSZikgvSSejJINMQ6b2Ks0kuzSbbEM22cZscktzKbeWV7kPe509no6eeDp44unoSYBLAO5l7nQO7Yybgxvu9u64ObjhZu+Gq73r2ZeLvQsudi5nE0h9v7ybWnGZmfSCUlLzjaQXlJKbkklJcgrmtDTIzsIuLxungjw8SwvwNhbhZSwkoNxAUDVXU2WOzhgcndEZihs9dpsmFCHEarQrDV8hRDLak1v2AFLK94GfgBuBOMAAPFixLlcIMQfYXbGrV6SUNVXuK4rSSErNpSQVJZFUlERyUTIpxSmkFqeSUpxCWkkaJaaSf2zj7eSNr7Mvvs6+dGzTER8nH3ycffB28sbLyQsvJy/aOLbBy9GrymQQFRXFkMghTfQJG5bRZCEp18CpHAOncw2kpuVQcvIU1uQk7DNSaVOQhb8hjwBDLl0MeThazedsb9XpMbq3wezlDSHtsfPxRfq1xcG/LW5+Pri6g72uGL01B70pE1GUROHpw3iEtWn0z2bThCKlvPsC6yUwpZp1K4AVjRGXoijnklKSacgkviCekwUnSchP4FThKRILE8kwZJxT1s3ejWC3YELdQxkYMJBA10ACXAPwd/XH38Wfts5tsdfb2+iTNJ3s4jJi04uIzyomPquE9MQUTPHxuKSdpl1RBiFFmXQuzuYKY8E525lcPbD4B2DXpSfOocG4tw/BJSQYO/8A7Pz9sPPxQVjLITsWMo9B9nHIPgTZcRCdAJayv3ems4c2oZjtPMBkaPTP3NxveSmK0sTKLGXE5cVxLPcYx3KPEZsXS1x+HIXlhWfLuDu4E+YZRmRgJO3c29Heoz2h7qGEuIfg6ehpw+ibnsli5URGMTFphcSkFnIsJZ/CuAR8004SXpBCWEEa1xem0qbs71tOFhdXRLsOuEQOxq1zJxw7dMChXTvs27VD7+Z27gGKMiDtAKRvgKOHIeMw5CbAmbojoQfvjuDbGTqPAO9O4NMJvMLAIwh0eg5FRTHEu2OjnwuVUBSlFbNYLcQXxHMo6xCHsw9zJOcIcXlxmKV2m8XV3pXObTpzfYfrCW8TTnib8LO3qC6VOomGJKUkIbuE/afzOZCUR3RKIaknU+icmUDXvNN0yz/NLfnJOJq0qwRp74Bdp064Db4Op65dcOzcGYdO4dj5ta36/JUVQcJWSNkDyXsgdT8Upf293qsD+PeEnndA28vAr5uWQOwcmuYEXIBKKIrSihjNRqKzo9mbsZd9GfuIzo6m2KT95ezh4EEPnx6M7zmebt7d6ObdjWD3YHSi9T7OXG62cig5n92JeexOzGXvqTwccjK5PDuOvnknGZV3Cp/8ilt+dnY4duuGy8g7cOrZE6cePXDsGIawr+H2XkkOJG6D039qr/Tov688fMIh7BoI7ANBfbRE4uTR+B/6IqiEoigtmMlq4nD2YXam7WRX2i4OZR3CZDUhEHT26syNYTfSx68Pvdv2pp17u1Z51VGZxSo5nFLAH3HZ7EzIYXdiLjpDCX2y4hhamMDDmbF45GoJRN+mDc79++PS736c+/bDqUd3dI6ONR+g3ACntkP8r3Dyd+32FYCdM4REwODpEBoJwf3AxbuRP23DUwlFUVqY9JJ0tqVsY3vKdnal7aLYVIxA0M2nG/d2u5cBAQPo49cHD4fm/dduU8kuLiPqeBZRxzPZHpdNnsFEUHEWN5fEMykjhraJxxBWKzoXF1wiI3G9YgIukZE4du6MqE1j1NyTELsRYjfBqR1apbmdk5Y4hr2oXYUE9YUW8KCCSiiKcomTUnI87zi/nf6N35J+42juUQACXAO4vsP1XBl0JZGBka2usrwmCVnFbDySzpaYDPYn5SOtkr6WXKYXx9Ijbg+OyVo7PscuXXCbOBHXwVfj0qcPwqEWdRVSareujq6Fo+sg65i23LcLDHwIOg2D9leCvXMjfkLbUAlFUS5BUkpicmPYlLiJLae2kFSUhEDQx68PT/R7gmtDrqVTm06t/hZWZScyilh3KI2Nh9OIzdDqjYa4l7Oo/DCdorcjTp8CnQ6XiAjcH7gXt2HDcAipQxeBmccg+is4/A3knQShg/ZXQb9x0HWU9iRWC6cSiqJcQk4VnuLHhB/ZcHIDiYWJ2Ak7IgMjmdBzAkNDh+Lj7GPrEJuV1PxSvj+QwtoDqRxLL0In4KpgF/7rlUznvVFYDh0AwGXAADwmPIj7yBHY+dThHBZnwaEv4OAayIjWkkjYtXD1E9B1NLhdsPurFkUlFEVp5orLi9mUuIkf4n9gf+Z+BIIBAQMY32M8I9qPULeyzmM0WfgpOo1v9iWzIz4HKaFfuza83tuRiIO/Ur5iE9JgQB8WhvdTT+J5003YBwbW/gBWC8Rtgb0fw4mfwWqG4P5ww+vQ41/g5tdon625UwlFUZohKSWHsw/zVexXbEzcSKm5lI6eHZnWfxqjw0bj7+pv6xCbnWPphazedZrv9qdQaDQT6u3Mf68N4+b8o+i/f4/SAwcod3LC48YbaTNmDM59+9TtlmBxFuz7BPZ+AgWnwdUPBv0H+t4Hbbs23ge7hKiEoijNSJmljA0nN/D50c85mnsUZztnbgy7kds7304v316qTuQ8JouVTUfSWbnjFH8l5uJgp+OGngGM7eFD591byJs/h/L0dBzat8f/2Zl4/utf6D3q+HRb+mHYuVSrH7GUabe0rp8LXW9sEU9mNSSVUBSlGSi0FLJo3yK+jv2avLI8wtuE80LkC4zuOBo3B7cL76CVKTFJFv8Wx8o/E8koLKOdtwvP39iN2zu5Yfl6NXlvrCKrqAiXgQMJmD0Lt2uuqd0jvmdIqbUT+eMtSPgN7F2g3/0Q+YjWxYlSJZVQFMWGEgoS+OTIJ6xNXosl2cK1oddyX7f7GBgwUF2NVCE5z8DybSdZvctAmeU4V4f78trtvbi6rT35H60g64nVSIMB9+uuw2fyZJx79qjbAaSE4xtg20JI2Qtu/jBiNvQfD85eF9hYUQlFUWzgWO4xPjz0IVtObcFB78Agt0HMHDmT9h7tbR1asxSXWcz7W+P5fn8KQsBAfzueH3MFXT305H7yMSdXfITVYMDjxhvxfeRhHDvX8SpCSq3hYdRrWkeMXh3gprfg8nvA3qlRPlNLpBKKojShmJwYlh5YSlRyFG72bkzqNYl7u91L9K5olUyqEJ9VzKJfTrD2YCqOdjruG9Seydd0JHbvDgK3/kj8u+9hyc3FfeQI2j7+eN0TCUDidtgyC5J3Q5v2cOsS6H0X6NXXY12pM6YoTeBE3gmWHFjCltNbcHdwZ0qfKdzT7R7V/Uk1TucYePuXWL7fn4KjnZ7J13TkocEd8XVzpPiP7fjMm0d6ahouERH4vb8U5969636QjBjYMhtObAL3ILh5EfS5R1W0XwSVUBSlEaWXpPPu/ndZF78OF3sXHr38Ue7vfj/uDu62Dq1Zyi4u471f41i16xQ6IZh4dRgPX9sJXzdHTKmpJD/7GkWbt4CvL8GL3sF95Mi61zWV5MBv82DvR+DoDiNehsiHW2RXKE1NJRRFaQRF5UUsi17GqphVAIzrMY6JPSfSxqnxh2G9FBlNFpb9nsD7W+Mxmq38OyKUJ0Z0xt/DCWkykbN8OVmLl4CUtJ02jcMdw+g1cmTdDmIxw57/acmkrBgGPARDZl6Svfo2V7YeU34U8A6gB5ZLKeeft/4tYGjFrAvgJ6VsU7HOAkRXrDstpbylaaJWlOpZrBa+OfENiw8sJs+Yx00db2Jq36kEuQXZOrRmSUrJ2oOpLNhwjNQCI9f38OfpUZfRqa32qLQxJobU51+g7OhR3IYPJ+C5Z7EPDoaoqLodKHkPrJ8G6Yeg41AY9Zo2OJXSoGyWUIQQemAxMBJIBnYLIdZKKWPOlJFSTqtU/jGgb6VdlEop+zRVvIpyIXvS9/DaX68RmxdLP79+LB2xlO4+3W0dVrN1OKWAWWuPsPdUHj2DPXjzrj4M6qj1o2UtLyd7yRJyli1H7+1F8LuL8KjrFQmAsVCrJ9mzAtwD4M6PofttoB7JbhS2vEIZCMRJKRMAhBBrgFuBmGrK3w3MaqLYFKXWskuzeWPPG6xPWE+gayALr13Ide2vU+1IqpFvKGfhz8f5fNdpvFwceP2O3ozpH4JOp50v4/FYUp9+mrLjx/H817/wn/kMes969FcW+zOsf0IbQjfyERj6XLMf8fBSZ8uEEgwkVZpPBiKrKiiEaA+EAb9WWuwkhNgDmIH5UsrvGytQRamKxWrhy9gvWbRvEWWWMh7q9RAP9X4IZztVuVsVKSU/HEhlzvoY8gzlPHBFB6aN7IKns/ZUlbRayV25kqw330Ln7k7IkiW4Dxt6gb1WoTQPNjyj9QLcthv8+1MI6d/An0apipBS2ubAQowBRkkpJ1XM3w9ESimnVlH2GSBESvlYpWXBUsoUIURHtEQzXEoZX8W2k4HJAP7+/v3XrFlTr3iLi4txc1NdYNRWSz9fqeWprM5ZTWJ5Il2duvJv73/jZ1//XmZb+vnKKLGyMqaMIzlWOnrqGN/DgXYe+rPrRVERnh9/jOORGMp696LgvvuQNfS5Vd358so9wGXHFmFvyud0uzGcan8nUqceA77Yn6+hQ4fulVJGXLCglNImL+AKYFOl+WeBZ6spux+4soZ9fQyMudAx+/fvL+vrt99+q/e2rVFLPV9l5jK5aN8i2eeTPnLw6sFybdxaabVaL3q/LfV8mS1W+cHWONnl+Z9kz5c2ypU7Tkqz5dzzVbJ7t4wdfI082qu3zP3881qdz3+cr3KDlD9Ol3KWh5TvRkiZvLcBP8Wl72J/voA9shbf67a85bUb6CyECANSgLHAPecXEkJcBngBf1Za5gUYpJRlQghf4Crg9SaJWmm1Dmcf5sXtLxKXH8ctnW5hRsQM9RhwDY6nF/H01wc5mFzAiG7+zL2tJwGef3djIqUkd8UKMt94E/vQEDq8vxqn7vV4iCHzGHw9ATKPQOSjMGKWalNiIzZLKFJKsxBiKrAJ7bHhFVLKI0KIV9Cy4dqKomOBNRVZ8oxuwAdCCCugQ6tDqa4yX1EuSrmlnCUHlvDRkY9o69yWJcOXMDhksK3DarYsVsmybQm8+XMsbk52LLq7Lzf3DjznIQWrwUDaCy9Q+NMG3K+/nsB5c9HX9ZaMlLBvpVZf4uAK934NnevxJJjSYGzaDkVK+RPw03nLXjpvfnYV2+0AejVqcIoCHM89znN/PEdsXiy3d76d6RHTVSv3GpzMLuGpLw+w73Q+o3oEMO9fPfFxczynTHlyMslTplJ24gR+05/Ce+LEOj8Rp7MY4buHtYr3jkPgXx+Cuxp0zNZUS3lFqYLFauGjIx+x+MBi2ji2YfHwxVwTco2tw2q2pJSs/iuJOetjsNcL3hnbh1suD/pHojDs20fylKlIq5XQDz/E7eqr6n6w7BP03zsDDEkw9HkYPB3qMtaJ0mhUQlGU86SXpPPstmfZk7GHke1H8tKgl1RdSQ1yist45ptothzN4OpwXxbeefk5dSVnFKxbT9pzz2EfFEToB+/j0KFD3Q92dD189zD2Ugf3fwed6vFYsdJoVEJRlEo2JW7i5T9fxmK1MOeqOdza6VbVQLEGf5zIZtqXBygwmHhhdDcmXBV2toHiGVJKspcuJXvRu7hERBD87iLsvOo4WJXVClsXwNb5ENyfvaH/4QqVTJodlVAUBSg1l/L67tf5OvZrevn2Yv7g+bTzaGfrsJotk8XKGz/H8sHv8XRq68bKCQPpFvjPdiPSYiF97lzyV6/B45abCZw7F52DQ90OVlYE3z4Mx3+EPvfC6Dcp276zgT6J0pBUQlFavbi8OGb8PoO4/Dgm9pzIlL5TsFeN4aqVlGvgsdX7OZCUz90D2/HSTd1xdtD/o5y1rIzU6TMo2rwZn0kTafvUU3W/2stPgtVjIfMojFqgdTOvrhibLZVQlFbtuxPf8equV3Gxd+GDER9wZfCVtg6pWfv5SDrTvzqIlLD4nn6M7h1YZTlrSQlJ/5mCYdcu/GY+g8/48XU/WPJeLZmYjXDvVxA+/OKCVxqdSihKq1RqLmXeznn8EP8DAwMGsuCaBfg6+9o6rGar3GxlwcZj/O+Pk/QK9uS9e/rS3se1yrKWwkKSJj9MaXQ0Qa8vwPOWeowsEbMWvn0I3Pxh3Drwu+wiP4HSFFRCUVqdxIJEpkVNIz4/nkcuf4RHej+CXvfPWzaKJq2glCmr9rHvdD4PXNGe50d3w9Gu6vNlzssjaeIkjCdOEPz2W/Xrcn7n+7BxJoQMgLtXg6tK9JcKlVCUVmXzqc28uP1F7HX2vD/ifXWL6wK2x2Xz+Or9GE0W3r27LzdfXv1AYebcXE6Pf5DyU6cIXfwebtfUsd2O1QqbX4Q/34PLboI7lqsuVC4xKqEorYLZauadfe/w8ZGP6eXbizeufYNAt6rv/yvao75Lt8azcNNxOrV1Y+l9/Qn3q75rFHNeHqcfnKAlkw/ex3XQoLod0FwOP/wHor+CgZNh1HxQV42XHJVQlBYv15jLjK0z+Cv9L+7qehdPD3gaB30dH11tRYqMJqZ/dZBNRzK4qXcgC+7ojatj9V8VZ5NJYiKhS5fUPZmUl8CXD0DcFhg+C66epp7kukSphKK0aIezDzMtahp5xjzmXjWXW8NvtXVIzVpcZjGTP93DqRwDL4zuxsSrw2p81NdSVETSxEmUJyQQsmQJrlfW8RaiIRc+vwtS9sDNi6D/uIv8BIotqYSitFjfx33PnD/n4Ovsy8obVqrx3S/g5yPpPPnlQRztdKyaFHl2fPfqWEtLSXrkUYyxsVqdSV375SrOhJW3Qc4JuPMT6F6Pp8GUZkUlFKXFMVlNLNy9kM+PfU5kQCT/d+3/4eVUx64+WhGrVfLOLyd455cT9A7x5P37+hPUpubKcFleTvJjj1O6bx/Bb76B27XX1u2g+Umw8lZtvPd7v9J6DFYueSqhKC1KrjGXp6KeYk/GHh7o/gDT+k/DTqd+zKtTZDQx7YuDbDmawZj+Icy9rSdO9jVXhkuLhZSnn6Hkjz8ImPMKHjfcULeD5ibAJ7eAsUDr4LFdHetclGZL/aYpLcax3GM8/uvj5BpzefXqV7m50822DqlZS8wu4aGVe0jILmHWzd0Zf2WHC3aNIqUk49XXKNq4Eb+nn8brzjvrdtDsE/DJzWAu0xosBvW5iE+gNDcqoSgtwsbEjbz4x4t4OnryyahP6OHbw9YhNWu/x2Yx9fN96HSCTycM5Mrw2jUezFm2nLxVq/B+8EF8JjxYt4NmHdeSidUC438Ef1Wn1dKohKJc0qzSynv732NZ9DL6tO3DW0PfUl2o1EBKyf/+OMmrPx2li787yx6IINTbpVbb5n//PVlvvonH6NH4zZhetwNnHtWSCUJLJqorlRZJJRTlklViKuHZbc/yW9Jv3N75dp6PfF61L6lBmdnC898d5uu9yVzfw583/92nxvYllZXs3EnaCy/iMmgQga+9iqjLCImZx+Djm0BnB+PXg2/nen4Cpbmz6biZQohRQojjQog4IcTMKtaPF0JkCSEOVLwmVVo3TghxouKlHl5vZZKKkrjvp/v4Pfl3Zg6cyewrZqtkUoPMQiNjP9zJ13uT+e/wziy9t3+tk0lZQgLJj/8Xx7AOhLy7qG7jmWTFalcmKpm0Cja7QhFC6IHFwEggGdgthFgrpYw5r+gXUsqp523rDcwCIgAJ7K3YNq8JQlds7K+0v3hy65Na9yAjlnJF0BW2DqlZi04u4KGVeygoNbH03n7c0Kv2Xc6Yc3NJmvwwwt6ekKXvo3d3r/2Bs0/AJzdp0+PWqWTSCtjyCmUgECelTJBSlgNrgNo2Y74e2CylzK1IIpuBUY0Up9KMrDm2hsmbJ+Pj5MPq0atVMrmAtQdTGfP+DvQ6wdePXlGnZGItLyd5ylTMWVmELlmMQ0hw7Q+cm6BdmUirdmXStks9olcuNbasQwkGkirNJwORVZS7QwhxDRALTJNSJlWzbZU/7UKIycBkAH9/f6KiouoVbHFxcb23bY0a+nxZpIWvc7/mj+I/6OHcg3Ee40jYl0ACCQ12DFtq6PNllZJvT5hYn2Cii5eOqX0EWbH7iYqt5Q6kxGPlpzjv30/+Q5NIz8uDWsbnaMyk7/7n0FuMHOgzj5IjaUBafT9KldTvY9001flq7pXy64DVUsoyIcTDwCfAsLrsQEr5IfAhQEREhBwyZEi9AomKiqK+27ZGDXm+8ox5PBn1JHuK9/Bgzwf5b9//trjxSxryfBWXmZn2xQE2J2RwV0Qoc27riYNd3W5G5Hz8MZl//onvlCl0e2zqhTc4ozAVPnoCRDlM/JEBgZfXMfraUb+PddNU58uWCSUFCK00H1Kx7CwpZU6l2eXA65W2HXLetlENHqFicyfyTvDYr4+RZchSjRVr4XSOgYdW7iEuq5jZN3dnXC0aK56v+I/tZL7+f7iPHIHvlP/UYcMsrTuVkmx44HtopGSiNF+2rEPZDXQWQoQJIRyAscDaygWEEJVv+N4CHK2Y3gRcJ4TwEkJ4AddVLFNakF9O/cK9P91LuaWcj0d9rJLJBeyIz+aWxX+QXmjk4wcHMP6qmnsKrkr5qVOkPPkkjuHhBM2fX/vHg0vz4bN/aX103fMFhETU4xMolzqbXaFIKc1CiKloiUAPrJBSHhFCvALskVKuBR4XQtwCmIFcYHzFtrlCiDloSQngFSllbpN/CKVRWKWVDw59wJIDS+jl24u3h76Nn4ufrcNqtqSUfLrzFC+viyHM15VlD0QQ5lv1eO81sRoMJE99DCEEIUsWo3Ot5T7KimHVnVp7k3vWQIc69jqstBg2rUORUv4E/HTespcqTT8LPFvNtiuAFY0aoNLkDCYDL2x/gc2nNnNLp1t46YqXcNQ72jqsZqvMbGHWD0dYszuJYZf58fbYPng42dd5P1JK0l58ibK4OEKXLcMhJKR2G5qMsOYebTyTOz+B8BF1PrbScjT3SnmlFUkqSuK/v/2X+Px4pkdM54HuD9T5lk1rkllk5NHP9rH3VB5ThnbiyZFd0evqd77yVq6k8McfaTttWu3HNbGY4ZuJcHIr3Pa+Gs9EUQlFaR52pu1k+tbpZxsrXhlUx5H/WpkDSfk88uleCkpNLL6nH6N71759yfkMe/eSUVEJ7zP5odptZLXC2sfg2HoYtQD63F3v4ysth0ooik1JKVkZs5I3975JR8+OLBq6iFCP0Atv2Ip9uSeJF747jJ+HI988eiXdgzzqvS9zTg4p057EISSEwNdeq90VoZTw8/Nw8HMY8iwMeqTex1daFpVQFJsxmAzM3jGbDYkbGNl+JHOvmouLfe16vm2Nys1W5v4Yw8o/T3FVuA/v3d0PL9f6918mLRZSZ8zAUlBA6LIP0bu51W7D3xfCziUQ+Shc+0y9j6+0PCqhKDZxuvA0T0Q9QVxeHE/0e4IJPSeo+pIaZBQa+c8qrb7kocFhPDPqMuz0F/fUf/aSpZTs+JPAeXNx6tq1dhvt/h/8Nhd6j4XrXwX1f6ZUohKK0uS2Jm3l2W3PotPpWDpiKVcFq8dMa7I7MZf/rNpHSZmZd+/uy82XB130Pkt27CB7yRI8b7sNz9tvr91Gh7+FH5+CLqPg1vegLl3YK62CSihKk7FYLSw5uIQPD31IN+9uvDX0LYLd6tDhYCtzZjCs1zYco523C6smRdLFvw69/VbDnJ1NytPP4NCpIwGzXqrdlWH8r/DtZGh3Bdz5Mejr/miy0vKphKI0iVxjLs/8/gw703Zya6dbeWHQCzjZOdk6rGaryGjimW8O8VN0Otd192fhvy+vV/uS80mrldSnn8FaVES7Ff9D5+x84Y1S9sKa+6BtV7h7NdjXYhulVVIJRWl0BzIP8NTWp8g35vPylS9ze+da3mJppY6mFTJl1T5O5Rp47sbLeGhwxwarX8pZ/j9Kduwg4JWXcepSiy7ls2LhszHg6gv3fQPObRokDqVlUglFaTRWaeWTI5/wzr53CHQN5LMbP6ObTzdbh9VsSSn5YncSs9YewcPZnlWTIhnU0adSAK4VAAAgAElEQVTB9m/Yv5+sd97B48YbaHPnnRfeoCAFPrtdG23x/u/APaDBYlFaJpVQlEZRbClm6i9T2ZayjZHtR/LylS/j7nDx9/9bqlKz5MkvD/Ld/hSuDvflrbv60Na94bqcsRQVkTp9BvaBgQS8/PKFr3gMuVoyMRZoA2T5dGqwWJSWSyUUpcHtTt/NgrQFlMgSnot8jrFdx6pHgmsQnVzA7B2lZJUamDaiC1OHhde7C5XqpL8yB1N6Ou0/+/TCw/iWl8Dnd0HuSe02l+qGXqkllVCUBmO2mll6cCnLDi3D186Xz274jO4+3W0dVrNltUpWbD/Jgo3HcLOD1Q8NIrIBb3GdUbB2LYXr1uH7+GO49O1bc2GLCb4a/3dnj2GDGzwepeVSCUVpEElFSTy37TkOZB3gtvDbuLrsapVMapBZaOSprw6y7UQ2I7r5c1tgUaMkk/KkJNJffgXniP74PvxwzYWtVvhhCpz4GW5+R3X2qNSZapmkXBQpJd+d+I4xa8cQnx/PgsELmHPVHBx1qsv56mw6ks71b//O7sRc5t7Wk2UP9MfNoeFvCUqzmdSnnwGdjuAFCxD6GoZNlhJ+fgEOfQHDXoD+4xs8HqXlU1coSr3llOYwZ+ccfjn9CwMCBjDvqnkEutW/19uWrtBo4pV1MXy9N5mewR68fVdfwv1q2X9WPeQsW0bp/v0ELVyIffAFGpBufxt2LobIR2Dw9EaLSWnZVEJR6mXzqc3M+XMOxaZinur/FA/0eACdUBe81dkel82Mrw6SXmhk6tBwHh/eGQe7xjtfpdGHyVq8BI/Ro/G8aXTNhfd9CltmQ88xcP1rqn8upd5UQlHqJM+Yx2t/vcaGkxvo7tOdeVfNI9wr3NZhNVtFRhOvbTjG57tO07GtK988eiV923k16jGtBgOpM2Zg5+tLwEsv1lz46DpY9zh0Gg63LVX9cykXxaYJRQgxCngHbUz55VLK+eetfxKYhDamfBYwQUp5qmKdBYiuKHpaSqlqEBuRlJINJzcw/6/5FJmKmNJnChN7TcRep/p0qs5vxzN57ttoMgqNTLo6jOnXd8XJvoZ6jAaSufANyk+dot1HH6H39Ky+4Mlt8PVECOoHd30KdvXvCl9RwIYJRQihBxYDI4FkYLcQYq2UMqZSsf1AhJTSIIR4FHgduKtiXamUsk+TBt1KpRWnMW/XPLYmb6WXby9evvJlOnt1tnVYzVZWURlzf4zhhwOphPu58fWjV9Kvka9Kzij+Yzt5n3+O97hxuA6KrL5g6gFYfTd4dYB7vwIH1yaJT2nZbHmFMhCIk1ImAAgh1gC3AmcTipTyt0rldwL3NWmErZzJamJVzCqWHFyClJLpEdO5r9t96HWN/1f2pchqlXy1N4lXfzqGodzM48PCmTIsHEe7pjlfloIC0p5/HodOnWj75LTqC2af0FrBO7fRulRx8W6S+JSWT0gpbXNgIcYAo6SUkyrm7wcipZRTqyn/HpAupZxbMW8GDqDdDpsvpfy+mu0mA5MB/P39+69Zs6Ze8RYXF+NW2xHtWoB4Yzxf5n5JqimVns49GeM9Bh+72reTaG3n61ShhU9jyonLt9LVS8e4Ho4EudW+PqIhzpfHio9w2rOH3Geexty+fZVlHI1Z9N0/E53VxP6+r1HqcmkOH9Dafr4u1sWer6FDh+6VUkZcqNwlUSkvhLgPiACurbS4vZQyRQjREfhVCBEtpYw/f1sp5YfAhwARERFyyJAh9YohKiqK+m57Kck0ZPLm3jf5MeNH/F38efuqtxnWblidu05pLecr31DOGz/HsmrXKbxcHPi/MZcxpn9Ik5+vwo2bSPnrL3wfm0r3ceOqLlSSDR/dAJTDhPVEXsJdqrSWn6+G0lTny5YJJQUIrTQfUrHsHEKIEcDzwLVSyrIzy6WUKRXvCUKIKKAv8I+EotSO0Wzks6OfsezQMkxWEw/1eohJvSapMd6rYbJYWbXzFG9tOUGR0cQDV3Rg2sgueDo3/UMK5uxs0mfPxqlnT3wnT666kLFAu82Vfxru+1b1z6U0ClsmlN1AZyFEGFoiGQvcU7mAEKIv8AHarbHMSsu9AIOUskwI4QtchVZhr9SRVVrZeHIjb+97m7SSNIaEDmFGxAzaebSzdWjNkpSSX49l8upPR4nPKuGqcB9eGN2dboEeNosnbfZsrAYDQQvmI+yrSGjlBq2zx4wjMHY1dFBDLiuNw2YJRUppFkJMBTahPTa8Qkp5RAjxCrBHSrkW+D/ADfiq4hbCmceDuwEfCCGsaN3HzD/v6TDlAqSU/Jn2J+/se4eYnBi6eXdj3tXzGBAwwNahNVt7T+WxYMMx/krMJczXleUPRDC8m59Ne1IuXLeO4i2/4DdjBo6dquhi3lwOX94Pp3fCmP9Bl+uaPkil1bBpHYqU8ifgp/OWvVRpekQ12+0AejVudC3XgcwDvLf/PXal7yLINYh5V8/jpo43qZbu1TiSWsDbW06wOSYDXzdH5t7Wk7sGhGKvt+35MmVkkD5nLs79+uE9vop6E4sZvpkAcVvg5kXQ846mD1JpVS6JSnmlYRzKOsSSA0vYnrodbydvZg6cyZ1d7sRBrxq0VeVYeiHvbDnBhsPpuDvZ8dTILky4OgxXR9v/2kgpSXvxRaTJRNCr8/7Z8aPVAt8/qrWEv/416F9NRb2iNCDb/2YojUpKyZ6MPSyPXs6O1B14OXoxrf80xnYdqyrcq3EwKZ/3fotjc0wGbo52PD4snImDO9qkwr06Bd9+S8nv2/B/7jkcOnQ4d6WUsH4aRH8Jw16EK/5jkxiV1kcllBbKYrWwNXkrKw6v4GDWQbydvHmi3xPcfdndKpFUQUrJthPZfPh7An/EZePhZMfjwzvz4JUd8HJtXldwprQ0Ml6bj8uAAXjdd++5K6WEDc/Avk+0XoOvUT0HK01HJZQWxmAy8EP8D3wW8xmni04T7BbM85HPc1v4bTjZOdk6vGbHaLKw7mAq//vjJMfSi/Bzd+SZUZdx36B2uDs1nyuSM7RbXS8hrVYCX52HqNyZ45kxTf76AAZN0cY1UZQmpBJKC5FYkMgXx7/gh7gfKDIV0du3N4/3e5zh7YZjp1P/zedLyS9l1c5TrNmdRG5JOV393Vl45+XccnlQo3Yrf7Hyv/6akj/+wP/FF3AIrdSMS0qtC/o/34OBk+H6eaobeqXJqW+aS1iZpYxfT//KNye+YVfaLux0dlzX/jruvuxu+vipfjPPZ7ZY+e14Fqv/Ok3Uca1Z04hu/oy/sgNXdPKx6eO/tWFKTSVz/gJcBg7E6+67/14hJfw2Txskq/+DcMPrKpkoNqESyiVGSklMTgzrEtaxPmE9BWUFBLkGMbXPVO7ocge+zr62DrHZOZFRxNf7kvl+fwoZhWX4uTvynyHh3DUglFDvS6M+SUpJ2kuzkFKee6tLSvh1LmxbCH3vh9FvqmSi2IxKKJeIpKIkNiVuYl38OhIKErDX2TM0dCh3dLmDQYGDVBuS86QXGFl/KJW1B1M5lFyAnU4wpGtbXrk1lGGX+dm8DUldFXzzjXar66UXcQgJ0RZKCb+8An+8Cf0egJveUQNkKTalEkozllSYxC+nf2FT4iYO5xwGoE/bPrx0xUtc1/46PB1rGDypFUorKGXj4XQ2HE5nd2IuUkKvYE9evKk7t/YJwtfN0dYh1ospLY2MM7e6xo7VFkoJW2bB9neg/3gY/ZZKJorNqYTSjFillZicGLYmb+WX079wIu8EAN19uvNk/ye5vsP1BLkF2TjK5kNKybH0In45msHmo5kcTMoHoIu/G9NGdOHmy4MI8720B446e6vLYiFw3lztVpfVChtnak9zRUyAG99QyURpFlRCsbE8Yx670naxPXU725K3kWPMQSDo69eXGREzGNZuGCHuIbYOs9koNJrYEZfN1thsfo/NIiW/FIDeIZ7MuL4ro3oG0Kltyxkno+Db7yjZtg3/55/XnuqyWmD9E7BvpfZosHqaS2lGVEJpYiWmEvZm7GVP+h52pe/iaM5RJBJ3B3euCrqKa0Ku4ergq/FyapohY5u70nIL+07n8Wd8DjviszmYXIDFKnFztOOKTj48NiycYZf54efR8trYmNLTyZg/H5eICLzuvQcsJvhhChz6Qmu0OOwFlUyUZkUllEYkpSTDkMHBrIMcyDzAvsx9HM89jkVasNfZ07ttb6b0mcKVQVfS3ae7GloXyCwysv90PvtO5fFXYi6HUwowWSR6neDyEE8evbYTgzv70q+91yVXsV4XUkrSZs1CmkzaU12WcvhqPMRu0LpTUS3glWZIJZQGlG/MJyYnhiM5RziSc4TorGgyS7X2Dk56J3q37c3EXhMZEDCAy9tejrOds40jtq2CUhNHUguITi7gUEoBB5PySc7TbmE56HX0DvFk0uCODOzgzYAwb9yaQaeMTaXg+x8o2fq71leXXxtYNQYS/4AbF8LAh2wdnqJUqfX8hjYgk8VEYmEi8fnxxObFcjzvOMdzj5NhyDhbJtQ9lAGBA+jl24vevr25zOcy7HXNryuPpmC2WEnMMRCbUcTx9CKOpRdyJLXwbPIACPFypneIJ+Ov7EDfdm3oEeSJk33rvGIzZWSS8dprOPfvj9dtI+GTmyE9Gm5fBr3vtHV4ilItlVBq4dfTv7Iubx3f/fodJwtPklSYhFmaAdALPWGeYUQERNDVqyvdfbrTzacbHg62GcHPVqSU5BlMnMwuITG7hKjYcr5I3ktcZjGJOSWYLBLQbvmH+bjSJ7QN90S2o3ugB71D2uDdzDpgtBUpJemzZiHLygiaMRnx0SgoSoexn0PXUbYOT1FqpBJKLXwf9z2/F/5Oe9Gejp4dGd5uOOFtwglvE06YZ1irGE9ESkluSTlpBUZS80tJyS8lOa+U5DwDp3NLSco1UFxmPlteJ6C9TxGd2roxvJs/4X5uXBbgTrifW6u98qiNwrVrKY6Kwu/R+3D4+UGQVhi3DkLVSJpK86cSSi3MuWoOe3bsYfjQ4bYOpcEZTRZyS8rJKS4nu6SM7KIysorLyCoqI7OojMxCI+mFRjIKyyg3W8/Z1tleT7CXM+28XYgM8ybU24UOPi6E+bqSEL2bEcOG2OQzXapMmZmkz3sV525heBe+B26+cP+34NvZ1qEpSq3YNKEIIUYB76CNKb9cSjn/vPWOwEqgP5AD3CWlTKxY9ywwEbAAj0spNzVWnJ6OnuhF8/yrWkpJqclCcZmZYqOZ4jIzRUYzhaUmioxmCkpNZ1/5pSbyDeXkG0zkGcrJKymnpNxS5X7dHO1o6+6Iv4cj/dt54e/hRKCnEwGezgS1cSK4jTPerg7Vdqh4WqceZ60TKUmf/TKytITAjvEIvx5wzxfgHmDryBSl1myWUIQQemAxMBJIBnYLIdZKKWMqFZsI5Ekpw4UQY4EFwF1CiO7AWKAHEARsEUJ0kVJW/e3YRKSUWKwSk0VSbrFiqniVm7X3MrM2XW7+e9potlBm0t6NJitGk4Uyk4VSkzZfarJQWm7BUG7GUG6peGnTxWVmSsrMWGXNcekEeDrb4+lsTxsXB3zcHAj3c8Pb1QFvVwe8XBzwdXPAx80RXzcH2ro74uKgLl6bktPOnRT/+it+fQpwjLgO7lgGDpd2K3+l9bngt4YQYjvwgpTytwY+9kAgTkqZUHGcNcCtQOWEciswu2L6a+A9of1JfCuwRkpZBpwUQsRV7O/PBo4RgFk/HOaXaAOOe6POJgyz1YrZIjFZrJitUpu2WpEX+HKvLWd7Pc4OepzsdDg76HFxsMPZQY+vmwMuji642OtxdbTDzdGu4l2Pu5M9bo52uDvZ4eFsf/bdzcEOnbpiaLZMp0/QZs2nOPmW4/3AOBg1F1SbJOUSVJs/Qx8GXhZCvICWWBrqSzsYSKo0nwxEVldGSmkWQhQAPhXLd563bXBVBxFCTAYmA/j7+xMVFVXnQAuzyvFxtOKgN6K3A70Q6HWgF9pf/3Y6gV7o0Ov02Fcs1wuBnQ5tXqeVsddRsUybtteBvf7vaYdK0/+8lWSueFWjrOIFlFa8Mqov3eiKi4vrda5bG+eSJPzffhXMktKxo9nqcj38vs3WYTV76uerbprqfF0woUgpDwN3CCH6Aa9UfNE9L6U82NjBNQQp5YfAhwARERFyyJAhdd7HkCEQFRVFfbZtrdT5qoXjG8j/v1mkJTliGj2Y3lMX0c3WMV0i1M9X3TTV+apL3xVxwBwgBdjbAMdOASqNYUpIxbIqywgh7ABPtMr52myrKM2T1QK/zsW04h4ydjvhfHkPckfffeHtFKWZu2BCEUL8JoRIRksizwG5wPgGOPZuoLMQIkwI4YBWyb72vDJrgXEV02OAX6WUsmL5WCGEoxAiDOgM/NUAMSlK4yrJhlV3Irf+H2nHeiCFI0ELFqru55UWoTZ1KE8BR6WUpRcsWQcVdSJTgU1ojw2vkFIeEUK8AuyRUq4F/gd8WlHpnouWdKgo9yVaBb4ZmGLrJ7wU5YISt8M3E8GQQ4H7BEqOb9S6pe/QARITbR2doly02tSh7Gusg0spfwJ+Om/ZS5WmjUCVnRdJKecB8xorNkVpMFYLbHsTol4FrzBMNy0mY9IzuERGat3SK0oLoRobKEpjyk+C7x6GU9uh5xjk6DdJffQJkJLAefO0ERgVpYVQCUVRGkv017D+SZAWuHUJ9LmHvFWfY9i5k4CXX8YhpMon3RXlkqUSiqI0tJJs+PEpiPkeQiPhXx+AdxhlJ0+SuXAhrtcMps2/VTf0SsujEoqiNKQj32nJpKxIG1nxqidAb4c0m0mdORPh6EjgnLnV9oGmKJcylVAUpSEUpMCGp+HYegjqq93i8u9+dnXO8uUYDx4i6I2F2Pv72TBQRWk8KqEoysWwWmD3cvjlFW16xGy44jHQ//2rZYyJIeu9xXjceAOeo0fbLFRFaWwqoShKfZ3eCT/NgPRD0Gk4jH4DvMPOKWI1GkmZ8TR2Xl74v/iijQJVlKahEoqi1FVhGmyZDYfWgEcwjFkBPW7Xxjc+T+Ybb1IeH0/o8uXYeXk1fayK0oRUQlGU2iorhh2LYMe7YDXD4Ke0VzXjlhRv307ep5/idd99uF19VRMHqyhNTyUURbkQcznsXwlRC6AkU7saGf7SP25vnbNJXh5pzz6HQ6dO+E1/qgmDVRTbUQlFUapjMUP0lxA1H/JPQbsr4O7VEBJR42ZSStJfmoU5L48OS5egc3JqooAVxbZUQlGU81lMcOgLrf+t3HgI7AOj34Tw4VXWk5wv/6uvKNq8Gb8Z03Hu0aMJAlaU5kElFEU5o7wE9q/S6kkKkiCgN9z1GVx2U60SCUBZwkkyXpuPyxWD8H7wwUYOWFGaF5VQFKUoHf76EHb/D4z5WncpN70F4SNqnUgAZHk5qdOno3N0JGj+fNXxo9LqqISitE5Swuk/4a9lcHSt1iix201ao8R2kfXaZeZbb2OMiSHkvXex9/dv4IAVpflTCUVpXQy5cOhL2LcSMo+AkycMfBgGTASfTvXebfHWreR+9BFe99yN+4gRDRiwolw6VEJRWj6LCeJ/0xoiHl0PljKtv62b34Fed1bbjqS2TBmZpM58FseuXfF75pkGClpRLj0qoSgtk9UKSbvgyLdw+FswZIOzN/QfD/3uh4BeDXIYabGQ+vTTWI1Ggt96E52jY4PsV1EuRTZJKEIIb+ALoAOQCPxbSpl3Xpk+wFLAA7AA86SUX1Ss+xi4FiioKD5eSnmgKWJXmjGLGU7vgGM/QswPUJQGekfoOgp6j9Uq2e0cGvSQ2e+/j2HXLgJffRXHjh0bdN+Kcqmx1RXKTOAXKeV8IcTMivnz7xUYgAeklCeEEEHAXiHEJillfsX6GVLKr5swZqU5MuRC/K9w4meI3aQ9paV3hM4joce/oMv14OjeKIcu+fNPst9bjMctN+P5r9sa5RiKcimxVUK5FRhSMf0JEMV5CUVKGVtpOlUIkQm0BfJRWi9zGSTvhoStkPAbpOwFadVuZ3W9AS4bDZ2GXXS9yIWYMjJJmT4Dh04dCZw9Ww2YVQWTyURycjJGo7HB9+3p6cnRo0cbfL8tVW3Pl5OTEyEhIdjb29frOLZKKP5SyrSK6XSgxmcshRADAQcgvtLieUKIl4BfgJlSyrJGiVSxrXKDljRO7YBT27VkYjKA0EFQP7jmae1qJKgv6PRNEpI0m0l96imsBgPtP/kYnYtLkxz3UpOcnIy7uzsdOnRo8IRbVFSEu3vjXHm2RLU5X1JKcnJySE5OJiys+n7qaiKklPXa8II7FmILEFDFqueBT6SUbSqVzZNSVtm3txAiEO0KZpyUcmelZeloSeZDIF5K+Uo1208GJgP4+/v3X7NmTb0+T3FxMW5ubvXatjWq1/mSVpxLU/EojMO9KBbPgmO4liSikxYkgmK3DhR49iDPqzcFnj0w29vm/8Pt2+9w/flnCh4cjzGyfm1WztcSf748PT3p1KlTo1y9WSwW9Pqm+QOiJajt+ZJSEh8fT0FBwTnLhw4duldKWXMndjTiFYqUstqH8YUQGUKIQCllWkVyyKymnAfwI/D8mWRSse8zVzdlQoiPgOk1xPEhWtIhIiJCDhkypM6fBSAqKor6btsaXfB8mYyQfRwyjkDaIUiP1gaqKivU1tu7QHB/6HsbhAxEtIvE3dkLdyCkKT5ANQo3/UzKzz/T5q676NaAjwi3xJ+vo0eP4uHh0Sj7VlcodVOX8+Xk5ETfvn3rdRxb3fJaC4wD5le8/3B+ASGEA/AdsPL8yvdKyUgAtwGHGz9kpV6MBZATD9kntASSHQuZx7ROF6VVK2PvAv49tDYhwf2021e+Xc8ZRrc5KEtIIO3ZZ3G6vDf+zz9n63AUpdmx1W/sfOBLIcRE4BTwbwAhRATwiJRyUsWyawAfIcT4iu3OPB68SgjRFhDAAeCRJo5fOcNihsIUrTPF/NOQdwryEumbeAB250BJ1t9lhR68O0LbrtoTWP7dwa+H1kK9ieo/6stSXELy1McQTk6EvPMOOoeGffxYUVoCmyQUKWUOMLyK5XuASRXTnwGfVbP9sEYNUNEaBhrzoTij4pWptesoTKt4T4GCFChO//tKAwABniFYdZ7QaRT4hGsJw6ezlkwauB1IU5BWK2nPzqQ8MZF2K1ZgH1BV1aDSXH3wwQfs37+f999/H5PJxIQJE9DpdCxfvrzeTzOdsXHjRv773/9isViYNGkSM2fOtFksABMmTGD9+vX4+flx+HDT37hpXvcUlIYnJZhKtboJYwGU5muJojQfSnOhNE9ry2HI0VqTG3K1q4qSbJCWf+7PwQ3cA8AjCDoOAc9g8AyFNqHg2U57t3PkYAuqE8hevISizVvwm/kMroMaphJeaTrR0dH07t2bwsJCbr/9diIjI5k3b95F79disTBlyhQ2b95MSEgIAwYM4JZbbqF79+5NHssZ48ePZ+rUqTzwwAMNts+6UAmlObBatPYVZmOld6OWCEylYK54NxnBVKJNl5doj8+WG6C8WHuVVXovK4SyIu1lNdVwcKF1kOjqCy4+0KadVo/h2lZ7ufmBm7/2cg9otEaCzVXhpp/JXrwYz9tuw3vcOFuHo9TDoUOHuPrqqxkyZAgPP/wwDz/8cIPs96+//iI8PJyOFT0kjB07lh9++OFsQrnxxhtZvnw5QUFBjR7LGddccw2JiYkNus+6UAmlNg6uISxhM5h+1b78rWatw0GruWLedO4yS/l50xXz5jLt3VKudVBorni3musXl9BrVwwOrtrL0U2b9wwBJw/ty9/RXUsYjh7au3MbcPLSpl28tfdmXn9hK8bjx0mtqIQPeFk1XrxURUdH89hjj7FixQpuvvnmWm0zePBgioqK/rF84cKFjKjoTTolJYXQ0NCz60JCQti1a9fZ+Z9++qnJYmkuVEKpjSPfEZq0GVLtQWenfZHr7bTpMy+9PejsK5bba/N6B+0LXWev1R3oq3jZOYKdk7bezqnSvBPYO1e8u2jT9s5/Tzu4aturL7lGYc7KIunRR9G7uRHy7ruq08eL9PK6I8SkFjbY/iwWC71CvZh1c81DLCclJeHm5kbnzp1JS0s7Z11ERASDBg3i6NGjLFq0iB6Vhmvetm1bg8Vam1gA+vfvz4ABAwCYOHHi2enGiKWxqIRSG/d8we8tqE5AqZnVaCRpylQsefm0//RT7P38bB2SUk/R0dFcfvnlLFu2jEGDBjFgwAD69u1LUlISAwcO5L333uOtt94iOTn5nIRSm6uC4OBgkpKSzq5LTk4mODi4zrGAlmwGDBjA+++//4/t1BWKolyipNVK6sxnMUZHE/LuIpx71vwXsFI7F7qSqKvaNtQ7dOgQvXr1IjAwkOXLl3PXXXexe/du9u7dS2xsLI888ghpaWlMmzbtnO1qc1UwYMAATpw4wcmTJwkODmbNmjV8/vnnZ9cPHz6clStXnk0y1cXi6enJ3r17OXr0KI888ggBAQHMnj27TrE0F2rQa0WpJOvtdyjauBG//2/v3qOjqq8Fjn93MgkhgOEdMEN4CCIgUCRRaXsrtNJo9IryyJVSFIMibe1dvWor6q1dVatoL+CyWINCC6tF0Yu9rcWCIgX7uC0I6DUYQCKBJCQEEh4hgSSTye/+MZOYhEwyM5yZM8nsz1qzkvM7j9nZzJrNOXPmtx9+WDsvdgG5ubmMH+/pfTN9+nSysrLIzs5mz549LFu2jJycHBwOB9XV1QEf2+FwsHLlSjIyMhgzZgxZWVlNZzkNDQ3k5+fTt2/fDmMB2LNnDy+88AI5OTktikmg5s6dy5QpUzh48CBOp5M1a9YEfaxg6BmKUl6nXnuNildeoXdWFn2z77E7HGWB9evXt1h++umnAc8dWOXl5cTExDBmzBh69AhudurMzEwyMzMvGs/Ly2PWrFl07969w1jAU1BKS0txOBxkZGRwxx13BBXP66+/HtR+VtGColVgCL8AABR8SURBVBRQuXUrZU89Tc+pUxn0xI/1jq4urq07sKx09dVXs3z5cr+3D3U84aKXvFTUO7/3I0oe/iEJE8aTsnwZ4tD/ZykVDC0oKqrVHPyMosWLiRs0iCEvv6y9TZS6BFpQVNSqO3qUwnsXEtO9O0PWrMHR7ANUpVTg9NxeRSVXWRmF2QvBVU/q+rXEO31/f0Ap5R8tKCrq1JeXU5i9EPeZM6SuW0e3K66wOySlugQtKCqq1J86ReE99+AqKSH1lVX6xUWlLKQFRUWN+tOnKVxwD3VFxQzJySHRO1eSUsoaWlBUVKg/fZrC7IXUHTnCkJyXta+JUiGgBUV1ea4TJyjMzsZVVIzzpZfo8eUv2x2SUl2SLbcNi0hfEdkqIoe8P/v42M4tIh97H283Gx8uIjtFJF9E3hCRztdXVoWFq7SUwvl34SopZciqVfT8l6/aHZIKs1WrVrF48WIAXC4X8+fP5+6778blaq/xnH+ys7MZOHAgV199te2xBBOP1ez6HsoSYJsxZhSwzbvclgvGmC95H7c1G38OWGGMGQmcBhaGNlzVGdUdOcLRed+mvqKC1NWr9TJXlGredvfmm28mNTWVdevWWdLDfcGCBWzZsiUiYgkmHqvZVVBmAOu8v68Dbvd3R/FMsvR1YGMw+6vocCF3H0e+NY+GCxdIXbuWxGsm2R2Sssknn3xC3759mTp1KnPmzLG0h/vXvva1FjMKN5eZmUlJSUnYYukonnCw6zOUZGNMY8uy40Cyj+0SRGQ3UA8sNcb8HugHnDHGNPbNLQb0W2mqSdXf/86x7/87sX36MGT1q3QbPtzukNTmJXA817LDdXfXQ8okuHlph9va1XZXWwBbSETeBwa1serx5gvGGCMixsdhhhpjjonICODPIpILnA0wjkXAIoDk5GR27NgRyO5Nqqqqgt43GtmVr4R//IPLfrue+kGDOPP9Byg+ehSOHg17HIHqiq+vpKSkpjfCbq46Ytz1HewRAAN1rjpq23ijba64uJgePXpwxRVXUFBQ0OKN+YYbbiAtLY3PPvuM559/njFjxjSta2/239Zv7lVVVTQ0NLT5pu9vLOA5u7jmmmsAmD9/PpMnTw44Fl/xuN3uDuNrVFNTE/xr0RgT9gdwEBjs/X0wcNCPfdYCswEBygGHd3wK8K4/zzt58mQTrO3btwe9bzQKd74a3G5TtmKFyRt9lTmyYIGpP3s2rM9/qbri6ysvLy9kx66srPRru3feecfccsstpqSkxKSmppq9e/caY4wpLCw03/nOd4wxxixfvtxs2bKlxX5f/epXzcSJEy96bN269aLnKCgoMOPGjQs6lsZ47r///jb3CyQWX/H4my9j2v53A3YbP95j7brk9TZwN7DU+/MPrTfw3vl13hhTKyL9ga8AzxtjjIhsx1NcNvjaX0WPhpoaSh97nMo//Ynec2Yz6IknEIs+5FSdWyhbAHdEWwCHz1JguogcAm70LiMiaSKy2rvNGGC3iPwfsB3PZyh53nWPAA+KSD6ez1TC2+dSRQxXSQlHvzWPys2bGfDQgwx68kktJqpJKFsAg++Wu9oCOIyMMRXAN9oY3w3c6/39f4HxPvY/DFwbyhhV5KveuYtjP/gBxuXC+cuX6DVtmt0hqQgT6hbAvlruagtgpToJYwynfvVrTixfTvzQoThXrqTbCL2TS/lPWwCHhhYU1am4z5yh5NHHqNq+nV7TpzP42WeI7dnT7rCUUmhBUZ2Ip/f7w7hOniT5scfoM//beL7nqpSKBFpQVMQzLhflL79Mec4q4i6/nGHrf0v3CRPsDksp1YoWFBXRag8XULJkCTWffELSHXeQ/PhjeolLqQilBUVFJFNfz6m1azn54i+Q7t1JeeEFLrspw+6wlFLt0IKiIk7NwYOU/uePqcnNpdf0Gxn0xBM4BgywOyylVAe0oKiI4a6qpnzlSk795jfEJiWRsmI5vW66ST94V6qT0IKibGeM4dyWLZQtfY76sjJ6Z2Ux8MH/ILZ3b7tDU0oFQAuKstWF3H2UPfssF/bupdtVV5HywgoSJ2nvEqU6I7vm8lJRrq6oiGM/+hFH5syh7uhRBj31JMPf2qjFRFkulG13t2zZwujRoxk5ciRLl3bcm8WuFsBbt24NKM5gaUFRYVV/8iTHn3yKzzNv4dy779Hvvvu44t0t9JkzB4mNtTs81QWFqu2u2+3me9/7Hps3byYvL4/XX3+dvLy8dvexowWw2+3moYceCijOYGlBUWHhKivj+M+eIf/G6Zx+8016z5rJFe+9x8CHHtTvlaiQClXb3V27djFy5EhGjBhBfHw8d955J3/4wxedNCKlBfCuXbsYMWKEzzitpJ+hqJCqLSjg1Np1nP3d7zDGkDTjNvovWkT80KF2h6bC6Lldz3Hg1AHLjud2uxk3YByPXPtIh9uGqu3usWPHGDJkSNM6p9PJzp07m5YjpQXwsWPHcDqdPuO0khYUZTljDHGHDlG0cSNV2/6MxMWRNHMm/e67j3hnit3hqShSVFREz549GTVqFKWlpS3WpaWlcf3117N//35efPFFxo0b17QuFE2t2osFYPLkyaSnpwOwcOHCpt87U4MtLSjKMg0XLnB20yZOr3+NvgcOcCEpiX6L76fvvHk4+ve3OzxlI3/OJAJx7tw5evXq1eF2ubm5TJw4kVdffZXrr7+e9PR0Jk2aRFFREddeey0rV65kxYoVFBcXtygo/pwVpKSkUFRU1LSuuLi4qTtjILGAp9ikp6eTk5Nz0X6XeoaSkpJCcXGx33FeCi0o6pIYY6j5NI8zb22kctM7NJw7R7fRo6mcN4/0hx4kJjHR7hBVFAtlC+D09HQOHTpEQUEBKSkpbNiwgddee61pfaS0AE5PT+fw4cM+47SSFhQVFFdZGZWb3uHsH/9I7YEDSLdu9Mr4Jn2ysug+eTKFH3ygxUTZLjc3l1tuuQVo2XZ37NixLFu2jIkTJzJr1iyqq6sD7trocDhYuXIlGRkZuN1usrOzm85yfLUAbiuWt956q6kF8KRLvG1+7ty57Nixg/LycpxOJz/96U9ZuHAhP//5z9uM02q2FBQR6Qu8AQwDjgBZxpjTrbaZBqxoNnQVcKcx5vcisha4ATjrXbfAGPNxiMOOevXl5Zx7/30qt7zL+Z07wRgSJkwg+Ykfk3TrrcRedpndISrVQqhbAGdmZpKZmXnReKS1AM7IyGD27NlBHTMQdp2hLAG2GWOWisgS73KLi6zGmO3Al6CpAOUD7zXb5IfGmI1hijdq1RUWUrV9O+e2/Znzu3dDQwPxw4bR/7vfJelfbyV+2DC7Q1QqYNoCODTsKigzgKne39cBO2hVUFqZDWw2xpwPbViqobaW87t3U/3Xv1H1t79Sl/85AN1GjaT/4sX0ysig25WjdMJGpdRF7CooycaYxvvmjgPJHWx/J9C63P9MRJ4AtgFLjDG1FscYFYzLRU1eHtX/3Mn5nf/k/N6PMDU1SFwcielp9MnKoue0acQ3u99eKaXaIsaY0BxY5H1gUBurHgfWGWN6N9v2tDGmj4/jDAY+AS43xriajR0H4oFXgM+NMU/62H8RsAggOTl58oYNG4L6e6qqqujZBb7RLdXVxB05Stzhw8Tl5xNfUIDU1QHgSknBdeWV1I4dS92Vo6Bbt6Cfp6vkK1y6Yr6SkpIYOXJkSI7tdruJ1al6/BZIvvLz8zl79myLsWnTpu0xxqR1tG/IzlCMMT5vkBaRMhEZbIwp9RaHE+0cKgv4n8Zi4j1249lNrYj8Gni4nThewVN0SEtLM1OnTg3gr/jCjh07CHZfu7jPnaP2wAEufPopNXl51OTuo66gwLMyJoZuV40mMSuLxLTJJKan4+jXz7Ln7oz5slNXzNf+/fv9+q5IMPz9HoryCCRfCQkJQd9tZtclr7eBu4Gl3p/tTSwzF3i0+UCzYiTA7cC+UAXaGTTU1lJXUEBt/ufUfp5P7WeHqD1wANexY03bOAYOJGHcOJJmzKD7xAkkjB+vc2gppSxlV0FZCrwpIguBo3jOQhCRNGCxMeZe7/IwYAjwQav914vIAECAj4HF4QnbPg21tbhKSnAVFVF3tJC6wkLqjhyhrqDAUzgaL13GxhI/dCjdJ06gd1YWCVeNJmHsWG2hq5QKOVsKijGmAvhGG+O7gXubLR8BLpojwBjz9VDGF27G7aa+ooL6shPUlx3HVVZGfWkprpJSXKWluI4do/5Ey6uCMYmJxA0bSvcJE0iaMYP4EcPpNnIU8cOHERMfb88fopSKavpN+RAwDQ00VFbiPnMG95kz1J85g/vUadynT1F/6hTuilOeAlJRTv3Jk7grTkFDQ4tjSFwcjsGDiRs0iB5f+QpxzhTinU7ihqQSnzqE2H799NZdpVRE0YLih/MffUS3Dz/kdNkJGqqraaiqoqG6CndVFQ3nqnCfq6Sh8hzuykrclZU0VFZ+cQmqFYmPJ7Z/Pxx9+xE3YCDdx43DMWCA55GcjGNgMnHJAz0FI0bb1SilOg8tKH4oz8mh9wd/4XizsZgePYjp2ZOYXj2JvSyJ2P79iB8+nNikJO/jMmJ79/7i0bcvsX36EtMjUc8slAqjVatW8dFHH5GTk4PL5SI7O5uYmBhWr159yZ0Ss7Oz2bRpEwMHDmTfvo7vDQplLO3Fs3XrVh599FHcbjf33nsvS5YsueTnaosWFD8Mevxxdk6dyvXTpnkKSWKitqtVqpNo3nZ35syZXHfddZZ1SlywYAEPPPAAd911l+2x+IqnsQXwtm3bcDqdpKenc9tttzF27FjLnreRXlPxQ3xqKm7v5xmxvXppMVGqEwll2922Wu420hbASikVAsefeYba/da1AK53u6m+ehyDHnusw23taLsL2gJYKaW6FG0BHF5aUJRSIefPmUQgIqEFcKC0BbBSSnVioWwB3JFobAGsH8orpbqs3Nxcxo8fD7Rsu7tnzx6WLVtGTk4ODoeD6urqoI4/d+5cpkyZwsGDB3E6naxZswbw3QK4rViAphbAOTk5LYqJFfE4HI6mFsBjxowhKyura7UAVkqpcAh1C2BfLXe1BbBSSkUJbQEcGnrJSymllCW0oCillLKEFhSllFKW0IKilAoZ42PWbRWZLvXfSwuKUiokEhISqKio0KLSSRhjqKioICEhIehj6F1eSqmQcDqdFBcXc/LkScuPXVNTc0lvfNHG33wlJCS0mPcrUFpQlFIhERcXx/Dhw0Ny7B07djRNW6I6Fq582XLJS0TmiMinItIgImntbHeTiBwUkXwRWdJsfLiI7PSOvyEi2kRdKaVsZtdnKPuAmcBffG0gIrHAS8DNwFhgrog0doR5DlhhjBkJnAYWhjZcpZRSHbGloBhj9htjDnaw2bVAvjHmsDGmDtgAzBBP/9yvAxu9260Dbg9dtEoppfwRyZ+hpABFzZaLgeuAfsAZY0x9s3GfczGLyCJgkXexSkQaC1kScLbV5q3Hmi/3B8oD/Bv81VYsVu3T3na+1vmTm7bGNF+BjUVyvvzdz6p8tTUebflqb32gr6fWy5ear6F+bWWMCckDeB/Ppa3WjxnNttkBpPnYfzawutnyfGClNzH5zcaHAPuCiO+VjsaaLwO7Q5iri2Kxap/2tvO1zp/caL66dr783c+qfHWUn2jIV6A5i5R8NX+E7AzFGBNcF5ovHMNTLBo5vWMVQG8RcRjPWUrjeKD+6MdYW9uEQjDP4+8+7W3na50/uWlrTPMV2Fgk58vf/azKV1vj0Zav9tYH83oKV76aiLd62UJEdgAPG2N2t7HOAXwGfANPwfgQ+JYx5lMR+W/gLWPMBhHJAT4xxvwyxLHuNsb4vCNNtaT5CozmKzCar8CEK1923TZ8h4gUA1OAd0TkXe/45SLyJwDv2ccDwLvAfuBNY8yn3kM8AjwoIvl4PlNZE4awXwnDc3Qlmq/AaL4Co/kKTFjyZesZilJKqa5D5/JSSillCS0oSimlLKEFRSmllCW0oFhARG4XkVe984p90+54Ip2IjBCRNSKyseOto5OI9BCRdd7X1Ty744l0+poKTKjes6K+oIjIr0TkhIjsazXe5sSUbTHG/N4Ycx+wGPi3UMZrN4vyddgYE3XzrwWYu5nARu/r6rawBxsBAslXtL6mmgswXyF5z4r6ggKsBW5qPuBrYkoRGS8im1o9Bjbb9T+9+3Vla7EuX9FmLX7mDs8XdhunHnKHMcZIshb/86WCy5el71mRPJdXWBhj/iIiw1oNN01MCSAiG/BMGfMscGvrY3gnrFwKbDbG7A1txPayIl/RKpDc4Zmjzgl8TJT+xy/AfOWFN7rIE0i+RGQ/IXjPisoXqh/ampjS5wSUwPeBG4HZIrI4lIFFqIDyJSL9vDMcTBKRR0MdXITzlbvfAbNE5GVsmEIjgrWZL31N+eTr9RWS96yoP0OxgjHmReBFu+PoLIwxFXiu3SofjDHVwD12x9FZ6GsqMKF6z9IzlLb5mphStU3zFTzNXWA0X4EJa760oLTtQ2CUt9VwPHAn8LbNMUUyzVfwNHeB0XwFJqz5ivqCIiKvA/8ARotIsYgs7GBiyqim+Qqe5i4wmq/AREK+dHJIpZRSloj6MxSllFLW0IKilFLKElpQlFJKWUILilJKKUtoQVFKKWUJLShKKaUsoQVFKaWUJbSgKKWUsoQWFKVsICLbRWS69/enReQXdsek1KXS2YaVssdPgCe9DccmEaVdGVXXolOvKGUTEfkA6AlMNcacazb+BvChMea/bAtOqSDoJS+lbCAi44HBQF2rYjID2ASMtys2pYKlBUWpMBORwcB6PK1rq0TkJu94AjDHGPMbIMnGEJUKihYUpcJIRBLxtPd9yBizH3gKz+cpAD8Eenpb2Y4Tke42halUUPQzFKUigIikAj8xxiz0Lv8E2GKM2WlvZEr5TwuKUkopS+glL6WUUpbQgqKUUsoSWlCUUkpZQguKUkopS2hBUUopZQktKEoppSyhBUUppZQltKAopZSyhBYUpZRSlvh/9WoXgyOirkoAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "KK_E = [1,10]\n",
    "KK_B = [0,1]\n",
    "Legend = []\n",
    "x_A = np.logspace(-2,2,100)\n",
    "x_B = 1\n",
    "for K_E in KK_E:\n",
    "    for K_B in KK_B:\n",
    "        Legend.append('$K_B = '+str(K_B)+'; K_E = '+str(K_E)+'$')\n",
    "        v = ECR_flow(x_A,x_B,K_E=K_E,K_B=K_B)\n",
    "        plt.semilogx(x_A,v)\n",
    "plt.legend(Legend)\n",
    "plt.ylim((-1,1))\n",
    "plt.xlabel('$x_A$')\n",
    "plt.ylabel('$v$')\n",
    "plt.grid()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note:\n",
    "    \n",
    "    1. When $K_B=0$, the flow $v$  depends on $x_A$ but not $x_B$\n",
    "    2. As $K_m = \\frac{K_c}{K_e}$ the two cases correspond to $K_m=1$ and $K_m=0.1$ respectively. A smaller $K_m$ leads to a steeper curve. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Simulated flows\n",
    "The theoretical values can be confirmed by simulation. However, the simulation is dynamic whereas the theoretical results correspond to the steady-state situation. This dichotomy is resolved by:\n",
    "\n",
    "1. Starting x_A at a constant value to allow the transient effects of initial conditions to die out.\n",
    "2. Then varying $x_A$ slowly:\n",
    "\\begin{equation}\n",
    "x_A = \n",
    "   \\begin{cases}\n",
    "   x_{min} & \\text{if $t<t_0$}\\\\\n",
    "   x_{min} + (x_{max}-x_{min})\\frac{t-t_0}{t_{max}-t_0} & \\text{otherwise}\n",
    "   \\end{cases}\n",
    "\\end{equation}\n",
    "3. It is not possible to set $K_B$ to zero as $\\ln K_B x_B$ is not defined in this case. Hence a small value is used instead."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.01 + (100.0-0.01)*np.heaviside(t-100.0,1)*((t-100.0)/99900.0)\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZgAAAEPCAYAAAB/WNKuAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzs3Xd8k1XbwPHfSZruQVvoLhRo2bsIZWkRB+IAFcVHRVBwPIIiiqKgKCAIbl9FAXnYKogTFUFRCip7bwq00L33TpPz/pGCBVratE0b2vPlk09y7yuHNFfu+z5DSClRFEVRlLqmaegAFEVRlMZJJRhFURTFIlSCURRFUSxCJRhFURTFIlSCURRFUSxCJRhFURTFIlSCURRFUSxCJRhFURTFIlSCURRFUSzCpqEDaEjNmzeXQUFBNdo2Pz8fJyenug2oEVPlZR5VXuZR5WWe2pTXvn370qSULaqzbpNOMEFBQezdu7dG20ZERBAeHl63ATViqrzMo8rLPKq8zFOb8hJCnK/uuuoSmaIoimIRKsEoiqIoFqESjKIoimIRKsEoiqIoFmFVCUYIsVQIkSKEOFrJciGE+D8hxBkhxGEhRK9yy8YIIU6XPcbUX9SKoihKRawqwQDLgaFXWX4bEFL2eAL4DEAI4QG8DvQF+gCvCyHcLRqpoiiKclVWVU1ZSrlNCBF0lVWGAyulaRjOnUKIZkIIXyAc+F1KmQEghPgdU6L6yrIRK4qiWAejNGKQBqSUGKTh4rTRaMSIEaM0IqXEKI1klmbWS0xWlWCqwR+ILTcdVzavsvlXEEI8gensB29vbyIiImoUSF5eXo23bYpUeZlHlZd5alteRmlEL/UXH6Wy9JLpyx8XlhukgVJZigGD6csc45XzpPHi61JKL37xX5h3YbuLSQCJkbJkUPbvwmuJKUFceF1+HXPo0OEeYfmLPNdagqk1KeViYDFA7969ZU0bG6mGXeZR5WUeVV5VM0ojuSW55JTksGX7FjzaeZBTkkNuSS4F+gLTo7SAfH3+xedCfeGV80oLKTWW1joeG2GDjcb00Gl0F1/baC+dZ6uxxUaUTWtt0AkdWo0WjdBcfGiFtsLX5acrna/RIhBXTAs0FOmNFOslZ8+eq5fP17WWYOKBwHLTAWXz4jFdJis/P6LeolIUpdaKDcWkF6abHkX/PmcUZZBemE5WcRY5JTnkFOdcTCQS+e8Okq7cp73WHkedI442jjjpnHDUOeJq64qPk8/FeQ42DtjZ2GGnreRRyTJbrS06jQ6dVoeNsEEIUX+FBegNRtLyiknOKSYlu5Ds1AzyE5MpSUunJD0DmZmCyM5Cm5uNXV4OjgU5eJfk41acTw9phPunWzzGay3BrAcmCiHWYLqhny2lTBRCbALmlruxfwvwSkMFqSjKpfJK8kjKTyKpIInE/ETT6/wkkvOTSSpIIq0wjXx9foXbOuuc8bD3oJl9MzzsPWjt1hpXW9eLDzc7N2IiY+gf2h83WzdcbF0uJg6tRlvP77RuFJcaSMgqIjajgLjMQlITUymKicWQkABpqdhkpuGUk4FnYTaehdm0KMom0KCvcF9F9k6UOLticG2GdPNC08yNHH0eRqMRjcay9bysKsEIIb7CdCbSXAgRh6lmmA5ASrkQ2AAMA84ABcCjZcsyhBCzgT1lu5p14Ya/oiiWJ6UkvSid8znnicmJITY3lvM554nNjSUuN45cfe4l62uEhhYOLfBx8qGDRweaOzTH094TTwdPPO098bD3wNPB9GxvY1/l8SPiIgj1DrXU27OI7EI9Z1PziErN51xqHulxSeijoxAJcTilJ+OTn45vfjpt8zPoqS+4ZFuDVkuJmycGj+Zo23XB1tcHez9fnP19cfRujtbFCa3IxkafhMiNgYxoyDwHmbshJ55SjT2aejjjsqoEI6X8TxXLJTChkmVLgaWWiEtRFBMpJUn5SZzOOs2ZrDOcyTzDmawznMs5R2Fp4cX1bIQN/i7+BLoE0r1Fd/yc/fBx8jE9HH1o4dgCG41Vff1YTGpuMccTczidnMvZ1HySz8VjOHuGZilxtMxNpmVuMoNyk3Et+TeJGDVaSlt4o2kTgGOrMFzbtMKuVUt0fv7ofLzRenggNBrQF0LycUg9CWmnIPVP2H0Sss6DLHfj38kL3IOgVX9wD+J0chEdpRGEZc/wmsb/sKIoZis1lnI26yzH0o9xLO0YkZmRnMk6Q54+7+I6Xo5ehDQLIdQ7lECXQFq5tqKlS0t8nH3QaXQNGH39MxolUWn5HE/M4XhCDicSskmMisU9LorgrDhCsuK4Myce98Kcf7dxdkHbpi0uHa7DITgY27ZtsG0VhM7HG2Fz2ddzUQ4kHYHTWyHxkOmRFgnSYFqutQXPYPDtDt3uh+btoEV7cG8Nds6X7Co5IoKO9XD5UCUYRVEASMhL4EDKAY6mHeVY+jFOZpy8eFbirHOmnXs7bm9zOyHNQghxD6Fts7a42bk1cNQNJ7tAz4HYTPbHZHEgJpPD59LxSomhc3o0XTLO8XTWOVwLTMlEajTYBAXhHBaOfadO2LVvj11wMFpPz4orB0gJaachZifE7jI90iL/Xe7iCz7doOOd4NsNWnQ0naForesr3bqiURSlXkgpic6JZl/yPvYn72df8j4S8xMBU82rDh4duDfkXjo370xnz860cm2FRlhbxx/1KzW3mB1R6ew4m8aec5lEJWXTLiuW0JRIRubG8GJqNLqSYgBsAgJwvOkGHLp2w75zZ+w7tEfj6Fj5zqWE5GMQtQXObzcllIJ00zIHdwjsazor8e1hSiwu3vXwjmtPJRhFaSLSCtPYnrCdf+L/YWfiTjKKTPVgPO096eXdizGdx9DLqxch7iFN5v7I1eQU6dl5Np3tZ9PZfjaNyOQ8vPPT6Z9xhmdyo2kdewJdYT4IgV3HDjgOvg/H0F449OqFzrsaCSAnAc5ugagI0yM/xTTfoy20G2pKKi3DwDMELFzby1LUp0hRGqlSYykHUw7yT8I//BP/DycyTgDgYe9BP79+XOd9HaHeobRybVXvbTis1bm0fP44mcIfJ5LZHZ1BqcFIp/xE7s87TY/zh3BKjAHAxtcX5ztuw2nAAJzCwtA2a1b1zqU03UM5tQFO/mx6DeDUAtqEQ5vB0OYGcAuw2PurbyrBKEojUmIoYWfiTjaf38yW2C1kFWdhI2zo7tWdSb0mMcBvAO092jf5y10XSCk5GJvFr0eT2HwimajUfJCSW0nmw8yTtD65F01KEmg0OIaG4jL2QZwGDcK2devqJWWjEWK2w/H1cOpXyI4BhOns5KaZEDwEvDpfs2coVVEJRlGucSWGEv6K+4tN5zexLW4b+fp8nHXO3BB4A0NaDqGfbz+cbZ2r3lETIaXkWEIOPx1O4JfDicRlFqLTCm5z1/OK4TCB+7YiE+IRtrY49e+Py6SJOA8ejI2HR3UPYDo7ObIOjn4LOfFgY286Q7nhJdPlL+cWln2TVkIlGEW5BkkpOZBygJ+ifmLTuU3kluTibufO0KChDGk5hL6+fbHV2jZ0mFYlNqOAb/bFsf5QAtFp+dhoBOEtnXlDE0fQ7j/RHzkMQuDYLwzXZyfictPNaJ2dqn+A3GQ4+AUcXmtql6KxgeCb4OZZ0P42sDVjX42ESjCKcg1JzEvkuzPf8dPZn4jPi8fBxoEbW97IHW3uIMw3TN2cv0yR3sCmY0l8vTeWf86kIwT0a+PJM221hB78k6IlP2PMy0MTEoLXi1NwveOO6t2gv8BohOgI2LvMdG/FWAqBYXD7e9DpbnDytNh7uxaoT6OiWDmD0cA/Cf/w9amv+Sv+L6SU9PXty9M9nmZIyyE46ZreL+OqnEnJZeWO8/xwIJ6colIC3B2YfGMwdxWeRfPNpxTs2UOhTofLbUNxf+A/OPTsYV5Fh8Is2Lcc9i0zdcHi4AF9n4LQsdA8xELv6tqjEoyiWKmMogy+jfyWbyK/ISE/AU97T8Z1GcfIdiPxc/Zr6PCsjtEoiYhMYdk/5/jrdBq2Wg23dfVhVHdvOh7dTsaH71EUFYXO3x+vKS/gds891b+vckHmOdi5EPavBH0+tBoIN75mavBoY2eR93UtUwlGUazMuexzrElfw5RvplBsKKavT1+e7/08NwbeiE7btLpfqY784lJ+O6fn9T0RnE8vwNvVjim3tGNUtxZofvqBjCeXkZSSgl3Hjvi99y6ut956ZTcsVUk4AH9/CCfWg9BAl5HQfyL4dLXMm2okVIJRFCtxMOUgy44uY0vsFrRoGR4ynEc6PUKbZm0aOjSrlF2gZ/n2cyzbHk1WgZ7QVk5MuaU9t7bzIO/bb0l7fSGG1DQcw8LwnTsXpwH9zW/vk3AAIuZB5Eawc4P+z0LfJ8FVnUFWh0owitKApJTsTNzJZ4c+40DKAVxtXRnfdTxBmUHc1f+uhg7PKqXkFvG/v6NZveM8+SUGburoRZhrDuPu7EP2jz9yfvICShMScezdmxYffohjaA268U84ABHzIfJXsG8GN74KfZ4Ee9e6f0ONmEowitIALk8sXo5evNznZe4OvhtHnWOtxpdvrNLzivks4iyrdp5HbzByRzc//hvelo6+rmxfuJDouz+i+PRp7Lt2xXfW7JqdsWREweaZcPwHlVjqgEowilLPdifuZsHBBexP2Y+XoxfT+07nnpB7VLuVSuQW6fn8r2j+91cUhXoD9/QKYOLgYIKaO1ESE0PsxGm4b/4Do78//h9+iMutt5ifWAoyYNu7sHsxaHVww8vQ72mwb7q9RdcFlWAUpZ6czjzN+/ve5+/4v/Fy9GJa32ncE3IPdlpV+6giRXoDK3ec49OIs2QV6BnW1Yfnb25PsJczxqIiUj74kIylS0GnI3fEcHrPnInGzsyyNOhNSWXrfCjOhZ4Pw+Dp4OJjkffU1KgEoygWllKQwoKDC/jhzA842TjxfOjzPNjxQZVYKiGlZMORJOZuOEF8ViHXt2vBi7e0p2uA6Wwif+cuEl+fgf58DK533YnXC1P458Rx85PL+e3w8/OQegLa3gi3vAnenS3wjpoulWAUxUIK9AUsPbqUlcdXojfqebDDgzzZ7Uma2Vej590m6mh8NrN+Ps7u6Aw6+Ljw5fi+9A9uDoAhO5vkd94h+5tv0QUG0nLp/3Dq39+04Ynj1T9Ifhr8PsPUrYtbIDzwFXQYZoF3o1hVghFCDAU+ArTAEinlvMuWfwAMLpt0BLyklM3KlhmAsv6viZFSqio4SoOQUrLp/Cbe2fMOKQUpDA0ayrO9niXQJbChQ7NaaXnFvLvpFGv3xuLuaMucu7vwwHUt0WpM91Lytm0jYfp0DBmZeI4fR/MJE9A4OJh3ECnhwGr47VUoyYOBk+H6F5tkH2H1xWoSjBBCCywAbgbigD1CiPVSyos/TaSUk8ut/wzQs9wuCqWUPeorXkWpSFR2FHN3zWVX4i46eHTgvRveo4eX+lhWxmiUrN0by1sbTlBQYmDcgNY8MyQENwdTg1JjYSEp77xL5pdfYhcSQstFi7Dv1Mn8A2XHw0/PwpnN0GoA3P4+eHWo43ejXM5qEgzQBzgjpYwCEEKsAYYDlZ37/gd4vZ5iU5SrKtAXsPDwQlYdX4WD1oFpfadxf7v70Wq0DR2a1YpMzmXad0fYez6Tvq09mHN3V4K9/h1WoPDYMRJefImSqCg8xo6lxeTnzL/PIqXpUtjGV0wdUQ57F3qPa7Tjr1gbIaVs6BgAEEKMBIZKKceXTY8G+kopJ1awbitgJxAgpTSUzSsFDgKlwDwp5Q+VHOcJ4AkAb2/v0DVr1tQo3ry8PJyd1Rgb1dWYy+tY4THWpq8l05BJmFMYd7nfhYvWpVb7bMzlVWKQrD+r59doPfY28EB7Wwb62/xbtVhKHLZsweXb7zC6uJAzdgwlHa5+tlFRedkWZ9D+1Cd4Zuwjy60LJzs8Q5GDqh0Gtft8DR48eJ+Usne1VpZSWsUDGInpvsuF6dHAJ5WsOxX4+LJ5/mXPbYBzQNuqjhkaGiprasuWLTXetilqjOWVUZghX972suyyvIsc/v1wuT95f53tuzGWl5RS7opKl9e//adsNfVn+fzagzItt+iS5aW5uTJ20nPyePsOMua/T8vSzMxq7feK8or8Tcr5baSc7S3lzoVSGgx19A4ah9p8voC9sprf69Z0iSweKH8XNKBsXkUeACaUnyGljC97jhJCRGC6P3O27sNUmjpZdhP/rV1vkVOcw5PdnuSJbk+ohpJXUaQ38M6mUyz9J5pAd8dLaoddXOdUJPGTJlESG4vXlBfwGDfO/AaTpcWmlvg7F4B3Fxi5FFq0r8N3opjDmhLMHiBECNEaU2J5AHjw8pWEEB0Ad2BHuXnuQIGUslgI0RwYALxdL1ErTUpKQQpv7nyTLbFb6OTZicU3L6a9h/oCu5r9MZlM+foQUWn5jA5rxcu3dcDJ7tKvnpxffyXhlWloXJxpuWwpTn36mH+gtDPw7WOQeAj6PAE3zwadfR29C6UmrCbBSClLhRATgU2YqikvlVIeE0LMwnRKtr5s1QeANWWnahd0BBYJIYyABtM9GDMqxitK1Tae28jsHbMpNhTzfOjzjO40Wo0geRVFegMfbj7N4m1n8XVz4IvxfRlw2VmLlJK0BZ+S9sknOPTqRcBHH2LTwvzx6punbofFC0zdvKh2LVbDqv46pJQbgA2XzZtx2fQbFWy3HVADMygWkVOSw1u73uLnqJ/p4tmFuYPm0tqtdUOHZdUik3N59qsDnEzK5T99Apk2rCMu9peOZWMsKiJx2jRyNvyK2/Dh+MyehcbWzMuMhlL4cxZdjn0EAdfBfSvAzb8O34lSG1aVYBTF2uxJ2sO0v6eRWpDK092fZny38eg0atCvykgp+WJXDLN/Po6LvQ3Lxl7H4A5eV6ynT0khbsJEio4erfn9lvx0+OZRiN5KvN9Q/MeuVKNKWhmVYBSlAsWGYj7e/zErj6+kpWtLVt62km4tujV0WFYtM7+Eqd8e5rfjydzQrgXv3tedFi5XfuGXnDtHzLjxlGZmEvDJx7gMGWL+wRIOwNrRkJcCwxdwOjsAf5VcrI5KMIpymajsKF7c+iKRmZGMaj+K50Ofx1Hn2NBhWbUdZ9OZvPYg6fnFvHp7Rx4b0BqN5sozksIjR4l98kmQklYrVuDQtYv5Bzv+I3z3JDg1h3GbwK8nqPFzrJJKMIpSzvqz63lz55vYa+1ZMGQB1wdc39AhWbVSg5GP/jjNJ1vO0NrTiSVjBtDFv+IxVPL++Yf4Z55F6+5O4JLPsWtt5n0sKeHvD+CPmRDQBx74EpzNrxCg1B+VYBQFU1cvc3bNYf3Z9YR6hzJ/0Hy8nbwbOiyrlppbzLNfHWBHVDr3hQbwxl2dr6h+fEHOr78S/9JU7Nq0IXDxYnTeV96XuarSEvh5MhxcDV1GwvAFqgryNUAlGKXJO5Vxihe3vci57HM81f0pnuz2pKp+XIU95zKY8MV+sgv1vHtfd0aGBlS6bvb69SS8/AoOvXoS+OmnaF3NHH64IMN0v+X836aRJsNfBnMrBCgNQv0VKU2WlJJ1keuYv3s+rnaufH7L5/T17dvQYVk1KSVL/opm3saTBLo7sOKxPnT0rTxhZH37HYmvvopj374EfroAjaOZ97KyYmH1PZB5Du75HLrdX7s3oNQrlWCUJimvJI83drzBpnOb6Ofbj7mD5tLcoXnVGzZhOUV6Xlx3iE3Hkhna2Ye37+uGq33lVbYzv/6apBmv49S/PwELPjF//JaUE7DqHijJh9E/QNCAWr4Dpb6pBKM0OWcyzzA5YjKxubFM6jWJx7o8hkao7tuv5kRiDv9dvY+4zEJevb0j4wa2vmq7lcyvviJp5iycrh9EwMcfm9/Nfswu+PJ+sLGHRzeATw1qmykNTiUYpUn5NfpXXt/+Oo42jiy5ZQm9farX63hT9uPBeKZ+exg3Bx1rngijd5DHVdfP+vY7kmbOwnnwYPw/+tD81vmnNsK6seDqB6O/B/dWNQ9eaVAqwShNgt6o5/2977P6xGp6tOjBe+Hv4eVoZk2mJsZglLy96SSLtkbRJ8iDBQ/1qrDhZHk5GzeS+NprOPXvX7PkcvBL+HEi+HaDh74xtXVRrlkqwSiNXmpBKlO2TmF/yn4e6vgQL4S+gE6runu5muwCPc+sOcC2yFRGh7XitTs6YWtz9cuIeVu3Ev/iSzj07EnAJx+bn1x2fw4bpkCbcBi1GuxqN2ib0vBUglEatf3J+3lh6wvk6/OZP2g+w9qoXnarEpmcyxMr9xKfVchb93TlP31aVrlN/q7dxD07CfuQEAIXfmZ+bbHtn8Bv06H9MLhvuepTrJFQCUZplKSUfHHiC97b+x7+Lv4svnkxIe4hDR2W1dt0LInn1x7EwdaGrx6v+n4LQOHRY8T997/oAgMI/N8StC5mnnlsewf+fBM6jYB7l5i63FcaBZVglEanQF/A69tfZ+O5jQwOHMycgXNwsVWXW67GaJR89MdpPvrjNN0D3Fg0ujc+blW3lC+JiyP2qafQNHOj5f/+h427e/UPKiVsmWNKMN0eMLXO16qvpMZE/W8qjcq57HNMjphMVHaUqoJcTXnFpTy/9iC/HU/m3l4BzLm7C/Y6bZXblWZmEvv4E8iSElqtWI7O24yudaSE31+D7R9Dr0fgjg9BU/UxlWuLSjBKoxERG8Erf72CjcaGhTctpJ9fv4YOyeqdS8vn8ZV7iUrLZ8YdnXh0QFC1xmUxFhURN2Ei+rg4Wi5bil3bttU/qJSw8RXY9Rlc9zjc9jZo1I+AxkglGOWaZ5RGPjv0GQsPLaSTZyc+CP8AP2e/hg7L6m2NTOWZL/ej0QhWPtbniuGMKyONRhJemkrh/v34f/gBjr3NaEt04cxl12cQ9jTcOlf1K9aIWdXPBiHEUCHEKSHEGSHEyxUsHyuESBVCHCx7jC+3bIwQ4nTZY0z9Rq40lJySHCb+MZGFhxYyvO1wVgxdoZJLFUz9iUXx6LLd+DVz4KeJA6udXABS5r9N7m+/4fXyVFyHDjXnwPDHLNNlseseV8mlCbCaMxghhBZYANwMxAF7hBDrpZTHL1t1rZRy4mXbegCvA70BCewr2zazHkJXGsjpzNNM2jKJxLxEpvedzqj2o8wfdreJKdIbmP79Ub7dH8fQzj68d3/3SrvYr0jmunVkrFiB++jReI4da97BI+bB3+9D6FjTZTH1f9XoWU2CAfoAZ6SUUQBCiDXAcODyBFORW4HfpZQZZdv+DgwFvrJQrEoD2xi9kRnbZ+Ckc2LZ0GX08OrR0CFZvZScIp5YtY+DsVk8d1MIz94YUuGok5Up2LuXpFmzcRo4EO+pL5l38G3vwNZ50ONhuP0Ddc+libCmBOMPxJabjgMq6jv9XiHE9UAkMFlKGVvJtv4VHUQI8QTwBIC3tzcRNRxqNS8vr8bbNkV1VV4GaWB91nr+zPmT1natGecxjqzjWUQcr/2+rUldf76isg383/5iCkolE3rY0cMmgW3bEqq9vSYtHc958zB6uBN9z91E/f13tbcNjPmOtlErSPIO56TbPbBtW03ewlWpv0fz1Fd5WVOCqY6fgK+klMVCiCeBFcCN5uxASrkYWAzQu3dvGR4eXqNAIiIiqOm2TVFdlFdGUQYvbX2JXTm7GNV+FFOvm9pou3ypy8/Xjwfjmb/5MM2d7fnykd508jNvwC9DXj7nH3wQvUZD0PLl5g11vONTiFoBXe7F5+7F+FionYv6ezRPfZWXNSWYeCCw3HRA2byLpJTp5SaXAG+X2zb8sm0j6jxCpcEcSzvG5IjJpBemM3vAbEYEj2jokKyewSh5Z9MpFm49S5/WHnz2UC88nc3rgkUajSRMnUrx2bMELl5kXnLZvxI2vQId74S7F6lGlE2QNV0I3QOECCFaCyFsgQeA9eVXEEL4lpu8CzhR9noTcIsQwl0I4Q7cUjZPaQS+P/09j/z6CAArh61UyaUacor0PL5yLwu3nuWhvi1ZPa6v2ckFIO2TBeT98QfeU6fiPMCMAb9O/AQ/TYK2N8K9S1X3L02U1fykkFKWCiEmYkoMWmCplPKYEGIWsFdKuR54VghxF1AKZABjy7bNEELMxpSkAGZduOGvXLv0Bj3z98xn7am19PXpy9s3vI2HfdV9YzV10Wn5jF+xh/PpBcwe0YXRYTUbTyVv61bSPv0UtxEjcB/9sBkBbINvHgP/ULh/FdiY2auy0mhYTYIBkFJuADZcNm9GudevAK9Usu1SYKlFA1TqTUpBCi9EvMDB1IOM7TyWSb0mYaOxqo+rVfrrdCoTvtiPViNYNa4v/dp61mg/JXFxxL80FbuOHfF54/XqV/+O3w9f/Qc82sKDX4Odc42OrzQO6i9WsTo7E3cyddtUCksLeef6dxja2ozGfE2UlJKl/5xjzi/HaeftwueP9CbQw8wu88sYi4qIe/ZZkJKA//sIjX3VnV4CkHYavhgJDh4w+jtwVGebTZ1KMIrVMEojiw8v5tODn9LarTVLb11K22Zm9HHVRBWWGJj2/RG+PxDPLZ28+WBUD7MaT5YnpSRp1myKj58gYOFn2AYGVr0RQHYcrBwBCHjkB9Nwx0qTpxKMYhUyizJ55a9X+CfhH4a1Hsbr/V7HUVezX+BNSWxGAU+u2seJpBxeuLkdEwYHm9V48nJZ69aR/d13NH/6v7hUtxprfjqsuhuKc2Dsz+CpfhQoJirBKA3uYMpBpmydQkZRBq+FvcZ97e5TXb5Uw9bIVJ796oDp8tiY6xjcwatW+ys8eozk2W/iNGAAzSdMqN5Gxbmmy2KZ502XxXy71yoGpXFRCUZpMOVHnfR28mbVsFV09uzc0GFZPSkln0ac5d3fTtHe24VFo0Np5elUq30a8vKIf/55tM2b4/fuOwhtNcZmKS2GtQ9D4iEYtRqCBtYqBqXxUQlGaRDZxdm8sf0NNsdsZnDgYGYPmI2bnVtDh2X1cov0TFl3iE3Hkrmrux/z7u2Ko23t/oyllCS9/gb6+HharVpZvVEpjQb47nGIioARn0GHYbWKQWmcVIJR6t3+5P1M/WsqaQVpvBD6AmM6j1GXxKrhTEoeT67ay7n0Al67oxOPVXNwsKpkf/cdOb/8Qos838GhAAAgAElEQVTnnsOxV6+qN5ASfp4Mx3+EW+ZAjwdrHYPSOKkEo9Qbg9HA4iOLWXhoIf7O/qwatoouzbs0dFjXhB8OxDPt+yM46LSsrkX7lssVnzlD0uw3cewXhufj46veAExjuuxfAYNegP4Tq15fabJUglHqRWZpJuN+G8e+5H3c2eZOpodNx0lXu/sGTUGR3sDMn47x1e5Y+gR58H//6YmPWzXbpVTBWFRE/OTn0Tg64jd/fvXuu2z/+N8xXW58rU7iUBovlWAUi/v9/O/MS5wHGpg7cC53tr2zoUO6JpxNzWPCF/s5mZTL0+Ftef7mdtho6677wOS35lF8+jSBn3+OzqsaNdAOfAG/vQqdhsPt76sBw5QqqQSjWEx2cTZv7X6LX6J+oaVtSz67/TNaurZs6LCuCTsTSpnw59/Y2mhY9uh1DG5fuyrIl8vZuJGstWvxHD8O50HVqP118hdY/wy0CYd7PgdNNc52lCZPJRjFIrbFbeON7W+QWZTJ092fpl1mO5VcqqGwxMCsn4/z1eFierdy5+MHe+Lr5lCnx9AnJJD42gzsu3ejxaRJVW9w7m9Y9yj49YBRX4CN+b0yK02TSjBKncoryeOdve/w3envCG4WzCdDPqGTZyc12mA1HI3PZtKaA5xNzWdYax0fjQ9DV4eXxKBsfJdXpoHBgP877yB0VXSjn3gIvnwA3IPgwXWq80rFLCrBKHVmZ+JOZvwzg+SCZMZ1GcfTPZ7GVqu6aq+KwShZvC2K938/haeTHV+M74s+7midJxeAjOUrKNi1C985b2LbsoozyrQzsOoecGgGo78Hp7qpuaY0HSrBKLVWoC/gg30fsObUGoJcg1h520q6t1BdhlRHfFYhz689yK7oDIZ19WHu3V1p5mhLRFzdH6vo1ClSP/gA55uG4HbPPVdfOSfB1L8Y0pRc3PzrPiCl0VMJRqmVLTFbmLt7Lkn5STzc8WGe7fUsDjZ1e8+gMZJSsv5QAq/+cBSjUfLOyG6MDA2wWINTY3ExCVNeRNPMDd/Zs69+nIIMU3IpzISxP0HzEIvEpDR+KsEoNZKcn8y83fPYHLOZ4GbBrLptFT28ejR0WNeElJwiXv3hKL8dT6ZXy2Z8MKpHrfsSq0rqBx+aqiQvXnT1rmCK8+CL+yAjGh7+Fvx6WjQupXFTCUYxS6mxlLWn1vLJgU/QG/VM6jWJMZ3GoFNjrldJSsm3++OZ9dMxikuNTBvWgccGtK7Tti0Vyd+xg4zly3F/8EGcr7++8hUvdF6ZsN801HHrQRaNS2n8rCrBCCGGAh8BWmCJlHLeZcufB8YDpUAq8JiU8nzZMgNwpGzVGCnlXfUWeBOxO3E3b+1+izNZZwjzDeO1sNdU1eNqSsgqZNr3R4g4lUrvVu68PbIbbVpYvkaWITubhJdfwbZNG7xenFL5ihc7r9wCwz+FjndYPDal8bOaBCOE0AILgJuBOGCPEGK9lPJ4udUOAL2llAVCiP8CbwOjypYVSinVNRoLSMhL4L297/Hb+d/wd/bnw/APubHljaqDymowGCVf7jrP/I2nMBglr9/ZiTH9gmo1KFh1SSlJmjmT0vR0ghYsQONQyb0xKeHn50ydV946F3o+ZPHYlKbBahIM0Ac4I6WMAhBCrAGGAxcTjJRyS7n1dwIP12uETUxeSR4rjq9g+dHlAEzoMYGxncdib1M3fWE1dodis3j1h6Mcic9mQLAnb93djZae9TdKZ87Pv5Cz4VdaPPccDl2uMs7O5tdh/0oYNAX6VXOgMUWpBiGlbOgYABBCjASGSinHl02PBvpKKSvsrlUI8QmQJKV8s2y6FDiI6fLZPCnlD5Vs9wTwBIC3t3fomjVrahRvXl4ezs6Ns9GZXur5J/cfNmVvIs+YR0/HnoxwH4GHjUeN99mYy+ty+XrJt5ElbIktxdVO8J8OtvT10Zp1xlfb8tJkZ+M5cxalPj5kTnkBNBXf5wmM+Za2USuJ97uN0yFPXrP9izWlz1ddqE15DR48eJ+Usnd11rWmM5hqE0I8DPQGbig3u5WUMl4I0Qb4UwhxREp59vJtpZSLgcUAvXv3luHVHXf8MhEREdR0W2tllEZ+ifqFBQcXEJ8XT1+fvjwX+lyddKnfGMvrcgaj5Nt9cczfeJLMglLGDghi8s3tcLU3vwJEbcpLSknc0xPINxgIXvAJdq1bV7zivuUQtRK63Iv/PUvwryQJXQuawuerLtVXeVlTgokHAstNB5TNu4QQ4iZgOnCDlLL4wnwpZXzZc5QQIgLoCVyRYJQrSSnZGreVTw58wqnMU3Tw6MCimxbRz6+fus9STX+fTuPNX45zMimXXi2bsXJcHzr7NcwIndk//kjeli14v/Jy5cnl2Pfw03MQfDOMWFjpGY6i1IY1JZg9QIgQojWmxPIAcMlQeUKInsAiTJfSUsrNdwcKpJTFQojmwABMFQCUqyg1lvLbud9YcnQJpzNPE+AcwPxB8xnaeigaob5wquN0ci5zN5xgy6lUAj0c+OTBntze1bfBErM+OZnkOXNx6B2K++jRFa905g/49nEI7Av3rwQb1Z2PYhlWk2CklKVCiInAJkzVlJdKKY8JIWYBe6WU64F3AGdgXdkf8IXqyB2BRUIII6DBdA/meIUHUigxlPDj2R9ZdnQZsbmxtHFrw9yBcxnaeig6jWrPUh1J2UX835+nWbsnFkdbLdOGdWBM/yDsbBquG3spJYmvvYYsLcVvzhxERWclMTtNbV1adIAH14Jt/VU6UJoeq0kwAFLKDcCGy+bNKPf6pkq22w50tWx0174CfQHrItex8thKUgpT6OLZhRcGv8DgwMHqjKWaUnKL+CziLF/sikFKyeiwVjw7JAQPp4Y/C8j+7jvyt/2F96uvYtuq1ZUrxO+D1SPB1c/USt+hWf0HqTQpVpVgFMuIyo5i3al1rD+7npySHPr69OXNgW8S5hum7rFUU3peMYu2RbFyxzn0BsnIXgFMvDGYQA/rOAPQJySQ/NY8HPv0wf3B/1y5QtIRU8/Ijh7wyHpw8a7/IJUmRyWYRqrEUMLm85tZF7mOvcl7sdHYMKTlEEZ3Gq16OjZDbEYBS/6KYu3eWEpKjYzo4c+zQ0IIam7ZvsPMIaUk8dXXkEYjvnMruDSWchJWDgdbZxjzk+oZWak3KsE0MjE5MXwT+Q0/nPmBzOJMApwDeK7Xc4wIHoGngxrPo7qOJ+SwaNtZfj6ciEbAiB7+PHlDW4K9rK+tRdbX68jfvh2fN17HNiDg0oXpZ2HlXaDRwZj14F7BpTNFsRCVYBoBvUHPltgtrItcx87EnWiFlsGBg7mv3X2E+YWp+yvVVGow8sfJFFbtOM/fZ9JwstXy2IAgHhvYus6HLa4rJXHxpMyfj1P/fjQbNerShZnnYcVdYCyFsRvAs23DBKk0WWYnGCHEW8BrZbW+NICzlDKn7kNTrkZv0LMzcSe/nf+NP2L+ILckF18nXyb2mMjdIXfj5ejV0CFeM1Jzi1m7J4Yvd8WQkF2Er5s9L97anof7tsLN0Xpr1UmjkcRXXwUhrhzjJTseVtwJJXkw9mfw6tBwgSpNVk3OYJyllKUAUkqjEOIDYFzdhqVUpKKk4qxzZnDgYIa2HsoAvwFoNQ1XTfZaYjRKdkals2ZPLL8eTURvkAwMbs6MOztzU0cvi3ehXxcy16yhYOdOfGbPQudf7r5Kdhwsv8M0YNgjP4KPqmCpNIyaJJjLOy/LrYtAlIpdSCqbzm3iz9g/L0kqtwbdSj+/fmrcezOcScnju/1x/HAgnoTsIlzsbHiobyseDmtllfdXKlMSG0vKO+/iNHAgzUaO/HdBVkxZcskyDXXs36vhglSavCoTjBBCK6U0lJu1XwjxLvAhYADUtZg6ppJK3UrOKWLj0SS+OxDPodgsNAKub9eCl4d15JZO3tjrrq2zPmk0kvjKNISNDb5vlrs0lnkOlt8JxdnwyA8quSgNrjpnMEuARy9MSCmXCyHCgFcBO+BNC8XWpBToC9iRsIM/Y/9kS+wWlVRqKS6zgI1Hk/j1aBL7YzKREjr4uDB9WEeG9/DDy/XaHXIgc/UXFOzdi+/cueh8fEwzM6JMyaUkz9TOxU8NjaQ0vOokGK0QYtZlLep3CiGOY+rORXXJUkMpBSlExEYQERvBrsRdlBhLcLF1ITwgXCUVMxmMkkNxWWw9lcqfJ1M4Ep8NQEdfVybf1I7buvgQ4u3SwFHWXsn586S8/z5ON1yP290jTDPTz5oui5UWmdq5+HZr2CAVpUx1EsxjwK9CiPFSyiUAQojuwFrga0sG19hIKYnMjGRL7BYiYiM4ln4MAH9nf+5vfz+DAwfT07un6g+smlJyi9gWmcbWyFT+Op1KVoEeIaBHYDNevq0DQzv7WFWDyNqSRiMJ06cjdDp8Z80yXRpLOQErR5RVRf4ZvK8ysJii1LMqE0xZdeR7gAghRALgj6m7/PFSys2WDvBapzfo2Zu89+KZSkJ+AgJB1xZdmdRrEuEB4bRt1lZ12VINKTlF7D6Xwe7oDHZFZXAq2VS/pLmzHUM6eBPevgUDg5vjbgX9gllC5hdfUrh3n+nSmLc3xO2DL+4FrV1ZVeSODR2iolyiOjf5PwUOY+r+fhlwBNNIk8kWju2alV2czd/xfxMRG8Hf8X+Tp8/DXmtPmF8YT3Z/kusDrqe5Q/OGDtOq6Q1GIpNzORyXzcGYLHafyyA6LR8AJ1stoUEeDO/pxw3tWtDRx7VexrhvSCUxMZdeGju7BdY8BM4tYPQP4FHJuC+K0oCqc4nsIKaeirsBOqAd8LkQ4ghwREpZszGHG5nY3NiLZyn7kvdhkAY87D24JegWwgPCCfMLw8HGOluDN7RSg5Fz6QUcic/iUGw2h+OyOJaQQ3GpEYBmjjp6t/LgwT4t6dPag85+rtdEO5W6Io1GEqdNN9UamzULceIn+HYceAabqiK7+DR0iIpSoepcIltcfloIEcC/Ced2oEkmGKM0cij10MWkcibrDADBzYJ5tMujhAeG07V5V9VNSzkGoyQ2o4DI5NyyRx6RyblEpeVTUpZMHG21dPFzY3RYK7oFNqN7gBstPRyb9CXEzC++NNUamzMHXfxv8NOz4B8KD35t6h1ZUaxUdS6R9QN2SiklgJQyDogDfrVwbFZJb9Tz5s43+T3ud3JjctEKLaHeobx03UuEB4QT6BpY9U4aMaNRkpJbzPn0fGIyCojJKOB8egGHogtJ2rzx4lkJgH8zB9p5O3NDuxaEeLvQ1d+NYC9ntI38cpc5Ll4aGzQIN69zsH4mtL0RRq0G28ZTgUFpnKpziewR4FMhxClgI7BRSplk2bCsl06j43zOeULsQ7g/9H4G+g/Eza5hxl6vb1JKcotLSc4uIjG7iKScIpLKPcdkFBCbUXBJEtFqBH7N7HGzFdzcrSXtfFxo5+1CsJczznaqr9WrueTS2PUaxB8zocu9MOIzsLFr6PAUpUrVuUT2XwAhRAfgNmC5EMIN2IIp4fxzWUv/Rm/ZrcvYunUr4W3CGzqUWisuNZBVoCcjv4TM/BIyCkrILNCbXueXkFlQQmpu8cUkUlBy5X+1h5Mt3q72tGnuxOD2LWjp6UQrD0daejji7+6ATqshIiKC8PBODfAOr12ZX35lujR2hx+6s1/BwMlw4wyoaChkRbFC5vyE/B+mXpQ/EEI4AIOB+4D3gd51EYwQYijwEaAFlkgp51223A5YCYQC6cAoKeW5smWvYOp00wA8K6XcVBcxVRKnpXZtllKDkUK9gUK9gYJiA3nFpeQWlZJbpCevuLTcdCl5xXpyi0rJKtCTWVByMaHkV5AwLnCxt8HDyRZPJ1s6+rgS3s4LHzc7vF3t8XVzwMfVHi9Xu2uuq5VrgTY1lZR338UpSIeb8wG440Po/WjVGyqKFTEnwTwJzBRCTAemSyk3ABvqKhAhhBZYANyM6R7PHiHE+st6ChgHZEopg4UQDwDzgVFCiE7AA0BnwA/YLIRoZ6kzK73BSFGpJKughBKDkZJSI3qDRF/2usRgRF82r8RgoKT032V6g+lRYpCXTl/YzmCkSF+WOErKHnoDRWWJpPy03nB5v6MV02kFLvY6nO1scHfU4eFkS3ALZ9ydbHF31OHuZIuHo63p2cmWZo463B1t0TWhmlrWRBqNuC/7HGEowjc0F/HQWgi5uaHDUhSzVTvBSCmPAvcKIXoBs8p+xU+XUh6qo1j6AGeklFEAQog1wHCgfIIZDrxR9vob4BNhCmQ4sEZKWQxECyHOlO1vRx3FdomwuX+Qnl8Cm3+vk/0JAbZaDbZaDTobDQ46LfY6DQ62Whx0WlzsbfByscPRVouDrRZ7nWm+g+7faUdb7cUk4mJvejjb2eBsb4OdjTrDuJZkvvcy2qhYvK+X6CZuUF2/KNesmtxlPQPMxtSFzL4a7qMi/kBsuek4oG9l65T1MJANeJbN33nZthUOPC6EeAJ4AsDb25uIiAizA701EAqKJI72dthoKHsI07O4clqrAV3ZtFaYXmsvbCe4Sq0pQ9njKi6sUgTkQDGmR5rZ78qy8vLyalTWTYo00urAUoqX78E2QMfxe2ZScioDTkU0dGRWT32+zFNf5VXt5CCE2AKEAIWYziqOA2MtE5bllLXrWQzQu3dvGR4ebvY+wqHsprX52zZVqryqUJSN/PZxYr7fi7BxJOGpGQwaem9DR3XNUJ8v89RXeZlz9vECcEJKWWihWOKB8o1IAsrmVbROnBDCBnDDdLO/OtsqinVKOQlfjyZzZyIFqS74zn6DxOaqKyHl2lftu7hSyv0WTC4Ae4AQIURrIYQtppv26y9bZz0wpuz1SODPsgag64EHhBB2QojWmM60dlswVkWpGwe+gM8HU5KcScrR5jgNHIjbSHXmojQOVtPSreyeykRgE6ZqykullMeEELOAvVLK9ZiqSq8qu4mfgSkJUbbe15gu25UCE5pa2xzlGlOSD79MgUNfIlsOJOF3Z4RNFL6zZ1lNNXhFqS2rSTAAFVV9vmygsyJMbW8q2nYOMMeiASpKXUg5AevGQuopuGEqGVHeFB54F9+33kLn69vQ0SlKnbGqBKMojZrRCLsXw+bXwc4FRn9PsTGA1EkjcR4yBLcRwxs6QkWpUyrBKEp9yI6HH5+GqAgIuQXu+hhp70nCqAfQODvjO2umujSmNDoqwSiKJUkJR9aZ7rcYS01dvoSOBSFI+7+PKTp+HP+P/w8bT8+GjlRR6pxKMIpiKdnx8OtLcPJnCOxr6gXZsy0AhYcPk7ZoEW7Dh+N6s+oGRmmcVIJRlLpmNMDepbB5pums5aaZ0P8Z0Ji67DEWFZEw9WVsWrTAe/q0Bg5WUSxHJRhFqUvJx2D9sxC/F9oMhjs+AI/Wl6yS8v77lERH03Lp/9C6ujZQoIpieSrBKEpdKMiArfNh9+fg0Azu+Ry63mfqybSc/J07yVy5CveHHsKpf/8GClZR6odKMIpSG4ZS2LcMtsyBomzoNQaGzABHjytXzc0lYdo0bIOC8JryQgMEqyj1SyUYRakJKSFqC2ycBqknIGgQDJ0HPl0q3SRp5ixKk1MI+mI1GgeHegxWURqGSjCKYq6YXfDnbDj3F7gHwajV0OGOKy6HlZe9fj05P/9M82efwaFHj/qLVVEakEowilJdCQfhzzfhzO/g5AW3vW1q02Jjd9XNSmJjSZo5C4fQUJo/+WT9xFqOXq8nLi6OoqKiej92fXFzc+PEiRMNHcY1ozrlZW9vT0BAADqdrsbHUQlGUaoStxf+/sDUnsXB3VTtuM/jYOtU5aaytJSEF18CjQb/t+cjtPU/umhcXBwuLi4EBQU12t4CcnNzcXFxaegwrhlVlZeUkvT0dOLi4mjdunWl61VFJRhFqYiUcPp3+OcjOP832DeD8Fcg7L9g71bt3aR9+hmFBw/i99676PwrHGTV4oqKihp1clHqnhACT09PUlNTa7UflWAUpTx9IRz7HrZ/AinHwDUAbn0Lej0Cds5m7apg717SFi7EbcQI3G6/3UIBV49KLoq56uIzoxKMogBkRJla3x9YDYWZ0KIjjFgIXUeC1vxr0IacHOJfegldQADer75qgYAVxfqpBKM0XaXFELkR9q+CM5tBaKDjHXDdeFO14xr+gpNSkvTGTFOV5C+/QOtc9b0aRWmMVIJRmhYpIW4PHPoKjn5rahzp4gs3TIXQMeDqV+tDZK1bR86GDbR47jkcunevg6AV5dqkaegAFMXipISko7BlLnzcC/53Mxz8CkJuhYe/g8nHYPArdZJcik6eJPnNOTgNGIDnE4/XQfCNx6JFi3jqqacAU9Xp0aNHM2bMGPR6fa33/fTTT+Pl5UWXLpU3dK2tjRs30r59e4KDg5k3b94ly7Kyshg5ciQdOnSgY8eO7Nixo9L9WLIcAB577DGLl0V1WUWCEUJ4CCF+F0KcLnt2r2CdHkKIHUKIY0KIw0KIUeWWLRdCRAshDpY9VEu2pk5KiNsHv88wJZWFA2Dr2+DqD8M/hRdPw72fQ/CQi70c15YhL4/4Sc+hbdYMv7fnIzRW8edlNY4cOUK3bt3Iycnhtttuo2XLlqxYsaJW7SwueOihh9i4cWMdRFkxg8HAhAkT+PXXXzl+/DhfffUVx48fv7h80qRJDB06lJMnT3Lo0CE6duxY6b4sWQ4AY8eOtWhZmMNa/gJeBv6QUoYAf5RNX64AeERK2RkYCnwohGhWbvmLUsoeZY+Dlg9ZsToFGabLXj9MgPc7wpIbYccCcG9tGuhrSiSM/Rl6PmQasrgOSSlJmjGDkthY/N97Vw0gVoHDhw/j4eFBeHg49913H3PmzKmzfQ8YMAAPjyv7fwOIjo5m+PDh9O7dmz59+nDq1Cmz9797926Cg4Np06YNtra2PPDAA/z4448AZGdns23bNsaNGweAra0tzZqZvpqGDRtGQkLCJfuyZDkAXH/99ZWWRX2zlnsww4HwstcrgAhgavkVpJSR5V4nCCFSgBZAVv2EqFid0mJIOABnt5hu0ifsB2k0tVlpO9h0Caz9UFPjSAvLWruWnA2/0mLyZByvu87ix6upmT8d43hCTp3us5OfK6/f2bnK9Y4cOcIzzzzD0qVLufPOO6u170GDBpGbm3vF/HfffZebbrqpyu31ej3jx49n8eLFtG3blg0bNjBv3jyWLVtWreNfEB8fT2Bg4MXpgIAAdu3aBZgSWIsWLXj00Uc5dOgQoaGhfPTRRzg5ObFhw4Yr9lWTcoDal0VDsJYE4y2lTCx7nQR4X21lIUQfwBY4W272HCHEDMrOgKSUxZVs+wTwBIC3tzcRERE1CjgvL6/G2zZFdVFe2tICXHNO4ZZ9nGZZx3HJjURrLEGiIcc1hMyW95Ph0Ysc12AQWtNPj12H6iT+q7GJicXj7bcp6dyJoyHBUAefi7r8fLm5uV38YtKX6DEYDHWy3wv0JfoKv/jKi4uLw8nJibZt2xIdHX3J+jfccAO9e/cmMjKSt99++5LLSxV9QV9Qfh8Gg4G8vDyMRuMl87///nuOHj3K3XffDUBpaSn9+/e/ZJ277rqL5OTkK/Y/Y8YMbi9rv1RYWIhe/+/7LCwspKSkhNzcXLKzs9m/fz/z5s3juuuu46WXXmLmzJm89tprZpUDwOrVqzl48CBSSpycnJg1a5bZZQFUWBblGQyGKv/PwNRItzafw3pLMEKIzYBPBYuml5+QUkohhLzKfnyBVcAYKaWxbPYrmBKTLbAY09nPrIq2l1IuLluH3r17y/DwcPPeSJmIiAhqum1TZHZ5lRZD8lHTGUr8AdNz6gnTGYrQgG936Pg4tOqHaDUAN0cP3IAgC8VfGUNODtFz30J6etJuyRJs3OvmbKkuP18nTpy42C3Im/c2zO3J6OhoevToweeff05YWBiDBg2iZ8+exMbGEhYWxqeffsoHH3xAZmbmJV2YVPdXe25uLs7Ozmg0mku2j4yMZO7cuRcvX1Vky5YtVcYfHBzM6tWrL+47PT2d1q1b4+LiQvv27QkICODGG28E4MEHH2TevHkVdsVSWTkA7Nixg9OnT7No0SIASkpKsLW1NbssgArLorzqdq1jb29/Mb6aqLcEI6Ws9BxOCJEshPCVUiaWJZCUStZzBX4Bpkspd5bb94Wzn2IhxDJgSh2GrliSlJCbCCknIPWk6TnpiGlkSGNZrRpHT/DrCR1uh5ZhENinzu+h1IQ0Gkl4aSr6hARarVxRZ8mlMTp8+DBdu3bF19eXJUuWMGrUKPbs2cO+ffuIjIzkqaeeIjExkcmTJ1+y3V9//VWr4/r6+rJp0yYeffRRNBoNR44coUuXLma3Ur/uuus4ffo00dHR+Pv7s2bNGr788ksAfHx8CAwM5NSpU7Rv354//viDTp06ATBkyBBWrlyJf1k3QZWVg5ubG8uXL2f69H9/b5dPLlD7smgI1nKJbD0wBphX9vzj5SsIIWyB74GVUspvLlt2ITkJYARw1PIhK+YQxlJIPwuZ0ZB2xnQ2knLS9FyU/e+Kjp7g3Rn6TTAlFb+e0KxljRs9WlLagk/Ji/j/9u49Osr6XPT498n9JuEqCRlEbkq4CG0jtAUpKhSNZxUr6JaNe3sAj1Cr9RxthVJbdDdsYa0iSlCRctFTreCWVqn1sjdoZNd6wOACwkVKkCAhoCSQkERym/mdP2YIkzBJZoZ5551kns9as2be6zzzMMyT931nfk8BfX/9OCnf/rbd4US0oqKi5tNNU6ZM4a677mLOnDkMHz6c5cuXM3r0aKZPn05tbS2pqYH/MHX27Nl8/PHHlJeX43A4ePLJJ5k7dy5z5szhww8/JDs7m+TkZEaOHMkrr7wS8P7j4uJYtWoVU6dOxel0MmfOHEaMuHjdKT8/n1mzZtHQ0MCgQYPYsGEDLpeL4uLiFhfc28rD5s2bqaurIy7u4key0+kkNojBUWfOnElBQcElubCDGNPm2ajwBSHSC4PdPrMAABd1SURBVHgduAo4BtxljDkjIjnAfGPMfSJyD7AB2O+16f80xuwWkQ9wX/AXYLdnm5qOnjcnJ8cUFhYGFbOeImvF2eg+Eqk6AedOQFWpu5icLYEzRzGVxxFcF9dP7uEejuXKYS3v0/rY9hICUf3BB5Q+8FPSf/xjMv99ScjH+gr1KbL2vjZrp9zcXAYMGEBMTAw9evQgLy8vqP1E4mjK+/btY/369Tz99NN+rb9//37y8vLo06cP1dXVrFixovnbaKHmb758vXdEZJcxJsef54mIIxhjTAVws4/5hcB9nsevAD7/9DDG3GRpgNHKGKivhtrT8E0F1JZ7HpdD9VfuQnKuzH1f8zXQ6o+V5J7QcyA4rudY+jiuHjPJ3aCr1xBIuzIij0r8Uf/FUcoeW0DSiBFkPLFYB5K8DO1duO7sRo4c6XdxARgxYgSvvfaahRGFX0QUGGUhl9NdJOqrof7cxcd1VV7zPcvOn/UqIhXue2eD7/0mdnP/8r1blvuUVrcsSM/yzHO475O6Na9eUlDA1d+aFJ7XbCFnTS2lDz2ExMfjyF9JTGL7zcaUimZaYMLNGHA1uT+4nY3um6vRM+2Z31TnvjWev3jv/djXfUOtjyJyDhpr/QhK3AUjOR1S+7jH5soYBam9IaW3e15qb6/p3hAffT3ljctF2YIFNJSUcNW6dcT3u/yhZZTqyrTABOPPP2HU8X/A8WfbLhKuxovLmpc3Xvxm1OWSWIhPgfgkiEt2f+AndXM3w0rv7/6WVWI397zEKy5O+5ofnwo6rEmHTj/9NDXbttF30SJSvzvO7nCUinhaYIJxrpSEhrNQLxCb4O7JHpvmfhwT576PjW/5ODYeYuJbPfaxXkycu1jEJbV/H0SPEhW8yjfeoGLtOrrPvJse/3KP3eEo1SlogQnGvX9hl36LLGrU7tjJySeeJPX73ydj0SK9qK+Un/S8iFLtaCgpofRnPyNhwACynlmBhGjEW6WigRYYpdrQdPYsx+f/BBGh/+oXiO3WreONlFLN9BSZUj64zp+ndP5PaCwr46r160jwGklXKeUfLTBKtWKamjjxfx7h/N69ZD3zDCk5fv1oWSnVip4iU8qLMYaTixdTU1BAxm9+TbepP7Q7pC4jklomR1PbYjtpgVHKy+mVK6na/Cd6P/ATesycaXc4XUoktUyOprbFdtJTZEp5nPnDK1S8sJrud95J74cesjsca7y70N0OIZQyRsGtSztcbe/evUyYMIFJkyYxb9485s2bF7IQxo8fT0VFhc9lubm5rF27ln5eIy9YGQu42xaXlJSEdJ+dkRYYpYCzr7/OV0uWcMWUyWQs/o3+1sUCdrRMBt8DakZT22I7aYFRUa9qyxZOLX6C1Ik30G/5ciSuC/+38ONIwwrHjx8nLS2NoUOHcvLkyRbLcnJy+O53v8vBgwdZuXJliz4rVjTZai8WgA0bNvDZZ581ty1etmyZpfF0ZV34f5JSHTv33vuULfwlKePG4Vi5kphWXQRVaBQVFTF69OjmVsHXX399c8vksWPHsmrVKlasWEFpaWmLAmPFEUNbsYC7bXFRURH5+fmAu22xNz2CCYwWGBW1qrdu5cTPf07ymDH0f/45YpKS7A6py7KrZTJo22I76bfIVFQ69847lD78v0kaMZz+a14kJiXF7pC6tKKiIkaNGgW0bBW8a9culi9fzurVq4mLi6O21p/2EpeaPXs23/ve9zh06BAOh4N169YBtNm22FcsgM+2xcGYOXOmz3iijR7BqKhT+ec3OfmrX5Hy7W/jWL2a2LTAe8CrwLz66qstpi+0Rs7NzaW8vJyYmBiys7NJTQ3u32LDhg0+WwAfOHCA6dOnk5x8sX9RW7EAPPbYYzz66KOX3ba4q3WmDJYWGBVVzm7cxKknniD1+9/D8dxzxCRHX+O0SGJ1y2RtW2yviDhFJiI9ReS/ROSw575HG+s5RWS357bFa/5AEdkhIsUisklE9EqtusSZl1/m1BNPkPaDH+B44QUtLkpZLCIKDLAQ2GaMGQps80z7ct4YM8Zz+5HX/GXACmPMEOAsMNfacFVnYlwuvv7d7/jqqaVcMWUKjvyVxCQm2h2WUl1epBSYacDLnscvA7f7u6G4fxF3E/BGMNurrs00NFC2YGFzN8qsZ1Yg+lVkpcJCjDF2x4CIVBpjunseC3D2wnSr9ZqA3UATsNQY86aI9Ab+n+foBRHpD7xrjPE5ypyI3A/cD9C3b9/vbNy4MaiYa2pqSEtLC2rbaGRHvqSujvQX15B48CA1P/oRtbfeAp3kF/qhzFd6ejpDhgwJyb4ildPpJDY21u4wOg1/81VcXExVVVWLeTfeeOMuY4xfQ4yH7SK/iGwFMnws+pX3hDHGiEhbVW+AMeaEiAwCPhCRIqCqjXV9MsasAdYA5OTkmGDbHhdoy+SAhDtfDaUnKH3gAeqPHCFzyRK6T78jbM8dCqHM18GDB31+w6orqa6u7vKvMZT8zVdSUlLzj1CDEbYCY4xp82euIvKViGQaY06KSCbwdRv7OOG5/0JECoBvAZuB7iISZ4xpAhzAiZC/ANVpfPPpp5T+7GGM00n/F18kbcJ4u0NSKipFyjWYLcC9nsf3Am+1XkFEeohIoudxb2A8cMC4z/F9CMxob3sVHSrfeINjc+YSm57O1Zs2anFRykaRUmCWAlNE5DAw2TONiOSIyFrPOtlAoYjswV1QlhpjDniWLQAeEZFioBcQnT+bjWKmoYFTeUs4+fivSR03jqtf30TiwIF2h6VUVIuIAmOMqTDG3GyMGWqMmWyMOeOZX2iMuc/z+O/GmFHGmNGe+3Ve239hjBlrjBlijLnTGFNv12tR4ddQWkrJrHs4+8or9Lz3XvqvfoHYbt3sDku1EkkdLYPx3nvvce211zJkyBCWLm05KnVlZSUzZsxg2LBhZGdn88knn7S5H6u7abYXp7dwdN2MiAKjVLCqt27l6B3TaSgpISt/JX1/ubBrD7ffiUVSR8tAOZ1OfvrTn/Luu+9y4MABXnvtNQ4cONC8/OGHH+aWW27h888/Z8+ePWRnZ7e5Lyvz0FGc3sLRdVMLjOqUTEMDXz21lNIHHyLhqqsY+KfNdJsyxe6wVDv27t1Lz549mTRpEnfeeSdLliwJ2b7Hjx/fYkBLb0ePHmXatGnk5OQwduxYDh06FPD+d+7cyZAhQxg0aBAJCQncfffdvPWW+1JvVVUV27dvZ+5c9++7ExISmscvy83NpaysrMW+rMxDe3G2NnHixDZzFir6p57qdOoPH6ZswULqDhygxz33cOVjv9A+Ln5atnMZn5/5PKT7HNZzGAvGLuhwPTs6WjY2NnLfffexZs0aBg8ezDvvvMPSpUvZsGGDX89/wYkTJ+jfv3/ztMPhYMeOHYC7gPXp04fZs2ezZ88evvOd7/Dss8+Smpoa9m6a7cVpBy0wqtMwTidnXv6/nH7mGWJSU8nKX6lHLZ2EXR0t33zzTfbv38/06dMBaGpq4oYbbmixzuTJkzl16tQl2y5ZsoRp06Z1+BxNTU189tln5OfnM27cOB5++GGWLl3Kb3/720vWjbZumlpgVKfQUFJC2eOPc75wF2k330zmk08Q17u33WF1Ov4caVjBro6We/bsYcmSJc2nr3zZunVrh/vJysri+PHjzdOlpaXNDcwcDgcOh4Nx48YBMGPGjDYvrlvdTbO9OO2gBUZFNFdDAxVr11Kx+kUkMZHMp54i/fZpSCcZ8kW52dXRMjMzk/fff5/Zs2cTExNDUVERI0eODPj9c/3113P48GGOHj1KVlYWGzdu5I9//CMAGRkZ9O/fn0OHDnHttdeybds2hg8fDoS/m2Z7cdpBL/KriPVNYSFHf3wH5SvzSbv5Jgb99W26//h2LS6dkF0dLefMmYPL5SI7O5sxY8awbNmyoN4/cXFxrFq1iqlTp5Kdnc1dd93V4kgrPz+fWbNmcd1117F7924WLVpkSzfNjuL0/tJBOLpuRsRgl3bJyckxhYWFQW2rY5EFJpB8NZ46xekVK6h6awvx/fqRsfg3pP3gB9YGGGFCPRZZe1+btVNubi4DBgwgJiaGHj16tOguGYhIHIts3759rF+/3u+GZ/v37ycvL++yu2n6w998+XrviEjkDXapVEdctbVUrFtPxfr14HLR6/776T1/HjEpKXaHpixidUdLO2k3TS0wKgKYpiaq3trC6Weeoen0abrl5tLnkUdIcNh3cVIpdfm0wCjbmKYmzv31r5Q//wINx46RNPo6slY+S8plDA+ulIocWmBU2JmmJqrefpuKF1bTcOwYicOGkZW/kismT9YL+Ep1IVpgVNi46uqo+stfqFi7lsZjX5KYnY1jVT5pN92ExOgXGpXqarTAKMs1njhB2p/+TPGChTirqkgcno3juVXuwqJHLEp1WVpglCWMMXyz81POvvIHqrd9QAqQMmUKPe+ZRXJOjhYWpaKAFhgVUk1nznDu7bep/I83qD98mNju3ek1dy6fD7ya4XfcYXd4Sqkw0gKjLpurvp6a7dupevMtaj76CJqaSBo+nMwleXS77TZikpI4UFBgd5hKqTDTAqOC4mpooPZvH3PuvXep2fYBrtpaYnv3pue//ivpt08j6Zpr7A5RKWWziPjqjoj0FJH/EpHDnvsePta5UUR2e93qROR2z7KXROSo17Ix4X8VXV/T2bNUbdnCiUce5fD3x1P6wAPUfLSdK269hf6//z1DCz6k72O/0OKifIqklsnR1LbYThFRYICFwDZjzFBgm2e6BWPMh8aYMcaYMcBNwDfAf3qt8osLy40xu8MSdRdnnE7O79tP+eoXKZn5zxweP4GyxxZQu2MHV0z9If1/v4Zr/vbf9MvLI+2GCdqqWLUrklomR1PbYjtFyifCNGCS5/HLQAHQXuOKGcC7xphvrA0ruhiXi/pDh/hm505qd+zkm8JCXOfOAZA0YgS9588n7cZJJI0Yob9bUQHbu3cvEyZMYNKkScybN4958+aFbN/jx4+noqLC57Lc3FzWrl1Lv379whKLd9tioLlt8YUh/L1NnDiRkpKSkD13pImUAtPXGHOhvdspoG8H698NtB5FbomI/AbPEZAxpt7XhiJyP3A/QN++fSkI8uJzTU1N0NtGCqmtJb7kGPHHSogvKSH+yBfEeIZLb+rTh4ZRo2i85hoahl2LKz3dvVFFBWzfHvBzdYV8hVMo85Went7cqKpy+dM0/OMfIdnvBQnXXEP3Rx/pcL2ioiIefPBBnn/+eW699VafzbNamzp1KjU1NZfMz8vL48Ybb2yedjqd1NTU4HK5Ltnvpk2bAFrMDyYWf+MpLi4mIyOjeZ+9evWisLCwzedoK24rOZ1Ov56vrq7ust6HYSswIrIVyPCx6FfeE8YYIyJt9hAQkUxgFPC+1+xf4i5MCcAa3Ec//+Zre2PMGs865OTkmGCHRO9sw/U7q6qoLy6mbt8+zu8t4nxREY1ffuleKELCoEEkT51K6rixpIwdS3xmZkifv7Ply26hHq7/wtDstQnxuGJjQ7LfC+IT4jsc+t27VXBlZWWL9dtrmfz3v//drxiqq6tJS0sjJibmsmKB9tsW+xNPcnIy8fEXc5KcnExCQkKbcfkbdyj5O1x/UlJSc8fNYIStwBhj2uxvKiJfiUimMeakp4B83c6u7gL+bIxpvhrndfRTLyIbgJ+HJOhOqOnsWRqOHKG++Aj1xcXUHymmofgITadPN68Tl5FB8qhRdJ8xg+TrRpE0YgSxEdZLQ1kjY9EiW57XrpbJgcQCXbNtsZ0i5RTZFuBeYKnn/q121p2J+4ilmVdxEuB2YJ9VgdrNNDbS+NXXNJadoOnkSRrLymgsO0nDl19Sf+QIzvLy5nUlJYXEwYNJHT+exCGDSRg8mKThI4jve6WNr0BFI7taJoO2LbZTpBSYpcDrIjIXOIb7KAURyQHmG2Pu80xfDfQHPmq1/asi0gcQYDcwPzxhh5YxBld1NY2nTrUoHo1lZTR6ppu+/hpcrhbbxfbsSbzDQdrEiSQOHkzi0CEkDh5MXGamXoxXEaGoqIjbbrsNaNkqePjw4SxfvpzRo0czffp0amtrSU1NDXj/s2fP5uOPP6a8vByHw8GTTz7J3Llz22xb7CuWzZs3+2xbHBvgKUXvtsVOp5M5c+Y0H5W1/sLBzJkzKSgouCTurkJbJlvQMtkYg6umBmdVFc6zle77qkqcle7HrqoqnJVVzdPOKs/jc+egdR/uuDjiMzKI79eP+MxM4rP6EZeZ6ZnuR3y/TGKSkoJ6DeGk12ACoy2TA9PWNYVIbltsJ22ZHMHOvv46qYWFfPXJJ+5CUdWqWFRVXVoovMSkpBDbvTsx3dOJTU8nMTOD2PR0YtO7E5ueTlzfK93FI6sfcb17IyG+KKtUpLC6ZbK2LbaXFpgglOevIu30aSo7KBSx3bsT61kW290zr1s3pNV5XaWU6oq0wARh4Ja3+NuuXUwK8lssSikVDfQKcBDievQAHRZFKaXapQVGqSgQzV/mUcEJxXtGC4xSXVxSUhIVFRVaZJTfjDFUVFSQdJnfUNXzPEp1cQ6Hg9LSUk57jebQ1dTV1V32h2E08SdfSUlJOByOy3oeLTBKdXHx8fEMHDjQ7jAsVVBQcFljZkWbcOVLT5EppZSyhBYYpZRSltACo5RSyhJRPRaZiJzGPbgmQDpQ5bW4o+neQDnWaP1codymvfXaWuZrfkfzNF+Bz/Oe1nxpviI1XwOMMX38WtMYozd3kV0T4HRhuGIJ5TbtrdfWMl/zO5qn+Qp8Xqv8ab40X50yX943PUV20V8CnLZSMM/l7zbtrdfWMl/zO5qn+Qp8XrhypvkKjOYrSFF9iuxyiEih8XPIaqX5CpTmKzCar8CEK196BBO8NXYH0MlovgKj+QqM5iswYcmXHsEopZSyhB7BKKWUsoQWGKWUUpbQAqOUUsoSWmCUUkpZQguMBUTkdhH5vYhsEpEf2h1PpBORQSKyTkTesDuWSCUiqSLysud9NcvueCKdvqcCY9VnlhaYVkRkvYh8LSL7Ws2/RUQOiUixiCxsbx/GmDeNMf8LmA/8k5Xx2i1E+frCGDPX2kgjT4C5uwN4w/O++lHYg40AgeQrWt9T3gLMlyWfWVpgLvUScIv3DBGJBZ4DbgWGAzNFZLiIjBKRt1vdrvTa9HHPdl3ZS4QuX9HmJfzMHeAAjntWc4YxxkjyEv7nSwWXr5B+ZmnDsVaMMdtF5OpWs8cCxcaYLwBEZCMwzRjzFPA/Wu9DRARYCrxrjPnM2ojtFYp8RatAcgeU4i4yu4nSPwwDzNeB8EYXeQLJl4gcxILPrKh8owYhi4t/PYL7P3tWO+s/BEwGZojIfCsDi1AB5UtEeonIauBbIvJLq4OLcG3l7k/AdBF5ARvGlIpgPvOl76k2tfX+suQzS49gLGCMWQmstDuOzsIYU4H73K9qgzGmFphtdxydhb6nAmPVZ5YewfjnBNDfa9rhmad803wFT3MXGM1XYMKaLy0w/vkUGCoiA0UkAbgb2GJzTJFM8xU8zV1gNF+BCWu+tMC0IiKvAZ8A14pIqYjMNcY0AQ8C7wMHgdeNMfvtjDNSaL6Cp7kLjOYrMJGQLx1NWSmllCX0CEYppZQltMAopZSyhBYYpZRSltACo5RSyhJaYJRSSllCC4xSSilLaIFRSillCS0wSimlLKEFRqkIICIfisgUz+M8Ecm3OyalLpeOpqxUZFgM/JunAdu3iNKulapr0aFilIoQIvIRkAZMMsZUe83fBHxqjPmdbcEpFQQ9RaZUBBCRUUAm0NCquEwD3gZG2RWbUsHSAqOUzUQkE3gVd6vfGhG5xTM/CbjTGPMHIN3GEJUKihYYpWwkIim42yE/aow5CPwW9/UYgF8AaZ7WvyNEJNmmMJUKil6DUSoCichVwGJjzFzP9GLgPWPMDnsjU8p/WmCUUkpZQk+RKaWUsoQWGKWUUpbQAqOUUsoSWmCUUkpZQguMUkopS2iBUUopZQktMEoppSyhBUYppZQl/j/1DEug+SXwDQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "##Time\n",
    "t_max = int(1e5)\n",
    "t = np.linspace(0,t_max,10000)\n",
    "t_0 = 1e2\n",
    "t_1 = t_max-t_0\n",
    "i_max = len(t)\n",
    "i_0 = int(i_max*t_0/t_max)\n",
    "i_1 = i_max-i_0\n",
    "\n",
    "## States\n",
    "x_min = 1e-2\n",
    "x_max = 1e2\n",
    "A_chemo = '{3} + ({0}-{3})*np.heaviside(t-{1},1)*((t-{1})/{2})'.format(x_max,t_0,t_1,x_min)\n",
    "#A_chemo = '1+t/1e5'\n",
    "print(A_chemo)\n",
    "X_chemo =  {'A':A_chemo}\n",
    "X0 = np.ones(s['n_X'])\n",
    "X0[s['spec_index']['E']] = 0.5\n",
    "X0[s['spec_index']['C']] = 0.5\n",
    "\n",
    "## Simulate\n",
    "Legend = []\n",
    "for K_C in [1,0.1]:\n",
    "    for K_B in [1E-6,1]:\n",
    "        Legend.append('$K_B = '+str(K_B)+'; K_C = '+str(K_C)+'$')\n",
    "        parameter = {'K_B':K_B,'K_C':K_C}\n",
    "        dat = st.sim(s,sc=sc,t=t,parameter=parameter,X0=X0,X_chemo=X_chemo,quiet=quiet)\n",
    "        v = dat['V'][:,s['reac_index']['r1']]\n",
    "        x_A = dat['X'][:,s['spec_index']['A']]\n",
    "        plt.semilogx(x_A[-i_1:],v[-i_1:]/K_C)\n",
    "plt.legend(Legend)\n",
    "#plt.ylim((-1,1))\n",
    "plt.xlabel('$x_A$')\n",
    "plt.ylabel('$v/K_C$')\n",
    "plt.grid()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Slowly-varying $x_A$\n",
    "This section explicitly shows the approach used to generate steady-state values using a slowly-varying signal.\n",
    "The first graph below shows the exogenous variable $x_A$: it has a constant value of 0.01  for $t<100$ and then rises as a slow ramp. The flow $v$ has a transient which settles well before $t=100$. In the above graphs, $v$ is plotted against $x_A$ for $t>100$ so that the initial transient does not appear."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZQAAAEMCAYAAADj8ECOAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzt3Xd0VPed///nWx0kUdSFBIjeMQjFYJvYwg3ci4hjO971fteOv8lJNsmm2msnm9j57jr7O1uyv+Rk15v1L/7ll42bcIxbCMbI3Rgkem8CRr0A6mVm3r8/5grLRMAgz2ja+3HOHObeuTN6v8VIL9175/O5oqoYY4wxn1VcqAswxhgTHSxQjDHGBIQFijHGmICwQDHGGBMQFijGGGMCwgLFGGNMQFigGGOMCQgLFGOMMQFhgWKMMSYgEkJdwEjKysrSoqKiYT23s7OT1NTUwBYUAazv2GJ9xxZ/+66srGxW1ewLbqiqYXkDVgH7gUPAw+fZrgxQoORCr7lkyRIdro0bNw77uZHM+o4t1nds8bdvYIv68Xs7LA95iUg88EvgBmAucI+IzB1iu3Tgm8Cmka3QGGPM2cIyUIBLgUOqekRV+4BngduG2O4J4GdAz0gWZ4wxkaTP7R2RrxOugVIAnBi07HLWnSEixcBEVX1tJAszxphI4PZ4eWtfA1/7XRWf/6e36HV7gv41I/KkvIjEAf8C/JUf2z4EPASQm5tLRUXFsL5mR0fHsJ8byazv2GJ9R74T7V7eq+nnw1oPbX1KWiIszU9g/VvvkJYkn9o20H2Ha6DUABMHLRc66wakA/OBChEByAPWisitqrpl8Aup6lPAUwAlJSVaWlo6rIIqKioY7nMjmfUdW6zvyNTc0cvL22opr3Sxp66ThDjh6tk5lC0pZMWsHJIShj4YFei+wzVQNgMzRGQKviC5G7h34EFVPQ1kDSyLSAXw3bPDxBhjolWv28Nbexspr3JRsb8Jt1dZUDCWH98yl1sXFZCRmjTiNYVloKiqW0S+DqwD4oGnVXW3iDyO7+Nra0NboTHGjDxVZbvrNOWVLl7ZUcuprn5y0pN5YPkUypYUMjM3PaT1hWWgAKjq68DrZ6370Tm2LR2JmowxJhTqT/ewZquL8koXh5s6SU6I4/p5eZQVF7B8ehYJ8eHx+aqwDRRjjIll3X0e1u2up7zKxXuHmlGFksnj+cc7p3LjgnzGjkoMdYl/xgLFGGPChNerbK5upbzKxes76+nodVM4fhR/c/UMyooLmJwZ3tPDWKAYY0yIHWvpZE1VDWu2ujjR2k1qUjw3LMhn9ZJCLi3KIC5OLvwiYcACxRhjQqCtp5/Xd9RRXuVic/VJROCKaVl8+7qZrJyXx+ikyPv1HHkVG2NMhPJ4lfcONVNe6WLd7np63V6mZqfyvZWzuGNxARPGjQp1iZ+JBYoxxgTZwYZ2Xqxy8YetNTS09TJ2VCJfKCmkrLiQRRPH4QzQjngWKMYYEwQnO/tYu72W8ioXO1yniY8TVszK5se3FHL1nBySE+JDXWLAWaAYY0yA9Hu8bNznG73+1r5G+j3K3Pwx/PDmudy2aAJZacmhLjGoLFCMMeYzUFV217bxYqWLtdtrae3sIystmfsvK6JsSSFz8seEusQRY4FijDHD0NjWwx+21VBeWcP+hnaS4uO4bm4uZUsKuHJGdtiMXh9JFijGGOOnnn4P6/c0UF7l4p0DTXgVFk8ax09vn88tCycwdnT4jV4fSRYoxhhzHqpK5bGTlFe5eHVHHe09biaMTeGrpdO4s7iQadlpoS4xbFigGGPMEE60dvHS1hrWVLmobuliVGI8N8zPo2xJIZdNzYyY0esjyQLFGGMcHb1u3thZx4uVLjYdbQXgsqmZfG3FdG5YkE9asv3KPB/77hhjYprXq+xp8bD2uW28saue7n4PRZmj+c51M7mjuIDC8aNDXWLEsEAxxsSkw00drKly8VJVDbWne0hPaeD2xQWsXlJA8aTxUTN6fSRZoBhjYsbprn7W7qhlTZWLrcdPESdw5cxsbitSvrl6BSmJ0Td6fSRZoBhjoprb4+Wdg02UV9awfk8DfR4vs3LT+bsbZ3P7ogJyxqRQUVFhYRIAFijGmKi0p7aN8ioXL2+robmjj4zUJO5dOonVSwqZN2GMHdIKAgsUY0zUaO7o5eVttZRXuthT10ZivHD17BzKigspnZVDUkLsjV4fSRYoxpiI1uv2sGFvI+WVLioONOHxKgsLx/KTW+dx6yUTGJ+aFOoSY4YFijEm4qgq212nebHyBK9sr+N0dz+5Y5J58PNTWF1cyIzc9FCXGJMsUIwxEaPudLfv2utVLg43dZKcEMfKeb7R68unZxFvo9dDygLFGBPWuvrc/HFXPWuqanj/cDOqcGlRBl/+/FRuXJjPmJTYnpAxnFigGGPCjterbDraSnmVizd21tHZ52Fixii+cfUM7iwuYHJmaqhLNEOwQDHGhI3q5k7WVLlYs7UG18lu0pITuGlhPmXFhXyuKMMmZAxzFijGmJBq6+nntR11lFe62HLsJCKwfHoW371+Fivn5TEqyQYcRgoLFGPMiPN4lXcPNlFeVcOfdtfT6/YyLTuV76+axR2LC8gfOyrUJZphsEAxxoyYAw3tlFe6eGlrDY3tvYwdlchdJRNZvaSQhYVjbfR6hLNAMcYEVWtnH2u31VBeVcPOmtMkxAmls7IpKy7k6jk5JCfYIa1oYYFijAm4PreXjft9o9c37m+k36PMzR/DD2+ey22LJpCVlhzqEk0QWKAYYwJCVdlV45uQce32Wlo7+8hKS+b+y4ooW1LInPwxoS7RBJkFijHmM2lo6zlz7fUDDR0kJcRx3dxcyooLuHJGNgnxNiFjrLBAMcZctJ5+D+t211NeVcN7B5vwKhRPGsf/uWM+Ny+YwNjRNno9FlmgGGP8oqpsOXaS8koXr+2oo73XTcG4UXxtxXTuWFzA1Oy0UJdoQmzEAkVEkoGFqrp5pL6mMeazO9Ha5ZuQcauLYy1djE6K54b5+ZQtKWDZlEwbvW7OCGqgiMi3gQXAOGAu8P8BFijGhLmOXjev7/SNXt90tBURuGxqJt+4egar5ueRmmwHN8yfC/a74lLgFVX9nYj8SlWf8OdJIrIK+DkQD/xaVZ886/FvAw8CbqAJ+GtVPRbY0o2JLR6vsrvZw8vPbeOPu+rp7vcwJSuV714/kzuKCykYZ6PXzfkFNVBU9W4RuVVEfgv4NT2oiMQDvwSuA1zAZhFZq6p7Bm22FShR1S4R+SrwT8AXA1y+MTHhcFPHmdHrdad7SE9p4I7iAsqKCymeNM5Grxu/BX2/VVXXisg64Psi8l+q+uULPOVS4JCqHgEQkWeB24AzgaKqGwdt/xFwX4DLNiaqnerq45XttbxYVcP2E6eIjxOunJHFHUXKN1avICXRRq+bixeUQBGRnwPfUlUFUNVe4AkRmejH0wuAE4OWXcDS82z/APDGcGs1Jlb0e7y8vb+J8ioXG/Y20ufxMjsvnUdvnMNtiyeQk55CRUWFhYkZNnF+5wf2RUV+ClwC3K2qnSKyEviRql7hx3NXA6tU9UFn+S+Apar69SG2vQ/4OnCVE1pDvd5DwEMAubm5S5599tlh9dTR0UFaWux9LNL6jnzH2jy8X+Pmozo3bX2QngTL8hNYXpDApPS4Tx3Siqa+L4b1fX4rVqyoVNWSC20XlD0UVX1MRO4FKkSkD+gAHvbz6TXA4D2ZQmfdp4jItcCjnCdMnFqeAp4CKCkp0dLSUj/L+LSKigqG+9xIZn1Hpqb2Xl7eVsOLlS721XeSGC9cM9t37fXSWdkknmP0eqT3PVzWd2AE65DXNcCXgU4gH9+nsPb7+fTNwAwRmYIvSO4G7j3r9RcD/4lvT6YxYIUbE8F6+j1s2NtIeZWLtw804fEqlxSO5fHb5nHLwgmMT00KdYkmygXrpPyjwA9V9T0RWQA8JyLfVtW3LvREVXWLyNeBdfg+Nvy0qu4WkceBLaq6Fvi/gDTgBWd3/biq3hqkXowJW6pK1fFTlFe5eHV7LW09bvLGpPDlz09l9ZICpuekh7pEE0OCdcjr6kH3d4rIDUA5cLmfz38deP2sdT8adP/aAJVqTESqOdXNS1Uu1lTVcKS5k5TEOFbN8x3SunxaFvE2et2EwIgMd1XVOucwmDFmmLr63Lyxs57yKhcfHmlBFS6dksFXrprGDQvySE+xCRlNaI3Y/Amq2j1SX8uYaOH1Kh8dbaG8soY3dtXR1edhUsZovnnNDO5cXMikzNGhLtGYM2xCHmPCUHVzJ+XOIa2aU92kJSdwy8IJlC0p5HNF4230uglLFijGhInT3f28tqOO8ioXlcdOEidwxfQsvr9qFtfPzWNUkg04NOHNAsWYEHJ7vLx7qJnyShd/2tNAn9vL9Jw0frBqNncsLiBvbEqoSzTGbxYoxoTA/vp2yqt8EzI2tfcybnQi93xuImVLCllQMNYOaZmIZIFizAhp6ehl7fZayqtc7KppIyFOKJ2Vw+olBayYnUNygh3SMpHNAsWYIOpze3lrXwMvVtZQsb8Rt1eZXzCGv79lLrdeMoHMtORQl2hMwFigGBNgqsoO12nKq1ys3V7Lqa5+stOT+evlUygrLmRWno1eN9HJAsWYAKk/3cNLW2sor3JxqLGD5IQ4rp+Xx53FBXx+ehYJ55iQ0ZhoYYFizGfQ3efhT3vqebHSxfuHmvEqlEwezz/euYAbF+QzdpSNXjexwwLFmIukqmyuPkl5pYvXdtbR0eumYNwovr5iOncWF1KU5dfVro2JOhYoxvjpRGvXmdHrx1u7GJ0Uz40L8ikrLmTplAzibEJGE+MsUIw5j/aefl7fWUd5ZQ0fV7ciApdPy+Rb185g5bw8UpPtR8iYAfbTYMxZPF5lV7OHPzy7lT/urqen38vUrFS+t3IWty8uoGDcqFCXaExYskAxxnGoscM3er2qhvq2HsakNFJWXEjZkkIWTxxno9eNuQALFBPTTnX18cr2Wl6sqmH7iVPExwlXzczmzqnKN8pWkJJoo9eN8ZcFiok5/R4vb+9vorzKxYa9jfR5vMzOS+exm+Zw66IJ5KSnUFFRYWFizEWyQDExQVXZXdvmG72+rZaWzj4yU5O4b9lkypYUMG/C2FCXaEzEs0AxUa2xvYeXt/omZNxX305SfBzXzMmhrLiQq2Zlk2ij140JGAsUE3V6+j2s39NAeZWLdw404VVYNHEcT9w+n1sW5jNudFKoSzQmKlmgmKigqlQdP0V5lYtXttfS3uMmf2wKX7lqGncWFzI9Jy3UJRoT9SxQTESrOdXNS1UuyqtqONrcyajEeFbNz6OsuJDLpmUSb6PXjRkxFigm4nT2uvnjrnrKq1x8eKQFVVg6JYOvlk7jxgX5pNnodWNCwn7yTETwepWPjrZQXlnDG7vq6OrzMDlzNN+6ZiZ3FhcwMWN0qEs0JuZZoJiwdrS5k/JK37XXa051k56cwK2XTKBsSSElk8fb6HVjwogFigk7p7v7eXVHLeWVLqqOnyJOYPmMbL6/ahYr5+XZgENjwpQFigkLbo+Xdw42UV5Vw/o9DfS5vczISePhG2Zzx+ICcsekhLpEY8wFWKCYkNpb10Z5pYs/bKuluaOX8aMTuffSSZQVFzK/YIwd0jImgligmBHX3NHLy9t8h7T21LWRGC+smJVD2ZJCVszKISnBRq8bE4ksUMyI6HV7eGtvI+VVLir2N+H2KgsKxvLjW+Zy66ICMlJt9Loxkc4CxQSNqrLddZryShev7KjlVFc/OenJPLB8CmVLCpmZmx7qEo0xAWSBYgKu7nQ3L22tobzSxeGmTpIT4rh+Xh5lxQUsn55Fgk3IaExUskAxAdHd52Hdbt/o9fcONaMKJZPH8493TuWmhfmMSUkMdYnGmCCzQDHD5vUqm6tbKa9y8frOejp63RSMG8XfrJjOncWFFGWlhrpEY8wIskAxF+14SxflVS7WbHVxorWb1KR4blyQT9mSQi4tyiDOJmQ0JiZZoBi/tPf08/rOOsora/i4uhURuGJaFt++biYr5+UxOsneSsbEurD9LSAiq4CfA/HAr1X1ybMeTwb+X2AJ0AJ8UVWrg1FLd5+HHrfS2esOxsuHLQV2Nbt56dmtrNtdT0+/l6nZqXxv5SzuWFzAhHGjQl2iMSaMhGWgiEg88EvgOsAFbBaRtaq6Z9BmDwAnVXW6iNwN/Az4YjDqufUX73GwsQveXBeMlw97Y0c1sXpJIWXFhSyaOM5GrxtjhhSWgQJcChxS1SMAIvIscBswOFBuA37s3H8R+IWIiKpqoIv58pVTqdy5j2nTpgb6pcNeW+1R/mb1CpITbEJGY8z5hWugFAAnBi27gKXn2kZV3SJyGsgEmgNdzF0lE8npOEzpldMC/dJhr6LihIWJMcYv4RooASMiDwEPAeTm5lJRUTGs1+no6Bj2cyOZ9R1brO/YEui+wzVQaoCJg5YLnXVDbeMSkQRgLL6T85+iqk8BTwGUlJRoaWnpsAqqqKhguM+NZNZ3bLG+Y0ug+w7XQNkMzBCRKfiC427g3rO2WQvcD3wIrAbeutD5k8rKymYROTbMmrIIwuG0CGB9xxbrO7b42/dkf14sLAPFOSfydWAdvo8NP62qu0XkcWCLqq4F/hv4rYgcAlrxhc6FXjd7uDWJyBZVLRnu8yOV9R1brO/YEui+wzJQAFT1deD1s9b9aND9HuALI12XMcaYodm0r8YYYwLCAsV/T4W6gBCxvmOL9R1bAtq3BGEcoDHGmBhkeyjGGGMCwgLFGGNMQFigXICIrBKR/SJySEQeDnU9gSYiT4tIo4jsGrQuQ0TWi8hB59/xznoRkX93vhc7RKQ4dJUPn4hMFJGNIrJHRHaLyDed9VHdN4CIpIjIxyKy3en9J876KSKyyenxORFJctYnO8uHnMeLQln/ZyEi8SKyVURedZajvmcAEakWkZ0isk1EtjjrgvJet0A5j0GzHt8AzAXuEZG5oa0q4H4DrDpr3cPABlWdAWxwlsH3fZjh3B4CfjVCNQaaG/iOqs4FlgFfc/5fo71vgF7galW9BFgErBKRZfhm6/5XVZ0OnMQ3mzcMmtUb+Fdnu0j1TWDvoOVY6HnAClVdNGjMSXDe66pqt3PcgMuAdYOWHwEeCXVdQeizCNg1aHk/kO/czwf2O/f/E7hnqO0i+Qa8jO9SCbHW92igCt/Eq81AgrP+zPse3+Diy5z7Cc52Eurah9FrofOL82rgVUCivedBvVcDWWetC8p73fZQzm+oWY8LQlTLSMpV1Trnfj2Q69yPuu+HczhjMbCJGOnbOfSzDWgE1gOHgVOqOnAFucH9fWpWb2BgVu9I82/A9wGvs5xJ9Pc8QIE/iUilM1kuBOm9HrYj5U14UFUVkaj8bLmIpAHlwLdUtW3whcOiuW9V9QCLRGQc8BIwO8QlBZWI3Aw0qmqliJSGup4QWK6qNSKSA6wXkX2DHwzke932UM7Pn1mPo1GDiOQDOP82Ouuj5vshIon4wuR3qrrGWR31fQ+mqqeAjfgO94xzZu2GT/d3pvfzzeod5q4AbhWRauBZfIe9fk5093yGqtY4/zbi+wPiUoL0XrdAOb8zsx47nwC5G98sx9FuYCZnnH9fHrT+L51PgiwDTg/abY4Y4tsV+W9gr6r+y6CHorpvABHJdvZMEJFR+M4d7cUXLKudzc7ufeB74tes3uFGVR9R1UJVLcL3M/yWqn6JKO55gIikikj6wH3gemAXwXqvh/qEUbjfgBuBA/iOMz8a6nqC0N/vgTqgH9/x0gfwHS/eABwE3gQynG0F36feDgM7gZJQ1z/MnpfjO668A9jm3G6M9r6dXhYCW53edwE/ctZPBT4GDgEvAMnO+hRn+ZDz+NRQ9/AZ+y8FXo2Vnp0etzu33QO/w4L1XrepV4wxxgSEHfIyxhgTEBYoxhhjAsICxRhjTEDE1DiUrKwsLSoqGtZzOzs7SU1NDWxBEcD6ji3Wd2zxt+/Kyspm9eMS6jEVKEVFRWzZsmVYz62oqKC0tDSwBUUA6zu2WN+xxd++ReSYP69nh7yMMcYEhAWKMcZEsbrT3fxx18iMw42pQ17GGBPtet0eKqtP8vaBJir2N7G/oR0R2PLotWSmJQf1a1ugGGNMhDvR2kXFgSbe3t/IB4db6OrzkBgvXDolg7Ils7lqZg4ZqUlBryMkgSIiGcBz+K7DUQ3cpaonh9jufuAxZ/GnqvrMWY+vxTctwvygFmyMMWGk1+3h46OtVOxvomJ/I4ebOgGYmDGKsuJCrpqZzWXTMklNHtlf8aHaQxm4WtiT4rus7sPADwZv4ITO3wMl+OZdqhSRtQPBIyJ3Ah0jW7YxxoSG62TXmQB5/1AL3f0ekhLiWDolg3uXTqZ0VjZTs1IZfBmGkRaqQLkN3yRtAM8AFZwVKMBKYL2qtgKIyHp8l6r9vXMdi2/ju0Tl8yNQrzHGjKh+j5ct1Sep2N/Ixv2NHGjw/f08MWMUXygppHRWNsumZjI6KXzOXISqknNdLWyw81057Angn4GuoFVojDEjrKm9l437G9m4r5H3DjbT3us+cy7krpKJrJidE/K9kPMJWqCIyJtA3hAPPTp4QfXirhYmIouAaar6t87lWy+0/UP49mTIzc2loqLC3y/1KR0dHcN+biSzvmOL9T2yvKpUt3nZ3uhhR5OHo22+KxSPSxYWZ8dzSXYyczPjGZXQA57jnNh9/FN/ZX9Wge47aIGiqtee6zERaRCRfFWtO+tqYYPV8MlhMfBdOawC39XlSpyrryUAOSJSoaqlDEFVnwKeAigpKdHhjoa1kbSxxfqOLSPZd0evm/cONvPWvgbe2tdEc0cvIrB44jjKluWwYnYOc/PHjMheSKD7DtUhr4GrhT3Jp68WNtg64B9EZLyzfD3wiHNO5VcAzh7Kq+cKE2OMCQcnWrvYsLeBDfsa2XSklT6Pl/SUBK6amc3Vs3MonTUyH+sNtlAFypPA8yLyAHAMuAtAREqAr6jqg6raKiJP4LsML8DjAyfojTEmnHm8yrYTJ3lzbyMb9jacOaE+NTuV+y+fzNWzcykpGk9ifHRNVhKSQFHVFuCaIdZvAR4ctPw08PR5XqcasDEoxpiQ6+x18+7BJtbv8X0qq7Wzj/g44dKiDB67aSLXzMllSlZ0z2gcPp83M8aYCFN3upsNext5c28DHxxqoc/jZeyoREpnZXPNnFyumpnN2FGJoS5zxFigGGOMn1SVPXVtvLnHFyI7a04DMDlzNH952WSunZtLyeTxJETZoSx/WaAYY8x59Hu8fHy0lfV7Gli/p4GaU91nPpX1g1WzuW5uDtOy08J2bMhIskDxg6ri9vo9VMYYE+E6et28c6CJP+2u5619jbT1uElOiOPzM7L4xjXTuXp2LtnpwZ25NxJZoPhh1b+9y/i4Xq69OtSVGGOCpam9lw17G/jTngbeO9RMn9tLRmoS18/L47q5uXx+RlZYTXMSjuy744fU5HjaO20PxZhoc6ylk3W763nho24OrXsTVd9cWX+5bDLXz8tjyeTxxMfZoSx/WaD4ISM1mf1/Nrm+MSbSDJxUX7e7gT/trmdffTsAk9Lj+NY1M7l+Xi6z89LtfMgwWaD4ITM1ifY+20MxJhJ5vcrWEyd5Y2c96/bUc6LVd1L9c5MzeOymOaycl8fhHR9TWjoj1KVGPAsUP2Sm+QJFVe0vF2MiQL/Hy6Yjrfxxdx3rdjfQ1N5LYrxwxfQsvlY6nWvn5pI16HK4h0NYazSxQPFDRmoSHoW2HndMDVIyJpL0uj28d7CZN3bV8+beBk519TMqMZ4Vs7NZOS+PFbNzGJNiP7/BZIHih8w036RtLR29FijGhJHuPg9vH2jkjV31bNjbSEevm/SUBK6dk8uq+XlcNTOblMT4UJcZMyxQ/JCZ6ts1bu3sY2p2iIsxJsZ19rrZuL+R13fWsXFfE939HsaPTuSmBfncsCCPy6dlkZQQmyPVQ80CxQ8D00q3dPaFuBJjYlNHr5sNext4fWcdFfub6HV7yUpLpmxJATfMz2fplIyYne4knFig+OGTQ14WKMaMlPaefjbsbeS1nXW8faCJPreXnPRk7v7cRG5ckE9JUYaNEQkzFih+GNhDae3sDXElxkS3gT2RV3d8EiJ5Y1L40tJJ3LggnyWTxhNnIRK2LFD8kJwQz6gEO+RlTDB09rrZsK+R13bUsnH/p0Pk5oX5LJ5oIRIpLFD8lJ4kdsjLmADp7vOwcX8jr+6o5a19jfT0e8kdk8yXlk7ipgX5FNueSESyQPHTmCSh1fZQjBm2XreHt/c38cqOOjbsbaCrz0NWWjJfLJnITQsnUDLZQiTSWaD4KT1J7JCXMRep3+Plg8MtvLK9lnW762nvcZORmsTtiwu4eWE+S6dk2on1KGKB4qf0JKHmtJ2UN+ZCvF7l4+pWXtleyxu76mnt7CM9OYGV8/O49ZIJXD4t0z7iG6UsUPw0Jkk42dVn83kZMwRVZWfNadZuq+XVHXXUt/UwKjGea+fmcsvCfK6alU1ygo1Yj3YWKH5KTxL6PWrzeRkzyKHGDtZuq2Ht9lqqW7pIjBeumpnN3900h2vn5NgFqWKM/W/7KT3Jt1di83mZWFd/uoe122t4eVstu2vbEIHLpmbylaumsWp+HuNGJ4W6RBMiFih+SncyxObzMrHodHc/b+ys4+VttXx0tAVVuKRwLI/dNIdbLplA7piUUJdowoAFip/GJDt7KPZJLxMjet0eNu5r5KWtNWzc10Sfx8uUrFS+ec0MbltUwJSs1FCXaMKMBYqfPjnkZYFiopdXlY+OtPDythpe21FHW4+brLRk7ls2mdsXT2BBwVj7UIo5JwsUPw0Eis3nZaLRwYZ21myt4fmPumnp+YjRSfGsmpfH7YsL7GO+xm8WKH5KjBPSkxPskJeJGk3tvazdXstLW13sqmkjPk6YlxHHj25fwHVzc+0TWuai2TvmImSkJdkhLxPRevo9/GlPA2uqXLx7sBmPV5lfMIYf3jyXWy+ZwO7KDyldVBDqMk2EskC5CBmpSTafl4k4qsrm6pOsqXLx2o462nvd5I9N4aErp3Ln4gJm5KaHukQTJSxH7pJNAAAQ2ElEQVRQLkJmajI1p7pDXYYxfjne0kV5lYs1W12caO1mdFI8N8zPp6y4gGVTM20iRhNwIQkUEckAngOKgGrgLlU9OcR29wOPOYs/VdVnnPVJwC+AUsALPKqq5cGuOzM1iR2uU8H+MsYMW0evm9d31PFilYuPj7YiAldMy+Jvr53Jynl5pCbb35AmeEL17noY2KCqT4rIw87yDwZv4ITO3wMlgAKVIrLWCZ5HgUZVnSkicUDGSBSdkZZk83mZsOP1+j7q+2Klizd21dPd72FqVirfWzmLOxYXMGHcqFCXaGLEBQNFRN4HHlPVjQH8urfh27sAeAao4KxAAVYC61W11aljPbAK+D3w18BsAFX1As0BrO2cMlOTbD4vEzZOtHbxQqWL8koXNae6SU9J4PbFBaxeUkjxpHH2R48Zcf7sofxv4Cci8hi+YPkwAF83V1XrnPv1QO4Q2xQAJwYtu4ACERnnLD8hIqXAYeDrqtoQgLrOKzPNN0eRzedlQqWrz80bO+t5ofIEHx3xHdJaPj2LH9wwm+vn5pKSaDP6mtARVfVvQ5Fi4HFn8VFV3X6B7d8E8oZ46FHgGVUdN2jbk6o6/qznfxdIUdWfOss/BLqB3wBNwBdU9UUR+TawWFX/4hx1PAQ8BJCbm7vk2WefvWCvQ+no6OBodwr/XNnLo0tTmDE+Nn5wOzo6SEtLC3UZIy6c+lZVDp/y8k6Nm4/r3PR4IGe0sLwggSsmJJA5KnCDDsOp75FkfZ/fihUrKlW15ELbXcw5lEPAE/gON1Ve6Lmqeu25HhORBhHJV9U6EckHGofYrIZPDosBFOI7NNYCdAFrnPUvAA+cp46ngKcASkpKtLS09FybnldFRQUrFi/mnyvfY9LMeZTOGyoro09FRQXD/Z5FsnDou7mjlzVVLp7f4uJQYxejk+K5eVEhd5VM5HNF44NySCsc+g4F6zsw/DmHshGYgW/vYI9z+6vP+HXXAvcDTzr/vjzENuuAfxCRgT2X64FHVFVF5BV8YfMWcI1TU9B9csjLxqKY4HB7vLx9oInnt5xgw95G3F5lyeTx/KxsATctnECafUrLhDF/3p3fAfaqaiAHYDwJPC8iDwDHgLsARKQE+IqqPqiqrSLyBLDZec7jAyfo8Z3A/62I/Bu+w1//K4C1nVNGqi9QbD4vE2jHW7p4fssJXqg8QUNbL1lpSfz18incVVLI9BwbeGgiwwUDRVWrAv1FVbUF357F2eu3AA8OWn4aeHqI7Y4BVwa6rgtJToi3+bxMwAxMg/Lc5uO8f6iFOIHSWTn85NaJXDMnh0SbkNFEGNt/vkg2n5f5rA40tPP7j4/z0tYaTnX1Uzh+FN+5biarSwrJH2tjRkzkskC5SDaflxmO7j4Pr+2s4/cfH6fy2EkS44Xr5+Vxz+cmcfk0mwbFRAcLlItk83mZi7G3ru3M3kh7j5upWak8euMc7iwuIDMtOdTlGRNQFigXyebzMhfS3efh1R21/G7TcbadOEVSQhw3zs/j7ksnsXRKho1gN1HLAuUi2Xxe5lwONLTzP5uOU17l8u2NZKfy2E1zKCsuZLzzCUFjopkFykWy+bzMYD39Hv64q57fbTrG5uqTJMXHsWp+Hvcutb0RE3ssUC6SzedlAI61dPI/m47zQqWL1s4+ijJH83c3zmb1kolnxisZE2ssUC5SRqrvRGprZx9Ts0NcjBlRbo+XDfsa+d2m47xzoIn4OOG6Obnct2yyfVLLGCxQLlqm89enDW6MHU3tvTy3+Tj/s+k4tad7yBuTwt9eO5Mvfm4ieWNTQl2eMWHDAuUi2XxesUFV+fhoK7/96Bh/3FVHv0f5/IwsfnTLPK6dk0OCjWI35s9YoFwkm88runX2uvnDthp+9X43rnUfkp6SwF8sK+JLyyYxLTv2pjc35mJYoFwkm88rOh1t7uS3Hx7jhcoTtPe4mZQex5N3zufWRRMYnWQ/Jsb4w35ShsHm84oOHq9Ssb+RZz48xjsHmkiMF26Yn8/9l0+m7ch2Vlw6KdQlGhNRLFCGwebzimynu/t5YcsJnvmwmhOt3eSOSebb183k7ksnkpPuO8lecdQ+sWXMxbJAGQabzysyHWxo5zcfVLOmqobufg+fKxrPD1bNZuW8PJsq3pgAsEAZBpvPK3J4vMrGfY385oNq3jvUTFJCHLddMoH7Ly9ifsHYUJdnTFSxQBkGm88r/LX39PPCFhe/+aCa461d5I9N4XsrZ3HPpZNsJLsxQWKBMgw2n1f4qm7u5DcfVPNipYuOXjclkwcOa+Xa2BFjgswCZRhsPq/woqp8eLiFp98/yoZ9jSTECTcvnMD/uqKIhYXjQl2eMTHDAmUYbD6v8NDr9rB2Wy1Pv1/N3ro2MlOT+JsV07lv2WRyxtiUKMaMNAuUYRiYz6vZxqKERHNHL7/76Di//egYzR29zMpN55/KFnLrogmkJMaHujxjYpYFyjAMHPKysSgj61BjO79+9yhrttbQ5/ZSOiubB5dP5YrpmfbhCGPCgAXKMNh8XiNHVfngcAv/9e4RKvY3kZwQR1lxIQ8sL2J6TnqoyzPGDGKBMgwD83nZIa/g6XN7eXVHLb9+9yh76trISkvi29fN5EtLJ5GZlhzq8owxQ7BAGaaMNJt+JRjaevp59uPjPP1eNfVtPczISeNnZQu4bVGBnR8xJsxZoAyTzecVWHWnu/l/3q/m95uO097r5rKpmfzjnQu4ama2XQnRmAhhgTJMmanJuE52hbqMiLevvo2n3jnC2m21eFW5aeEEHvr8VBYU2rQoxkQaC5Rhsvm8hk9V2XS0lf94+zAV+5sYlRjPfcsm88DyKUzMGB3q8owxw2SBMkwD51BsPi//ebzK+j31/MfbR9h24hSZqUl89/qZ3LdsMuNG2/xaxkQ6C5RhykxNwu1V2rrdjB1t06+cT6/bw0tVNTz1zhGONHcyKWM0T9w+ny8sKbQT7cZEEQuUYTozn1dnrwXKOXT0uvmfTcf49btHaWzvZX7BGH5x72JumJ9PvJ1oNybqWKAMk83ndW4tHb385oNqnvmgmrYeN1dMz+Sf77qE5dOz7PCgMVEsJIEiIhnAc0ARUA3cpaonh9jufuAxZ/GnqvqMs/4e4O8ABWqB+1S1OfiVf8Lm8/pzNae6+a93jvDs5uP0ur2snJvHV0qnsWiizfhrTCwI1R7Kw8AGVX1SRB52ln8weAMndP4eKMEXHJUishZoB34OzFXVZhH5J+DrwI9HsH6bz2uQI00d/KriMC9trQHg9sUFfOWqqTY1ijExJlSBchtQ6tx/BqjgrEABVgLrVbUVQETWA6uAFwEBUkWkBRgDHAp+yZ9m83nBnto2fllxiNd31pEUH8d9yybz5SunUjBuVKhLM8aEQKgCJVdV65z79UDuENsUACcGLbuAAlXtF5GvAjuBTuAg8LVgFjuUWJ7Pa+vxk/zirUNs2NdIWnIC//vKaTywfArZ6TbHljGxLGiBIiJvAnlDPPTo4AVVVRHRi3jdROCrwGLgCPB/A48APz3H9g8BDwHk5uZSUVHh75f6lI6Ojj977qh4D/uOnqCiomlYrxkJBvd94KSHlw/1sbvFS2oi3DE9kWsnJ5KaWM/uyvrQFhpgQ/1/xwLrO7YEuu+gBYqqXnuux0SkQUTyVbVORPKBxiE2q+GTw2IAhfgOjS1yXv+w81rP4zsHc646ngKeAigpKdHS0tJzbXpeFRUVnP3cgj3vk5CUQGnp0mG9ZiTYuHEjSRPn8+8bDvLRkVay0pJ45Iap3LdsMqnJ0fshwaH+v2OB9R1bAt13qH4jrAXuB550/n15iG3WAf8gIuOd5evx7YmkAHNFJFtVm4DrgL3BL/nPZaYms7v2NL9+90govnzQeVV58eMeDqzbRE56Mj+8eS73XjqJUUk2GNEY8+dCFShPAs+LyAPAMeAuABEpAb6iqg+qaquIPAFsdp7z+KAT9D8B3hGRfuf5fzXSDQDMzU/nzb0N/PS1kOTZiMhIER6/bR53lUy0Ue3GmPMKSaCoagtwzRDrtwAPDlp+Gnh6iO3+A/iPYNboj7+9biYPXjk11GUExLmGG27+8D2uvqxoJEsxxkSo6D0IPgJEhDEp0T3tSpyNbDfG+Cku1AUYY4yJDhYoxhhjAkJU/R4CEvFEpAnfSfzhyAJGdL6wMGF9xxbrO7b42/dkVb3gNLgxFSifhYhsUdWSUNcx0qzv2GJ9x5ZA922HvIwxxgSEBYoxxpiAsEDx31OhLiBErO/YYn3HloD2bedQjDHGBITtoRhjjAkIC5QLEJFVIrJfRA45V5eMKiLytIg0isiuQesyRGS9iBx0/h3vrBcR+Xfne7FDRIpDV/nwichEEdkoIntEZLeIfNNZH9V9A4hIioh8LCLbnd5/4qyfIiKbnB6fE5EkZ32ys3zIebwolPV/FiISLyJbReRVZznqewYQkWoR2Ski20Rki7MuKO91C5TzEJF44JfADcBc4B4RmRvaqgLuN/iuhDnYwCWaZwAb+OTyADcAM5zbQ8CvRqjGQHMD31HVucAy4GvO/2u09w3QC1ytqpfguxTEKhFZBvwM+FdVnQ6cBB5wtn8AOOms/1dnu0j1TT49M3ks9DxghaouGvQR4eC811XVbue4AZcB6wYtPwI8Euq6gtBnEbBr0PJ+IN+5nw/sd+7/J3DPUNtF8g3f5ROui8G+RwNVwFJ8g9sSnPVn3vf4LiNxmXM/wdlOQl37MHotdH5xXg28im8+1KjueVDv1UDWWeuC8l63PZTzG/IyxCGqZSSd6xLNUff9cA5nLAY2ESN9O4d+tuG7sN164DBwSlXdziaD+zvTu/P4aSBzZCsOiH8Dvg94neVMor/nAQr8SUQqnSvYQpDe6zbbsDkv1Yu7RHMkEZE0oBz4lqq2yaCZlaO5b1X1AItEZBzwEjA7xCUFlYjcDDSqaqWIlIa6nhBYrqo1IpIDrBeRfYMfDOR73fZQzq8GmDhoudBZF+0anEszc9YlmqPm+yEiifjC5HequsZZHfV9D6aqp4CN+A73jBORgT8wB/d3pnfn8bFAywiX+lldAdwqItXAs/gOe/2c6O75DFWtcf5txPcHxKUE6b1ugXJ+m4EZzqdBkoC78V2+ONoNXKIZPn2J5rXAXzqfBFkGnB602xwxxLcr8t/AXlX9l0EPRXXfACKS7eyZICKj+OQS2huB1c5mZ/c+8D1ZDbylzsH1SKGqj6hqoaoW4fsZfktVv0QU9zxARFJFJH3gPr5Lqe8iWO/1UJ8wCvcbcCNwAN9x5kdDXU8Q+vs9UAf04zte+gC+48UbgIPAm0CGs63g+9TbYWAnUBLq+ofZ83J8x5V3ANuc243R3rfTy0Jgq9P7LuBHzvqpwMfAIeAFINlZn+IsH3IenxrqHj5j/6XAq7HSs9Pjdue2e+B3WLDe6zZS3hhjTEDYIS9jjDEBYYFijDEmICxQjDHGBIQFijHGmICwQDHGGBMQFijGGGMCwgLFGGNMQFigGBNiIlIoIl8MdR3GfFYWKMaE3jVAxF60y5gBNlLemBASkeX45lE6BbQDd6rqkdBWZczwWKAYE2Ii8kfgu6q664IbGxPG7JCXMaE3C9h3wa2MCXMWKMaEkIhk4Zsi3H3BjY0JcxYoxoRWEVAb6iKMCQQLFGNCax+QJSK7ROTyUBdjzGdhJ+WNMcYEhO2hGGOMCQgLFGOMMQFhgWKMMSYgLFCMMcYEhAWKMcaYgLBAMcYYExAWKMYYYwLCAsUYY0xA/P/h43wBvh+nGgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.subplot(2,1,1)\n",
    "plt.plot(t[:5*i_0],x_A[:5*i_0])\n",
    "plt.xlabel('$t$')\n",
    "plt.ylabel('$x_A$')\n",
    "plt.grid()\n",
    "plt.subplot(2,1,2)\n",
    "plt.plot(t[:5*i_0],v[:5*i_0])\n",
    "plt.xlabel('$t$')\n",
    "plt.ylabel('$v$')\n",
    "plt.grid()\n",
    "plt.show()"
   ]
  }
 ],
 "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.6.8"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}