doc/example/conversion_reaction.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Conversion reaction\n",
    "==================="
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import amici\n",
    "import amici.plotting\n",
    "import pypesto"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Compile AMICI model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "libSBML Warning (SBML unit consistency): In situations where a mathematical expression refers to a compartment, species or parameter, it is necessary to know the units of the object to establish unit consistency. In models where the units of an object have not been declared, libSBML does not yet have the functionality to accurately verify the consistency of the units in mathematical expressions referring to that object. \n",
      " The units of the <parameter> 'R1_k2' cannot be fully checked. Unit consistency reported as either no errors or further unit errors related to this object may not be accurate.\n",
      "\n",
      "libSBML Warning (SBML unit consistency): If the attribute 'units' on a given Parameter object has not been set, then the unit of measurement associated with that parameter's value is undefined.\n",
      "Reference: L3V1 Section 4.7.3\n",
      " The <parameter> with id 'R1_k2' does not have a 'units' attribute.\n",
      "\n",
      "libSBML Warning (SBML unit consistency): In situations where a mathematical expression refers to a compartment, species or parameter, it is necessary to know the units of the object to establish unit consistency. In models where the units of an object have not been declared, libSBML does not yet have the functionality to accurately verify the consistency of the units in mathematical expressions referring to that object. \n",
      " The units of the <parameter> 'R1_k1' cannot be fully checked. Unit consistency reported as either no errors or further unit errors related to this object may not be accurate.\n",
      "\n",
      "libSBML Warning (SBML unit consistency): If the attribute 'units' on a given Parameter object has not been set, then the unit of measurement associated with that parameter's value is undefined.\n",
      "Reference: L3V1 Section 4.7.3\n",
      " The <parameter> with id 'R1_k1' does not have a 'units' attribute.\n",
      "\n",
      "libSBML Warning (SBML unit consistency): In situations where a mathematical expression contains literal numbers or parameters whose units have not been declared, it is not possible to verify accurately the consistency of the units in the expression. \n",
      " The units of the <kineticLaw> <math> expression 'compartment * (R1_k1 * A - R1_k2 * B)' cannot be fully checked. Unit consistency reported as either no errors or further unit errors related to this object may not be accurate.\n",
      "\n",
      "libSBML Warning (Modeling practice): As a principle of best modeling practice, the units of a <parameter> should be declared rather than be left undefined. Doing so improves the ability of software to check the consistency of units and helps make it easier to detect potential errors in models.\n",
      " The <parameter> with the id 'R1_k2' does not have a 'units' attribute.\n",
      "\n",
      "libSBML Warning (Modeling practice): As a principle of best modeling practice, the units of a <parameter> should be declared rather than be left undefined. Doing so improves the ability of software to check the consistency of units and helps make it easier to detect potential errors in models.\n",
      " The <parameter> with the id 'R1_k1' does not have a 'units' attribute.\n",
      "\n"
     ]
    }
   ],
   "source": [
    "import importlib\n",
    "import os\n",
    "import sys\n",
    "import numpy as np\n",
    "\n",
    "# sbml file we want to import\n",
    "sbml_file = 'conversion_reaction/model_conversion_reaction.xml'\n",
    "# name of the model that will also be the name of the python module\n",
    "model_name = 'model_conversion_reaction'\n",
    "# directory to which the generated model code is written\n",
    "model_output_dir = 'tmp/' + model_name\n",
    "\n",
    "# import sbml model, complile and generate amici module\n",
    "sbml_importer = amici.SbmlImporter(sbml_file)\n",
    "sbml_importer.sbml2amici(model_name,\n",
    "                         model_output_dir,\n",
    "                         verbose=False)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Load AMICI model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYkAAAEZCAYAAABiu9n+AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4wLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvqOYd8AAAIABJREFUeJzt3XmcHVWZ//HPt/csnYUskH1hTQBBaAmb7GjQYRlEkEXZJIID4ziO/GBcUEfHZZTBcRgwLIKIRFnUjKLgCAgmISYhbAkgIRDSSSAhe0K6k+5+fn9UdXLT9E36Jn3v7eX7flGve+vUqaqnmnQ9fWo5RxGBmZlZa0qKHYCZmXVcThJmZpaVk4SZmWXlJGFmZlk5SZiZWVZOEmZmlpWThFmeSZon6YQOEMetkr5S7Disc3GSsA5F0rGSpktaK2mVpGmSPpAuu0TSX3LY1mhJIalsF2PZrfWbRcSBEfHE7mxD0l2SvrmbcVwZEf+2O9uw7me3/vGbtSdJfYDfAlcBvwQqgA8C9cWMa0cklUVEQ7Hj2BlJpRHRWOw4rPNxS8I6kv0AIuK+iGiMiE0R8WhEPC9pHHArcJSkDZLWAEj6qKS5ktZJWizpaxnbezL9XJOuc1S6zmWSXpK0WtIjkkZliec966etmWmS/lPSSuBrkvaW9JiklZLekXSvpH7NG5H0hqRT0u8lkq6T9Fpa/5eS9sio29ySWpMezyWSJgEXAtemcfxvWnecpCfSuvMknZGxnbsk3SLpYUkbgRNbtkYk/Z2kZ9P1p0t6X8ay/ydpiaT1kl6RdHIO/x+tK4kIT546xAT0AVYCdwOnAf1bLL8E+EuLshOAg0n+4Hkf8DZwVrpsNBBAWUb9M4EFwDiSlvSXgelZ4mlt/UuABuCadP0ewD7AqUAlMIgkudyUsc4bwCnp988BTwPD0/o/Bu5Ll40C1gPnA+XAAODQdNldwDcztlmeHse/krS4TkrX3T+j/lrgmPRnU5W5DeD9wHJgAlAKXJzGWQnsDywGhmb8HPYu9r8PT8WZ3JKwDiMi1gHHkpyYbwNWSJoqac8drPNERLwQEU0R8TxwH3D8DnZzJfDtiHgpkstE/w4cuoPWRGuWRsSPIqIhktbOgoj4Y0TUR8QK4MYdxHAl8KWIqI2IeuBrwDnpfY8LgP+LpCW1JSJWRsSzWbZzJNAb+E5EbI6Ix0gu1Z2fUec3ETEt/dnUtVh/EvDjiJgZSavtbpLLekcCjSTJYryk8oh4IyJey+HnY12Ik4R1KOnJ+5KIGA4cBAwFbspWX9IESY9LWiFpLclJeOAOdjEK+GF6iWUNsAoQMCyHMBe3iGFPSVPSyzPrgJ/tIIZRwK8y9v8SyUl5T2AE0NaT8VBgcUQ0ZZQtanEci8luFPCF5jjSWEaQtB4WAP9EksCWp8c2tI1xWRfjJGEdVkS8THKJ5KDmolaq/RyYCoyIiL4k9y20g/qLgc9ERL+MqUdETG8thGyhtZj/97Ts4IjoA1yUEUNr+z+txf6rImJJumzvNu5zKTBCUubv8EhgSRvib47jWy3i6BkR9wFExM8j4liSZBLAd3ewLevCnCSsw5B0gKQvSBqezo8guXzydFrlbWC4pIqM1aqBVRFRJ+kIkks2zVYATcDYjLJbgeslHZjuo6+kj2cJqbX1W1MNbADWShoGfHEHdW8FvtV8eUvSIElnpsvuBU6RdK6kMkkDJB2aLnu7RRwzgXdJbmaXp+9hnA5M2UmszW4DrkxbYpLUK30IoFrS/pJOklQJ1AGb0p+DdUNOEtaRrCe5kTozfSLnaeBF4Avp8seAecBbkt5Jyz4LfEPSeuCrJI/OAhAR7wLfAqall1SOjIhfkfxVPCW9NPQiyU3y92ht/Sxxfx04jORG8e+Ah3ZwjD8kafk8msb8dHrMRMSbwEfS410FPAsckq53B8k9gjWSfh0Rm0mSwmnAO8D/AJ9KW187FRGzgSuA/wZWk9wEvyRdXAl8J93uW8Bg4Pq2bNe6HkV40CGzfJL0JnBRRDy508pmHYxbEmZ5JGkQyWOxbxQ5FLNd4iRhlidKuhN5FfhReinJrNPx5SYzM8vKLQkzM8uq03fwN3DgwBg9enSxwzAz6zTmzJnzTkQMakvdTp8kRo8ezezZs4sdhplZpyFpUVvr+nKTmZll5SRhZmZZOUmYmVlWnf6ehJlZPm3ZsoXa2lrq6lr2tt7xVVVVMXz4cMrLy3d5G04SZmY7UFtbS3V1NaNHj0bK1rlvxxMRrFy5ktraWsaMGbPL2ynY5SZJd0paLunFLMsl6b8kLZD0vKTDChWbmVk2dXV1DBgwoFMlCABJDBgwYLdbQIW8J3EXMHEHy08D9k2nScAtBYjJzGynOluCaNYecRcsSaQ9YK7aQZUzgZ9G4mmgn6Qh+Yrn7ulvMOO1lfnavJlZl9CRnm4axvbDLdaSZUhJSZMkzZY0e8WKFTnvqG5LI/fOXMTFP/krf5z/9q5Fa2bWDXSkJNFmETE5ImoiombQoDa9Wb6dqvJSfjHpKMbtVc2VP5vDr+bW5iFKM7POryMliSUkA7E3G8724/W2q/69Krj3iiOZMGYPPv+L57h7+hv52pWZ2W478cQT+eMf/wjAl7/8Za655pqC7LcjPQI7Fbha0hSS4RzXRsSyfO6wd2UZd17yAa65by43TJ3H2k1buOakfTrtTSozy6+v/+885i9d167bHD+0DzecfuDO9/31r/PVr36V5cuXM3fuXKZOndqucWRTsCQh6T7gBGCgpFrgBqAcICJuBR4mGd93AckA75cWIq6q8lJuufAwrn3weW78499Y8+4WvvzRcZSUOFGYWcdx3HHHERHceOONPPHEE5SWlrJx40Y++9nPUlFRwQknnMCFF17Y7vstWJKIiPN3sjyAfyhQONspKy3h++ccQp+qcu6c9jrr6rbwnbMPpqy0I12NM7Nia8tf/PnywgsvsGzZMgYMGEB1dTUADz30EOeccw6nn3465513Xl6ShM+CqZISccPp4/n8KfvxwJxaPnvvM9RtaSx2WGZmLFu2jAsvvJDf/OY39O7dmz/84Q9A8jb4iBHJrdzS0tK87NtJIoMkPnfKvtxw+ngenf82l989iw31DcUOy8y6sXfffZezzz6bH/zgB4wbN46vfOUrfP3rXwdg+PDh1NYmT2c2NTXlZf+dfozrmpqayMegQw/OqeXaB5/noGF9ueuSD9C/V0W778PMOr6XXnqJcePGFTuMVm3cuJGrr76aqqoqjj322FYvN7UWv6Q5EVHTln10pKebOpSPHT6c6qoyrr5vLudNnsE9l09gzz5VxQ7LzGyrXr168ZOf/CSv+/Dlph340IF7cdelH2DJ6k187JbpLFq5sdghmZkVlJPEThy990B+fsWRbKxv4JxbZ/DyW+37jLSZWUfmJNEGh4zoxy8/cxQlgvN+/DTPvLm62CGZmRWEk0Qb7btnNQ9ceTT9epZz4W0zeerV3DsWNDPrbJwkcjBij57cf+VRjBrQk8vumsXvX8hrryFmZkXnJJGjwdVV/GLSUbxveD/+4efP8MtZi3e+kplZJ+UksQv69iznnsuP4Jh9BnLtg89z25MLix2SmVleOEnsop4VZdx+cQ0fPXgI33r4Jb7/yCt09hcTzcxa8st0u6GyrJT/Ov/9VFeV8d+PL2Dtpi18/YwD3YOsmXUZbknsptIS8e2zD+Yzx43lnqcX8flfPsuWxvz0oWJm3ZcHHerEJHH9R8bRt2c53/vDK6yva+DmCw6jR0V+emU0syL5/XXw1gvtu829DobTvrPTasUadMgtiXb02RP24ZtnHcTjryzn4jv/yrq6LcUOycy6iMxBh6ZMmUJpaSkLFy7k8ssv55xzzsnbft2SaGcXHTmKPj3K+edfPMsFtz3N3ZcewYDelcUOy8zaQxv+4s+X1gYdGjt2LHfccUdek4RbEnlwxiFDue1TNbz69gY+/uMZLFmzqdghmVknlm3QoUJwksiTEw8YzD2XT2DFuno+fst0XluxodghmVkntKNBhwrBSSKPjhizB/dNOpL6hibOvXUGLy5ZW+yQzKyT6dmzJzNmzODUU08FknsTM2bMAGDlypVceeWVzJ07l29/+9t52b+TRJ4dNKwv9195FJVlJZw/+Wn++vqqYodkZl3EgAEDuPXWW3nttde4/vrr87IPJ4kCGDuoNw9cdTSD+lTyyTtm8tjLbxc7JDOzNnGSKJCh/Xpw/2eOYt89e3PZXbO54qezmbPIrQoz69icJApoQO9Kpkw6in88aR9mvbGKj90yg4/fOp3/m/82TU3u98nMOh4niQLrXVnGP39of6ZfdxI3nD6epWvq+PRPZ/Phm57k/tmL2dzgLj3MOprO2nlne8TtJFEkPSvKuPSYMTzxxRO46bxDKS0RX3zgeY773uNMfvI11vttbbMOoaqqipUrV3a6RBERrFy5kqqqqt3ajjrbgbdUU1MTs2fPLnYYuy0iePLVd7j1ideYsXAl1VVlXHTkKC49ZjSDq3fvf7KZ7botW7ZQW1tLXV1dsUPJWVVVFcOHD6e8vHy7cklzIqKmLdtwkuiAnlu8hslPLuT3Ly6jrKSEjx0+jCs+OJaxg3oXOzQz6wKcJLqIN97ZyG1PLeT+ObVsaWziw+P34jPHj+X9I/sXOzQz68ScJLqYFevruXv6G/x0xhusq2vgiDF7cNXxe3PC/oOQPMCRmeXGSaKL2lDfwC9mLeaOpxaydG0d++9ZzaTjxnLGoUMpL/UzCGbWNrkkiYKeWSRNlPSKpAWSrmtl+UhJj0uaK+l5SR8pZHwdXe/KMi4/dgx/vvZEbjz3EAC+cP9zHP+9x7n9qYVsqG8ocoRm1tUUrCUhqRT4G3AqUAvMAs6PiPkZdSYDcyPiFknjgYcjYvSOttudWhItRQRPvLKCW//8GjNfX0WfqjI+ddRoLj56NIOqPYaFmbUul5ZEIQcdOgJYEBELASRNAc4E5mfUCaBP+r0vsLSA8XU6kjjxgMGceMBg5r65mh//eSE3P7GAyU8t5OOHD+eKD45l9MBexQ7TzDqxQrYkzgEmRsSn0/lPAhMi4uqMOkOAR4H+QC/glIiY08q2JgGTAEaOHHn4okWLCnAEncPCFRu47amFPDhnCVuamjjtoL34zHF7c8iIfsUOzcw6iA57T6INzgfuiojhwEeAeyS9J8aImBwRNRFRM2jQoIIH2ZGNHdSbb5/9Pv5y3YlcdfzePPXqO5x58zTOn/w0T7yy3H1EmVlOCnm5aQkwImN+eFqW6XJgIkBEzJBUBQwElhckwi5kcHUV1048gM+euA/3zXyTO/7yOpf8ZBb9epZz2Mj+HD4qmQ4Z3o8eFaXFDtfMOqhCJolZwL6SxpAkh08AF7So8yZwMnCXpHFAFbCigDF2Ob0ry7jiuLFcfPRofv/iMqYvWMnsRat47OUk75aViAOH9eXwkf2pGd2fmlH9GdzH3YCYWaKg70mkj7TeBJQCd0bEtyR9A5gdEVPTJ5puA3qT3MS+NiIe3dE2u/PTTbtj9cbNPPPmamYvWs2cRat5bvEa6tMeaEfs0YPDR/bn8NF7UDOqP/vtWU1piV/aM+sq/DKd5WxzQxPzlq5lTpo0Zi9azYr19QBUV5Zx6Mh+1Izag8NH9efQkf3oXVnIRqiZtScnCdttEcHiVZuYvWjV1sTxytvriYASwbghfagZ1Z/DRvWnZvQeDOvXo9ghm1kbOUlYXqzdtIVnF69hzhurmL1oNc8uXsO7mxsBGNK3KkkYo/pTM2oPxg2ppsxdhZh1SB31ZTrr5Pr2KOf4/QZx/H7JY8cNjU28/NZ6ZqdJ45lFq/nd88sA6FFeyqEj+lEzOnmK6n3D+9G/Z7k7JDTrZNySsHa1dM2mrQlj9qJVzF+6juZXM3qUlzKkXxVD+/ZgSN8qhvTrwdD0c1i/Kob07UEv3+swyzu3JKxohvbrwRn9enDGIUMB2FjfwHOL1zB/2TqWra1j2dpNLFlTx9/+toIVG+pp+TdKn6oyhvbblkSGNX/v24Oh/arYq28VlWV+r8OsUJwkLK96VZZx9D4DOXqfge9ZtrmhibfX1W1NHkvXbPtcumYTzy5ew+p33zvW98DeFdsSSZo8mj+H9uvB4OoqP7Jr1k6cJKxoKspKGLFHT0bs0TNrnU2bG1m2dhPL1iaJI7M1snDFRqYtWPmeLtJLS8Se1ZUMSRNJv57l9K4sp7qqbOvUu7Kc3pVlGWXJfEWZb7abZXKSsA6tR0UpYwf13uH43uvqtrBsTR1L125i2ZrmJJJ8n7d0Hes2bWF9XQObG5t2ur/KspLtkkZzIuldVUafqvLt5quryqmubP7eXDep45aMdRVOEtbp9akqp89e5ey/V/UO69U3NLK+roENdQ2sr2tgff2Wrd831Dewvm4L6+sbMupsYUN9A2+uejepn863pY/EXhWl9Kgoo7KshIqykq2fFaXp59ay0q1llVnrlW79XlHaol5mefm2bZWViNLmSaLESct2kZOEdRuVZaVU9i5lYO9dH5ApInh3c+O2pJKRZDbUNbAuTSTr6xp4d3MjmxuaqG9IPjc3NiWfDU1JyyajrL6hic0NjclnY9N7bujvLglKlSSNspIkabRMJKWlyqhTsrXO1rrKqJ8xZdYpkZCgRKJEIERJSTL2SUlaLprn0zrN5Zl19N51Skq2bVuw3b4ktm63+XvzgSs9frH9suYyMtfNXJaxDulyIOv2mtdja61t89u2u622WqyzrW7GFrOs3zz/wX0H5v2xcicJsxxIoldlGb0qy9gzTx0hRgQNTbE1odSnn5sbG7d+31q2XaJp3LqsKd1GU9P2n43p1NAU2etE0NjYWp0mGpuC+oZGGgMam5poaNxWpzGtHwER0BTb5psiOa6kLPm+rQ4EkVGHrevZjr3+7fyP8OwkYdbBSKK8VJSXltCrm49C2zJpNG2XXIIAoolt39PPZN0k+aT/bZvf+n1bsmoug6TOdsszlrHdsm3b2279jO2w3bq8t26LOrFdndiuLDOGQnKSMLMOSxKlglLye0nFsvPzfmZmlpWThJmZZeUkYWZmWTlJmJlZVjknCUm9JLmHNTOzbmCnSUJSiaQLJP1O0nLgZWCZpPmS/kPSPvkP08zMiqEtLYnHgb2B64G9ImJERAwGjgWeBr4r6aI8xmhmZkXSlvckTomI9/TXHBGrgAeBByWVt3tkZmZWdDttSbSWIHaljpmZdT47bUlIWk/r74ELiIjo0+5RmZlZh7DTJBERO+5/2czMuiy/J2FmZlm1uYM/STXAl4BR6XrNl5vel6fYzMysyHLpBfZe4IvAC8DOx4E0M7NOL5cksSIipuYtEjMz63BySRI3SLod+BNQ31wYEQ+1e1RmZtYh5JIkLgUOAMrZdrkpACcJM7MuKpck8YGI2H9XdyRpIvBDoBS4PSK+00qdc4GvkSSf5yLigl3dn5mZ7b5cksR0SeMjYn6uO0l7jb0ZOBWoBWZJmpq5LUn7kvQPdUxErJY0ONf9mJlZ+8olSRwJPCvpdZJ7Erk8AnsEsCAiFgJImgKcCWQmnCuAmyNiNcmGl+cQm5mZ5UFbuuU4iqS314m7sZ9hwOKM+VpgQos6+6X7m0ZySeprEfGHLDFNAiYBjBw5cjfCMjOzHWnLG9efAuYA3wVOBOojYlHz1I6xlAH7AicA5wO3SerXWsWImBwRNRFRM2jQoHYMwczMMrWl76arACQdAJwG3CWpL8k4E38ApkVE4042swQYkTE/PC3LVAvMTHuUfV3S30iSxqy2HIiZmbW/NvfdFBEvR8R/RsRE4CTgL8DHgZltWH0WsK+kMZIqgE8ALV/M+zVJKwJJA0kuPy1sa3xmZtb+crlxvVVEbAIeTqe21G+QdDXwCMn9hjsjYp6kbwCz0ze5HwE+JGk+0Ah8MSJW7kp8ZmbWPhTR2lARGRW2H09CLRYXfTyJmpqamD17djFDMDPrVCTNiYiattT1eBJmZpZVTpebJB0CfDCdfTIinm//kMzMrKNo841rSZ8j6S58cDrdK+mafAVmZmbFl0tL4nJgQkRsBJD0XWAG8KN8BGZmZsWXy/ClInnqqFkj772RbWZmXUguLYmfADMl/SqdPwu4o/1DMjOzjqLNSSIibpT0Z+CYtOjSiJibn7DMzKwjyOnppoiYQ9KPk5mZdQNtThKSaoAvAaPS9XLpKtzMzDqhXFoS9wJfBF5g2/ClZmbWheWSJFakfSyZmVk3kUuSuEHS7cCfSEamAyAiHmr3qMzMrEPIJUlcChwAlLPtclMAThJmZl1ULkniAxGxf94iMTOzDieXN66nSxqft0jMzKzDyaUlcSTwnKSFJPck/AismVkXl0uSmNhK2Y5HLDIzs04tlyQxiO1fpmvmloSZWRfll+nMzCwrv0xnZmZZ+WU6MzPLyi/TmZlZVn6ZzszMsvLLdGZmllWuL9M9K+l1/DKdmVm3sLsv05mZWRe20yQhSZFYtLM67RuamZkVW1vuSTwu6RpJIzMLJVVIOknS3cDF+QnPzMyKqS2XmyYClwH3SRoDrAGqgFLgUeCmiJibvxDNzKxYdpokIqIO+B/gfySVAwOBTRGxJt/BmZlZceVy45qI2AIsy1MsZmbWweTynsRukzRR0iuSFki6bgf1PiYpJNUUMj4zM9tewZKEpFLgZuA0YDxwfmsv50mqBj4HzCxUbGZm1rqck4SkXukJP1dHAAsiYmFEbAamAGe2Uu/fgO8CdbuwDzMza0c7TRKSSiRdIOl3kpYDLwPLJM2X9B+S9mnjvoYBizPma9OyzH0dBoyIiN/tJKZJkmZLmr1ixYo27t7MzHLVpvckgL2B64G9ImJERAwGjgWeBr4r6aLdDURSCXAj8IWd1Y2IyRFRExE1gwYN2t1dm5lZFm15uumUiNgiaXREbB2RLiJWAQ8CD6aPxu7MEmBExvzwtKxZNXAQ8IQkgL2AqZLOiIjZbdi+mZm1s522JNLHXqGVcSMkHdmizo7MAvaVNEZSBfAJYOtIdxGxNiIGRsToiBhN0kpxgjAzK6K23JM4V9J3gGpJ49LLQs0mt3VHEdEAXA08ArwE/DIi5kn6hqQzcg3czMzyry2Xm6aRdMPxaZJ7BvtLWgMsBTblsrOIeBh4uEXZV7PUPSGXbZuZWftrS7ccS4CfSnotIqYBSBoAjCZ50snMzLqoXLoKn9ZcFhErgZUt6+QpRjMzKxJ3FW5mZlntalfhPUgSjLsKNzPrwtxVuJmZZZVzV+GSrgLKJD0LPBsRf8tPaGZmVmw5JQlIHlmVtCdwKPD3kvaJiCvaPzQzMyu2NicJST8E/il90ultkpfiHslbZGZmVnS5dBW+nqQvpV4Akj4sadpO1jEzs06szS2JiPiypAtIOuDbDGwAso4uZ2ZmnV8ul5tOBq4ANgJDgMsi4pV8BWZmZsWXy+WmLwFfSftUOgf4haST8hKVmZl1CLlcbjop4/sLkk4jGU/i6HwEZmZmxZfzGNfNImIZcHI7xmJmZh3MLicJgIjIqatwMzPrXHJ+mc7MrFuKSCZafm6tsK1e1rKMZW0p29n2e/SHZLjnvHGSMOssmpqgcTM01kNDOjVuTj8zy7ZAU0MyRWP6vSl72db5xmTaOt9a2c7qNJKcPJsypvRkmlmWlzotTuA0f2Q7uWeedHdSp6O6If9d6DlJmO2qpibYvB7q1kH9+nRal0x162DLuxkn7npoaHmCzyzbDA117y3LrN/UlqHk20lJGag0+SwpTaas82mZSrZ9Nk9o2/eSUlB58pdvtjpbl+2oTgmIVparlU9aL9tpnbZsR0kcZPwln1k3W9l2f/m3pWwn288zJwnrfiKSE3jzib1u3baTe33mST/jhL91PmP55vVt32dJOZRVQmlFxmcVlFVAaWVSVlkNZYMyypo/q1opS6dWyyoyTuaZJ/uSbSf1lmWZCaD5JG2Gk4R1FQ2bYf1SWLsE1i2BtYuT7+uXwaY17z3hR+POt1neKzlxV/VJPiv7QPWQdD6jbLs6fbfNl/fcdiIv2a1nRMyKxknCOr6mRtjwdpoAapPPtbXbvq9bAhuW857rxz36Q59hyWe/Uduf7Cur05N539ZP9hXVUOpfDzP/FlhxRcC7K9OT/pKMRJCRANYvS26QZqronSSAvsNgzwOh7/Bk6jMs/RwKFb2Kc0xmXYiThOXfumXw9ovbJ4K1i5Pv65YmN2wzlVZsO9mPOiZJBM3zzYmgqq+vm5sVgJOEta+mRnh7HiyemUxvzoS1b25brhKoHpqc+IccCgd8FPoMz0gEI6DXQCcAsw7CScJ2T/16qJ2VJIPFM6F29ranfqqHwIgJcNRnk4TQbwT03svX+s06Ef+2WttFJJeJ3pwJi59OPpfPS15kUgkMPhAOOQ9GHAkjjoB+I90iMOvknCQsu8Yt8NYL6WWjp5PP9cuSZRW9YXgNHHctjJwAw2qSJ4PMrEtxkrBtNq2GxbO23U9YMid56Qyg78jkJvLII5NLSIPH+7KRWTfg3/LuKgJWLdz+BvOKl5JlKoW9DobDPpUkhBETkhvLZtbtOEl0J8ueh4VPbEsMG1ck5ZV9k3sIB30svXR0uN8xMDOggElC0kTgh0ApcHtEfKfF8n8GPg00ACtIxtBeVKj4uqwIePVRmPZDWDQtKes/BvY5JUkMI46EQQe42wgza1VBkoSkUuBm4FSgFpglaWpEzM+oNheoiYh3JV0FfA84rxDxdUkNm+HFB2H6f8Hy+cm7CB/+dtJaqN6z2NGZWSdRqJbEEcCCiFgIIGkKcCawNUlExOMZ9Z8GLipQbF1L/Xp45qcw4+bkjebBB8LfT4aDzobS8mJHZ2adTKGSxDBgccZ8LTBhB/UvB36fbaGkScAkgJEjR7ZHfJ3fhuUw81aYdTvUrYXRH4TTf5hcVvK7Cma2izrcjWtJFwE1wPHZ6kTEZGAyQE1NTQcfOirPVr4G038Ez/48GaVs3OlwzOeSdxjMzHZToZLEEmBExvzwtGw7kk4BvgQcHxH1BYqtc1oyJ7kZPX9q0iHeoefD0f8IA/YudmRm1oUUKknMAvaVNIYkOXwCuCCzgqT3Az8GJkbE8gLF1blEwII/wbSb4I04KhV0AAAJvUlEQVSnkkdXj/08TLjSN6PNLC8KkiQiokHS1cAjJI/A3hkR8yR9A5gdEVOB/wB6A/cruYb+ZkScUYj4OrzGLfDiQ0nLYfm8pBfVD30LDr84GSDHzCxPCnZPIiIeBh5uUfbVjO+nFCqWTqN+A8y9J3lSae3i5H2Gs26Bg85JxjU2M8uzDnfj2oANK+Cvk5Opbg2MPBo+8n3Y90N+6c3MCspJoiNZtTBpNcz9GTTUJwPyHPO55M1oM7MicJLoCJbOTZ9U+g2UlMH7zkueVBq0X7EjM7NuzkmiWCLgtceS5PD6n6GyT5IYJlwJfYYUOzozM8BJovCaGmHer5LHWN96IRnO89RvwOGXQFXfYkdnZrYdJ4lC2lIHD14OL/8WBu4HZ/w3vO9cKKssdmRmZq1ykiiU+vUw5QJ4/Un48L/DhKv8pJKZdXhOEoXw7ir42cdg2XNw1q1JFxpmZp2Ak0S+rVsK9/w9rHodzrsneazVzKyTcJLIp1UL4adnJi2Jix6AMccVOyIzs5w4SeTLWy/Cz85O+l26eGoybrSZWSfjO6f5sPivcNdHQKVw6e+dIMys03KSaG+vPZZcYuqxB1z2Bxh8QLEjMjPbZU4S7Wn+b+Dec2GPsXDZI9B/VLEjMjPbLU4S7eWZe+D+S2DYYXDJbz0IkJl1CU4S7WH6j2Dq1TD2RPjkr6BH/2JHZGbWLvx00+6IgMf+DZ76AYw/C86+zYMBmVmX4iSxq5qa4OF/gdl3wGGfgr+7CUpKix2VmVm7cpLYFY1b4NdXwQv3J917n/oNSMblNjPrUpwkcrX53eQG9auPwMk3wAf/udgRmZnljZNELurWws8/AW/OgL/7T6i5rNgRmZnllZNEW218J+mob/l8OOcOOOhjxY7IzCzvnCTaYs3iJEGsrYXzp8C+pxY7IjOzgnCS2Jl3XoWfngX16+CTD8Goo4sdkZlZwThJ7MjSZ5PBgiB5i3rIIcWNx8yswPzGdTZvTIO7T4fyHkk/TE4QZtYNOUm05m+PJmNBVO+V9OQ6cJ9iR2RmVhROEi298ABMOR8G7Z+MBdF3eLEjMjMrGieJTLNuhwc/DSMmwMW/hV4Dix2RmVlR+cY1JB31PfWDpLO+/SbCx+9K7kWYmXVzThIR8MevJN19H3wunPU/UFpe7KjMzDqEgl5ukjRR0iuSFki6rpXllZJ+kS6fKWl0XgNqaoSp1yQJ4gNXwN//2AnCzCxDwZKEpFLgZuA0YDxwvqTxLapdDqyOiH2A/wS+m7eAGurhgUth7j1w3BfhI/8BJb5FY2aWqZBnxSOABRGxMCI2A1OAM1vUORO4O/3+AHCylIc+uBvq4b5PJGNSf+hbcNKX3dW3mVkrCnlPYhiwOGO+FpiQrU5ENEhaCwwA3smsJGkSMAlg5MiRuUdSWgF77A0Hng2HfTL39c3MuolOeeM6IiYDkwFqamoi5w1I8NHvt3dYZmZdTiEvNy0BRmTMD0/LWq0jqQzoC6wsSHRmZvYehUwSs4B9JY2RVAF8Apjaos5U4OL0+znAYxGRe0vBzMzaRcEuN6X3GK4GHgFKgTsjYp6kbwCzI2IqcAdwj6QFwCqSRGJmZkVS0HsSEfEw8HCLsq9mfK8DPl7ImMzMLDu/GGBmZlk5SZiZWVZOEmZmlpWThJmZZaXO/oSppBXAol1cfSAt3ubuBnzM3UN3O+budrywe8c8KiIGtaVip08Su0PS7IioKXYcheRj7h662zF3t+OFwh2zLzeZmVlWThJmZpZVd08Sk4sdQBH4mLuH7nbM3e14oUDH3K3vSZiZ2Y5195aEmZntgJOEmZll1S2ThKSJkl6RtEDSdcWOJ98kjZD0uKT5kuZJ+lyxYyoUSaWS5kr6bbFjKQRJ/SQ9IOllSS9JOqrYMeWbpM+n/65flHSfpKpix9TeJN0pabmkFzPK9pD0R0mvpp/987HvbpckJJUCNwOnAeOB8yWNL25UedcAfCEixgNHAv/QDY652eeAl4odRAH9EPhDRBwAHEIXP3ZJw4B/BGoi4iCSYQi64hADdwETW5RdB/wpIvYF/pTOt7tulySAI4AFEbEwIjYDU4AzixxTXkXEsoh4Jv2+nuTEMay4UeWfpOHAR4Hbix1LIUjqCxxHMi4LEbE5ItYUN6qCKAN6pKNZ9gSWFjmedhcRT5KMsZPpTODu9PvdwFn52Hd3TBLDgMUZ87V0gxNmM0mjgfcDM4sbSUHcBFwLNBU7kAIZA6wAfpJeYrtdUq9iB5VPEbEE+D7wJrAMWBsRjxY3qoLZMyKWpd/fAvbMx066Y5LotiT1Bh4E/iki1hU7nnyS9HfA8oiYU+xYCqgMOAy4JSLeD2wkT5cgOor0OvyZJAlyKNBL0kXFjarw0mGe8/I+Q3dMEkuAERnzw9OyLk1SOUmCuDciHip2PAVwDHCGpDdILimeJOlnxQ0p72qB2ohobiU+QJI0urJTgNcjYkVEbAEeAo4uckyF8rakIQDp5/J87KQ7JolZwL6SxkiqILnJNbXIMeWVJJFcp34pIm4sdjyFEBHXR8TwiBhN8v/4sYjo0n9hRsRbwGJJ+6dFJwPzixhSIbwJHCmpZ/rv/GS6+M36DFOBi9PvFwO/ycdOCjrGdUcQEQ2SrgYeIXkS4s6ImFfksPLtGOCTwAuSnk3L/jUdc9y6lmuAe9M/gBYClxY5nryKiJmSHgCeIXmKby5dsIsOSfcBJwADJdUCNwDfAX4p6XKS4RLOzcu+3S2HmZll0x0vN5mZWRs5SZiZWVZOEmZmlpWThJmZZeUkYWZmWTlJmJlZVk4SZmaWlZOE2S6QNFzSeVmW9ZD057Rb+taWV0h6Mu211KxDc5Iw2zUnk71fpMuAhyKisbWFaRf1fwJaTTJmHYmThFmOJB0L3AicI+lZSWNbVLmQtB8dSb0k/U7Sc+nIac2J4ddpPbMOzc1dsxxFxF8kzQL+JSJezFyW9pk0NiLeSIsmAksj4qPp8r5p+YvABwoUstkuc0vCbNfsD7zcSvlAIHM0uBeAUyV9V9IHI2ItQHoparOk6vyHarbrnCTMciRpIMkIaA2tLN4EVDXPRMTfSO5dvAB8U9JXM+pWAnX5jNVsd/lyk1nuRpNlHOWIWC2pVFJVRNRJGgqsioifSVoDfBpA0gDgnXSgHLMOyy0Js9y9TNKv/4uSWhsF7VHg2PT7wcBf03E8bgC+mZafCPwu75Ga7SaPJ2HWziQdBnw+Ij65gzoPAdell6PMOiy3JMzaWUQ8Azy+o5fpgF87QVhn4JaEmZll5ZaEmZll5SRhZmZZOUmYmVlWThJmZpaVk4SZmWXlJGFmZln9f5uPV4rZWeJHAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# load amici module (the usual starting point later for the analysis)\n",
    "sys.path.insert(0, os.path.abspath(model_output_dir))\n",
    "model_module = importlib.import_module(model_name)\n",
    "model = model_module.getModel()\n",
    "model.requireSensitivitiesForAllParameters()\n",
    "model.setTimepoints(amici.DoubleVector(np.linspace(0, 10, 11)))\n",
    "model.setParameterScale(amici.ParameterScaling_log10)\n",
    "model.setParameters(amici.DoubleVector([-0.3,-0.7]))\n",
    "solver = model.getSolver()\n",
    "solver.setSensitivityMethod(amici.SensitivityMethod_forward)\n",
    "solver.setSensitivityOrder(amici.SensitivityOrder_first)\n",
    "\n",
    "# how to run amici now:\n",
    "rdata = amici.runAmiciSimulation(model, solver, None)\n",
    "amici.plotting.plotStateTrajectories(rdata)\n",
    "edata = amici.ExpData(rdata, 0.2, 0.0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Optimize"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create objective function from amici model\n",
    "# pesto.AmiciObjective is derived from pesto.Objective, \n",
    "# the general pesto objective function class\n",
    "objective = pypesto.AmiciObjective(model, solver, [edata], 1)\n",
    "\n",
    "# create optimizer object which contains all information for doing the optimization\n",
    "optimizer = pypesto.ScipyOptimizer(method='ls_trf')\n",
    "\n",
    "#optimizer.solver = 'bfgs|meigo'\n",
    "# if select meigo -> also set default values in solver_options\n",
    "#optimizer.options = {'maxiter': 1000, 'disp': True} # = pesto.default_options_meigo()\n",
    "#optimizer.startpoints = []\n",
    "#optimizer.startpoint_method = 'lhs|uniform|something|function'\n",
    "#optimizer.n_starts = 100\n",
    "\n",
    "# see PestoOptions.m for more required options here\n",
    "# returns OptimizationResult, see parameters.MS for what to return\n",
    "# list of final optim results foreach multistart, times, hess, grad, \n",
    "# flags, meta information (which optimizer -> optimizer.get_repr())\n",
    "\n",
    "# create problem object containing all information on the problem to be solved\n",
    "problem = pypesto.Problem(objective=objective, \n",
    "                          lb=[-2,-2], ub=[2,2])\n",
    "\n",
    "# maybe lb, ub = inf\n",
    "# other constraints: kwargs, class pesto.Constraints\n",
    "# constraints on pams, states, esp. pesto.AmiciConstraints (e.g. pam1 + pam2<= const)\n",
    "# if optimizer cannot handle -> error\n",
    "# maybe also scaling / transformation of parameters encoded here\n",
    "\n",
    "# do the optimization\n",
    "result = pypesto.minimize(problem=problem, \n",
    "                          optimizer=optimizer, \n",
    "                          n_starts=10)\n",
    "# optimize is a function since it does not need an internal memory,\n",
    "# just takes input and returns output in the form of a Result object\n",
    "# 'result' parameter: e.g. some results from somewhere -> pick best start points"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Visualize"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.axes._subplots.AxesSubplot at 0x7fb9a6218748>"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAEWCAYAAABliCz2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4wLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvqOYd8AAAIABJREFUeJzt3Xl8VNX9//HXB0iAhJ0ERCAGZRMRUCPuCyCKS2211gXbr0pbqrZ1q7bVb3+tbb+ltVqpVqu1StWvIu7Vr3UBcUOtyiICgsgqYU/YIZD18/vj3tghJjBZZm6SeT8fj3nkzpl77/lMCJ85c+6555i7IyIizV+LqAMQEZHkUMIXEUkRSvgiIilCCV9EJEUo4YuIpAglfBGRFKGEL1IDM/sfMys0s/Vx7Puwmf1PuH2qma2uY52Xm9m7dTlWZH+U8CVyZnazmb1SpWxJDWUXx3G+W83ssXrGlAP8BBjk7gfU51yJYmZvmdn3oo5Dmg4lfGkM3gGON7OWAGbWA0gDjqhS1jfcN6HMrBWQA2xy942Jrk8kWZTwpTGYSZDgh4XPTwLeBBZXKVvm7msBzOwuM8s3s+1mNtvMTgrLxwC3ABeZ2U4z+yQs72hmD5nZOjNbE3bXVH6YXG5m75nZRDPbBLwFTAMODM/xcLjf02a23sy2mdk7ZnZYXd6smbmZXWNmy8Muo9vNrNr/i2Z2vJnNDOucaWbHh+W/C38n94Qx3lOXWCS1KOFL5Ny9BPgQODksOhmYAbxbpSy2dT+T4MOgCzAZeNrM2rj7q8AE4El3b+fuQ8P9HwbKCL4lHAGcDsR2hxwDLAe6A6OBM4G14TkuD/d5BegHdAPmAI/X422fB+QBRwJfB8ZV3cHMugD/Au4GugJ3Av8ys67u/t8Ev6MfhTH+qB6xSIpodAnfzCaZ2UYzW9BA53vVzLaa2UtVyvuY2YdmttTMnjSz9FqcM8fMpprZIjNbaGa51exzQ/jaPDObbmYHxbx2m5ktCB8XxZSPNLM5YfkjYddC5UXAbWY2N3z8MixvY2YfmdknZvapmf065lwPheXzzOwZM2sXR1zlMXW8WM17utvMdsb7e6qlt/lPcj+JIJnNqFL2duXO7v6Yu29y9zJ3/xPQGhhQ3YnNrDtwFnCdu+8Ku2kmArHXA9a6+1/C8+2u7jzuPsndd7h7MXArMNTMOtbx/d7m7pvdfRXwZ+CSavY5G1ji7v8bxvUE8BnwtTrWKSmu0SV8gpbYmAY83+3Ad6opvw2Y6O59gS3Ad6vuEF4Uy63m2EeB2939UGA4UF0/78dAnrsPAZ4B/hie82yCVt0wglbljWbWIfxK/whwsbsPBr4ALos53wx3HxY+fhOWFQMjw1bsMGCMmR0bvna9uw8N618FVLYAq40rtDumjnOr/C7ygM7VvM+G8g5wYtiqzXb3JcD7BH37XYDBxLTwzezG8AN3m5ltBToCWTWc+yCCLqN14Yf/VuBvBC31Svn7Cs7MWprZH8xsmZltB1aGL9VU5/7E1vcFcGA1+xwYvkaVfXvWsU5JcY0u4bv7O8Dm2DIzOyRsqc82sxlmNrAW55sO7KhyPgNGEiQ8CBLtN+I5n5kNAlq5+7Tw/Dvdvaiaet+MKf8A6BVuDwLeCVtsu4B5BB9wXYESd/883G8a8M39vDd398oWd1r48PC17THvtW1MeU1x7es9tyT44Pzp/vath38TJO3vA+/Bl+9hbVi21t1XhPGcFMZyIdDZ3TsB2wALz1V1Cth8gg/HLHfvFD46uHtsH/z+po0dS9D1cloYZ25YbjUdsB+9Y7ZzCN5nVWsJPqyosu+acFtT3UqtNLqEX4MHgB+7+1HAjcBf63m+rsBWdy8Ln68m/lZTf2CrmT1nZh+HF9xa7ueY7xL0/wJ8QtASzzCzLGAEwX/+QqBV2JIGuIC9k8JxYRfNK7EXC8OW51yCbxnT3P3DmNf+AawHBgJ/2U9cAG3MbJaZfWBmsR+APwJedPd1+3mfdRZ2o8wCbiDoyqn0blgW23/fnqA/voDgd/ZLoEPM6xuA3MoLoWHcU4E/VX6bChsRp9QixPYEHxqbgAyC6wT1cZOZdTaz3sC1wJPV7PMy0N/MxppZq7D7bxBQ2T25ATi4nnFICmn0CT/sez6e4KLcXIKv4j3C186P6QuPfbxWj/quqOzHJrio9nL4/Plwl1YE/ck3AkcT/Ie7fB/n+3Z4ntsB3H0qwX/k94EnCFq25R4sTHAxMNHMPiL4VlIenmYOcFDYdfMX4J+V53f3cncfRtBSH25mg2Neu4KgW2AR8OW1guriCh3k7nkErdk/h0nxQOBbVP+B0dDeJuhmib3xaEZYFpvwXwNeBT4n6OLYw95dJE+HPzeZ2Zxw+7+AdGAhQRfeM4R/R3F6NKxrTXiOD2pxbHVeAGYDcwkuzD5UdQd33wScQ3A/wCaCbzXnuHthuMtdwAVmtsXM7q5nPJIK3L3RPQi+Li8ItzsA6+p5vlOBl2KeG2GLOnx+HPBaNce9BeRWKTsWeDvm+XeAe2uo9zSCZNttH7FNBs6qpvx04KkajllJ0D1RtfyXwI3VlJ9c5f3HE9fDBN8yzib4lrAyfFQAS6P+G2nKD4KumL5Rx6FH6j0afQvfg37cFWb2LQj6pM1s6H4O2985nWCc9wVh0WUELa54zAQ6mVl2+HwkQYtvL2Z2BMG3kXM95uadsAuma7g9BBhC0N2AmXULf7YGfgbcHz4/IOyLx8yGE3wz22Rm2WbWKSxvSzCc8LPwd9Q3LDfgXILRHfuKq3NYL2FX0wnAQnf/l7sf4O657p4LFHlwoVtEmphWUQdQlZk9QdAiz7JgPpJfAZcC95nZLwguTE4h6AuP53wzCPqw24Xn+667v0aQUKdYMP/Jx1Tzlbo67l5uZjcC08NkOhv4e1jXb4BZ7v4iQVdJO4KuKIBVHox8SQNmhGXbgW/7f64l3GRm5xAk9Pvc/Y2w/ALgKjMrA3YTjORxC+4+fSS8htCC4BvBS5UjfsysA8G3mU+Aq8Jz1RTXocDfzKwiPNcf3P0rH2Qi0nRZ0NgVEZHmrtF36YiISMNoVF06WVlZnpubG3UYIiJNxuzZswvdPXv/ezayhJ+bm8usWbOiDkNEpMkws6p3Y9dIXToiIilCCV9EJEUo4YuIpAglfBGRFKGELyKSIpTwRUQiMn3yDC7NvYrTW17IpblXMX3yjP0fVA+NalimiEiqmD55BhPH309xUQkAG1cVMnH8/QCMGntSQupUC19EJAKTbpn8ZbKvVFxUwqRbJiesTiV8EZEIFORvqlV5Q1DCFxGJQHbvrrUqbwhK+CIiEfiv316Mp+19GbV1RjrjJoxNWJ1K+CIiEdjYL5f1ZxxHhx6dMTO65WRx/QNXJuyCLWiUjohI0u0uKecv05dw6Fl5PPXCjwkXI0o4tfBFRJLs4fdXsnFHMT8dMzBpyR6U8EVEkmrb7lLuf3sZIwZkc3Rul6TWrYQvIpJED7yzjG27S7nxjAFJr1sJX0QkSTbu2MOkd1fytaEHctiBHZNevxK+iEiS3PvGUkrKK7hhdP9I6lfCFxFJgvzNRUz+aBUX5vWmT1ZmJDEo4YuIJMHE1z+nhRnXjuoXWQxK+CIiCbZ4/Q6e/3gNlx2fywEd20QWhxK+iEiC3TF1Me3SW3HVKYdEGkdCE76ZXW9mn5rZAjN7wsyi+2gTEYnAnFVbmLZwA98/+WA6Z6ZHGkvCEr6Z9QSuAfLcfTDQErg4UfWJiDQ27s7try6ma2Y6407sE3U4Ce/SaQW0NbNWQAawNsH1iYg0Gu8uLeTfyzfxwxF9adc6+qnLEpbw3X0NcAewClgHbHP3qVX3M7PxZjbLzGYVFBQkKhwRkaRyd25/bTE9O7Xl0mNzog4HSGyXTmfg60Af4EAg08y+XXU/d3/A3fPcPS87OztR4YiIJNWrC9Yzb/U2rj2tH61btYw6HCCxXTqnASvcvcDdS4HngOMTWJ+ISKNQVl7BHVMXc0h2Jucf0TPqcL6UyIS/CjjWzDIsmP9zFLAogfWJiDQKz328hmUFu7jx9AG0atl4Rr8nsg//Q+AZYA4wP6zrgUTVJyLSGBSXlXPX60sY0qsjYwYfEHU4e0noZWN3/xXwq0TWISLSmDz+wSrWbN3NH755eFIXN4lH4/muISLSxO0sLuPeN5dy3MFdObFvVtThfIUSvohIA5n07go27SrhpjEDGl3rHpTwRUQaxJZdJfz9neWMHtSdI3M6Rx1OtZTwRUQawH1vL2NnSRk3np78pQvjpYQvIlJP67ft4ZH3V3LesJ4MOKB91OHUSAlfRKSe7pq+hAp3ro9o6cJ4KeGLiNTDisJdPDUrn0uG59C7S0bU4eyTEr6ISD3cOe1z0lu24Ecj+0Ydyn4p4YuI1NGna7fxf5+s5YoTcunWvvGv76SELyJSR3e8tpgObVrxg5OjXbowXkr4IiJ1MHPlZt5cXMCVpx5Cx4y0qMOJixK+iEgtuTt/fPUzstu35orjo1+6MF5K+CIitfTW4gJmrtzCNSP70ja9cSxuEg8lfBGRWqiocP742mJ6d2nLRUc3jqUL46WELyJSCy/NX8eiddu5YXR/0ls1rRTatKIVEYlQaXkFd05dzIDu7Tl3aONZujBeSvgiInF6etZqVm4q4sYzBtCyReOb/nh/lPBFROKwp7Scu6Z/zpE5nTjt0G5Rh1MnSvgiInF49N8r2bC9mJvOGNgoFzeJhxK+iMh+bN9Tyl/fWsZJ/bI47pCuUYdTZ0r4IiL78eA7y9laVMpPzxgYdSj1ooQvIrIPhTuLefDdFZx1+AEc3qtj1OHUixK+iMg+3PvmUvaUlnPD6Ma7dGG8lPBFRGqweksRj3+wiguO6kXfbu2iDqfelPBFRGpw1+tLALj2tMa9dGG8lPBFRKqxdOMOnp2zmm8fexA9O7WNOpwGkbCEb2YDzGxuzGO7mV2XqPpERBrSn6Z+Ttu0lvxwRNNY3CQerRJ1YndfDAwDMLOWwBrg+UTVJyLSUOat3sorC9Zzzah+dG3XOupwGkyyunRGAcvc/Ysk1SciUme3v7aYzhlpfP+kprO4STySlfAvBp5IUl0iInX2/rJCZiwp5OpT+9K+TdNYujBeCU/4ZpYOnAs8XcPr481slpnNKigoSHQ4IiI1CpYuXMwBHdrwneMOijqcBpeMFv6ZwBx331Ddi+7+gLvnuXtednZ2EsIREanetIUbmJu/lWtP60ebtKazdGG8kpHwL0HdOSLSyJVXOHdMXUyfrEy+dVSvqMNJiIQmfDPLBEYDzyWyHhGR+nph7ho+37CTG0b3p1XL5nmLUsKGZQK4+y6g6c4lKiIpoaSsgomvf86gHh04+/AeUYeTMM3zY0xEpBamzFxF/ubd3DRmAC2a4NKF8VLCF5GUVlRSxt3TlzI8twun9m/eA0eU8EUkpf3jvZUU7izmp2MGNNmlC+OV0D58EZHGaPrkGUy6ZTIF+YWUdWjHMReeSl5ul6jDSji18EUkpUyfPIOJ4+9n46pC3KHltp3seGwq0yfPiDq0hFPCF5GUMumWyRQXlexVVrq7hEm3TI4oouSJO+GbWUYiAxERSYaC/E21Km9O9pvwzex4M1sIfBY+H2pmf014ZCIiDay4rJz0rA7Vvpbdu/nfMhRPC38icAawCcDdPwFOTmRQIiINbfH6HXz9nvdYecwQWrTeexbM1hnpjJswNqLIkieuLh13z69SVJ6AWEREGlxFhTPp3RV87Z53KdxZzF/+eDE/fegquuVkYWZ0y8ni+geuZNTYk6IONeHiGZaZb2bHA25macC1wKLEhiUiUn8btu/hxqc/YcaSQk47tBt/+OYQstq1hoHdUiLBVxVPwr8SuAvoSbBM4VTgh4kMSkSkvl6Zv46bn5/PntJyfnfeYMYOz2n2N1btz34TvrsXApcmIRYRkXrbWVzGr1/8lKdnr2ZIr45MvGgYh2S3izqsRmG/Cd/M/gF41XJ3H5eQiERE6mj2F1u4/sm5rN5SxI9G9OXa0/qR1kynOq6LeLp0XorZbgOcB6xNTDgiIrVXWl7BX95Yyj1vLOHATm158gfHcXQKTJVQW/F06Twb+9zMngDeTVhEIiK1sLJwF9c9OZe5+Vs5/8ie3HruYXRoZouPN5S6TJ7WD+jW0IGIiNSGu/PkzHx+89JC0lq24J6xR3DOkAOjDqtRi6cPfwdBH76FP9cDP0twXCIiNdq8q4SfPzuPqQs3cPwhXfnThUPp0bFt1GE1evF06bRPRiAiIvF4a/FGbnpmHtuKSvnvsw7luyf2adarVDWkGhO+mR25rwPdfU7DhyMiUr09peX84ZXPePj9lfTv3o5HrhjOoAOrnxdHqrevFv6f9vGaAyMbOBYRkWp9unYb102Zy5KNO7nihFx+NmYgbdJaRh1Wk1Njwnf3EckMRESkqooK5+8zlnPH1MV0zkjn0XHDObmZrzubSHGN0jGzwcAggnH4ALj7o4kKSkRk7dbd3PDUXD5YvpkzDuvO788fQpfM9KjDatLiGaXzK+BUgoT/MnAmwTh8JXwRSYgXP1nLL56fT1mF88dvDuFbeb1Sfh6chhBPC/8CYCjwsbtfYWbdgccSG5aIpKLte0r51Quf8vzHazgipxN/vmgYB3XNjDqsZiOehL/b3SvMrMzMOgAbgd4JjktEUsxHKzZz/ZNzWb99D9ed1o8fjehLK82D06DiSfizzKwT8HdgNrAT+HdCoxKRlFFSVsGfX/+c+95eRk6XDJ6+8jiOzOkcdVjNUjw3Xl0dbt5vZq8CHdx9XjwnDz8oHgQGEwzlHOfu+rAQSWHTJ89g0i2TKcjfROeeXdg+Io9FPXpwUV5v/t/XBtGudV1mfJF4xLOI+YtmNtbMMt19ZbzJPnQX8Kq7DyS4DqCVskRS2PTJM5g4/n42rirE3dm8ehOlU6ZzbecKbrtgiJJ9gsXTQfYn4ERgoZk9Y2YXmFmb/R1kZh0JFjt/CMDdS9x9a72iFZEmbdItkykuKtmrzErLmHXfyxFFlFr2m/Dd/e2wW+dg4G/AhQQXbvenD1AA/MPMPjazB83sK5fbzWy8mc0ys1kFBQW1DF9Emor8zUVszC+s9rWC/E1JjiY1xXUJ3MzaAt8kWN/2aOCROA5rBRwJ3OfuRwC7gJ9X3cndH3D3PHfPy87WHXQizc2arbu5+bn5jLjjLcraVz/EMrt31yRHlZriufHqKWA48CpwD/C2u1fEce7VwGp3/zB8/gzVJHwRaZ7Wb9vDvW8uZcrMVRjG2GNyGHrwd/jHdZP26tZpnZHOuAljI4w0dcRzheQh4BJ3L6/Nid19vZnlm9kAd18MjAIW1iVIEWk6Nu7Yw31vLePxD1dRUeFceHRvfjiiLz07tQUG0ykj/ctROtm9uzJuwlhGjT0p6rBTgrl/ZX3yhju52TCCYZnpwHLgCnffUtP+eXl5PmvWrITFIyKJU7izmL+9vYz//eALSsudC47sxY9G9qV3l4yoQ2vWzGy2u+fFs29Cx0C5+1wgrkBEpGnasquEv72znEfeX0lxWTnfOKIn14zsR26WpkRobDToVUTqZFtRKQ++u5xJ766gqLScc4ceyDWj+nFIdruoQ5MaxDs9ck/goNj93f2dRAUlIo3X9j2lTHp3BQ/NWMGO4jLOPrwH153Wj37dtRpqYxfPKJ3bgIsILrhWXrh1QAlfJIXsLC7jkfdX8sA7y9m2u5QzDuvOdaf159AeWmawqYinhf8NYIC7Fyc6GBFpfIpKynj031/wt7eXsaWolNMO7cZ1p/VncM+OUYcmtRRPwl8OpAFK+CIpZHdJOY9/+AX3v72Mwp0lnNI/m+tH92dY705RhyZ1FE/CLwLmmtl0YpK+u1+TsKhEJDJ7Sst54qNV/PWtZRTsKObEvllcP7ofRx3UJerQpJ7iSfgvhg8RacaKy8p5atZq7n1jKeu37+GYPl2455IjOOZgTXvQXMQzH/4jZpYO9A+LFrt7aWLDEpFkKS2v4JnZq7nnjaWs2bqbow7qzJ0XDuW4Q7pqHdlmJp5ROqcSTJa2EjCgt5ldpmGZIk1L7MIj2b27ctn/XMy2gQdz9xtLyN+8m2G9O/H78w/npH5ZSvTN1H6nVjCz2cDYcD4czKw/8IS7H9XQwWhqBZHEqFx4JHbSMk9rxfozjqPP6CO5YXR/Th2QrUTfBDX01ApplckewN0/N7O0OkcnIklX08IjA2Yv4JkXfqxEnyLiXcT8QeCx8PmlgJrhIk1IQQ0Lj+xYv1XJPoXEswDKVQR32V4TPhaGZSLSBCzduBPvVP20B1p4JLXEM0qnGLgzfIhIE/LC3DXc/Nx8Oo7Io/PL71G6RwuPpLIaW/jhSleY2Xwzm1f1kbwQRaS2isvK+cU/53PtlLkcdmAHnn9oHD958Eq65QQjcLrlZHH9A1dq4ZEUU+MoHTPr4e7rzOyg6l539y8aOhiN0hGpv/zNRVz9+Bzmr9nGD04+mBvPGEBay7iWr5YmqEFG6bj7unDzanf/WZUKbgN+9tWjRCRK0xZu4CdPzQXg7/+Vx+hB3SOOSBqTeD72R1dTdmZDByIidVdaXsGElxfx/UdncVDXTP51zUlK9vIVNbbwzewq4GrgkCp99u2B9xMdmIjEZ/22Pfz4iTnMXLmF7xx7EL8451Bat2oZdVjSCO1rlM5k4BXg98DPY8p3uPvmhEYlInGZsaSA66bMZXdpOXddPIyvD+sZdUjSiO2rD38bsM3M7gI2u/sOADPrYGbHuPuHyQpSRPZWXuH85Y0l3DV9Cf26teOvlx5F325aS1b2LZ47be8Djox5vrOaMhFJksKdxVz/5FxmLCnk/CN78j/fGExGelzLU0uKi+evxDxm7Ka7V5iZ/rpEIjBz5WZ+NHkOW4tKue2bh3NhXm9NjSBxi2eUznIzu8bM0sLHtQTLHopIkrg7D7yzjIsf+IC2aS15/uoTuOjoHCV7qZV4WupXAncDvwAcmA6MT2RQIvIf24pK+cnTn/D6og2cOfgAbrtgCB3aaMJaqb145tLZCFychFhEpIr5q7dx9eTZrNu6h1+eM4grTshVq17qLJ4Vr7KB7wO5sfu7+7g4jl0J7ADKgbJ4b/8VSXXuzmMfruK3/7eQrHbpPHXlcRyZ0znqsKSJi6dL5wVgBvA6QeKurRHuXv1k3CLyFbuKy7j5ufm8+MlaTh2QzcQLh9E5Mz3qsKQZiCfhZ1SdS0dEEmPx+h1c/fhsVhTu4qYzBnDVKYfQooW6cKRhxDNK5yUzO6uO53dgqpnNNrNqL/Sa2Xgzm2VmswoKCupYjUjT9+zs1Xz93nfZtruMx753DD8c0VfJXhpUPIuY7wAygWKgFDDA3b3Dfk9u1tPd15hZN2Aa8GN3f6em/TU9sqSiPaXl3Prip0yZmc8xfbrwl0uOoFuHNlGHJU1Egy5i7u7Vr40WB3dfE/7caGbPA8OBGhO+SKpZWbiLqx+fw8J12/nhiEO4/rT+tNLc9ZIg8YzSObm68n211MPjMoEW7r4j3D4d+E2dohRphl6Zv46bnplHq5bGPy4/mhEDu0UdkjRz8Vy0vSlmuw1BK302MHI/x3UHng/HDLcCJrv7q3UJUqQ5KSmr4PevLOIf761kWO9O3DP2CHp1zog6LEkB8XTpfC32uZn1Bv4cx3HLgaF1D02keZg+eQaTbplMQf4muvTsws5RR7Og+wFcfnwut5x1KOmt1IUjyVGXv7TVwKENHYhIczR98gwmjr+fjasKcXc2rd7E7snT+HGncm499zAle0mqePrw/0IwvBKCD4hhwJxEBiXSXEy6ZTLFRSV7lbUoLWfO/a/Az8+NKCpJVfH04ceOkywDnnD39xIUj0izsjG/+pvMC/I3JTkSkX2vaZvj7qvc/ZFkBiTSHGzbXcofXvmM0vaZpG3f9ZXXs3t3jSAqSXX76kD8Z+WGmT2bhFhEmoVXF6xj9J1v8+TMVQz5/hm0zth7HpzWGemMmzA2ougkle2rSyf2nu6DEx2ISFO3YfsefvnCAl77dAODenTgwcvyGNKrE9OP6PXlKJ3s3l0ZN2Eso8aeFHW4koL2lfC9hm0RiVFR4Uz+aBW3vfIZJeUV/GzMQL53Uh/SwjtmR409SQleGoV9JfyhZradoKXfNtyGWsylI9LcLd24k5ufm8fMlVs4/pCuTDjvcHKzMqMOS6RaNSZ8d2+ZzEBEmpKSsgruf3sZ97yxlLbpLfnjBUP41lG9tBqVNGrxDMsUkRizv9jCzc/N4/MNOzlnSA9+9bXDyG7fOuqwRPZLCV8kTjuLy7j91c949IMvOKBDGx66LI9Rh3aPOiyRuCnhi8Rh+qIN/OKfC1i/fQ+XHZfLjWcMoF1r/feRpkV/sSL7ULCjmFv/71P+NW8d/bu3495Lj9di4tJkKeGLVMPdeXrWan738iJ2l5Tzk9H9+cEph2iyM2nSlPBFqlhZuItbnp/P+8s2MTy3CxPOP5y+3dpFHZZIvSnhi4RKyyv4+4zl3PX6EtJbtuB35w3mkqNztJC4NBtK+CLAvNVb+dmz81m0bjtjDjuAX3/9MLprIXFpZpTwJaUVlZRx59TPmfTeCrLateb+bx/FmMEHRB2WSEIo4UvKevvzAv77+fms3rKbscfk8LMxA+nYNi3qsEQSRglfUs7mXSX89qWFPP/xGg7OzuSpHxzH8D5dog5LJOGU8CVluDv/nLuG3/zfQnYWl3HNyL5cPaIvbdI0bZSkBiV8abamT57x5Tz0XXp2oWzMsczqks0ROZ34w/lDGHBA+6hDFEkq3UUizdL0yTOYOP5+Nq4qxN3ZtHoTWx55lcvbFPPMlccr2UtKUsKXZunBmx+nuKhkr7IWpeUseXgaLTWuXlKUunSk2Sgrr+CdJQU88VE+BfmbqC6tF+RvSnpcIo2FEr40efmbi3hqVj5Pz1rN+u17yGqXzgHZHSku2PaVfbN7d40gQpHGIeEJ38xaArOANe5+TqLrk9RQXFbOtIUbeHJmPu8uLQTglP7Z3HruIEYO7M6Mg9swcfz9e3XrtM5IZ9yEsVGFLBK5ZLTwrwUWAVoDV+pt6cYdTPkon+c+XsPmXSX07NSWa0f141t5venZqe2X+1UuGl4HHDbTAAAQnElEQVQ5Sie7d1fGTRirxcQlpSU04ZtZL+Bs4HfADYmsS5qvopIy/jVvHU/OzGfWF1to1cIYPag7Fw/P4cS+WTVehB019iQleJEYiW7h/xn4KVDjGDgzGw+MB8jJyUlwONJUuDsL1mzniZmreHHuWnYWl3Fwdia3nDWQ84/sRVY7rSErUlsJS/hmdg6w0d1nm9mpNe3n7g8ADwDk5eV5ouKRpmHb7lJemLuGKR/ls3Dddlq3asHZQ3pw8dE5HJ3bGTMNqRSpq0S28E8AzjWzs4A2QAcze8zdv53AOqUJcnc+WrGZJ2fm86/56yguq2BQjw789uuHce6wnprQTKSBJCzhu/vNwM0AYQv/RiV7iVWwo5jn5qzmyZn5LC/cRfvWrbjgqF5cMjyHwT07Rh2eSLOjcfiSVOUVzowlBTw5M59pCzdQVuEcnduZq0f05azDDyAjXX+SIomSlP9d7v4W8FYy6pLoxU5aVjkccuBZeTwd3hy1ZutuumSmc8UJuVx0dG/6dtO8NiLJYO6N5zppXl6ez5o1K+owpB4qJy2LveHJ0lux7vTj2HnYIZzYN4uLj85h9KDupLfSVE4i9WVms909L5599f1Zaq20vIItRSVs2VXK5l0lbCkqCX7uKuG16x/5yqRlXlJG7kfz+OvD36N3l4yIohYRJfxmpLqulP3deFRe4WzfXcrmoiBhbwoTd+XzzbtK/5PQw5879pTVeL6+BduqnbSsuGC7kr1IxJTwm4mqXSkbVxVyx/fu46Plm8g6eUhMAo9pkReVsrWohIoaevVat2pB18x0Omem0yUznZwuGXTJTKdzRjpdMtOC8oz/vN4pI41xU15i46rCr5xLk5aJRE8Jv5mYdMvkr3SllO0p5bU7/snKnWmktbQwUQcJe+ABHeicmbZXwv7y9TCRt02v/dJ/4yaM1aRlIo2UEn4zsTH/q61qgLQdu5h/6+m0a90qKXepatIykcZLCb+JK69w7n1zKaXtM0nbvusrr3frnUX7Nsm9U1WTlok0ThoX14St2bqbSx74gDunfc5B3x5F67bpe72urhQRiaUWfhP10ry13PzcfNxh4kVDOe+IXkw/PlddKSJSI9141cTsKi7j1hc/5enZqxnWuxN3X3wEOV013FEkVenGq2bqk/ytXDvlY77YXMSPR/blmlH9SGupXjkRiY8SfhNQXuH87Z1l3Dn1c7q1b82U7x/LMQdrXLuI1I4SfiO3bttubnjyE/69fBNnH96DCecdTscMzQ8vIrWnhN+IvbpgHT97dj6l5RX88YIhfOuoXlrxSUTqTAm/ESoqKeO3Ly3iiY9WMaRXR+66+Aj6ZGVGHZaINHFK+I3MgjXbuGbKx6wo3MWVpxzCDaP7axphEWkQSviNREWF89C7K/jja5/RJTOdx797DMf3zYo6LBFpRpTwG4GN2/fwk6c/YcaSQk4f1J3bvjmEzpnp+z9QRKQWlPAj9vrCDfz02XkUlZQx4bzDuWR4b12YFZGEUMKPyJ7Scia8vIhH//0Fg3p04O5LhmltVxFJKCX8CCxat51rp3zM5xt28r0T+3DTmAG0blX7uedFRGpDCT+J3J2H31/J71/5jA5t0nh03HBO7p8ddVgikiKU8JOkcGcxNz39CW8uLmDkwG788YIhZLVrHXVYIpJClPCT4K3FG7nx6Xls31PKb75+GN859iBdmBWRpFPCT6DisnJue2Uxk95bwYDu7Xn8e8cw4ABdmBWRaCjhJ8iSDTu4ZspcFq3bzuXH5/LzMwfSJk0XZkUkOkr4DWT65BlfrjaV0a0jXxwzlBZHDWDS5XmMHNg96vBERBK3pq2ZtTGzj8zsEzP71Mx+nai6ojZ98gwmjr+fjasKcXd2bdhK1qvv8YsDWyjZi0ijkchZuYqBke4+FBgGjDGzYxNYX9KUVzj5m4t45/MCHnl/JROvf4TiopK9dyop45nfPhNNgCIi1UhYl44Hi+XuDJ+mhY/Gs4Dufrg7m3eVsKJwF8sLd7GicBcrCsKfm3ZRUlbx5b59C7ZR3ZibgvxNyQtYRGQ/EtqHb2YtgdlAX+Bed/+wmn3GA+MBcnJyEhlOtYpKyoIkHpPQKxP8tt2lX+6X1tLI6ZJBn6x2nDIgmz5ZmfTJyuTg7Eyun/ISG1cVfuXc2b21DKGINB4JTfjuXg4MM7NOwPNmNtjdF1TZ5wHgAYC8vLxafwOIvVia3bsr4yaMZdTYk/bap6y8gtVbdrO8cCfLK1vp4WPdtj177dujYxsOzs7ka0N70CerHQeHib1X57a0qmHB8HETxjJx/P17deu0zkhn3ISxtX07IiIJk5RROu6+1czeBMYAC/a3f7wqL5ZWJtqNqwr50/fu5/2lhbQ4asCXrfVVm4ooq/jPZ0nHtmkcnJ3JcYd0DRN6O/pkZZKblUFGeu1/JZUfMPv74BERiZIFXe0JOLFZNlAaJvu2wFTgNnd/qaZj8vLyfNasWXHXcWnuVdV2pZR2yGTtjy6iT9egyyW2+6VPVjs6Z6TpTlcRaRbMbLa758WzbyJb+D2AR8J+/BbAU/tK9nVR00XRtB27+Ow3Y2jRQkldRKRSIkfpzAOOSNT5IbgoWl0Lv1vvLCV7EZEqmvTq2OMmjKV1xt5LAepiqYhI9Zr01Aq6WCoiEr+EXbSti9petBURSXW1uWjbpLt0REQkfkr4IiIpQglfRCRFKOGLiKQIJXwRkRTRqEbpmFkB8EUdD88CvnoXVvI1hjgaQwygOKpSHHtrDHE0hhigfnEc5O7Z8ezYqBJ+fZjZrHiHJjX3OBpDDIpDcTSFOBpDDMmMQ106IiIpQglfRCRFNKeE/0DUAYQaQxyNIQZQHFUpjr01hjgaQwyQpDiaTR++iIjsW3Nq4YuIyD4o4YuIpIgmn/DNbJKZbTSzBlsrtw4x9DazN81soZl9ambXRhRHGzP7yMw+CeP4dRRxxMTT0sw+NrMGXemsljGsNLP5ZjbXzCKZitXMOpnZM2b2mZktMrPjIohhQPg7qHxsN7Prkh1HGMv14d/nAjN7wszaRBTHtWEMnybzd1FdzjKzLmY2zcyWhD87J6LuJp/wgYcJFkePUhnwE3cfBBwL/NDMBkUQRzEw0t2HAsOAMWZ2bARxVLoWWBRh/ZVGuPuwCMdb3wW86u4DgaFE8Dtx98Xh72AYcBRQBDyf7DjMrCdwDZDn7oOBlsDFEcQxGPg+MJzg3+QcM+ubpOof5qs56+fAdHfvB0wPnze4Jp/w3f0dYHPEMaxz9znh9g6C/9A9I4jD3X1n+DQtfERyVd7MegFnAw9GUX9jYWYdgZOBhwDcvcTdt0YbFaOAZe5e17va66sV0NbMWgEZwNoIYjgU+NDdi9y9DHgbOD8ZFdeQs74OPBJuPwJ8IxF1N/mE39iYWS7BWr4fRlR/SzObC2wEprl7JHEAfwZ+ClREVH8lB6aa2WwzGx9B/X2AAuAfYffWg2aWGUEcsS4GnoiiYndfA9wBrALWAdvcfWoEoSwATjKzrmaWAZwF9I4gjkrd3X1duL0e6J6ISpTwG5CZtQOeBa5z9+1RxODu5eHX9l7A8PCra1KZ2TnARnefney6q3Giux8JnEnQ1XZykutvBRwJ3OfuRwC7SNDX9XiYWTpwLvB0RPV3JmjN9gEOBDLN7NvJjsPdFwG3AVOBV4G5QHmy46iOB2PlE/LNXAm/gZhZGkGyf9zdn4s6nrDb4E2iub5xAnCuma0EpgAjzeyxCOKobFHi7hsJ+qyHJzmE1cDqmG9azxB8AETlTGCOu2+IqP7TgBXuXuDupcBzwPFRBOLuD7n7Ue5+MrAF+DyKOEIbzKwHQPhzYyIqUcJvAGZmBH20i9z9zgjjyDazTuF2W2A08Fmy43D3m929l7vnEnQfvOHuSW/FmVmmmbWv3AZOJ/gqnzTuvh7IN7MBYdEoYGEyY6jiEiLqzgmtAo41s4zw/80oIrqwb2bdwp85BP33k6OII/QicFm4fRnwQiIqaZWIkyaTmT0BnApkmdlq4Ffu/lCSwzgB+A4wP+w/B7jF3V9Ochw9gEfMrCXBh/lT7h7ZkMhGoDvwfJBXaAVMdvdXI4jjx8DjYXfKcuCKCGKo/NAbDfwgivoB3P1DM3sGmEMwuu1jopve4Fkz6wqUAj9M1sX06nIW8AfgKTP7LsEU8RcmpG5NrSAikhrUpSMikiKU8EVEUoQSvohIilDCFxFJEUr4IiIpQglfGoyZ9TKzF8IZ/5aZ2V3hUMR4jn3YzC5IUFwrzSyrAc5zuZkdGPP8wdpMkmdmeWZ2d33jEKkrJXxpEOFNNM8B/wxn/OsPtAN+V82+9b7/oyHOUQeXE0wHAIC7f8/d476Jyt1nufs19Qlgf+87vAdDpFpK+NJQRgJ73P0fEMzpA1wPjAvvqrzczF40szeA6Ra4x8wWm9nrQLfKE5nZUWb2djjh2Wsxt5y/ZWZ/Due1vza8s/hZM5sZPk4I9+tqZlPDec4fBKy6gM3sEgvmyl9gZrfFlO80s4nh8dPDei4A8ghuoJprZm3DePJijrk9POZ1Mxsevr7czM4N9znVwrUBzOxl+8/c9NvM7LJw4rvbw/cyz8x+EHPcDDN7kWru0g3r/pOZfQIcF/uNJvxW8Va4fasFc7FXxlWvDx9pgtxdDz3q/SCY43xiNeUfA0MIWsergS5h+fnANIL50A8EtgIXEEzp/D6QHe53ETAp3H4L+GvMuScTTI4GkEMwtQXA3cAvw+2zCSaiyqoS14EEt/lnE9yF+wbwjfA1By4Nt38J3BNTf17MOb58Hh5zZrj9PMGkXGkEc63PDctPBV6qEsdRwDygIzAe+EVY3hqYRTDJ2KkEk671qeF378CFMc9XVr5fgg+pt8LtW8PfbWsgC9gEpEX9t6NH8h5NfmoFaVKmuXvlPOAnA0948E1gbdjyBxgADAamhVMitCSYRrfSkzHbpwGDwv0AOlgwY+nJhHObu/u/zGxLNbEcTZAICwDM7PHwuH8STOlcWc9jBF1V+1NCMOsiwHyg2N1LzWw+kFvdAWEr/H8JkvU2MzsdGBJzLaMj0C8890fuvqKGussJJu6Lx7/cvRgoNrONBNNPrI7zWGnilPCloSwkaKF/ycw6ELS8lxLMELkrjvMY8Km717QMYOw5WgDHuvueKvXGG3O84pl/pNTdK/erIFh9DHevqK7fPexrnwL8xt0rJ3Qz4Mfu/lqVfU9l37+7PeEHZ6Uy/tNdW3X5wOKY7XKUA1KK+vCloUwHMszsv+DLhPYn4GF3L6pm/3eAi8J+6x7AiLB8MZBt4bqvZpZmZofVUOdUgonJCPcdFnPusWHZmUB164N+BJxiZllhrJcQrHoEwf+Lyg+vscC74fYOoH0NsdTWH4B57j4lpuw14CoLptrGzPpb3RZLWUnQVQTwzXpFKc2KEr40iLB1ex7wLTNbQjC3+B7glhoOeR5YQvDN4FHg3+F5SgiS7W3hRci51Dxf+jVAXniBcyFwZVj+a+BkM/uUoGtnVTXxriNYiORN4BNgtrtXTkm7i2DxmAUEF6N/E5Y/DNxfedF237+R/boROD3mwu25BMtBLgTmhHX/jbq1wH8N3BVe3G4Ui3pI46DZMkWqMLOd7t4u6jhEGppa+CIiKUItfBGRFKEWvohIilDCFxFJEUr4IiIpQglfRCRFKOGLiKSI/w8JCO338REkdAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEWCAYAAACJ0YulAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4wLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvqOYd8AAAIABJREFUeJzt3XmYHFW9//H3JwQkgbAmIEtCUCSCgAiDLOo1sihBFPmx5CIgixIQFLiiIMKVIAgoqwhhE25kERM2MSwCAdlkM8HIDiJbAghJ2AITlpDv74+qGZowM90zmarTPf15PU8/091VU/WZfpLz7VPLOYoIzMzMAPqlDmBmZvXDRcHMzNq5KJiZWTsXBTMza+eiYGZm7VwUzMysnYuC1Q1JX5L0eOocHZE0UtKM1DnMiuaiYAtN0jOS5kp6s+JxRg2/F5LWaHsdEXdExIiCMo6XdGwR225kLna2oP6pA1if8Y2ImJw6RKOT1D8i5qXOUatGy2vVuadghZK0hqTbJL0uaZakCfn7t+er/DPvWYxe8Ftr3gP5iaQHJL0l6XxJK0q6XtIcSZMlLVux/mWS/pPv63ZJn8nfHwPsChya72tS/v7Kkq6QNFPS05IOrNjWgLx38aqkR4CNqvydIelASU/lf+eJkvrlyz4p6RZJs/Nll0haZoG/8zBJDwBvSeov6aeS/p3/nY9I2r5i/T0l/U3SqZJey/e5Wf7+dEkvS9qjYv2PSTpJ0nOSXpJ0dv73LQFcD6xc0cNbWVK/iv3PljRR0nL5tobnf+t3JT0H3CJpcUkX5+u+Junvklbs3r8UqxsR4YcfC/UAngG27GTZpcARZF9AFge+WLEsgDUqXo8EZiyw3XuAFYFVgJeB+4HP5du6BTiqYv29gUHAx4DTgGkVy8YDx1a87gdMBX4OLAZ8AngK+Fq+/ATgDmA5YCjwUGW2Dv7OAP6arz8MeAL4Xr5sDWCrPNcQ4HbgtAX+zmn5fgbk7+0ErJznHA28BayUL9sTmAfsBSwCHAs8B5yZ7+OrwBxgyXz9U4E/59kGAZOA4zv6zPP3Dso/91Xz7Z0DXJovG57/rRcCSwADgH3zbQ7M82wILJX636UfPfz/nDqAH43/yBu1N4HXKh775MsuBM4FVu3g92opCrtWvL4COKvi9Q+BP3WSaZl8+0vnrxcsChsDzy3wO4cD/5c/fwrYumLZmBqKQuX6+wM3d7Lut4B/LPB37l3lM54GbJc/3xP4V8WydfP9r1jx3mxgfUB5QflkxbJNgac7+szz9x4Ftqh4vRLwHtnh5rai8ImK5XsDdwHrpf636MfCP3xOwXrLt6LjcwqHAscA90l6FTg5Ii7oxnZfqng+t4PXSwJIWgT4Jdk37CHA/HydwcDrHWx3NbLDJq9VvLcIWe8Asm/p0yuWPVtD1gXXXznPtiLwG+BLZN/U+wGvdvG7SPoO8COyRhiyv3NwxSoLfg5EREefzRCyb/BTJbVvnuxv7cxqwFWS5le89z5Zj62jvBeR9XL+mB8Wuxg4IiLe62IfVqd8TsEKFRH/iYh9ImJlssMM4yqvOOpF3wa2A7YEluaDxrStJVxwOODpZN+Wl6l4DIqIbfLlL5I1dG2G1ZBhwfVfyJ8fl+9/3YhYCtitIleb9nySVgPOA34ALB8Ry5Advlrwd2oxi6xAfKbi71w6IpZccL8VpgOjFvhsFo+I5zvKGxHvRcTREbE2sBmwLfCdHmS1OuCiYIWStJOkVfOXr5I1Jm3fQF8iO5bfGwYB75AdNhlI1hBXWnBf9wFz8hO8AyQtImkdSW0nlCcCh0taNs//wxoy/CRffyjZcfkJFdneBF6XtArwkyrbWYLsc5oJIGkvYJ0a9v8RETGfrMCcKmmFfHurSPpavspLwPKSlq74tbOBX+bFCUlDJG3X2T4kfUXSunlv7Q2yQ03zO1vf6puLgvWWSfrwfQpX5e9vBNwr6U2yk50HRcRT+bKxwO/zK1Z2Xsj9X0h2yOZ54BGyE6WVzgfWzvf1p4h4n+wb7frA02TfqH9H1ssAODrf3tPAjWSHSKq5muzk9TTg2nyfbdvagOww1rXAlV1tJCIeAU4G7iZrtNcF/lbD/jtzGPAkcI+kN4DJwIh8X4+RXQzwVP7ZrEx2qOvPwI2S5pB9lht3sf2PA5eTFYRHgduo7fOyOqQIT7JjtrAkBfCpiHgydRazheGegpmZtXNRMDOzdj58ZGZm7dxTMDOzdg1389rgwYNj+PDhqWOYmTWUqVOnzoqIIdXWa7iiMHz4cKZMmZI6hplZQ5FUy135PnxkZmYfcFEwM7N2LgpmZtbORcHMzNq5KJiZWbvCioKkoZL+mk8l+LCkgzpYR5JOl/SksikXNygqD8DYsWOL3HxD2arfTh95mFl9K6MNK+yOZkkrkU0feL+kQWSjR34rHwGybZ1tyIYk3oZsFMbfRERXozHS0tISPb0kVRK+g5suC8BN8y8rMYmZdcfCtGGSpkZES7X1CuspRMSLEXF//nwO2ZC6qyyw2nbAhZG5B1gmLyZmZpZAKecUJA0nm2z93gUWrcKHp/WbwUcLB5LGSJoiacrMmTO7te+xY8ciibapCNue+1CSmTWCstuwwgfEk7Qk2aQbv4yIKxdYdg1wQkTcmb++GTgsIjo9PtTTw0eTJ09mq6228uEjfPjIrBEtbBuW/PBRHmJR4ArgkgULQu55Pjyv7ar5e73u2GOPLWKzZmalKKsNK/LqI5FNR/hoRJzSyWp/Br6TX4W0CfB6RLxYVKbVVlutqE03lM56A+4lmNW3MtqwInsKXwB2BzaXNC1/bCNpP0n75etcBzxFNn/secD+BebBo6t+4Kb5l3HstT/90Gszq29ltGGFjZKanydQlXUCOKCoDNa1jUdtmDqCmdUZ39FsZmbtGm4+hZ4655xzUkeoay88N5OVh1Wdf8PMEimrDWuansKIESMYMWJE6hj1Jz/Ad/zoU9PmMLMuldWGNU1RmDRpEpMmTUodo+4MWmYJAJ78x9OJk5hZV8pqw5qmKJx88smcfPLJqWPUnVH7bAnAvHfnJU5iZl0pqw1rmqJgHdvnhN1SRzCzOuKiYGZm7VwUrF1ra2vqCGaWmIuCtRt30PjUEcwssaa5T+Giiy5KHaFuLTZgUd6d+x63T7yHH59X6EgjZtZDZbVhTVMUhg4dWn2lJrXulz7N1BsfZO6cuamjmFknymrDmubw0YQJE5gwYULqGHXpxxe7d2BW78pqw5qmp3DWWWcBMHr06MRJ6s/gwYNTRzCzKspqw5qmp2BmZtW5KNiH3PD7W1JHMLOEXBQMgH79s38K/3fEHxMnMbOUXBQMgJU+sQIAr/zntcRJzCylpjnRfPnll6eOUNcOHLcPh215DDE/Ukcxsw6U1YY1TVHwFTZd22Dz9VJHMLMulNWGNc3ho/HjxzN+/PjUMczMeqSsNsxFwT7iuceeSx3BzBbgomDly6fmPG63M9LmMLNkXBSs3VLLLQnAsw+5p2DWrFwUrN32B40CYN677ydOYmapuChYu92O3Dl1BDNLrGkuSb3uuutSRzAz67Gy2rCmKQoDBw5MHaGhtLa2+jMzqyNl/X9smsNH48aNY9y4caljNIzT9zsvdQQzq1BWG9Y0RWHixIlMnDgxdYy6t9iARQH425+mJE5iZpXKasOapihYbTbYYh0A3n7z7cRJzCwFFwX7kMP/eHDqCGaWUGFFQdIFkl6W9FAny5eWNEnSPyU9LGmvorJY7Xxy2ay5FdlTGA9s3cXyA4BHIuKzwEjgZEmLFZjHzMyqKOyS1Ii4XdLwrlYBBkkSsCTwCjCvqDy33nprUZvus6455wa23fdrqWOYGeW1YSnPKZwBrAW8ADwIHBQR8ztaUdIYSVMkTZk5c2aZGZtS29ScF//CExOZNZuUReFrwDRgZWB94AxJS3W0YkScGxEtEdEyZMiQHu3spJNO4qSTTupx2Gay6pofBzw1p1k9KasNS1kU9gKujMyTwNPAp4va2TXXXMM111xT1Ob7lP85bz8AwjNzmtWNstqwlEXhOWALAEkrAiOApxLmsdw6m66VOoKZJVLYiWZJl5JdVTRY0gzgKGBRgIg4GzgGGC/pQbLpXQ6LiFlF5TEzs+qKvPpolyrLXwC+WtT+rXc8/fCzrP6Z1VLHMLOSNM0dzQMGDGDAgAGpYzSOfGrOE3b/bdocZgaU14Y1zdDZ119/feoIDWXpIUvx+stv8OxD01NHMTPKa8Oapqdg3bPD/2wLwPvzOrx1xMz6qKYpCscccwzHHHNM6hgNY5fDtk8dwcwqlNWGNU1RuPnmm7n55ptTxzAz65Gy2rCmKQrWc62trakjmFlJXBSsqlP2Pjt1BDMriYuCdepjAz8GwL3X3p84iZmVpWmKwvLLL8/yyy+fOkZDafnqugC8/dY7iZOYWVltmKLBRj1raWmJKVM8qXwZWltb2W7JPQC4af5lidOY2cKQNDUiWqqt1zQ9Bes+T81p1nyapigcfvjhHH744aljmJn1SFltWNMMc3H33XenjtDQLjtlEjv96BupY5g1rbLasKo9BUlbdvDeHsXEsXqzSD4158QTr06cxMzKUMvho59LOkvSEpJWlDQJ8FfGJjF07VUAeO3l1xMnMbMy1FIUvgz8m2w+5TuBP0TEjoWmsrrx4wsOyJ401kVqZtZDtZxTWBb4PFlhWBVYTZKiwa5lXXXVVVNHaEgjNvhk6ghmRnltWC1F4R7ghIi4QNIA4FfA34DNCk3Wyy6++OLUEczMeqysNqyWorBlRDwHEBFzgQMl/VexsawePX7/v91zMOvjajmnMF3SbpJ+DiBpGPB2sbF638EHH8zBBx+cOkZjyqfmPHHPM9PmMGtiZbVhtfQUxgHzgc2BXwBzgCuAjQrM1eumTZuWOkLDWnbFpXn1P68z/bHnU0cxa1pltWG19BQ2jogDyHsHEfEqsFihqayu7HLEDgDM99ScZn1eLUXhPUmLkF+UKGkIWc/BmsT2B4xKHcHMSlJLUTgduApYQdIvye5VOK7QVGZmlkTVcwoRcYmkqcAWZKccvxURjxaerJetueaaqSP0Ca2trR491SyBstqwTudTkLRcV78YEa8UkqgKz6eQxlb9dgJg0+1a+MVVhyVOY2bdVet8Cl31FKaSnUcQMAx4NX++DPAcsHov5LQGsfiSH+PtN99h6o0PpI5iZgXq9JxCRKweEZ8AJgPfiIjBEbE8sC1wY1kBe8uYMWMYM2ZM6hgNa5OvbwjAu3PfTZzErDmV1YbVcqJ5k4i4ru1FRFxPgw1xAfDEE0/wxBNPpI7RsP7n/H1TRzBramW1YbXcvPaCpCOBtoE3dgVeKC6S1SOfXDZrDrX0FHYBhpBdlnoVsEL+XpckXSDpZUkPdbHOSEnTJD0s6bZaQ5uZWTGqFoWIeCUiDoqIz+WPg2q88mg8sHVnCyUtQzaExjcj4jPATrWGtrQu/dVVqSOYWUGqHj6StCbwY2B45foRsXlXvxcRt0sa3sUq3waurBiB9eXqcXtu/fXXL3LzTaH/oosw7733ufyUa9jlsO1TxzFrKmW1YbWcU7gMOBv4HfB+L+57TWBRSbcCg4DfRMSFHa0oaQwwBmDYsGE92tlpp53Ws5TWbthnVuWpac/yxqw3UkcxazpltWG1nFOYFxFnRcR9ETG17dEL++4PbAh8Hfga8L95r+QjIuLciGiJiJYhQ4b0wq6tJ474w4HZk4aac8/MuqOWnsIkSfuTnWR+p+3NXrijeQYwOyLeAt6SdDvwWaCQa6522203wDOwLYxhn+5ZL83MFl5ZbVgtRWGP/OdPKt4L4BMLue+rgTMk9Scbintj4NSF3GanZsyYUdSmzcwKV1YbVsuAeD0azkLSpcBIYLCkGcBRwKL5Ns+OiEcl/QV4gGwo7t9FRKeXr1p9eejuR1ln07VSxzCzXtZpUZC0eUTcIun/dbQ8Iq7sasMRUfVehog4ETixakqrG5KICE757llc8MjpqeOYWS/rqqfwZeAW4BsdLAugy6JgfdNyH1+a2S++xvNPvpQ6ipkVoNOiEBFH5T/3Ki9OcTbddNPUEfqE3Y/eidPGnOepOc1KVlYb1ul8CvXK8ymk1za3wk3zL0ucxMxqVet8CrXcp2BmZk2iy6IgqZ+khhsmuyM77LADO+ywQ+oYfcqcOXNSRzBrGmW1YV0WhYiYD5xZeIoSzJ49m9mzZ6eO0aec8N+/SR3BrGmU1YbVcvjoZkk7SFLhaawhLL7k4gBMu/WRxEnMrLfVUhT2JRsU711Jb0iaI8kjojWx/9pxYwDenfte4iRm1ttqmU9hUET0i4hFI2Kp/PVSZYSz+nTAGXunjmBmBallPgWRTcG5ekQcI2kosFJE3Fd4ul60xRZbpI7QZ3hqTrPyldWGVb1PQdJZZGMTbR4Ra0laFrgxIjYqI+CCfJ9CffC9CmaNpTfvU9g4Ig4A3gaIiFfJRjU14+KjJ6SOYGa9qJai8J6kRcinVpE0hKzn0FBGjRrFqFGjUsfoM/ovtggAV/72L4mTmDWHstqwWorC6WQT7Kwg6ZfAncDxhaYqwNy5c5k7d27qGH3GJz47HIA5r76ZNohZkyirDavl6qNLgEPJCsGLwLciYmLRway+HXHZIdmTxho6y8yqqOXqo4siYnfgsQ7esya18jDPlW3WF9Vy+OgzlS/y8wsbFhPHzMxS6mrmtcOBnwED8juY24a5eBc4t4RsvWrbbbdNHaHPmjp5GhtuuX7qGGZ9WlltWC33KRwfEYeXkqYGvk+hfny1/87E/GClNVbkwifOSB3HzLrQm/cpHCFpN0n/m294qKTPL3RCa3jLr7wsAC89MzNxEjPrLbUUhTOBTYFv56/fpAGH0x45ciQjR45MHaNP2efX2bUGnprTrHhltWG+o9l6bPP//mLqCGbWy5rmjmYzM6uup3c0H1doKms4s2bNSh3BzHpBT+9o9tCY9iG/3sVXH5n1BVXvaM69BNyRrz9A0gYRcX9xsXrfzjvvnDpCnzRg0ADmzpnLQ397InUUsz6trDaslmEujgH2BP7NByPdBLB5cbF63/777586Qp+0xbe/yDXn3MR7b3tqTrMildWG1dJT2Bn4ZES8W3SYIrW2tgKeNay37XPyblxzzk2pY5j1eWW1YbUUhYeAZYCXC01SsG222QaAW2+9NW2QPsZF1qwcZbVhtRSF44F/SHoIeKftzYj4ZmGpzMwsiVqKwu+BXwEP4vsTrAvnHnohY379ndQxzGwh1HKfQmtEnB4Rf42I29oe1X5J0gWSXs57GF2tt5GkeZJ2rDm11ZX+i2XfLf5y/i2Jk5jZwqqlKNwh6XhJm0raoO1Rw++NB7buaoX8TulfATfWsD2rUyNaPgnAnNfeSpzEzBZWLYePPpf/3KTivaqXpEbE7ZKGV9n2D4ErgI1qyLFQ9txzz6J30bQO/cNB7DF8f0/NaVagstqwqkUhIr5SxI4lrQJsD3yFKkVB0hhgDMCwYcN6tD8XheJ4ak6z4tVNUQCQ9HWyaTkXb3svIn6xkPs+DTgsIuZL6nLFiDiXfLa3lpaWHn0fbRubZ/DgwT35dTOzpMpqw2q5o/lsYCDZN/rfATsC9/XCvluAP+YFYTCwjaR5EfGnXtj2R+y4Y3Ye2/cpFOve66ey8ShP4W3W28pqw2o50bxZRHwHeDUijiabcGfNhd1xRKweEcMjYjhwObB/UQXBiqd+WW/vjB9ekDiJmS2MWg4fvZ3/bJW0MjAbWKnaL0m6FBgJDJY0AzgKWBQgIs7uUVqrW0OGLs/Lz87ipWc9NadZI6ulKEyStAxwInA/2TUm51X7pYjYpdYQEbFnretafdrvlD34xQ4nE+/7EiSzRtZlUZDUD7g5Il4DrpB0DbB4RLxeSjprGF/afpPqK5lZ3euyKORXBp1Jfq9CRLxDxfhHjeT73/9+6ghmZj1WVhtWy+GjmyXtAFwZEQ17bGD06NGpIzSNWbNm+dJfs15WVhtWy9VH+wKXAe9IekPSHElvFJyr102fPp3p06enjtEUjtvp9NQRzPqcstqwWu5oHlR4ihLsvvvugO9TKNISywzgrdfm8ug9/0odxazPKasNq6WngKRlJX1e0n+1PQpNZQ1py92+DMC8d+YlTmJmPVW1KEj6HnA7cANwdP5zbLGxrBH94PTvpo5gZguplp7CQWQD1j2bD473OeC1QlOZmVkStRSFtyPibQBJH4uIx4ARxcYyM7MUarkkdUZ+R/OfgJskvQo8W2ys3nfIIYekjtBUzjjwfB9OMutFZbVh6s6tB5K+DCwN/CUi3i0sVRdaWlpiypQpKXZtNRg1YBfmvTOPgUsP4OpXL0wdx8xykqZGREu19To9fCRpcUkHSzpD0r6S+ufzM/85VUFYGI8//jiPP/546hh93lobrQFA6xtzEycx61vKasO6Onz0e+A94A5gFLA22UnnhrTvvvsCvk+haD+78iB2WeH7nprTrJeV1YZ1VRTWjoh1ASSdT+9MrGN9nIe3MGtsXV199F7bk4jw3UhmZk2gq6Lw2XysozckzQHWa+Sxj6x8t19xV+oIZtZNnRaFiFgkIpbKH4Mion/F86XKDGmNRYtkU3OefchFiZOYWXfVcp9Cn3DkkUemjtA0VlxtBf7z1EvMmjE7dRSzPqOsNqxpisKWW26ZOkLT+MFv9+LIr59AzPclSGa9paw2rKZRUvuCadOmMW3atNQxmsLGozZMHcGszymrDWuansLBBx8M+D4FM2tMZbVhTdNTsDReeG5m6ghm1g0uClaM7AIkfrXLaWlzmFm3uChYIQYtswQAT0x9KnESM+sOFwUrxKh9sisl5r3rm+HNGknTnGg+7rjjUkdoKvucsBsTf3116hhmfUZZbVjTFIXNNtssdQQzsx4rqw1rmsNHd911F3fd5bF4UmhtbU0dwazhldWGNU1P4Wc/+xng+xRSOPeQizn4rDGpY5g1tLLasKbpKVj5Fl18UQBu+cOdiZOYWa0KKwqSLpD0sqSHOlm+q6QHJD0o6S5Jny0qi6WxzhfWBGDuHE/NadYoiuwpjAe27mL508CX89ndjgHOLTCLJXDopT9IHcHMuqmwcwoRcbuk4V0srzxjcg+walFZLA1PzWnWeOrlRPN3ges7WyhpDDAGYNiwYT3awWmnebgFM2tcZbVhyYuCpK+QFYUvdrZORJxLfnippaWlR4P0r7/++j3KZ71j8iV3sOWuX0odw6xhldWGJb36SNJ6wO+A7SKi0Gm6Jk+ezOTJk4vchXWgX//sn9j5P704cRKzxlZWG5aspyBpGHAlsHtEPFH0/o499ljAM7CVbcXhQ3jxyZeY/eKrqaOYNbSy2rAiL0m9FLgbGCFphqTvStpP0n75Kj8HlgfGSZomaUpRWSydg8/Oblrz1JxmjaHIq492qbL8e8D3itq/1YcNNl8vdQQz6wbf0WxmZu1cFKw0M57+T+oIZlZF8ktSy3LOOeekjtC8BAQcP/pUzrzvV6nTmDWkstqwpikKI0aMSB2haQ1adgnmvPIWTz3wXOooZg2rrDasaQ4fTZo0iUmTJqWO0ZS+cUA2BJan5jTrubLaMEU01qWCLS0tMWVK969eHTlyJOD5FFLZqt9OANw0/7LEScwa08K2YZKmRkRLtfWapqdgZmbVuShYqTw1p1l9c1GwUv12/wtSRzCzLrgoWCkWG5BNzXnnVfcmTmJmXWmaS1Ivuuii1BGa2voj1+a+6//J23PeTh3FrCGV1YY1TU9h6NChDB06NHWMpnXEZT9KHcGsoZXVhjVNUZgwYQITJkxIHaNpDRw4MHUEs4ZWVhvWNIePzjrrLABGjx6dOImZWfeV1YY1TU/B6se1v7sxdQQz64SLgpWmbWrOi8ZenjiJmXXGRcFKs8oaKwLwyouvJU5iZp1xUbDS/Oj87wPQaONtmTWTpjnRfPnlPmSR2jqbrpU6glnDKqsNa5qiMHjw4NQRzMx6rKw2rGkOH40fP57x48enjmG55x7zhDtm3VFWG+aiYOVS9uOX3z49bQ6zBuOiYH3SUoOXAuDZh2ckTmJmHXFRsFLt+KNtAXj/vfcTJzGzjrgoWKl2OWz71BHMrAsuCmZm1q5pLkm97rrrUkewBbS2tnr0VLMaldWGNU1PYeDAgW6A6swpe5+dOoJZwyirDWuaojBu3DjGjRuXOoYBiw1YDIB7r70/cRKzxlFWG9Y0RWHixIlMnDgxdQwDNvraegC8/dY7iZOYNY6y2rCmKQpWPw69+IepI5hZJworCpIukPSypIc6WS5Jp0t6UtIDkjYoKkubZ555puhdWA22W3KP9udb9duJrfrtlDCNWeMoow0rsqcwHti6i+WjgE/ljzHAWQVmAeDZZ58tehdWRWcFwIXBrLoy2rDCikJE3A680sUq2wEXRuYeYBlJKxWVx8zMqkt5TmEVYHrF6xn5ex8haYykKZKmzJw5s1s7GTt2LJK47bbb2raFJMaOHduz1GZmJSq7DVORs2BJGg5cExHrdLDsGuCEiLgzf30zcFhETOlqmy0tLTFlSperdJXHs34l1tVhopvmX1ZiErPGszBtmKSpEdFSbb2UPYXngaEVr1fN3zMzs0RSFoU/A9/Jr0LaBHg9Il4scodHHXVUkZu3GnTWG3Avway6Mtqwwg4fSboUGAkMBl4CjgIWBYiIsyUJOIPsCqVWYK9qh45g4Q4fmZk1q1oPHxU2IF5E7FJleQAHFLV/MzPrPt/RbGZm7VwUzMysnYuCmZm1c1EwM7N2hd68VgRJM4GeDgAyGJjVi3GagT+z7vHn1T3+vLpnYT6v1SJiSLWVGq4oLAxJU2q5JMs+4M+se/x5dY8/r+4p4/Py4SMzM2vnomBmZu2arSicmzpAA/Jn1j3+vLrHn1f3FP55NdU5BTMz61qz9RTMzKwLLgpmZtau6YqCpBMlPSbpAUlXSVomdaZ6JmknSQ9Lmi/Jlw52QtLWkh6X9KSkn6bOU+8kXSDpZUkPpc5S7yQNlfRXSY/k/xcPKnJ/TVcUgJuAdSJiPeAJ4PDEeerdQ8D/A25PHaReSVoEOBMYBawN7CJp7bSp6t54smHzrbp5wCERsTawCXBAkf++mq4oRMRQwyWNAAAEOUlEQVSNETEvf3kP2Yxv1omIeDQiHk+do859HngyIp6KiHeBPwLbJc5U1yLiduCV1DkaQUS8GBH358/nAI/SyXz2vaHpisIC9gauTx3CGt4qwPSK1zMo8D+tNa983vvPAfcWtY/CJtlJSdJk4OMdLDoiIq7O1zmCrFt2SZnZ6lEtn5eZpSVpSeAK4OCIeKOo/fTJohARW3a1XNKewLbAFuEbNap+XlbV88DQiter5u+Z9QpJi5IVhEsi4soi99V0h48kbQ0cCnwzIlpT57E+4e/ApyStLmkx4L+BPyfOZH1EPp/9+cCjEXFK0ftruqIAnAEMAm6SNE3S2akD1TNJ20uaAWwKXCvphtSZ6k1+4cIPgBvITgJOjIiH06aqb5IuBe4GRkiaIem7qTPVsS8AuwOb523WNEnbFLUzD3NhZmbtmrGnYGZmnXBRMDOzdi4KZmbWzkXBzMzauSiYmVk7FwVrGJLezy/He0jSZZIGps4EIOlnCfY5XtKOZe/X+j4XBWskcyNi/YhYB3gX2K/WX8xHMi1Kt4tCwXnMesxFwRrVHcAaAJL+JGlqPtb8mLYVJL0p6WRJ/wQ2lfRzSX/Pexrn5neKIulWSadKmiLpUUkbSbpS0r8kHVuxvd0k3Zf3Vs6RtIikE4AB+XuXdLZeR3kqtvtpSfdVvB4u6cH8eYeZK0l6RtLg/HmLpFvz50vk8xbcJ+kfkjxyq1XlomANR1J/srkLHszf2jsiNgRagAMlLZ+/vwRwb0R8NiLuBM6IiI3ynsYAsvGv2rwbES3A2cDVwAHAOsCekpaXtBYwGvhCRKwPvA/sGhE/5YMezK6drddJHgAi4jFgMUmr52+NBibkz7vKXM0RwC0R8XngK8CJkpboxu9bE+qTA+JZnzVA0rT8+R1k48FAVgi2z58PBT4FzCZrkK+o+P2vSDoUGAgsBzwMTMqXtY1V9CDwcES8CCDpqXybXwQ2BP6ef1kfALzcQcYtulhvwTyVJpIVgxPyn6NryFzNV4FvSvpx/npxYBjZUBxmHXJRsEYyN//23U7SSGBLYNOIaM0PnSyeL347It7P11scGAe0RMR0SWMr1gN4J/85v+J52+v+gIDfR0S1mfq6Wq89TwcmAJdJuhKIiPhXDZnbzOODXn/lcgE7eJIk6w4fPrJGtzTwal4QPk02XWFH2hrLWfm49N29cudmYEdJKwBIWk7Savmy9/Khjaut16mI+DdZT+J/+eDQUa2ZnyHrnQDsUPH+DcAPK86dfK5aDjMXBWt0fwH6S3qU7NDLPR2tFBGvAeeRzTl9A9lw1zWLiEeAI4EbJT1ANtf3Svnic4EHJF1SZb1qJgC7kR1K6k7mo4HfSJpCVljaHAMsmmd7OH9t1iWPkmpmZu3cUzAzs3YuCmZm1s5FwczM2rkomJlZOxcFMzNr56JgZmbtXBTMzKzd/wcIzGNUMH9n1QAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# waterfall, parameter space, scatter plots, fits to data\n",
    "# different functions for different plotting types\n",
    "import pypesto.visualize\n",
    "\n",
    "pypesto.visualize.waterfall(result)\n",
    "pypesto.visualize.parameters(result)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Data storage"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# result = pypesto.storage.load('db_file.db')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Profiles"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# there are three main parts: optimize, profile, sample. the overall structure of profiles and sampling\n",
    "# will be similar to optimizer like above.\n",
    "# we intend to only have just one result object which can be reused everywhere, but the problem of how to \n",
    "# not have one huge class but\n",
    "# maybe simplified views on it for optimization, profiles and sampling is still to be solved\n",
    "\n",
    "# profiler = pypesto.Profiler()\n",
    "\n",
    "# result = pypesto.profile(problem, profiler, result=None)\n",
    "# possibly pass result object from optimization to get good parameter guesses"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Sampling"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# sampler = pypesto.Sampler()\n",
    "\n",
    "# result = pypesto.sample(problem, sampler, result=None)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "# open: how to parallelize. the idea is to use methods similar to those in pyabc for working on clusters.\n",
    "# one way would be to specify an additional 'engine' object passed to optimize(), profile(), sample(),\n",
    "# which in the default setting just does a for loop, but can also be customized."
   ]
  }
 ],
 "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.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}