Tutorial 1 -- The Basics.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "7167adf3",
   "metadata": {},
   "source": [
    "# Counting Subgraphs with PACT\n",
    "\n",
    "This notebook is intended as a quick practical introduction for the most basic way to count subgraphs using **PACT (PAttern Counting Toolkit)**.\n",
    "\n",
    "We will work through a running example in which we want to count the how often all 7 node trees occur as subgraphs in an example host graph. The folder `data/tutorial` should contain files describing all of our pattern graphs (the trees with 7 nodes), the full *spasm base* of our patterns, and a host graph file `schi_pombe.gml`. The pattern and spasm base graphs were created with nauty (their creation process is described in detail at the end of this tutorial).  Our example host graph is the PPI network of *Schizosaccharomyces pombe* (fission yeast), derived from the [BioGRID](https://thebiogrid.org/) dataset.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "866790a0",
   "metadata": {},
   "source": [
    "# Part 1: Pattern Precomputations\n",
    "\n",
    "Our overall intention is to compute the number of times a pattern $P$ occurrs as a subgraph in hostgraph $H$, in short we write $\\#Sub(P \\rightarrow H)$ for this number. We will do so using the following identity:\n",
    "\n",
    "$$\\#Sub(P \\rightarrow H) = \\sum_{F \\in Spasm(P)} \\mu(F,P) \\cdot \\#Hom(F \\rightarrow H)$$\n",
    "\n",
    "where\n",
    "  - $Spasm(P)$ is the set of all homomorphic loop-free images of $P$,\n",
    "  - $\\mu(F,P)$ is a coefficient depending on $P$ and $H$, and\n",
    "  -  $\\#Hom(F \\rightarrow H)$ is the number of homomorphisms from $F$ into $H$.\n",
    "  \n",
    "We will refer to the set $Spasm(P)$ together with the corresponding coefficients $\\mu$ as the *homomorphism base* of $P$. Computing the homomorphism base as well as $\\#Hom$ is computationally expensive, but more flexible and cheaper (for larger patterns) than other known methods of obtaining $\\#Sub(P \\rightarrow H)$. \n",
    "\n",
    "Note that the identity is not a heuristic or similar, the resulting count is the **exact** number of non-induced subgraph occurrences of $P$ in $H$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "ff036059",
   "metadata": {},
   "outputs": [],
   "source": [
    "%load_ext rich\n",
    "%load_ext autoreload\n",
    "%autoreload 2\n",
    "import networkx as nx\n",
    "import matplotlib.pyplot as plt\n",
    "from itertools import islice\n",
    "plt.rcParams[\"figure.figsize\"] = (2, 2)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0cc33c45",
   "metadata": {},
   "source": [
    "## Loading the Patterns\n",
    "\n",
    "For convenience we use the settings module for a unified interface to input files for different scenarios. Ideally, simply switching the imported settings allows us to seemlessly change the analysed patterns and spasm graphs.\n",
    "\n",
    "The individual graphs are stored as `GraphWrapper` objects, which wraps networkx graphs with various additional information that will be useful for our subgraph counting process, like a unique `id`. We do nothing else special with the patterns and just store them in a standard python contain (here a dictionary, indexed by the graphs internal id)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "71693ed9",
   "metadata": {},
   "outputs": [],
   "source": [
    "import settings.tutorial as settings\n",
    "from pact.graphwrapper import GraphWrapper\n",
    "\n",
    "patterns = dict()\n",
    "\n",
    "with open(settings.raw_pattern_filename(), 'r') as f:\n",
    "    lines = map(str.rstrip, f.readlines())\n",
    "    for line in lines:\n",
    "        G = GraphWrapper.from_g6str(line, **settings.PATTERN_GRAPHWRAPPER_PARAMS)\n",
    "        patterns[G.id] = G"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "b90daee9",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
      ],
      "text/plain": []
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAANwAAADcCAYAAAAbWs+BAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/H5lhTAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAVwUlEQVR4nO3de1RUVd8H8O8ZBkZB0gQxUwgR8DWWF0yUpanNq7HsKV1Cq6Umlcv3XeaFMsq017tm+rjSZUtZ78rM9ZhXNPNCF1t4IfMKqA9o9K7OgMKYmsSooKA4l/P+YeMtFIbZ+5x9Zn6ff0f22TJ+nTl7n/37SYqiKCCEqMKg9QQI8ScUOEJURIEjREUUOEJURIEjREUUOEJURIEjREUUOEJURIEjREUUOEJURIEjREUUOEJURIEjREVGrSdA+Kqtd6DcVovbDheCjAZEh4UgxERvu1boN++DLJevY1O+FXm/VcJ6pQ73n7+SAES1DYa5awTG9otCXPtQrabplyQ6D+c7zl+pw8ydZ3CotAoBBglO16PfWvfrA2PDsTi1OyLbBqs4U/9FgfMR2YVWzMspgcOlPDZoDwswSDAaJCwYkYDRSVEcZ0gACpxPyMqzYFmu7PU401LikWGOYzAj8ii0Sqlz2YVWJmEDgGW5MrYWWpmMRRpGgdOx81fqMC+nhOmYc3NKcP5KHdMxyT0UOB2bufMMHB7crzWFw6Vg5s4zTMck91DgdMpy+ToOlVZ5tEDSFE6XgkOlVSitvM50XHIHBU6nNuVbEWCQuIwdYJCw8Tjdy/FAgdOpvN8qmX+6uTldCvLkSi5j+zsKnA7dqHfAynlhw2qrQ229g+s1/BEFTocqbLXgvXmqACi31XK+iv+hwOnQbYfLp67jTyhwOhRkVOdtU+s6/oR+ozoUHRYCPuuT90h/XYewpcrxHDqTxVaIyYiotsGo4LhwEhUWTO8RB9x+o3Qmiy9z1whsyK/gsjUQYJBgjo9gPi7hcFqAzmSpw3L5Ol787Gdu4+/LHITYCPqPkDWm93DZhVYMXXEQR8/aAKDR/33drx89a8PQFQeRTU+qN1lc+1AMjA1n/rRJgEHCwNhwChsnzAKXlWfBRzvOoN7h8vhrjtOloN7hwkc7ziArz8JqSj5vcWp3GBkHzmiQsDi1O9MxyT1MAkdnsrQR2TYYC0YkMB1z4YgE+mrPkdeBozNZ2hqdFIVpKfFMxvowpStGUZkFrrwOHJ3J0l6GOQ7/TOsOk9Hg8T1dgEGCyWjA0rTumGKO5TRD4uZV4OhMljhGJ0VhX+Zg9I8JA4BGg+d+vX9MGPZlDqZPNpV4tS0wP6eE617QG/2ewXzG9yj+wL0Hun7fSbiC2wLSvfBJuLOpbY6PQHpyFK1GqsyrwA3+NI/r0w7PhAXj4DQzt/F9maIoCA0Nxaz5CzHi9f+ip3wE0ezfvJpnsugfiOcuXbqE2tpaJMTHIuHp1lpPh/yl2fdwdCZLbLJ8Z5smPp7NCiZho9mBozNZYpNlGQaDATExMVpPhdyn2YGjM1lis1gsiI6ORlBQkNZTIfdp9r9mOpMlNlmW6eukgJodOPeZLJ7oTFbzUeDE5NX3NXPXCK61EelMVvM4HA6UlZVR4ATkVeDG9oviWhsxPZmefmgOq9UKu92OuDjqhCMarwJHZ7LERFsC4vJ6CZD1mSxFUWCAQmeyvCDLMkwmEyIjI7WeCnmI14FjfSZLkiRU/fi/OPDt18zG9DeyLCM2NhYBAQFaT4U8hMkmF8szWZn/2QVpvZ7CuHHjkJmZCYeDym17ymKx0P2boJjtKrM6kzX1xf/Al19+iVWrVmHVqlUYNmwYbDYbq2n6BdoSEBfTxzhYncmSJAkZGRnYu3cvioqKkJSUhDNn6EBqU9y6dQsVFRUUOEExL5PndrcupVwJq62BupRNPJN17tw5jBw5EmVlZVi/fj3S0tJ4TNdn/Prrr0hISMDPP/+MgQMHaj0d8jBFBTdu2ZVfLlxTTlVcUX65cE25ccvu2c/fuKG89tprCgBl7ty5itPp5DRT/du5c6cCQLl06ZLWUyENUOW5qRCT0aszWSEhIdi6dSt69eqF2bNno7i4GBs2bEBoKO3TPUyWZYSGhqJ9+/ZaT4U0QDeP4kuShJkzZ2L37t04cOAAkpOTUVpaqvW0hONeMJEk3o+Wk+bQTeDchg8fjvz8fNjtdiQlJSE3N1frKQnFYrHQgonAdBc4AOjWrRsKCgqQnJyMl156CcuXL4fi4dpPbb0DJRer8W/rVZRcrPaZ9rqyLNMenMB0e/alTZs2+O677zBz5kxMmzYNRUVF+OKLL9CyZctH/oyvd/SpqanBH3/8QZ9wAuO2LaCmLVu2YPz48UhISMCuXbvQqVOnB173l44+J0+eRJ8+fVBQUICkpCStp0MaoMuvlA8bM2YMjhw5gsrKSvTp0wdHjhy5+5o/dfSxWO40QqGvlOLyicABQO/evXHixAnEx8fDbDZjzZo1ftfRR5ZltGvXDm3atNF6KuQRdHsP15CIiAjs27cP7733Ht7P2o6wfzzNZNxluTLatTIJXw6cnqEUn898wrkFBQXhfxYtQ/uX3/F45fJx9NDRhwInPp8LHHCno48iGZhu/ore0UdRFDqWowM+Fzh/7ehTVVWFa9eu0Sec4HwucJvyrVwriW08LuaqJdUx0QefC1zeb5VcK4nlyZVcxvaWe0sgNpaaKorMpwKnZkcf0ciyjMjIyMc+aUO051OB8+eOPrRCqQ8+FTh/7uhDgdMHnwqcv3b0cblcKC0tpcDpgFj/crzkrx19Lly4gJs3b9IenA74VOD8taMPbQnoh08FDvDPjj6yLMNoNCI6OlrrqZBG+Fzg/LGjj8ViQUxMDAIDA7WeCmmEzwXOHzv6UFkF/fC5wAH+19GHtgT0wycDx6OjT+UPWfhu61dMj/ywYLfbcfbsWQqcTvhk4AC2HX3eHxKLN5+PxeTJk/H222+jvr6eybgslJeXw+l00ldKnRBrfZuxDHMcwluZMC+nBA6X4tFiSoBBgtEgYeGIhDsnvYeuQmJiIiZNmoRffvkF33zzDTp06MBx9k1DWwL64rOfcG6sOvoAwPjx43Hw4EGUl5ejT58+yM/P5zfxJpJlGS1btkTHjh21ngppAp8ok9dUrDr6XLp0CWlpaTh16hRWr16NcePG8Z76I02ePBlHjhxBcXGxZnMgTedXgbtfbb0D5bZa3Ha4EGQ0IDosxKMnSOrr6zFlyhSsXbsW7777LpYtW6bJPtjQoUPRpk0bbN++XfVrE8/59D3c43jb0cdkMmHNmjXo3bs3pk6ditOnT2Pbtm1o164dw1k2TpZlpKenq3pN0nw+fw/HkyRJmDx5Mvbv34+SkhIkJSWhqKhItevX1dXh/PnztGCiIxQ4BgYNGoQTJ04gLCwM/fv3R3Z2tirXLSsrA0ArlHpCgWMkKioKhw4dQmpqKsaMGYMZM2bA6XRyvaZ7S4D24PTDb+/heAgODsbGjRvRu3dvTJ8+HcXFxdiyZQuefPJJLteTZRlt2rRBeHg4l/EJe/QJx5gkSfjggw/w448/oqCgAH379kVJSQmXa1G3U/2hwHHy4osvorCwEC1atEBycjJ27drF/BrU7VR/KHAcdenSBceOHUNKSgpSU1Mxf/58uFzsChDRsRz9ocBx1qpVK3z99df4+OOPsXDhQqSlpaGmpsbrca9evYo///yTPuF0hgKnAoPBgNmzZ2P37t04cOAAkpOT71ZKbi73z1Pg9IUCp6Lhw4ejoKAATqcTSUlJ2LNnj8dj1NY7UHKxGgeKzyIwojM6PhPDYaaEF799llJL1dXVGDt2LH744QcsXrwYM2bMeOxK492Hrn+rhPVKAw9dtw2GuWsExvaLQlx78UpAkHsocBpxOp2YN28ePvnkE4waNQpr165FSMiD9S7PX6nDzJ1ncKi0CgEG6bHn+dyvD4wNx+LU7ojkXC6QNA8FTmPbt2/HW2+9hbi4OOzatetuqbvsQqtXB2cXjEjAaMFbJPsjCpwATp8+jZEjR6Kmpgbbtm3Dr1IkluXKXo87LSUeGWbaNhAJBU4QNpsNo0aNQoEtEG2HZTAbd2la9wdOrRNtUeAEcu7PGgxZ/hOcYNef3GQ0YF/mYLqnEwRtCwhkbs7/QQowMn020uFSMHPnGWbjEe9Q4ARhuXwdh0qrmJdpd7oUHCqtQmnldabjkuahwAliU76VaxOSjcetXMYmnqHACSLvt0quTUjy5EouYxPPUOAEcKPeAeuVOq7XsNrqUFvv4HoN0jgKnAAqbLXgvVSsACi31XK+CmkMBU4Atx3szsiJcB3yaBQ4AQQZ1Xkb1LoOeTR6BwQQHRYC3lVJpL+uQ7RFgRNAiMmIKM5PgkSFBXtUyp3wQYEThLlrBNd9OHN8BJexiWcocIIY2y+K6z5cejI9wCwCCpwg4tqHYmBsOPNPuQAJGBgb/tj2W0Q9FDiBLE7tDiPDwCmKAof9NkZ0oP03UVDgBBLZNhgLRiQwG0+SJLSz5mHUK0MxZ84cOBz0pInWKHCCGZ0UhWkpbErffZjSFfmbV2DRokVYsmQJBg0ahPLyciZjk+ahA6iC8ramycIRCQ+c9D5+/DjGjBmDq1evYvXq1Rg1ahSPaZNGUOAExrpqV3V1NSZOnIjs7GyMHz8eK1eu/FulMMIXBU4HLJev41+Hy/BVbiEC23YA7nsuRcKdTW1zfATSk6MaXY1UFAVfffUVMjIy0LFjR2RnZyMxMZHvX4DcRYHTiaKiIiQmJiLv0FG0i3kWtx0uBBkNiA4LadYTJLIsY/To0SgpKcHSpUsxdepUanulAlo00Ql3t9Mez3ZFwtOtkRj1JBKebt3sx7Xi4+Nx7NgxZGRkIDMzEy+//DIqK+mQKm8UOJ2QZRlhYWFo27YtszFNJhOWL1+OPXv24OTJk+jRowf27t3LbHzydxQ4neDZfHHYsGEoLi5Gz549kZKSgunTp+P27dtcruXvKHA6wbv54lNPPYU9e/bg008/xWeffYYBAwagtLSU2/X8FQVOJ9z9vHkyGAyYNm0ajh49imvXriExMREbNmzgek1/Q4HTAZvNhitXrqjWfLFPnz44deoUXn31Vbz55ptIT09n0rWVUOB0wd3tVM1+3qGhoVi3bh02bdqEnJwcJCYmoqCgoNnjuRtJ/tt6FSUXq/22ghgdAdYB95ZAbGys6td+/fXXkZycjDFjxmDAgAFYtGgRPvzwQxgMjf9fTY0k/442vnVg9uzZWLduHX7//XfN5mC32zF37lwsXboUQ4YMwfr169GhQ4cG/yw1knw0+kqpAzy3BJoqMDAQS5Yswd69e1FSUoIePXrg+++//9ufyy60YuiKgzh61gYAjT547X796Fkbhq44iOxC3y7JToHTAd5bAp4YMmQIiouL0a9fP7zyyiuYOnUqbt26BQDIyrPgox1nUO9weVwuwulSUO9w4aMdZ5CVZ+ExdSFQ4ASnKIoqWwKeaNeuHb799lusXLkSn3/+OZKTk7EiJ59J11YAWJYrY6uPftJR4AR38eJF1NXVCRU44M5p8nfeeQcFBQW4GRCCzw5dABgWbJ+bU4LznPstaIECJzgttgQ80bNnTzz39qcwBAQCDMvZ+mojSQqc4GRZhsFgQExMjNZTaZDl8nUcPXcVisT2n5KvNpKkwAlOlmV07twZQUFBWk+lQdRI0jMUOMGJsCXwONRI0jMUOMGJtCXwMGok6TkKnMAcDgfKysqE/YSjRpKeo8AJrKKiAna7XdjAUSNJz1HgBObeEhA1cNRI0nO+8zfxQbIsw2QyITIyUuupNIgaSXqOAicwWZYRGxvbpKMwWqBGkp4T850kANQpq+AtaiTpGQqcwCwWi7BbAm7USNIzFDhB3bp1CxUVFcJ/wnFrJGmQfLKRJAVOUGVlZVAURfjAAXwaSbocdkw3i7lY5A0KnKBE3xK4H49GknU/r8M/BvfDwYMHmY0rAgqcoGRZRmhoKCIi9LFowLqR5IltK9GlSxeYzWbMmTMHdrudydhao8AJyr1CqaeONhnmOPwzrTtMRoPH93QBBgkmowFL07pjijkWkZGR2L9/Pz7++OO73VvPnj3LaebqocAJSg9bAg0ZnRSFfZmD0T8mDAAaDZ779f4xYdiXOfiBrq0BAQGYNWsWDh8+jMuXL6NXr17YtGkTv8mrgMrkCapDhw6YMGECFixYoPVUmu1uXUq5ElZbA3UpPWgkWV1djcmTJ2Pz5s144403kJWVhSeeeILr/HmgwAmopqYGrVu3xsaNGzF27Fitp8NEbb0D5bZarxtJbty4EZMmTUJERAQ2b96Mfv36cZgtP/SVUkB6WqFsqhCTkUkjyfT0dBQVFSE8PBzPP/88lixZAqfT2ayxtCi/7jsPqfkQ0QsHaa1Lly44fPgw5s+fj1mzZiE3NxcbNmxAp06dGv1Zrcuv01dKAS1cuBBZWVnUArgJfvrpJ6Snp6Ourg5r165Fampqg39OlPLr9JVSQHpdodTCCy+8gOLiYgwePBhpaWmYOHEi6uoeLPsgUvl1CpyAKHCeCQsLw44dO7B69WqsX78ezz33HIqLiwGIV36dAicYRVF0cUpANJIkYcKECTh58iRMJhP69u2L/17yL+HKr1PgBFNVVYVr167RJ1wzdevWDcePH8dbU95Hrq01wHCJgkX5dQqcYNzNFylwzdeiRQvcTBgJY6AJYPhoHIvy6xQ4wWjZ7dRXWC5fx6HSKrCu9cWi/DoFTjAWiwWRkZFo2bKl1lPRLZHLr1PgBEMrlN4Tufw6BU4wFDjviF5+nQInEJfLhdLSUgqcF0Qvv06BE8iFCxdw8+ZN2oPzgujl1ylwAqEtAe+JXn6dAicQWZZhNBoRHR2t9VR0S/Ty6xQ4gVgsFnTu3BmBgYFaT0W3RC+/ToETCK1QsiFy+XUKnEAocGyIXH6dAicIu92Oc+fOUeAYELn8OgVOEOXl5XA4HLQlwAjr8usAYDRIWJza3asxKHCCoC0BtliXXweAhSMSvC63QIEThCzLaNmyJTp27Kj1VHwG6/Lr9xepbS6q2qWh+2s1niy7jC5dnxW226leZZjjEN7KhHk5JXC4FI8WUwIMEowGCQtHJDAJG0BVu1T3uDJtUBQ8ExbCtUybvxKlahcFTiWivOH+jmX59eagwKkgu9Dq1VeaBSMSMJrRVxpyD6vy656gwHGWlWdhUjlqWko8Msy0ZaB3dIfOUXahVbgybURbFDhOzl+pw7ycEqZjsijTRrRFgeNk5s4zcDB+no9FmTaiLQocB+4ybawfoGVRpo1oiwLHgchl2oi2KHAciFymjWiLAseY6GXaiLYocIyJXqaNaIsCx5joZdqItihwjIlepo1oi941xkQv00a0RYFjTPQybURbFDgORC7TRrRFgeNA5DJtRFsUOA5ELtNGtEWB40TUMm1EWxQ4TkQt00a0RYHjSMQybURbVGJBBd7WNGFZpo1oiwKnEqraRQAKnOq0LtNGtEWB05AWZdqItihwhKiIVikJUREFjhAVUeAIUREFjhAVUeAIUREFjhAVUeAIUREFjhAVUeAIUREFjhAVUeAIUREFjhAV/T9HgXHIWcyUbAAAAABJRU5ErkJggg==",
      "text/plain": [
       "\u001b[1m<\u001b[0m\u001b[1;95mFigure\u001b[0m\u001b[39m size 20\u001b[0m\u001b[1;36m0x200\u001b[0m\u001b[39m with \u001b[0m\u001b[1;36m1\u001b[0m\u001b[39m Axes\u001b[0m\u001b[1m>\u001b[0m"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "ex_pat = list(patterns.values())[-2]\n",
    "nx.draw(ex_pat.graph)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "eab29362",
   "metadata": {},
   "source": [
    "## The Spasm Space\n",
    "\n",
    "In the basic usage of PACT for subgraph counting we also explicitly specify a *spasm space*. Formally this is a set of graphs $S$ such that $Spasm(P) \\subseteq S$ for each pattern $P$ that we want to count.\n",
    "\n",
    "Since our patterns are all the trees with 7 nodes, we know that the spasm space must contain all connected graphs with at most 7 nodes and 6 edges.\n",
    "\n",
    "The spasm space is stored in a special `SpasmSpace` object, which implements various techniques for optimising later computations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "098a4f5f",
   "metadata": {},
   "outputs": [],
   "source": [
    "from pact.spasmspace import SpasmSpace\n",
    "\n",
    "spasm_space = SpasmSpace()\n",
    "\n",
    "for sf in settings.raw_spasm_filenames():\n",
    "    with open(sf) as f:\n",
    "        lines = list(map(str.rstrip, f.readlines()))\n",
    "        spasm_space.add_from_g6lines(lines, settings.SPASM_GRAPHWRAPPER_PARAMS)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "dbbbb64a",
   "metadata": {},
   "source": [
    "We'll get an example graph `ex_graph` from the spasm space to illustrate some operations going forward.\n",
    "\n",
    "To get an interesting graph we can use the `iter_by_ev(x, y)` method which gives us an iterator over all graphs in the spasm space with *x* edgs and *y* vertices. We can then draw the graph using networkx by accessing the graph attribute, which stores the networkx graph inside our wrapped graph."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "58f515ec",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
      ],
      "text/plain": []
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAANwAAADcCAYAAAAbWs+BAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/H5lhTAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAjv0lEQVR4nO3deUBUVf8/8PedGRhkU0HAVEZRYJgBXDBETSoSfTKfh0JNVChT0xZMwzTNCtNHMbdEpZ/5GC6RSy7k11wSSR4wFdwFUUAeZQnRiUVFBgZm5v7+MEwFBpi5c+8w97z+qrkz53woPpzPnHvOuRRN0zQIgmCFgOsACIJPSMIRBItIwhEEi0jCEQSLSMIRBItIwhEEi0jCEQSLSMIRBItIwhEEi0jCEQSLSMIRBItIwhEEi0jCEQSLRGx0Uq1So6C8GnVqLSxFAvRytIGNmJWuCcKkGO23/sbdKuzIKEJKrgJFFUo8uQeIAiBxsEaQ1BnhARJ4uNgZKwyCMCkU0/vhiiuUWPhzFk7ml0EooKDRNt98w/VA9y6ICfWFq4M1k6EQhE5cVF6MJtzuc0VYdDAbai2tM9GeJRRQEAkoLA7xxgR/CVPhEEQjXFdejCVcXMoNrE7KM7iduSM9MTPIg4GICOJvplJ5MZJwu88VYUFiVqPX+/boiLF+PTCktyN6dO6ASmU9LhVVYs3xPNwqq262vRVjfBFGRjqCIaZUeRmccMUVSgSvTYVKrW107f9N8sPzPTvj8NVS5JRWwclOjMlDesLaUoTQjaeQd/dhk22KRQIkR71EvtMRBjO1ysvghHsrPgOnb5Y3+ZfDT9IZWSX3UK/5+1ovR2scm/0ijly9g6g9l5tsUyigMLS3IxKmBRgSGsFzzVVe1pZCvPdib/R37YR+PTqhk7Ul5u69gn0X/9DZHhOVl0E3vm/crcLJ/LJmh+mLRZVPJRsAFJQrkad4CHdn22bb1WhpnMwvQ76iypDwCB4rrlBi0cHsJq85WFti9nBP9HGyxfXS1v+ORR/MRnGF0qC4DEq4HRlFEAqoNn+ui60lKqvrdL5HKKDwY3qRvqERPLfw5yyomxkIFFUq+C9LxrCVKYg5er3Vbaq1NBb+3HjEbAuDEi4lV9GmL6EA8Eb/7niuYwf8knlb5/s0WhopeQpDwiN4qqXKq06jxZ8PVW1ul4nKS++Ee6hSo6iNw2sfJxssed0bFworsb+FehkAisqVqFap9Q2R4Cl9K6/WMLTy0jvhCsur0ZaxzclWjC2T/VFVq8YHOy6gNQMjDaCgvPnbBwTRFH0qr9YytPLSO+HqmrgN0Bw7sQjbpvjDvoMFJm89C0VV64fztvRDEPpUXm1lSOWld8JZilr3UbFIgO8nPw+3LjaYtv0c8hVN33sztB+CANpeeenDkMpL79/mXo42aKlKFlBA3MQB8JN0xoc7L+Ji0b029UH91Q9BtBZbFZG+/ei9NNpGLILEwRqFOobvL0bLMULeFcev3UWnDpZ4o3/3p64fuFyisw+JozXZN0e0CVsVkb79GPTbHCR1RkJGYbNfUOXP2QMARshdMELu0ui6roQTCigEeTobEh7BQw2VlzHLSkMqL4MSLjxAgm1nCpq9PmFzut5ta7Q0IgaTBcxE27Sm8gKAt4f0hL2VBVzsrQAAw2XO6Nrx0T9vP12AKh2TIoZUXgYlnIeLHQLduzS7llJftEaNTnUK2GgeAiC7wYm2aanyAoAZgb3Ro/Pfi+NH+TyHUT7PAQAOXCppNuEMrbyMultAXyKKxoNd86CquI2lS5figw8+gEhEvssRrXPjbhVGxKYZrf3kqBfh7qzfQGDwN0xXB2ssDvE2tJmnLAvti5wLpzBx4kTMnj0bgwYNwtmzZxntgzBfXSzV6FRbClrD7ColoYBCoHsXvZMNYOiYvAn+Eswd6clEU5g3UoowfwkcHBzw3Xff4fTp06BpGoMHD8b777+PyspKRvohzNPRo0fh4+ODon0rYCESMtq2SEAhJtTXoDYYm0OdGeSBr8f4QiwStHkdm1BAQSwSYMUYX0QGuT91bfDgwTh37hxiY2Oxc+dOSKVSbN++HQyffUS0cw8ePMD06dPx2muvwdfXF5lnUrA0tC+jfSwJ8TZ4U3S7OrWrtLQUc+bMwe7duxEYGIiNGzfC25vZcpZof1JSUjBlyhSUl5djzZo1mD59Oijq0R99pnZ8zxspbTQY6IPxhGvw+HSkPAWKyps4HcnRGkGezogYLGlzTZycnIzIyEjcvHkTUVFRiI6Ohq1t8xtaCfOkVCrx2WefYf369XjppZewdetWuLm5NXqfoWeaLAnxZuyMHaMl3JOMcf6fSqXCqlWrsGzZMjg5OWHdunV44403Hv9lI8zb6dOn8c4776C4uBhff/01PvroIwgEzX9DMqtTu7h08+ZNfPTRRzhy5AhGjx6NDRs2NPlXjjAPKpUKixYtwqpVq+Dv74/t27dDKpW2+vNPVl4FZdVP/YE2tPJqFdoMaLVaev/+/XSPHj1oKysreunSpXRtbS3XYREMu3DhAu3t7U1bWFjQy5cvp+vr6w1qz9q+M/356o30xcIK+mrJPfphrWHttYZZJFyDqqoqet68ebRIJKKlUimdnJzMdUgEA+rq6uivvvqKFolEdP/+/enMzEyD21QqlTQAOiEhgYEIW8+sEq5BVlYWPWzYMBoAPXHiRPr27dtch0ToKSsri/bz86OFQiEdHR1Nq1QqRtotKSmhAdCHDx9mpL3WMsvdnT4+PkhLS8PWrVtx/PhxeHl5YcOGDdBoNFyHRrSSRqPBypUrMXDgQNTU1CA9PR2LFy+GpaUlI+1XVFQAADp37sxIe63GanpzoLy8nJ4xYwYNgPbz86MzMjK4DoloQV5eHj1kyBCaoih67ty5dE1NDeN9pKWl0QDonJwcxtvWxSxHuCc5ODhg06ZNOHPmDLRaLQYPHowPPviALBEzQVqtFhs2bEC/fv2gUChw8uRJrFq1ClZWVoz31TDCOTg4MN62TqymN8fq6+vpdevW0XZ2drSTkxO9fft2WqvVch0WQdP0rVu36KCgIBoAHRkZST98+NCo/W3ZsoUGQNfV1Rm1n2eZ/Qj3JJFIhFmzZiEnJwfDhw/H5MmT8fLLLyM7u+kjsQnjo2ka33//PXx9fZGfn4/k5GTExcXBxsa4Z9lUVlbCzs4OFhYWRu3nWbxKuAbdunXDrl27kJSUhNLSUvTv3x/z589HdTU5A5NNt2/fxujRozF9+nSMHz8eWVlZGD58OCt9V1RUsF9OgqcJ12DEiBHIysrCokWLsH79esjlchw4cIDsRDAymqaxY8cOeHt749KlSzh06BDi4+PRsWNH1mKoqKhgf4YSPE84ABCLxfjiiy+QnZ0Nb29vhIaGIiQkBLdu3eI6NLOkUCgwbtw4REREYNSoUbh69SpGjx7NehyVlZVkhONS7969cfjwYezfvx+XL1+Gt7c3YmJioFK1/aEPRNMSExPh4+OD1NRU7NmzBzt37oSjoyMnsZCS0gRQFIUxY8bg+vXriIyMRHR0NPr164cTJ05wHVq7VllZiYiICIwdOxYvvPACsrOz8eabb3IaEykpTYitrS1WrVqFS5cuwcnJCcOHD0d4eDju3LnDdWjtTsORB4cOHcIPP/yAxMREuLg0PqOUbaSkNEG+vr5ITU3F1q1bkZSUBKlUiri4OLJErBWePfLg6tWreOutt0xmvyJXJSWvbnwboqys7KklYmfPnuU6JJN14sQJumfPnrSNjQ29adMmk1tcoNFoaIqi6P/85z+s901GuFZydHR8vERMo9EgICAAH374IVki9gSlUolZs2bhlVdeQa9evZCZmYkZM2aYzKjW4P79+6Bpmoxw7UV9fT0dGxtLlog94dSpU7SHhwdtZWVFx8bG0hqNhuuQmpWfn08DoE+cOMF632SE04NIJMLs2bORk5ODV155BZMnT0ZQUBCuXbvGdWisU6lUWLBgAQIDA+Hg4IDLly9j9uzZOs8X4RpnW3MAMsIxISkpifbw8KBFIhE9f/58oy+8NRXnz59n9MgDthw7dowGQBcWFrLeN0k4htTU1NBLliyhxWIxLZFI6AMHDnAdktHU1dXRixYtenzkwZUrV7gOqU127dpFA6CrqqpY75skHMPy8/PpUaNG0QDof/3rX/StW7e4DolRxjrygE3ffvstLRKJOPnebbqFdjvVp08fHD58GPv27cPFixchl8sRExODuro6rkMziEajwYoVK4x25AGbGm56czJ7ynqK88iDBw/oTz75hBYKhbSXlxcns2JMyM3NNfqRB2yaM2cO7eXlxUnfZIQzIjs7O6xevRqXLl2Co6MjXnnlFURERLSbJWINRx7079/f6EcesImrdZQAWdrFCl9fX6SlpWHLli349ddf4eXlhW+//dakl4gVFBQgODgYs2bNwtSpU3HlyhW88MILXIfFCK7WUQIk4VgjEAgwZcoU5ObmYvz48Zg5cyYCAgJw7tw5RtqvVqmRffs+LhVVIvv2fVTreEa1LjRHRx6wibN1lDCDZwu0V+np6Xj//feRmZmJ999/H8uWLWtzmfP4nPxcBYoqmnhCkYM1gqTOCA+QwMOl5XPyS0pKMH36dBw9ehRTp07FN998w+oubLb4+Phg+PDhWLduHet9kxGOI4MHD8b58+fxzTffICEhAV5eXkhISGjV8Q7FFUq8FZ+BEbFpSMgoROEzyQYANIDCCiUSMgoxIjYNb8VnoLhC2WR79F9HHvj4+HB25AGbSEnJUyKRCB9//DFycnIQFBSEt99+u8UlYrvPFSF4bSpO3ywHgBafd9Zw/fTNcgSvTcXuc0VPXVcoFBg7diznRx6wicuSkiScCejevTt2796NY8eOoaSkBP369cNnn33W6BSxuJQbWJCYBZVa26YHCwKPEk+l1mJBYhbiUm4AeHTkgbe3N9LS0jg/8oAtNTU1qK2tJbOUBDBy5EhkZWXhyy+/xNq1a+Ht7Y2DBw8CeDSytebRuZEvu6Ng+Wgcm/1is+9ZnZSH4PeiMXbsWAwbNswkjjxgC2cnLv+FJJyJsbKyQnR0NLKzsyGTyfD666/j1bHhiP6/qy1+tqu9FSKD+rQ4Q0nTNG7YD8C67380mSMP2NKwf5EkHPGUPn364MiRI9i3bx+u2faDqr7laf7PX5PhUtE9ZJXc1/k+iqIgsrDEObib3OZQY+N0aw5Iwpk0iqLQd9hICLp5gxIIdb53UC8HjPLpiiWHWrcnT0MDJ/PLkK+oYiLUdoOUlIROOzKKIBToHoUEFLA4xBs/nS9G7t3WJ5BQQOHH9KKW32hGGkpKMsIRTUrJVbQ4Ixke0BPdO3XAmuMtT6o8SaOlkZKnMCS8dqeiogJ2dnYQiUSc9E8SzoQ9VKlR1MzN6gadrC0wZ4Qn1p+4gYrqtm8BKipX6r0MrD3i8h4cQBLOpBWWVzdaQfKsuSOkuKesx/YzBXr1QQMoKOfPU4O4XGUCANyMq0Sr1Km1Oq/3crTGxEESLDl0DS52f2+ZEYsEEAkp9OjUAVUqNe7X1BvUjznhcmsOQBLOpFmKdBcgXe2tIBRQWBzijcUh3o2u/z7/FWw5davFmcuW+jEnXJeUJOFMWC9HG1BAs2Vl7t0qzEg43+j1T0ZIYSMWYsmhaygs1/0dkPqrH76orKyEm5sbZ/2ThDNhNmIRJA7WKGxm4qRSWY+ka3cbvT71BTcAlk1ee5bE0Ro2Yv78GnBdUvKnlmingqTOLd6H05dQQCHI09kobZsqUlISOoUHSLCtjTOQEzant+p9Gi2NiMESPaJqnzQaDe7fv09uCxDN83CxQ6B7F+ZHOa0G8s4U3J1b3gluLhoe4kFKSkKnmFBfiBhOOIrW4njMVISFhaGkpITRtk0V1+soAZJw7YKrg3WT0/6GWP6mH7Z9uwapqanw8vLC2rVroVab94oTrrfmACTh2o0J/hLMHenJSFvzRkoxwV+C8PBw5OTk4J133sHcuXMxcOBAnDp1ipE+TBHXW3MAknDtyswgD3w9xhdikaDN3+mEAgpikQArxvgiMsj98eudOnXChg0bcPbsWYjFYgwbNgzTpk1DWVkZ0+FzjpSURJtN8JcgOeolDO396OyRlhKv4frQ3o5IjnoJYf5Nz0oOHDgQZ86cwXfffYfExERIpVJs3rwZWq35LPuqrKyEhYUFp2dskoRrh1wdrJEwLQDHP34RbwX0RE9HazROOxo9Ha3xVkBPJEe9iIRpAXB1sNbZrlAoxHvvvYfc3FyEhIRgxowZGDp0KC5dumSsH4VVDTe9Od3lzskTDQjGPaytp6+W3KMvFlbQAa+Opd+cGGFwm2lpabSPjw8tEAjoWbNm0ffu3WMgUu5ERUVx9hCPBmSEMxM2YhG8u3XEAEln9O/ZBXnXsgxuMzAwEBcvXsTKlSsRHx8PLy8v7Nq1q1WH1ZoirrfmAKSkNEsymQw5OTmMPCzEwsICn3zyCXJycjBs2DBMmjQJwcHByMnJYSBSdnG9jhIgCWeWZDIZVCoVCgoKGGuzR48e2Lt3L44ePYrCwkL07dsXn3/+OZRK3bsRTAnX6ygBknBmSS6XA4DOI9P19eqrr+Lq1av4/PPPsWbNGsjlcvzyyy+M92MMpKQkjKJ79+6ws7PD9evXjdK+lZUVFi1ahKtXr8LLywshISF4/fXXGR1RjYGUlIRRUBQFmUxmtIRr4O7ujqNHj2Lfvn24cOHC4+eZq1Qqo/arD5qmSUlJGI9MJjNKSfksiqIwduxYXL9+HZGRkYiOjka/fv3w22+/Gb3vtqipqYFKpSIJRxhHwwjH1hS+nZ0dVq1ahcuXL8PJyQnBwcGYNGkSSktLWem/JVwfANuAJJyZksvlqKqqYn3rjY+PD9LS0rBt2zYkJyfDy8sL69ev53wngimsowRIwpktmUwGAEb/HtcUiqIwefJk5OTkYNKkSfj444/h7++P9PTW7UQ3BpJwhFG5ublBLBZzknANHBwcsHHjRqSnp0MgEGDIkCGYMWMGysvLWY+FlJSEUQmFQkilUlYmTloyaNAgnD17FnFxcdizZw+kUim2bNnC6k4EU9gLB5CEM2ts3BpoLaFQiMjISOTk5GDUqFGYNm0aAgMDkZmZyUr/FRUVsLe35+whHg1IwpkxuVxuMgnXoGvXrkhISEBKSgoqKyvh5+eHOXPmoKrKuM+pq6ys5Hx0A0jCmTWZTIY///zTJHdvv/zyy7h8+TJiYmKwadMmeHl5Yc+ePUa7jWEKN70BknBmjcuZytawtLTEp59+iuvXryMgIABhYWH4xz/+gby8tj3nrjVIwhFG5+HhAaFQaBITJ7pIJBIkJibi0KFDyM/Ph6+vL6Kjo1FTU2Nw29UqNbJv38dtlSXELn04fxYeRbfX3YREq0ilUowaNQqxsbFch9IqSqUSy5cvx8qVK9G9e3fExcXhtddea1MbN+5WYUdGEVJyFSiqUD71MBQKgMTBGkFSZ4QHSODhwu5BuCThzFxoaCiUSiWOHTvGdShtkpeXh8jISCQnJyM0NBSxsbGQSHQfy15cocTCn7NwMr8MQgGl81HNDdcD3bsgJtS3xfNemEJKSjPH1iJmpnl6eiIpKQk//fQT0tPTIZPJsHLlStTVNf1Y5d3nihC8NhWnbz66qd7Sc9Ebrp++WY7gtanYfa6I2R+gGSThzJxMJsMff/yBBw8ecB1Km1EUhfHjxyMnJwfvvfceFi5ciAEDBiA1NfWp98Wl3MCCxCyo1NoWE+1ZGi0NlVqLBYlZiEu5wWT4TSIlpZm7cOECnn/+eWRkZGDQoEFch2OQK1eu4IMPPsCZM2cQERGB1atXI6VIhQWJjQ9M8nC2xcfBnvDt3hFOtmLU1GtwQ1GF/6TdxG85imb7WDHGt9mzO5lARjgz5+XlBcB0bw20Rb9+/fD7778jPj4eR48ehez5Yfgi8UqT7+3euQNsxULsv/gHFh/KxoYTj0av+Mn+mOjv2mwf0QezUdzMAzCZQEY4HujVqxcmTJiAr7/+mutQGFNeXo6RMQdRJnQAJWzdci0BBRyaOQxikRDD16Y2+R6hgMLQ3o5ImBbAZLh/x2CUVgmT0l4nTnSpUFui3NK51ckGAFoauH2/FvYdmv+MRkvjZH4Z8hXGWWpGEo4HTGkRM1N2ZBS16oEmHSyE6GxtAYmDNaa94IaXPZ1w6n+6twcJBRR+TDfOrCV55DAPyOVyrFu3DrW1tbCysuI6HEak5CpaNSP5xWgZwgN6Ang0ev2afQfR/3dV52c0WhopeQp8BWafyQeQhOMFmUwGrVaLvLw89O3bl+twDPZQpUZRKyc2tpy6hSNZpXCxt8Lovs9BSFGwFLVc2BWVK1GtUsNGzGyKkJKSB0x9EXNbFZZXo7Uzff/7sxqn/leOxEslmLb9PKzFQnz/tn+Ln6MBFJRXGxRnU0jC8YCDgwNcXFzMZuKkTq3/TvGjWaXo79oJvbu0/Iw4Q/ppDkk4njCniZPWlITNsbIQAgDsrFouFQ3ppzkk4XhCLpebzQjXy9GmiQdQPs3RxrLRayIBhTF+PVBTp8ENxUOdn6f+6odpZNKEJ2QyGTZv3gy1Ws35uR6GshGLIHGwRqGOiZOYUF/YikU4W1CBO/dr4WQnxhv9u8Pd2Rb/PnwNyjrdj/KSOFozPmECkBGON+RyOerr63Hz5k2uQ2FEkNRZ5324Q5m3oaVpRARIsPQNH7w7zA2l92vw7g/nEP/7LZ1tCwUUgjydmQ4ZABnheKNhpvLatWvw9PTkOBrDhQdIsO1MQbPXf8ksxS+Z+h2zrtHSiBhsnAXMZITjia5du6Jjx45mM3Hi4WKHQPcurVpt0hZCAYVA9y5wdzbOTnCScDxBUZRJHptniJhQX4gYTjiRgEJMqC+jbT6JJByPmNsiZlcHaywOYXb51ZIQb6Met0ASjkfkcjlycnJYPWLc2Cb4SzDrZbdH/2LgTrN5I6VG3XwKkITjFZlMhurqahQXF3MdCmNomsaVnStQ9dt3sBQJ2vydTiigIBYJsGKMLyKD3I0U5d9IwvGIua2pBID4+Hj88MMPiP1oPH6b8zKG9nYEgBYTr+H60N6OSI56yegjWwOy45tHtFotbG1tsWzZMkRFRXEdjsEuX76MwYMHY/Lkydi0adPj1x+fS5mnQFF5E+dSOlojyNMZEYMlRpuNbA5JOJ7x8/PDwIEDsXnzZq5DMci9e/fw/PPPw97eHqdPn252n1+1So2C8mrUqbWwFAnQy9HGKCtIWovc+OYZc1jETNM0pkyZgrKyMiQlJencVGsjFsG7W0cWo9ONfIfjmYZFzO25sFm7di0OHDiA7du3o3fv3lyH0yYk4XhGJpOhsrISCkXzZzOaslOnTuHTTz/FvHnz8Prrr3MdTpuRhOMZuVwOoH3OVCoUCoSFhWHIkCFYtmwZ1+HohSQcz/Tp0wcikajdrTjRaDQIDw9HfX09fvrpJ1hYWHAdkl7IpAnPWFhYwMPDo92NcEuWLMGJEyeQlJSEbt26cR2O3kjC8VB72/197Ngx/Pvf/8aSJUswfPhwrsMxCCkpeag93RooLi5GeHg4Xn31VSxcuJDrcAxGEo6H5HI5SktLce/ePa5D0amurg7jx4+HtbU1EhISIBC0/1/X9v8TEG3WXtZUzp8/HxcuXMDevXvh6OjIdTiMIAnHQ1KpFBRFmXTC7du3D7GxsVizZg0CAozzJBsukITjoQ4dOsDNzc1kJ07y8vIwdepUjB8/HjNnzuQ6HEaRxcs89c9//hM0TePw4cNch/IUpVKJIUOGoLa2FufPn4edHbur+Y2NjHA8Zarnm8ycORM3btzAvn37zC7ZAJJwvCWTyVBQUACl0niP122rLVu2YOvWrdi4cSN8fY13kA+XSMLxlEwmA03TyM3N5ToUAMCVK1cQGRmJd999F5MnT+Y6HKMhCcdTpnRr4P79+xg3bhy8vLywfv16rsMxKrK0i6c6duyIbt26cT5TSdM0pk2bBoVCgV9//RUdOnTgNB5jIwnHY6YwcbJu3Trs378fiYmJ6NOnD6exsIGUlDzG9cGwp0+fxrx58/DJJ58gNDSUszjYRBKOx2QyGfLz81FfX89633/++SfCwsIQEBCA5cuXs94/V0jC8ZhcLodarUZ+fj6r/Wo0GkREREClUrXrzaT6IAnHY08+wopNS5cuxfHjx7Fz5050796d1b65RhKOx5ycnODo6MjqxMnx48exePFifPXVVwgODmatX1NB1lLyXGBgIFxdXbFz506j9/XHH39gwIABGDhwII4cOWIW+9vain8/MfEUtnZ/19fXIywsDFZWVvjxxx95mWwASTjea3iElUaj+yHzhlqwYAHOnj2LPXv2oEuXLkbty5SRhOM5mUyG2tpaFBYWGq2PxMREfPPNN1i9ejWGDBlitH7aA5JwPGfsNZX5+fmYMmUKxo0bh1mzZhmlj/aEJBzPubq6wtbW1ii3BmpqajBu3Di4uLggPj4eFMXs87jbI7KWkucoioKXl5dRRriPPvoIubm5yMjIgL29PePtt0ck4QijLGLetm0b4uPjsXXrVvTt25fRttszUlISjxcxM3VLNjMzEx9++CGmTp2Kd955h5E2zQVJOAIymQwPHjxAaWmpwW09ePAA48aNg4eHB+Li4hiIzryQkpJ4/Aira9euGfSgDJqm8e677+LOnTu4cOGC2W8m1QcZ4Qi4ubnB0tLS4O9xGzZswN69e7F161Z4eHgwFJ15IQlHQCQSQSqVGpRw6enpmDt3LqKiojB27FgGozMvZPEyAQAICwvD3bt38d///rfNny0rK4Ofnx969OiB1NRUXu1vaysywhEA9F/ErNVqERERgZqaGt5tJtUHSTgCwKOJE4VCgfLy8jZ9LiYmBklJSdixYwdcXV2NFJ35IAlHANBvTWVycjKio6MRHR2NkSNHGis0s0ISjgAAeHp6QiAQtDrhSkpKMGnSJAQHB+PLL780cnTmg9yHIwAAYrEYffr0Qea1XGTfvo86tRaWIgF6OdrARvz0r0l9fT0mTJgAS0tL7NixA0KhkKOo2x8yS0ngxt0q7Mgows7UK6izsAeeWNVPAZA4WCNI6ozwAAk8XOwwb948xMbGIjU1FUOHDuUu8HaIJByPFVcosfDnLJzML4NQQEGjbf5XoeG61F6LEyumY+Wi+YiKimIxWvNASkqe2n2uCIsOZkP9V5LpSrYnr+fco9HjvU3o+sIAo8dojsgIx0NxKTewOinP4HbmjvTEzCCyhKstyAjHM7vPFTWbbJZCAeaM8ETogO7o2MECOXceYHVSHn7PL2vy/auT8uBkK0aYv8SYIZsVcluAR4orlFh0MLvZ66vf7Itpw9xw4HIJFv+SDY0W2PqOP57v2bnZz0QfzEZxhek8RdXUkYTjkYU/Zz3+zvasfj06IqRfd6w8lovlR3Ow61wxJn2fjpJ7NfhslKzZNtVaGgt/zjJWyGaHJBxP3LhbhZP5Zc1OjozyeQ5qjRa7zhY9fk2l1mLPuWIM7NkZz3W0avJzGi2Nk/llyFdUGSVuc0MSjid2ZBRBKGj+1Czvbva4VVaNhyr1U69f/uMeAED+XPOHAAkFFH5ML2r2OvE3knA8kZKr0Dn172xnBUWVqtHrDa+52Dc9wgGPRrmUPIXhQfIASTgeeKhSo6iFiQ0rCwHqNNpGr6vqNY+v61JUrkT1M6Mj0RhJOB4oLK9GSzdba+u1sBQ2/nUQWwgfX9eFBlBQXq1nhPxBEo4H6tS6kwUAFFW1cLYTN3q94bW7D2oZ6YfvSMLxgKWo5f/N10ofwK2LDWyf2RnQ37XT4+tM9MN35L8QD/RytEFLp/ofvXoHIqEAEwf9vWrEUijAmwNdcamoEqX3dY9w1F/9ELqRpV08YCMWQeJgjUIdEyeXi+/hUOZtfPoPKRxtLFFYXo2xfj3Qo3MHzN+f2WIfEkfrRvvmiMbIfyGeCJI6IyGjUOetgU/2XkHJvRqM+Wst5fU7VZi2/RzOFlTobFsooBDk6cx0yGaJ7BbgiRt3qzAiNs1o7SdHvQh3ZzujtW8uyHc4nvBwsUOgexedq030IRRQCHTvQpKtlUjC8UhMqC9EDCecSEAhJtSX0TbNGUk4HnF1sMbiEG9G21wS4g1XB2tG2zRnJOF4ZoK/BHNHejLS1ryRUrL5tI3IpAlPPXmmSUvnmTxJKKAgElBYEuJNkk0PJOF4TJ9TuwLduyAm1JeUkXoiCUc8PpcyJU+BonLlUwudKTy6qR3k6YyIwRIyG2kgknDEU6pVahSUV+s8eZnQH0k4gmARmaUkCBaRhCMIFpGEIwgWkYQjCBaRhCMIFpGEIwgWkYQjCBaRhCMIFpGEIwgWkYQjCBaRhCMIFpGEIwgW/X9TahbmILJRXgAAAABJRU5ErkJggg==",
      "text/plain": [
       "\u001b[1m<\u001b[0m\u001b[1;95mFigure\u001b[0m\u001b[39m size 20\u001b[0m\u001b[1;36m0x200\u001b[0m\u001b[39m with \u001b[0m\u001b[1;36m1\u001b[0m\u001b[39m Axes\u001b[0m\u001b[1m>\u001b[0m"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "ex_graph = next(spasm_space.iter_by_ev(6,5))\n",
    "nx.draw(ex_graph.graph, with_labels=True, font_color='w')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "94fe069e",
   "metadata": {},
   "source": [
    "## Computing the Homomorphism Base of Patterns\n",
    "\n",
    "The homomorphism base of a graph $P$ are those graphs in the spasm space for which the coefficient $\\mu_F$ in the sum above is not zero. This can be quite costly to compute for larger patterns, but it is independent of the host graph and can therefore be pre-computed once.\n",
    "\n",
    "The `hombase_coeffs` function is the basic way to compute the homomorphism base. The function takes as a second parameter a spasm space that must contain all possible graphs in the spasm of the pattern $P$. The homomorphism base is returned as a dictionary with elements of the form $id(F): \\mu(F,P)$ for each graph $F$ in the spasm of $P$. The id's refer to the id of $F$ in the given spasm space.\n",
    "\n",
    "Since we will have to deal with very large numbers of homomorphisms, high precision arithmetic is important here. For this we use the [gmpy2](https://gmpy2.readthedocs.io/en/latest/index.html) package which provides python bindings for the [GNU Multiprecission](https://gmplib.org/) library. In particular, we use the `mpq` type for rational numbers as the coefficients are guaranteed to be rational."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "0ff80895",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
      ],
      "text/plain": []
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "\n",
       "\u001b[1m{\u001b[0m\n",
       "    \u001b[1;36m244471403881519891900181704987527009525\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-3\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m200805017286104562580457396622371148304\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m35\u001b[0m,\u001b[1;36m2\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m280376347674979585194257264981957219449\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m10\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m258662326575684880685635720424942876393\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-8\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m146011057545644161557291642472811544920\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-3\u001b[0m,\u001b[1;36m2\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m245413697315286802743807153654882190908\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-79\u001b[0m,\u001b[1;36m4\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m273134957091060411677327821996853583980\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-9\u001b[0m,\u001b[1;36m2\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m66719643163285929538505745610842253983\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-27\u001b[0m,\u001b[1;36m4\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m84273283320066270384244169320407146834\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m3\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m193101649776123763688982186209706020543\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m3\u001b[0m,\u001b[1;36m2\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m69590959277943715886386460662349719302\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m9\u001b[0m,\u001b[1;36m2\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m194004892045500727445598428862645766084\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m5\u001b[0m,\u001b[1;36m4\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m55146198916550377801394465856303775856\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m2\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m54695293067644195284398106656028363738\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m1\u001b[0m,\u001b[1;36m2\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m241758357588706690769845049306420308138\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m9\u001b[0m,\u001b[1;36m2\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m24883373538079352055065965295148335320\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m7\u001b[0m,\u001b[1;36m4\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m183865322710435675491774000165040442587\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m1\u001b[0m,\u001b[1;36m2\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m134672613886082743296041172130123887562\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-5\u001b[0m,\u001b[1;36m4\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m249250281349047097037132682571385289116\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-1\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m74311936876948642230131623326800756430\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-1\u001b[0m,\u001b[1;36m2\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m100146012587478022275270086846012758251\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-1\u001b[0m,\u001b[1;36m4\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m333447451475214955032543438464829452839\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-1\u001b[0m,\u001b[1;36m4\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m219831900909649418933876728263860555244\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-1\u001b[0m,\u001b[1;36m2\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m12314306155591027532348329433199392317\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m1\u001b[0m,\u001b[1;36m4\u001b[0m\u001b[1m)\u001b[0m\n",
       "\u001b[1m}\u001b[0m"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from pact.hombase import hombase_coeffs\n",
    "\n",
    "hombase_coeffs(ex_pat, spasm_space)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "af56ff19",
   "metadata": {},
   "source": [
    "The standard pattern for computing all homomorphism bases of our patterns is the following.\n",
    "\n",
    "The iteration here is wrapped inside a `rich.Progress` track as computation can take a while for large (or many) patterns and a progress indicator is helpful.\n",
    "\n",
    "Homomorphism bases are stored in the `hombase` attribute of a `GraphWrapper`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "05336354",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "63af2cae8d26462495cd1b9987c5d92a",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "Output()"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
      ],
      "text/plain": []
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\">\n",
       "</pre>\n"
      ],
      "text/plain": [
       "\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from pact.ui import default_progressbar\n",
    "\n",
    "with default_progressbar() as progress:\n",
    "    track = progress.track(patterns.values(), total=len(patterns))\n",
    "    for pattern in track:\n",
    "        pattern.hombase = hombase_coeffs(pattern, spasm_space)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9b6ab694",
   "metadata": {},
   "source": [
    "## Precompute Decompositions for the Spasm Space\n",
    "\n",
    "As outlines above, we will count the subgraph occurrences of our patterns by counting the number of homomorphisms from graphs in the spasm space. \n",
    "\n",
    "To count the homomorphisms effectively we use dynamic programming on tree decompositions. This is again independent of the host graph and can be done in another pre-computation phase (note that there is no dependence on the homomorphism base computation above).\n",
    "\n",
    "There are some technical details involved in this process that will be documented separately in the future. For now, this standard pattern will suffice to compute decompositions. (The `times` parameter scales the number of attempts to find a good decomposition. Better decompositions will lead to faster homomorphism counting.)\n",
    "\n",
    "Decompositions are stored in the `td` attribute of a `GraphWrapper`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "db358cda",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "53b543e1e7ba41eab16774cd1c19ad78",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "Output()"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
      ],
      "text/plain": []
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\">\n",
       "</pre>\n"
      ],
      "text/plain": [
       "\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from pact.balgowrapper import balgo_multitry_for_cheapest_decomp\n",
    "\n",
    "with default_progressbar() as progress:\n",
    "    track = progress.track(spasm_space.graphs_iter(), len(spasm_space))\n",
    "    for G in track:\n",
    "        decomp, badness = balgo_multitry_for_cheapest_decomp(G, times=2)\n",
    "        G.td = decomp\n",
    "        G.td_badness = badness"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "50ce4138",
   "metadata": {},
   "source": [
    "Decompositions support some basic pretty printing for inspection."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "87f0fb93",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
      ],
      "text/plain": []
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "\n",
       "B: \u001b[1m{\u001b[0m\u001b[1;36m0\u001b[0m, \u001b[1;36m3\u001b[0m, \u001b[1;36m4\u001b[0m\u001b[1m}\u001b[0m    λ: \u001b[1m[\u001b[0m\u001b[32m'E_5'\u001b[0m, \u001b[32m'E_0'\u001b[0m\u001b[1m]\u001b[0m\n",
       "  B: \u001b[1m{\u001b[0m\u001b[1;36m2\u001b[0m, \u001b[1;36m4\u001b[0m\u001b[1m}\u001b[0m     λ: \u001b[1m[\u001b[0m\u001b[32m'E_4'\u001b[0m\u001b[1m]\u001b[0m\n",
       "  B: \u001b[1m{\u001b[0m\u001b[1;36m0\u001b[0m, \u001b[1;36m4\u001b[0m\u001b[1m}\u001b[0m     λ: \u001b[1m[\u001b[0m\u001b[32m'E_1'\u001b[0m\u001b[1m]\u001b[0m\n",
       "  B: \u001b[1m{\u001b[0m\u001b[1;36m1\u001b[0m, \u001b[1;36m3\u001b[0m, \u001b[1;36m4\u001b[0m\u001b[1m}\u001b[0m  λ: \u001b[1m[\u001b[0m\u001b[32m'E_3'\u001b[0m, \u001b[32m'E_2'\u001b[0m\u001b[1m]\u001b[0m"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ex_graph.td"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "035be344",
   "metadata": {},
   "source": [
    "## From Decompositions to Homomorphism Counting Plans\n",
    "\n",
    "In a final preperatory phase we translate decompositions into execution plans. These plans abstractly describe the steps that should be taken to compute homomorphism count. \n",
    "\n",
    "Plans are expected to be stored in the `plan` attribute of `GraphWrapper` objects."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "bbd067f7",
   "metadata": {},
   "outputs": [],
   "source": [
    "from pact.planner import node_to_ops\n",
    "\n",
    "for G in spasm_space.graphs_iter():\n",
    "    G.plan = node_to_ops(G.td)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9513f2fc",
   "metadata": {},
   "source": [
    "Plans support human readable ouput. Technically, a plan is simply a `deque` of `code.Operation` objects. The plan is to be read in sequence but ultimately it is up to the execution engine how to execute it. A more detailed explanation of the individual operations will follow."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "eacfa951",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
      ],
      "text/plain": []
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "\n",
       "\u001b[1;35mdeque\u001b[0m\u001b[1m(\u001b[0m\u001b[1m[\u001b[0m\n",
       "    ⍴: _edge_base -> \u001b[1;35mE_5\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m3\u001b[0m, \u001b[1;36m4\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    ⍴: _edge_base -> \u001b[1;35mE_0\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m0\u001b[0m, \u001b[1;36m3\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    node$\u001b[1;36m0\u001b[0m = E_5 ⨝ E_0 \u001b[1m(\u001b[0mon \u001b[1m{\u001b[0m\u001b[1;36m3\u001b[0m\u001b[1m}\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    node$\u001b[1;36m0\u001b[0m = project node$\u001b[1;36m0\u001b[0m to \u001b[1m{\u001b[0m\u001b[1;36m0\u001b[0m, \u001b[1;36m3\u001b[0m, \u001b[1;36m4\u001b[0m\u001b[1m}\u001b[0m,\n",
       "    ⍴: _edge_base -> node$\u001b[1;35m1\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m2\u001b[0m, \u001b[1;36m4\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    ⍴: _edge_base -> node$\u001b[1;35m2\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m0\u001b[0m, \u001b[1;36m4\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    ⍴: _edge_base -> \u001b[1;35mE_3\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m1\u001b[0m, \u001b[1;36m4\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    ⍴: _edge_base -> \u001b[1;35mE_2\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m1\u001b[0m, \u001b[1;36m3\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    node$\u001b[1;36m3\u001b[0m = E_3 ⨝ E_2 \u001b[1m(\u001b[0mon \u001b[1m{\u001b[0m\u001b[1;36m1\u001b[0m\u001b[1m}\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    node$\u001b[1;36m3\u001b[0m = project node$\u001b[1;36m3\u001b[0m to \u001b[1m{\u001b[0m\u001b[1;36m1\u001b[0m, \u001b[1;36m3\u001b[0m, \u001b[1;36m4\u001b[0m\u001b[1m}\u001b[0m,\n",
       "    node$\u001b[1;36m0\u001b[0m = node$\u001b[1;36m0\u001b[0m ⋉ node$\u001b[1;36m2\u001b[0m \u001b[1m(\u001b[0mon \u001b[1m{\u001b[0m\u001b[1;36m0\u001b[0m, \u001b[1;36m4\u001b[0m\u001b[1m}\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    node$\u001b[1;36m1\u001b[0m = Count extensions of \u001b[1m{\u001b[0m\u001b[1;36m4\u001b[0m\u001b[1m}\u001b[0m in node$\u001b[1;36m1\u001b[0m,\n",
       "    node$\u001b[1;36m0\u001b[0m = Merge count from node$\u001b[1;35m1\u001b[0m\u001b[1m(\u001b[0m\u001b[1m{\u001b[0m\u001b[1;36m4\u001b[0m\u001b[1m}\u001b[0m\u001b[1m)\u001b[0m into node$\u001b[1;36m0\u001b[0m,\n",
       "    node$\u001b[1;36m3\u001b[0m = Count extensions of \u001b[1m{\u001b[0m\u001b[1;36m3\u001b[0m, \u001b[1;36m4\u001b[0m\u001b[1m}\u001b[0m in node$\u001b[1;36m3\u001b[0m,\n",
       "    node$\u001b[1;36m0\u001b[0m = Merge count from node$\u001b[1;35m3\u001b[0m\u001b[1m(\u001b[0m\u001b[1m{\u001b[0m\u001b[1;36m3\u001b[0m, \u001b[1;36m4\u001b[0m\u001b[1m}\u001b[0m\u001b[1m)\u001b[0m into node$\u001b[1;36m0\u001b[0m\n",
       "\u001b[1m]\u001b[0m\u001b[1m)\u001b[0m"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ex_graph.plan"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ad76458d",
   "metadata": {},
   "source": [
    "# Part 2: Execution on Host Graphs\n",
    "\n",
    "With everything now in place, all that is left is to actually compute the $\\#Hom(F \\rightarrow H)$ terms of our original formula for $\\#Sub(P \\rightarrow H)$ for specific host graphs $H$."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7a8df0f0",
   "metadata": {},
   "source": [
    "## Loading the Host Graph\n",
    "\n",
    "The required format for the host graph depends on which execution engine is used. In this tutorial we will use the *naive pandas execution engine*, which takes as input a `pandas.DataFrame` that represents the edge relation of our host graph. The dataframe must have attributes `s` and `t`, source and target, respectively. \n",
    "\n",
    "Note that an edge $\\{a,b\\}$ in an undirected graph must have two entries in the relation, one $(a,b)$ and one $(b,a)$. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "df694045",
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "\n",
    "host_nx = nx.read_gml('data/tutorial/schi_pombe.gml', label=None)\n",
    "\n",
    "# IMPORTANT: for undirected graphs we need both directions in our edge relation\n",
    "edges = [{'s': a, 't': b} for a,b in host_nx.edges()] + [{'s': b, 't': a} for a,b in host_nx.edges()]\n",
    "host_df = pd.DataFrame(edges).drop_duplicates()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d8c555dd",
   "metadata": {},
   "source": [
    "## Executing Counting Plans\n",
    "\n",
    "The plan of a graph G form a list of instructions which can be executed to count the number of homomorphisms from G into the host graph. In particular, for a `GraphWrapper` object with a plan, we can execute the plan on a given host using different execution backends.\n",
    "\n",
    "At the moment only one such backend exists: the naive pandas backend executes the instructions directly as given using pandas for the heavy lifting. The homomorphism count can be computed via this backend with the `naive_pandas_homcount(G, H)` function, which directly returns $\\#Hom(G \\rightarrow H)$. Note that as described above $H$ is passed as an edge relation in a `pandas.DataFrame`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "95b2aa56",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
      ],
      "text/plain": []
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "\u001b[1;36m454750\u001b[0m"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from pact.naive_exec import naive_pandas_homcount\n",
    "\n",
    "naive_pandas_homcount(ex_graph, host_df)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "08dcf356",
   "metadata": {},
   "source": [
    "## Counting a Pattern\n",
    "\n",
    "We finally arrive at counting the number of times a pattern occurs as a subgraph of the host. Let us recall our original formula for computing the subgraph count.\n",
    "$$\\#Sub(P \\rightarrow H) = \\sum_{F \\in Spasm(P)} \\mu(F,P) \\cdot \\#Hom(F \\rightarrow H)$$\n",
    "In the previous steps we have computed $Spasm(P)$ as well as all $\\mu$ coefficiens, as well as tree decomposition based execution plans for counting $\\#Hom(F \\rightarrow H)$ for all graphs that occur in the spasms of our patterns.\n",
    "\n",
    "What is left is to apply the plans to compute the values of $\\#Hom(F \\rightarrow H)$ for concrete host graphs $H$.\n",
    "The code below demonstrates how to gather the (non-zero) terms of the above sum for a pattern, and consequently compute the full subgraph count, for our example pattern `ex_pat`. \n",
    "\n",
    "Note that because our coefficients use the `gmpy2.mpq` type for precision, this produces `mpq` values which are dispayed as (numerator, denominator) pairs."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "f5312a41",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "#Sub(P -> H): 62638680\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
      ],
      "text/plain": []
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "\n",
       "\u001b[1m{\u001b[0m\n",
       "    \u001b[1;36m244471403881519891900181704987527009525\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-9246\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m200805017286104562580457396622371148304\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m73815\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m280376347674979585194257264981957219449\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m238400\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m258662326575684880685635720424942876393\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-175888\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m146011057545644161557291642472811544920\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-109845\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m245413697315286802743807153654882190908\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-1094229\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m273134957091060411677327821996853583980\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-1738197\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m66719643163285929538505745610842253983\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-2348271\u001b[0m,\u001b[1;36m2\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m84273283320066270384244169320407146834\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m1534746\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m193101649776123763688982186209706020543\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m2688111\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m69590959277943715886386460662349719302\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m3754314\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m194004892045500727445598428862645766084\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m2565735\u001b[0m,\u001b[1;36m2\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m55146198916550377801394465856303775856\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m909500\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m54695293067644195284398106656028363738\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m111938\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m241758357588706690769845049306420308138\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m12445155\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m24883373538079352055065965295148335320\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m2177574\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m183865322710435675491774000165040442587\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m4875256\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m134672613886082743296041172130123887562\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-60067845\u001b[0m,\u001b[1;36m2\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m249250281349047097037132682571385289116\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-10527076\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m74311936876948642230131623326800756430\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-10360537\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m100146012587478022275270086846012758251\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-12737789\u001b[0m,\u001b[1;36m2\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m333447451475214955032543438464829452839\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-22169881\u001b[0m,\u001b[1;36m2\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m219831900909649418933876728263860555244\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m-35095743\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m12314306155591027532348329433199392317\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m280639315\u001b[0m,\u001b[1;36m2\u001b[0m\u001b[1m)\u001b[0m\n",
       "\u001b[1m}\u001b[0m"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "hombase_pattern_terms = dict()\n",
    "\n",
    "for spasm_id, coeff in ex_pat.hombase.items():\n",
    "    spasm_graph = spasm_space[spasm_id]\n",
    "    spasm_count = naive_pandas_homcount(spasm_graph, host_df)\n",
    "    hombase_pattern_terms[spasm_id] = coeff * spasm_count\n",
    "\n",
    "print('#Sub(P -> H):', sum(hombase_pattern_terms.values()))\n",
    "\n",
    "# we also display the individual terms of the sum for each F in Spasm(P)\n",
    "hombase_pattern_terms"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "019f9894",
   "metadata": {},
   "source": [
    "## Counting all Patterns\n",
    "\n",
    "When we want to count the subgraph occurrences for a number of different patterns we can be a bit smarter.\n",
    "In settings such as our example, many graphs occur in the homomorphism base of multiple patterns. \n",
    "\n",
    "Rather than counting their homomorphism counts repeatedly for each pattern, we can instead count the homomorphisms once for every graph in the spasm space and compute the patterns counts all at once.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "2bec5099",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "280b07e150e44eb089528c1174a855c7",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "Output()"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
      ],
      "text/plain": []
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\">\n",
       "</pre>\n"
      ],
      "text/plain": [
       "\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
      ],
      "text/plain": []
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "\n",
       "\u001b[1m{\u001b[0m\n",
       "    \u001b[1;36m190384357088276965136693023116342970814\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m68828058\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m87968536505999792063251622645303437152\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m109656596\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m116130532718508524490079362148244371643\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m30382730\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m304019808730398462645962837222145370967\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m53045833\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m169753036403415102524309189092001454411\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m82915983\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m17902548692492397906272402987230910192\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m56009486\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m144812595679229883217623805773300660911\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m50062057\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m129632563928494135457881829008522067115\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m93357910\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m245024156124066838666528091004667761126\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m14232418\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m14689815793139536459496296402553498920\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m62638680\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m,\n",
       "    \u001b[1;36m98135613334211697414771868798897728850\u001b[0m: \u001b[1;35mmpq\u001b[0m\u001b[1m(\u001b[0m\u001b[1;36m9372935\u001b[0m,\u001b[1;36m1\u001b[0m\u001b[1m)\u001b[0m\n",
       "\u001b[1m}\u001b[0m"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pattern_counts = {pattern_id: 0 for pattern_id in patterns.keys()}\n",
    "\n",
    "with default_progressbar() as progress:\n",
    "    track = progress.track(spasm_space.graphs_iter(), total=len(spasm_space))\n",
    "    for spasm_graph in track:\n",
    "        hom_count = naive_pandas_homcount(spasm_graph, host_df)\n",
    "\n",
    "        for pattern in patterns.values():\n",
    "            if spasm_graph.id in pattern.hombase:\n",
    "                coeff = pattern.hombase[spasm_graph.id]\n",
    "                to_add = hom_count * coeff\n",
    "                pattern_counts[pattern.id] += to_add\n",
    "\n",
    "pattern_counts"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d99203d1",
   "metadata": {},
   "source": [
    "## Plot the Most Common Patterns\n",
    "\n",
    "From here on you can analyse our results with whichever tools you prefer. Here is an example of plotting the top counts and plotting the top 3 patterns."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "872519c7",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
      ],
      "text/plain": []
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>pattern</th>\n",
       "      <th>#Sub</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>87968536505999792063251622645303437152</td>\n",
       "      <td>109656596</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>129632563928494135457881829008522067115</td>\n",
       "      <td>93357910</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>169753036403415102524309189092001454411</td>\n",
       "      <td>82915983</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>190384357088276965136693023116342970814</td>\n",
       "      <td>68828058</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>14689815793139536459496296402553498920</td>\n",
       "      <td>62638680</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "\n",
       "                                   pattern       #Sub\n",
       "\u001b[1;36m1\u001b[0m   \u001b[1;36m87968536505999792063251622645303437152\u001b[0m  \u001b[1;36m109656596\u001b[0m\n",
       "\u001b[1;36m7\u001b[0m  \u001b[1;36m129632563928494135457881829008522067115\u001b[0m   \u001b[1;36m93357910\u001b[0m\n",
       "\u001b[1;36m4\u001b[0m  \u001b[1;36m169753036403415102524309189092001454411\u001b[0m   \u001b[1;36m82915983\u001b[0m\n",
       "\u001b[1;36m0\u001b[0m  \u001b[1;36m190384357088276965136693023116342970814\u001b[0m   \u001b[1;36m68828058\u001b[0m\n",
       "\u001b[1;36m9\u001b[0m   \u001b[1;36m14689815793139536459496296402553498920\u001b[0m   \u001b[1;36m62638680\u001b[0m"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "result_df = pd.DataFrame({'pattern': pid, '#Sub': int(c)} for pid, c in pattern_counts.items())\n",
    "\n",
    "# show the top patterns\n",
    "result_df.sort_values('#Sub', ascending=False).head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "id": "2e594ea6",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
      ],
      "text/plain": []
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABW0AAAErCAYAAACl2aHRAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/H5lhTAAAACXBIWXMAAA9hAAAPYQGoP6dpAACe9ElEQVR4nOzdd1hTZxsG8DsBmSIKihtxIFW0TsQJ4kIFcaHiXq2tW6vWrWito9pqK45aZxFXVRT3xIkijlbFhVUJTgQUERAMOd8ffkmJbEhIQu7fdXm1niTveYLJwznPec/zigRBEEBEREREREREREREWkGs6QCIiIiIiIiIiIiI6D8s2hIRERERERERERFpERZtiYiIiIiIiIiIiLQIi7ZEREREREREREREWoRFWyIiIiIiIiIiIiItwqItERERERERERERkRZh0ZaIiIiIiIiIiIhIi7BoS0RERERERERERKRFWLQlIiIiIiIiIiIi0iIs2hIRERERERERERFpERZtiYiIiIo4kUiUqz9nzpzRWCyLFy9W+74zs3r1amzevDnD9jt37sDX1xdPnjwp9JiIiIiIiAw1HQARERERqZe/v7/S3//880+cOHEiw/ZatWoVSjzt27fHoEGDlLY1aNCgUPb9udWrV6N06dIYMmSI0vY7d+5g3rx5aN26Nezs7DQSGxERERHpLxZtiYiIiIq4AQMGKP398uXLOHHiRIbthaVmzZoa27emJSYmwtzcXNNhEBEREZGWY3sEIiIiIkJiYiImTZqEypUrw9jYGA4ODli2bBkEQVB6nkgkwpgxYxAQEAAHBweYmJigUaNGOHfuXJ72l5ycjA8fPuTpNb6+vhCJRLh37x569+6NEiVKwNraGuPHj88w1qZNm9CmTRvY2NjA2NgYtWvXxpo1a5SeY2dnh/DwcJw9e1bRpqF169bYvHkzevXqBQBwc3PLtH3EkSNH0KpVK5ibm8PCwgIeHh4IDw9XGn/IkCEoXrw4/v33X3Tu3BkWFhbo378/gP9+jvv27UOdOnVgbGwMR0dHHD16VGmMhIQETJgwAXZ2djA2NoaNjQ3at2+P69ev5+lnR0RERES6hTNtiYiIiPScIAjw8vJCcHAwhg8fjvr16+PYsWOYMmUKnj17huXLlys9/+zZs9i5cyfGjRsHY2NjrF69Gh07dsSVK1dQp06dHPe3efNmrF69GoIgoFatWpg1axb69euX63h79+4NOzs7LFq0CJcvX8Zvv/2GN2/e4M8//1Q8Z82aNXB0dISXlxcMDQ1x4MABjBo1CjKZDKNHjwYArFixAmPHjkXx4sUxc+ZMAEDZsmVRvXp1jBs3Dr/99htmzJihaBsh/6+/vz8GDx4Md3d3LFmyBElJSVizZg1atmyJGzduKLVTkEqlcHd3R8uWLbFs2TKYmZkpHrtw4QL27t2LUaNGwcLCAr/99ht69uwJiUQCa2trAMC3336L3bt3Y8yYMahduzZiY2Nx4cIF3L17Fw0bNsz1z4yIiIiIdIxARERERHpl9OjRQvrDwH379gkAhAULFig9z9vbWxCJRMLDhw8V2wAIAISrV68qtkVGRgomJiZC9+7dc9x38+bNhRUrVgj79+8X1qxZI9SpU0cAIKxevTrH186dO1cAIHh5eSltHzVqlABA+OeffxTbkpKSMrze3d1dqFatmtI2R0dHwdXVNcNz//rrLwGAEBwcrLQ9ISFBKFmypPD1118rbX/58qVgaWmptH3w4MECAGHatGkZxgcgGBkZKf1s//nnHwGAsHLlSsU2S0tLYfTo0RleT0RERERFG9sjEBEREem5w4cPw8DAAOPGjVPaPmnSJAiCgCNHjihtb9asGRo1aqT4u62tLbp27Ypjx44hLS0t231dvHgR48ePh5eXF7799ltcu3YNderUwYwZM5CcnJyreOUzZeXGjh2reB9ypqamiv+Pj49HTEwMXF1d8ejRI8THx+dqP5k5ceIE3r59i759+yImJkbxx8DAAM7OzggODs7wmpEjR2Y6Vrt27VC9enXF37/88kuUKFECjx49UmwrWbIkQkND8fz583zHTERERES6h0VbIiIiIj0XGRmJChUqwMLCQmm7vB1AZGSk0nZ7e/sMY9SsWRNJSUl4/fp1nvZtZGSEMWPG4O3bt7h27VquXvP5/qtXrw6xWIwnT54otl28eBHt2rWDubk5SpYsiTJlymDGjBkAUKCibUREBACgTZs2KFOmjNKf48ePIzo6Wun5hoaGqFSpUqZj2draZthWqlQpvHnzRvH3n376Cbdv30blypXRpEkT+Pr6KhV1iYiIiKhoYk9bIiIiItKoypUrAwDi4uLy9XqRSKT093///Rdt27bFF198gV9++QWVK1eGkZERDh8+jOXLl0Mmk+U7Vvlr/f39Ua5cuQyPGxoqH14bGxtDLM58noSBgUGm24V0i7/17t0brVq1QmBgII4fP46lS5diyZIl2Lt3Lzp16pTft0FEREREWo5FWyIiIiI9V6VKFZw8eRIJCQlKs23v3buneDw9+WzT9B48eAAzMzOUKVMmz/uXzxzN7WsjIiJQtWpVxd8fPnwImUymWADswIEDSElJQVBQkNJs1sxaF3xe8M1pu7ydgY2NDdq1a5ereAuqfPnyGDVqFEaNGoXo6Gg0bNgQP/74I4u2REREREUY2yMQERER6bnOnTsjLS0Nfn5+StuXL18OkUiUoTh46dIlXL9+XfH3qKgo7N+/Hx06dMhy9iiATFsnJCQkYMWKFShdurRSn9zsrFq1SunvK1euBABFnPIY0s9YjY+Px6ZNmzKMZW5ujrdv32a6HUCGx9zd3VGiRAksXLgQHz9+zPC6vLaHyE5aWlqGVg42NjaoUKECUlJSVLYfIiIiItI+nGlLREREpOe6dOkCNzc3zJw5E0+ePEG9evVw/Phx7N+/HxMmTFBaLAsA6tSpA3d3d4wbNw7GxsZYvXo1AGDevHnZ7mfVqlXYt28funTpAltbW7x48QIbN26ERCKBv78/jIyMchXv48eP4eXlhY4dO+LSpUvYunUr+vXrh3r16gEAOnToACMjI3Tp0gXffPMN3r9/jz/++AM2NjZ48eKF0liNGjXCmjVrsGDBAtSoUQM2NjZo06YN6tevDwMDAyxZsgTx8fEwNjZGmzZtYGNjgzVr1mDgwIFo2LAhfHx8UKZMGUgkEhw6dAgtWrTIUPzOr4SEBFSqVAne3t6oV68eihcvjpMnTyIsLAw///yzSvZBRERERNqJRVsiIiIiPScWixEUFIQ5c+Zg586d2LRpE+zs7LB06VJMmjQpw/NdXV3RrFkzzJs3DxKJBLVr18bmzZvx5ZdfZrufFi1aICQkBOvXr0dsbCzMzc3RpEkTbNy4EW3atMl1vDt37sScOXMwbdo0GBoaYsyYMVi6dKnicQcHB+zevRuzZs3C5MmTUa5cOYwcORJlypTBsGHDlMaaM2cOIiMj8dNPPyEhIQGurq5o06YNypUrh7Vr12LRokUYPnw40tLSEBwcDBsbG/Tr1w8VKlTA4sWLsXTpUqSkpKBixYpo1aoVhg4dmuv3kRMzMzOMGjUKx48fx969eyGTyVCjRg2sXr0aI0eOVNl+iIiIiEj7iIT0940REREREWVDJBJh9OjRKptNmhe+vr6YN28eXr9+jdKlSxf6/omIiIiICgt72hIRERERERERERFpERZtiYiIiIiIiIiIiLQIi7ZEREREREREREREWoQ9bYmIiIiIiIiIiIi0CGfaEhEREREREREREWkRFm2JiIiIiIiIiIiItAiLtkRERERERERERERahEVbIiIiIiIiIiIiIi3Coi0RERERERERERGRFmHRloiIiIiIiIiIiEiLsGhLREREREREREREpEVYtCUiIiIiIiIiIiLSIizaEhEREREREREREWkRFm2JiIiIiIiIiIiItAiLtkRERERERERERERahEVbIiIiIiIiIiIiIi3Coi0RERERERERERGRFmHRloiIiIiIiIiIiEiLsGhLREREREREREREpEVYtCUiIiIiIiIiIiLSIizaEhEREREREREREWkRFm2JiIiIiIiIiIiItAiLtkRERERERERERERahEVbIiIiIiIiIiIiIi3Coi0RERERERERERGRFmHRloiIiIiIiIiIiEiLsGhLREREREREREREpEV0tmg7ZMgQFC9eXNNhkA7x9fWFSCTSdBhUBDD/UF5t3rwZIpEIT5480XQopOOYfyivmH9IVZh/KK+GDBkCOzs7TYdBOo65h/KqKOWefBdt5QeA8j+GhoaoWLEihgwZgmfPnqkyRpV78eIFpk2bBjc3N1hYWEAkEuHMmTMFHvf48eMYPnw46tSpAwMDg3x9SIKCgtCwYUOYmJjA1tYWc+fOhVQqzfC8EydOoGXLljAzM0OpUqXg7e2d4WD8zJkzSv9Gn//58ccfFc+VFzSz+vPy5cs8v5fPyWQy/PTTT6hatSpMTEzw5ZdfYvv27Zk+d9euXWjatClKliwJa2truLq64tChQznuIykpCb6+vir59yTtpcv559SpUxg2bBhq1qwJMzMzVKtWDV999RVevHhRoHE/fvyIefPmoVq1ajA2Nka1atWwYMGCTPPH55KTkxW5y9LSEsWLF0e9evXw66+/4uPHj0rPPXfuHLy8vFC5cmWYmJigXLly6NixIy5evJjp2KmpqVi4cCG++OILmJiYoGzZsvDw8MDTp08Vz/n83/PzP5cvXy7Qz0Zuw4YNqFWrFkxMTGBvb4+VK1dmeI6dnV2Wcdjb2+e4j4ULF2Lfvn0qiZe0ky7nn7x+f3NLJpNh7dq1qF+/PooXL46yZcuiU6dOCAkJyfG1UVFRmDdvHpo0aYJSpUqhdOnSaN26NU6ePJnp83Nz/BMbG4ulS5fCxcUFZcqUQcmSJdG0aVPs3Lkzw3jMP6RLdDn/fO7rr7+GSCSCp6dngcdKSEjA999/j6pVq8LY2BgVK1aEt7c3kpKScnztixcvMGLECFStWhWmpqaoXr06vvvuO8TGxmZ4bm7OT54/f44BAwbAwcEBFhYWKFmyJJo0aYItW7ZAEASl5w4ZMiTL77yJiUnBfij/l5KSgqlTp6JChQowNTWFs7MzTpw4ofScJ0+eZJsHv/7662z38fz5c/j6+uLvv/9WScykfXQ592T3e76gNY785p579+7h+++/R/369WFhYYHy5cvDw8MDV69ezfI1O3fuRLNmzWBubo6SJUuiefPmOH36dJbPv3DhguJ9xsTEKD3G3KNbDAs6wPz581G1alV8+PABly9fxubNm3HhwgXcvn1bZf/gqnb//n0sWbIE9vb2qFu3Li5duqSScbdt24adO3eiYcOGqFChQp5ff+TIEXTr1g2tW7fGypUrcevWLSxYsADR0dFYs2aN4nkHDx5E165d0bBhQyxevBjv3r3Dr7/+ipYtW+LGjRsoU6YMAKBWrVrw9/fPsB9/f38cP34cHTp0yPDYmjVrMr2KVbJkyTy/n8/NnDkTixcvxtdffw0nJyfs378f/fr1g0gkgo+Pj+J5K1euxLhx4+Dh4YHFixfjw4cP2Lx5Mzw9PbFnzx706NEjy30kJSVh3rx5AIDWrVsrPTZr1ixMmzatwO+DtIcu5p+pU6ciLi4OvXr1gr29PR49egQ/Pz8cPHgQf//9N8qVK5evcQcMGIC//voLw4YNQ+PGjXH58mXMnj0bEokE69aty/a1ycnJCA8PR+fOnWFnZwexWIyQkBBMnDgRoaGh2LZtm+K5Dx48gFgsxrfffoty5crhzZs32Lp1K1xcXHDo0CF07NhR8dyPHz/Cw8MDISEh+Prrr/Hll1/izZs3CA0NRXx8PCpVqqQUh/zf83M1atTI188kvd9//x3ffvstevbsie+++w7nz5/HuHHjkJSUhKlTpyqet2LFCrx//17ptZGRkZg1a1amOfNzCxcuhLe3N7p166a0feDAgfDx8YGxsXGB3wtpB13MP3n5/ubFlClT8Msvv2DAgAEYNWoU3r59i99//x2urq64ePEimjRpkuVr9+/fjyVLlqBbt24YPHgwpFIp/vzzT7Rv3x4bN27E0KFDFc/N7fHPpUuXMHPmTHTu3BmzZs2CoaEh9uzZAx8fH9y5c0dxnJAe8w/pEl3MP+ldvXoVmzdvVkms8fHxcHV1xdOnTzFixAjUqFEDr1+/xvnz55GSkgIzM7MsX/v+/Xs0a9YMiYmJGDVqFCpXrox//vkHfn5+CA4OxrVr1yAWf5rjlNvzk5iYGDx9+hTe3t6wtbXFx48fceLECQwZMgT379/HwoULlWIwNjbG+vXrM8RmYGBQ4J8N8Kk4s3v3bkyYMAH29vbYvHkzOnfujODgYLRs2RIAUKZMmUzPGY8ePYqAgIAc88/z588xb9482NnZoX79+kqP/fHHH5DJZCp5L6R5upx7Mvs9X5AaR0Fyz/r167Fhwwb07NkTo0aNQnx8PH7//Xc0bdoUR48eRbt27ZSe7+vri/nz58Pb2xtDhgzBx48fcfv27SwL5jKZDGPHjoW5uTkSExMzfQ5zjw4R8mnTpk0CACEsLExp+9SpUwUAws6dO/M7dK4MHjxYMDc3z9dr3717J8TGxgqCIAh//fWXAEAIDg4ucEzPnj0TUlNTBUEQBA8PD6FKlSp5en3t2rWFevXqCR8/flRsmzlzpiASiYS7d+8qPa9GjRpCSkqKYtvff/8tiMVi4bvvvstxPzVq1BDs7e2Vts2dO1cAILx+/TpPMctfm9N7ffr0qVCsWDFh9OjRim0ymUxo1aqVUKlSJUEqlSq229vbC05OToJMJlNsi4+PF4oXLy54eXllu5/Xr18LAIS5c+fm+X2Q7tDl/HP27FkhLS0twzYAwsyZM/M15pUrVwQAwuzZs5W2T5o0SRCJRMI///yTr3HHjBkjABBevHiR7fMSExOFsmXLCu7u7krblyxZIhQrVkwIDQ3N9vVZ/Xvmhvy12UlKShKsra0FDw8Ppe39+/cXzM3Nhbi4uGxf/8MPPwgAhIsXL+YYj7m5uTB48OAcn0e6S5fzT2ay+v7m1sePHwVTU1PB29tbafujR48EAMK4ceOyff3t27czHHt8+PBB+OKLL4RKlSopbc/t8c+jR4+EJ0+eKL1WJpMJbdq0EYyNjYX3798rtjP/kC4pCvlHJpMJzZo1E4YNGyZUqVIlw3cjr0aOHCmULFlSePToUZ5fGxAQIAAQDh48qLR9zpw5AgDh+vXrim0FOT8RBEHw9PQUzM3Nlc55CvLzHDx4sODq6prtc0JDQwUAwtKlSxXbkpOTherVqwvNmjXLcR9t27YVSpQoISQnJ2f7vLCwMAGAsGnTptyETjpIl3NPQX7PZ6cguefq1atCQkKC0raYmBihTJkyQosWLZS2X7p0SRCJRMIvv/yS6/HXrFkjWFtbC+PHj8+0xsPco1tU3tO2VatWAIB///1XsS01NRVz5sxBo0aNYGlpCXNzc7Rq1QrBwcFKr5VPj162bBnWrVuH6tWrw9jYGE5OTggLC8tx33///TfKlCmD1q1bZ5ipkJ6FhQWsrKxy9X5iYmJw7969XN1eU6FCBRQrVixX437uzp07uHPnDkaMGAFDw/8mQI8aNQqCIGD37t0AgLi4ONy5cwfdu3eHkZGR4nn16tVDrVq1sGPHjmz3c+XKFTx8+BD9+/fPV5z5tX//fnz8+BGjRo1SbBOJRBg5ciSePn2qNNv53bt3sLGxUeo/W6JECRQvXhympqZZ7uPJkyeKWTbz5s1TTKv39fUFkHlPW5FIhDFjxuCvv/5C7dq1YWpqimbNmuHWrVsAPs2OqVGjBkxMTNC6detM+8GFhoaiY8eOsLS0hJmZmWJmUXoJCQmYMGEC7OzsYGxsDBsbG7Rv3x7Xr1/P3Q+QckUX8o+Li4ti1kb6bVZWVrh7967S9tzmn/PnzwOA0ox1+d8FQcj0luDckLd4efv2bbbPMzMzQ5kyZZSeJ5PJ8Ouvv6J79+5o0qQJpFJprvKoOgQHByM2NlYp/wDA6NGjkZiYmGPrlW3btqFq1apo3rx5ts8TiURITEzEli1bFPlnyJAhADLvKWlnZwdPT0+cOXMGjRs3hqmpKerWrato77J3717UrVsXJiYmaNSoEW7cuJFhn/fu3YO3tzesrKxgYmKCxo0bIygoSOk58tYZ9vb2MDExgbW1NVq2bJnhFiUqGF3IP5nJ7PsL5D7/fPz4EcnJyShbtqzSdhsbG4jF4mx/bwOAo6MjSpcurbTN2NgYnTt3xtOnT5GQkAAgb8c/VatWRZUqVZTGFIlE6NatG1JSUvDo0aNsY1Il5h/mn8KgS/nH398ft2/fVmrT9rkXL17g3r17GVo0fe7t27fYtGmTor1BamoqUlJScoxB7t27dwCQIX+VL18eAJTyV37PT+Ts7OyQlJSE1NTUXMdXULt374aBgQFGjBih2GZiYoLhw4fj0qVLiIqKyvK1L168QHBwMHr06JHtDMozZ87AyckJADB06FBF/tm8eTOAjH0l03/eVq1ahWrVqsHMzAwdOnRAVFQUBEHADz/8gEqVKsHU1BRdu3ZFXFxchv0eOXIErVq1grm5OSwsLODh4YHw8HCl57x8+RJDhw5FpUqVYGxsjPLly6Nr167s761CupR7gE/n5GlpaVk+Xli5p1GjRhnubra2tkarVq0ynA+uWLEC5cqVw/jx4yEIQo7vNS4uDrNmzcL8+fNVcrd0fjD3qDb3qLxoKw+kVKlSim3v3r3D+vXr0bp1ayxZsgS+vr54/fo13N3dM+0/sW3bNixduhTffPMNFixYgCdPnqBHjx7ZfnnCwsLQpk0bNGjQAEeOHFFZo2o/Pz/UqlULV65cUcl4WZEfDDdu3Fhpe4UKFVCpUiXF4/JkkNnBgZmZGZ4/f55tb5aAgAAAyLJoGxcXh5iYGKU/ORVscuPGjRswNzdHrVq1lLbLb5lMfzLQunVrHD16FCtXrsSTJ09w7949jB49GvHx8Rg/fnyW+yhTpoyijUT37t3h7+8Pf3//bNspAJ8KXpMmTcLgwYPh6+uLu3fvwtPTE6tWrcJvv/2GUaNGYcqUKbh06RKGDRum9NrTp0/DxcUF7969w9y5c7Fw4UK8ffsWbdq0UfrMfPvtt1izZg169uyJ1atXY/LkyTA1Nc2QlKlgdDX/vH//Hu/fv89QuMht/skqL8hvy7l27Vqu4khNTUVMTAyioqIQGBiIZcuWoUqVKpneHvzu3TtFUWfGjBm4ffs22rZtq3j8zp07eP78Ob788kuMGDEC5ubmMDc3x5dffpnhoE0uPj4+Q/7JrKdcXmWVXxs1agSxWJxpMSL9a+/evYt+/frluB9/f38YGxujVatWivzzzTffZPuahw8fol+/fujSpQsWLVqEN2/eoEuXLggICMDEiRMxYMAAzJs3D//++y969+6tdJtPeHg4mjZtirt372LatGn4+eefYW5ujm7duiEwMFDxPF9fX8ybNw9ubm7w8/PDzJkzYWtry4tGKqZL+Sen7y+Q+/wj71O2efNmBAQEQCKR4ObNmxgyZAhKlSqldMCeFy9fvoSZmZkijxX0+Ec+JoAMuRZg/mH+0W26kn8SEhIwdepUzJgxI9t2UNOnT0etWrVy7JV54cIFfPjwATVq1IC3tzfMzMxgamqKFi1a5KrHofxC+vjx43H58mU8ffoUhw8fxo8//ohu3brhiy++UDw3r+cnycnJiImJwZMnT7BlyxZs2rQJzZo1yzSHfZ57YmJiFAXlgrhx4wZq1qyJEiVKKG2Xn39l9zPasWMHZDJZjhN9atWqhfnz5wMARowYocg/Li4u2b4uICAAq1evxtixYzFp0iScPXsWvXv3xqxZs3D06FFMnToVI0aMwIEDBzB58mSl1/r7+8PDwwPFixfHkiVLMHv2bNy5cwctW7ZUKor07NkTgYGBGDp0KFavXo1x48YhISEBEokk29go93Ql9wCAm5sbSpQoATMzM3h5eSEiIiLDcwor92Tl5cuXGY5RTp06BScnJ/z2228oU6aMogeun59fpmPMnj0b5cqVy/EYAGDu0Znck98puvJp5idPnhRev34tREVFCbt37xbKlCkjGBsbC1FRUYrnSqVSpVvZBEEQ3rx5I5QtW1YYNmyYYtvjx48FAIK1tbXS7WL79+8XAAgHDhxQbEs/pfvChQtCiRIlBA8PD+HDhw95eh85tUeQtw3Ia/uEvLZHWLp0qQBAkEgkGR5zcnISmjZtKgiCIKSlpQklS5YU2rZtq/ScmJgYwdzcXAAgXL16NdN9SKVSoWzZskKTJk0yPCZ/n5n9cXBwyDb23LRH8PDwEKpVq5Zhe2JiogBAmDZtmmLbq1evhLZt2yrFULp0aSEkJCTbfQhC9u0R5O8xPQCCsbGx8PjxY8W233//XQAglCtXTnj37p1i+/Tp0wUAiufKZDLB3t5ecHd3V7pVKikpSahatarQvn17xTZLS0ul1hBUMEUl/8jJb389deqU0vbc5p89e/YIAAR/f3+l7WvXrhUACHXq1MlVHNu3b1f63jVu3Fi4efNmps91d3dXPM/IyEj45ptvlG5h2bt3r+LnaW9vL2zatEnYtGmTYG9vLxgZGSm1bJD/e2b2x9jYONuYc3N78ujRowUDA4NMHytTpozg4+OT5WsnTZokABDu3LmT7T7ksro9WR5n+lxTpUoVAYBSbjt27JgAQDA1NRUiIyMV2+V5Kf1noW3btkLdunWVPncymUxo3ry5UgucevXqFfj2U/pPUcg/OX1/BSFvxz8RERFCw4YNlb671apVE+7du5frmD4fz8TERBg4cKBiW0GOfwRBEGJjYwUbGxuhVatWStuZf5h/dImu55/JkycLVatWVTw/q/YIgwcPzvCZzcwvv/yiiL1JkyZCQECAsHr1aqFs2bJCqVKlhOfPn+cY0/r164WSJUsqffcHDx6s1K5OEPJ+frJo0SKl57Zt2zbDeZ78fWb2J6eWNbm5RdnR0VFo06ZNhu3h4eECAGHt2rVZvrZRo0ZC+fLlM7T0ykx2tygPHjxY6TxR/nkrU6aM8PbtW8V2+XnW560C+/btKxgZGSk+MwkJCULJkiWFr7/+Wmk/L1++FCwtLRXb37x5k+H2bMo/Xc49O3fuFIYMGSJs2bJFCAwMFGbNmiWYmZkJpUuXzvI7WRi553Pnzp0TRCKRUru7uLg4xX6KFy8uLF26VNi5c6fQsWPHTL/D//zzj2BgYCAcO3ZMEISsW2Ay9+hW7inwQmSfN0m2s7PD1q1blRaYMTAwUDQ0lslkePv2LWQyGRo3bpzplfY+ffooXa2RT7vP7Ha24OBgdOnSBR06dMCOHTuUbplTBV9fX8Xt9eqUnJwMAJkuEmFiYqK44iEWi/HNN99gyZIlmD59OoYNG4Z3797h+++/V9xuIx/rc6dOncKrV68wY8aMLOPYs2dPhisi5ubmSn//fPXBpKQkyGSyDNstLCwU7yc5OTnL9/Z5zGZmZnBwcEClSpXg6emJhIQELF++HD169MD58+dVsihIem3btlWaOu/s7Azg0xUSCwuLDNsfPXoEOzs7/P3334iIiMCsWbMyzMZp27Yt/P39IZPJIBaLUbJkSYSGhuL58+f5WqSOMlcU8s+5c+cwb9489O7dG23atFF6LLf5p3PnzqhSpQomT54MMzMzNGrUCKGhoZg5cyYMDQ2zzAmfc3Nzw4kTJ/D27VucOnUK//zzT5bN6xcvXoxJkyYhKioKW7ZsQWpqKqRSqeJx+a07CQkJuHHjBipXrgwAaNOmDWrUqIGffvoJW7duVRpz1apVqFmzptK2z5vhv3nzRum2Jvl+Ps8/6WfoJScnZ/lvY2JikuXPRyaTYceOHWjQoEGGuwRUpXbt2mjWrJni7/I806ZNG9ja2mbY/ujRI7Ru3RpxcXE4ffo05s+fj4SEBMUt5ADg7u6OuXPn4tmzZ6hYsSJKliyJ8PBwRERE5GoFesodXc4/OX1/gbwd/1hYWMDR0RHNmjVD27Zt8fLlSyxevBjdunXD+fPnM53ZmpWkpCT06tULpqamWLx4sWJ7QY5/5DM23r59i5UrV2b6HOYf5h9doov558GDB/j111+xffv2HBfF27x5s+IW1+zIv4MikQinTp1SzLZr0KABmjVrhlWrVmHBggXZjlGxYkU0adJEcSx1/vx5/PbbbyhdujSWLVumeF5ez0/69u2Lxo0b4/Xr1zh48CBevXqV6ffdxMQEBw4cyLA9fd6UyWQZbtNNSUnBx48fM+QfS0tLRbu+vJx/pffgwQNcu3YNEydOzNDSS1V69eoFS0tLxd/leWbAgAFKrQKdnZ2xfft2PHv2DNWqVVMcp/bt21fpvRsYGMDZ2VlxN5epqSmMjIxw5swZDB8+XOmzTfmni7mnd+/e6N27t+Lv3bp1g7u7O1xcXPDjjz9i7dq1iscKM/ekFx0djX79+qFq1ar4/vvvM+wnNjYWO3bsQJ8+fQAA3t7eqFu3LhYsWKA0o3bcuHHo1KlTrhYvZe75RBdyT4GLtvKD3Pj4eGzcuBHnzp3L9B9oy5Yt+PnnnzP0CMlspd70B4rAf9Pt37x5o7T9w4cP8PDwQKNGjbBr1y6lH7Kukd8qk1kvlA8fPijdSjN//nzExMTgp59+UpzQdOjQAcOHD8fatWuzvD0gICAABgYGii97ZlxcXHI8uZL3jc1p+6ZNmxT91ExNTbN8b/LH5Xr16gVDQ0OlJNK1a1fY29tj5syZ+e7PmZXPP2/yL7G8yPT5dvnnUH5LxeDBg7McOz4+HqVKlcJPP/2EwYMHo3LlymjUqBE6d+6MQYMGoVq1aip7H/pI1/PPvXv30L17d9SpUyfT1Ttzy8TEBIcOHULv3r3Rs2dPAJ8uAP3000/48ccfc92uoWzZsoq+bt7e3li4cCHat2+PiIiIDLcxpl+hc8CAAWjYsKFilVDgv+90ixYtlL5Ltra2aNmyJUJCQjLsv0mTJhluIf5cgwYNEBkZmWH75/ln7ty5ioKTqalplj3kPs+v6Z09exbPnj3DxIkTs42pIPKbfx4+fAhBEDB79mzMnj0707Gjo6NRsWJFzJ8/H127dkXNmjVRp04ddOzYEQMHDsSXX36p6rejV3Q5/+T0/c0LqVSKdu3aoXXr1koF0Xbt2sHR0RFLly7FkiVLcjVWWloafHx8cOfOHRw5ciTDRc78Hv+MHTsWR48exZ9//ol69epl+hzmH+YfXaKL+Wf8+PFo3ry54jhFFeTfny5duih9/5s2bYqqVatmeqyR3sWLF+Hp6YnLly8rvv/dunVDiRIlMG/ePAwbNgy1a9cGkPfzkypVqih6a/ft2xcjRoxAu3btcP/+faXvvYGBQYZC2OckEkmm/2ZAxvwTHByM1q1bA8jb+Vd6ObXTU4WCnn99PtFBTj75yNjYGEuWLMGkSZNQtmxZNG3aFJ6enhg0aFC2rTkoe7qYezLTsmVLODs74+TJk/l6fUFzT3qJiYmKC0EXLlxQGk++n2LFisHb21uxXSwWo0+fPpg7dy4kEglsbW2xc+dOhISE4Pbt27naL3PPJ7qQewpc5Ux/kNutWze0bNkS/fr1w/379xUfuK1bt2LIkCHo1q0bpkyZAhsbGxgYGGDRokVKTavlPp/ZICcIgtLf5QtV7N+/H0ePHoWnp2dB347GyBvev3jxIsMH5sWLF4r+HwBgZGSE9evX48cff8SDBw9QtmxZ1KxZE/369YNYLM50JmpycjICAwPRrl27DM328+rzxSP+/PNPHD9+PMOsOUdHR6X3FxwcDEEQlBr4v3jxAgAUJ2aPHj3C0aNHsW7dOqWxrKys0LJlywwLfKlCVp+3nD6H8t5uS5cuVToBTk/+HejduzdatWqFwMBAHD9+XHESu3fvXnTq1KmA70B/6XL+iYqKQocOHWBpaYnDhw8rzerOD0dHR9y+fRt37tzBmzdvFAvrTZw4Ea6urvka09vbGzNnzsT+/fuz7YtkZGQELy8vLF68GMnJyTA1NVV8pzPLNzY2Ntn2ccxOQECA0tVZ+ffp87yU/oJI+fLlkZaWhujoaNjY2Ci2p6amIjY2NsvZ7wEBARCLxejbt2++Ys2NguafyZMnw93dPdPnyn8XuLi44N9//8X+/ftx/PhxrF+/HsuXL8fatWvx1VdfFfQt6C1dzj/pZfb9zYtz587h9u3b+OWXX5S229vbo1atWnn6vf3111/j4MGDCAgIyPSgOD/HP/PmzcPq1auxePFiDBw4ME/v7XPMP8w/2kLX8s/p06dx9OhR7N27V6nvn1QqRXJyMp48eQIrK6sMd/vlJKdjjc+LPp/7/fffUbZs2QwXbLy8vODr64uQkBDUrl1bJecn3t7e+OOPP3Du3LksvzdZKVeuXIY8s3TpUrx8+RI///yz0vb0F6bKly+faW/Oz8+/Prdt2zY4ODigUaNGeYozLwqaf/z9/TMtgKQv5E2YMAFdunTBvn37cOzYMcyePRuLFi3C6dOn0aBBg4K+Bb2ka7knO5UrV8b9+/fz9dqC5h651NRU9OjRAzdv3sSxY8dQp04dpcflC32WLFkyw89Jfkzx5s0b2NraYsqUKejVqxeMjIwUeVa+NlFUVBRSU1PzfMcvc4925B6VTk2Vfxnliw1MmzYNwKfV46pVq4a9e/cqFezmzp1boP2JRCIEBASga9eu6NWrF44cOaKo7usaedHv6tWrSgXa58+f4+nTp5ku5JF+VlxaWhrOnDkDZ2fnTGeaBAUFISEhQSVXLT6/InPhwgWYmJhke6Wmfv36WL9+Pe7evau4Yg0AoaGhiscB4NWrVwCQ6aqOHz9+zHD75ufSf77UrXr16gA+XVXJ6SoV8Cl5jRo1CqNGjUJ0dDQaNmyIH3/8kUVbFdGl/BMbG4sOHTogJSUFp06dUly0KSiRSKR0seTw4cOQyWS5+nxmRl6ciI+Pz9VzBUFAQkKCYhXyYsWKZfoL+/nz51nO2M9JixYtlP7+9OlTABnzUnrp82vnzp0V269evQqZTJbpRZeUlBTs2bMHrVu3ztMBTmHlIHlRqFixYrn697WyssLQoUMxdOhQvH//Hi4uLvD19WXRREV0Kf9k5vPvb14U9Pe23JQpU7Bp0yasWLEix0Jlbo9/Vq1aBV9fX0yYMAFTp07NVRzZYf75hPlHu+hC/pEvvpLZ4sDPnj1D1apVsXz5ckyYMCFPschP7LM61ki/kFhmXr16lWXuAqDIX6rIc3k5pvpcZudZW7duRUpKSo75Jzg4GO/evVMqiH9+/pVeaGgoHj58qFjgJzc0cf5lY2OTq/xTvXp1TJo0CZMmTUJERATq16+Pn3/+OcNkI8o7Xcg92Xn06FG+z0cKmnuAT0XAQYMG4dSpU9i1a1emk2zEYjHq16+PsLAwpKamKrWDeP78OYD/ZrxGRUVh27Zt2LZtW4ZxGjZsiHr16uV5kTTmnv9oMveovFFE69at0aRJE6xYsUIx/VletU5/tSQ0NBSXLl0q8P6MjIywd+9eODk5oUuXLjmucpxX8tWVk5KSVDbmx48fce/ePcWVBuDTLLkvvvgC69atUzogWLNmDUQikdJ0+MwsW7YML168wKRJkzJ9fNu2bTAzM0P37t1V8ybyqGvXrihWrBhWr16t2CYIAtauXYuKFSuiefPmAD7NzBCLxdi5c6fS5+Xp06c4f/58jlcm5D3k5FeV1KlRo0aoXr06li1bpug3k97r168BfDrA+/wAzcbGBhUqVMj0tgHKP13IP4mJiejcuTOePXuGw4cPZ9vjryD5Jzk5GbNnz0b58uWVCiBJSUm4d++eUi+emJiYDFezAShaNqSfgRIdHZ3heW/fvsWePXtQuXJlxVVfCwsLdO7cGSEhIbh3757iuXfv3kVISAjat2+f5/eUX23atIGVlRXWrFmjtH3NmjUwMzODh4dHhtccPnwYb9++zfOFLnNz80LJPzY2NmjdujV+//13pd8lcvL8AyBDz+3ixYujRo0azD8qpgv5J7ffXyD3+UfeB3bHjh1K269fv4779+8r/d7OLP8An2ZtLFu2DDNmzMh0FfbsZHX8s3PnTowbNw79+/fPMAu4MDH/MP8UBm3PP23atEFgYGCGP2XKlEHjxo0RGBiILl26KJ7/4sWLDLdVZ8bBwQH16tXD/v37lfLK8ePHERUVpXSsER8fj3v37ikdk9esWROvXr3CmTNnlMbdvn07ACjyV17OT9J//tPbsGEDRCIRGjZsmO17UiVvb2+kpaUpzRBOSUnBpk2b4OzsnOHuTgCKgk+/fv1yvR/5+ieFkX/c3d1RokQJLFy4MNPPh/znn5SUpPguyFWvXh0WFhbMPyqk7bkHyPw7efjwYVy7dg0dO3ZU2l5YuQf41Lpp586dWL16daYXtOT69OmDtLQ0bNmyRbHtw4cPCAgIQO3atRUXdzPLsfK2mH/++SeWL1+e7XtSJeYe1eYetTSBlU/N3rx5M7799lt4enpi79696N69Ozw8PPD48WOsXbsWtWvXzrTYlVempqY4ePAg2rRpg06dOuHs2bMZppZ/Tt4YOjw8HMCnac4XLlwAAMyaNUvxPD8/P8ybN0+pR0dWbt68iaCgIACf+n3Fx8cr9lOvXj3FwcizZ89Qq1YtDB48WKnR9dKlS+Hl5YUOHTrAx8cHt2/fhp+fH7766iulRSi2bt2KPXv2wMXFBcWLF8fJkyexa9cufPXVV5n2iYqLi8ORI0fQs2fPHHtb7t69O9PntG/fvkBtFSpVqoQJEyZg6dKl+PjxI5ycnLBv3z6cP39e0WsX+HSlaNiwYVi/fj3atm2LHj16ICEhAatXr0ZycjKmT5+e7X5MTU1Ru3Zt7Ny5EzVr1oSVlRXq1KmT4+chP8RiMdavX49OnTrB0dERQ4cORcWKFfHs2TMEBwejRIkSOHDgABISElCpUiV4e3ujXr16in+zsLCwDLcVUMFpe/7p378/rly5gmHDhuHu3bu4e/eu4rHixYujW7duir/nJf/07t0bFSpUQO3atfHu3Tts3LgRjx49wqFDh5RaL1y5cgVubm5KPRe3bt2KtWvXolu3bqhWrRoSEhJw7NgxnDhxAl26dFG6TblTp06oVKkSnJ2dYWNjA4lEgk2bNuH58+cZ+rktXLgQp06dQps2bTBu3DgAwG+//QYrK6tMF0Q8cuSIUoFXrnnz5gXq/2xqaooffvgBo0ePRq9eveDu7o7z589j69at+PHHH2FlZZXhNQEBATA2Ns5z771GjRrh5MmT+OWXX1ChQgVUrVpV0eBe1VatWoWWLVuibt26+Prrr1GtWjW8evUKly5dwtOnT/HPP/8A+LTYUOvWrdGoUSNYWVnh6tWr2L17N8aMGaOWuPSZtuefvHx/c5t/GjVqhPbt22PLli149+4dOnTogBcvXmDlypUwNTVVmjmXWf4JDAzE999/r2in8PkMhPTHH7k9/rly5QoGDRoEa2trtG3bVtEjTS6znML8kzfMP9pHm/OPra1thj6CwKdbSMuWLat07AMA06dPx5YtW/D48WOlxYIzs3z5crRv3x4tW7bEN998g/j4ePzyyy+oWbMmRo4cqXheYGAghg4dqrTmxpgxY7Bp0yZ06dIFY8eORZUqVXD27Fls374d7du3V3x/8nJ+8uOPP+LixYvo2LEjbG1tERcXhz179iAsLAxjx47N0MZFKpVmOfOqe/fuGRaEzgtnZ2f06tUL06dPR3R0NGrUqIEtW7bgyZMn2LBhQ4bnp6WlYefOnWjatKliVlluVK9eHSVLlsTatWthYWEBc3NzODs7Z9kLsyBKlCiBNWvWYODAgWjYsCF8fHxQpkwZSCQSHDp0CC1atICfnx8ePHiAtm3bonfv3qhduzYMDQ0RGBiIV69ewcfHR+Vx6TNtzj3Ap9/jDRo0QOPGjWFpaYnr169j48aNqFy5cobzkcLKPStWrMDq1avRrFkzmJmZZcgB6b/733zzDdavX4/Ro0fjwYMHsLW1hb+/PyIjI5V6bH+eRwEoZtZ26tQpw7pFzD15o9HcI+TTpk2bBABCWFhYhsfS0tKE6tWrC9WrVxekUqkgk8mEhQsXClWqVBGMjY2FBg0aCAcPHhQGDx4sVKlSRfG6x48fCwCEpUuXZhgTgDB37lzF3wcPHiyYm5srPScmJkaoXbu2UK5cOSEiIiLb+AFk+Se9uXPnCgCE4ODgXP9MMvszePDgDO8z/Ta5wMBAoX79+oKxsbFQqVIlYdasWUJqaqrSc0JDQwUXFxehVKlSgomJiVCvXj1h7dq1gkwmyzSutWvXCgCEoKCgLGOXv8+s/mT3/ufOnav075iVtLQ0xefAyMhIcHR0FLZu3ZrheR8/fhRWrlwp1K9fXyhevLhQvHhxwc3NTTh9+nSO+xAEQQgJCREaNWokGBkZKX1u5O8xPQDC6NGjlbZl9TkMDg4WAAh//fWX0vYbN24IPXr0EKytrQVjY2OhSpUqQu/evYVTp04JgiAIKSkpwpQpU4R69eoJFhYWgrm5uVCvXj1h9erVuXo/lJEu558qVapk+T37/HuUl/yzZMkS4YsvvhBMTEyEUqVKCV5eXsKNGzcyPE/+OU7/fsLCwoRevXoJtra2grGxsWBubi40bNhQ+OWXX4SPHz8qvd7Pz09o2bKlULp0acHQ0FAoU6aM0KVLF+HcuXOZxnXt2jWhXbt2grm5uWBhYSF07dpVePDggdJzssudAIRNmzZl+b7lr82NdevWCQ4ODoKRkZFQvXp1Yfny5Znmzfj4eMHExETo0aNHrsZN7969e4KLi4tgamqqlOflcT5+/Fjx3CpVqggeHh4ZxshLXvr333+FQYMGCeXKlROKFSsmVKxYUfD09BR2796teM6CBQuEJk2aCCVLlhRMTU2FL774Qvjxxx8z/G6h3NHl/JOX729e8k9SUpIwf/58oXbt2oKpqalgaWkpeHp6ZshBmeWfvBx/5Pb4Jy85hflHGfOPdtPl/JOZrD6HgwcPzvCZzc6JEyeEpk2bCiYmJoKVlZUwcOBA4cWLF0rPkf/sPv9O37t3T/D29hYqV64sFCtWTKhSpYowefJkITExUel5uT0/OX78uODp6SlUqFBBKFasmGBhYSG0aNFC2LRpU4bvvPx9ZvUnu/c/ePBgwdXVNcefTXJysjB58mShXLlygrGxseDk5CQcPXo00+cePXpUACD89ttvOY77uf379wu1a9cWDA0NlX7Ouf28ZXWeldVnPjg4WHB3dxcsLS0FExMToXr16sKQIUOEq1evCoLw6XM5evRo4YsvvhDMzc0FS0tLwdnZWdi1a1ee3xvpdu6ZOXOmUL9+fcHS0lIoVqyYYGtrK4wcOVJ4+fJlhucWVu7J63f/1atXwuDBgwUrKyvB2NhYcHZ2zvJ7nJ78GOv169eZvk/mHt3IPSJByOSeWCIiIiIiIiIiIiLSCJX3tCUiIiIiIiIiIiKi/GPRloiIiIiIiIiIiEiLsGhLREREREREREREpEVYtCUiIiIiIiIiIiLSIizaEhEREREREREREWkRFm2JiIiIiIiIiIiItAiLtkRERERERERERERahEVbIiIiIiIiIiIiIi3Coi0RERERERERERGRFmHRloiIiIiIiIiIiEiLsGhLREREREREREREpEVYtCUiIiIiIiIiIiLSIizaEhEREREREREREWkRQ00HQNorMUWKJ7GJSJXKYGQohp21OcyN+ZEhItI05mciygnzBBFpAnMPEekCXclV2hcRaVTEqwQEhEoQfD8akrgkCOkeEwGwtTKDm4MN+jvbwr6shabCJCLSO8zPRJQT5gki0gTmHiLSBbqYq0SCIAg5P42Kuqi4JMwIvIXzD2NgIBYhTZb1x0L+eKsapbGwe11UtjIrxEiJiPQL8zPRf3RlVkRhY54gIk1g7iEiXaDLuYpFW8KOMAnmBoVDKhOy/fB+zkAsgqFYhHlejvBxslVjhERE+on5mUg3Z0UUJuYJItIE5h4i0gW6nqtYtNVzfsERWHb8QYHHmdyhJsa42asgIiIiApifiXR5VkRhYZ4gIk1g7iEiXVAUcpVYI3slrbAjTKKSDzAALDv+ADvDJCoZi4hI3zE/k77bESZBu+VnEfIoFgBynBkhfzzkUSzaLT+LHXrwmWeeICJNYO4hIl1QVHIVi7Z6KiouCXODwlU65pygcETFJal0TCIifcP8TPrOLzgC0/beQopUlqfb2IBPxdsUqQzT9t6CX3CEmiLUPOYJItIE5h4i0gVFKVexaKunZgTegjSPJ0I5kcoEzAi8pdIxiYj0DfMz6bOiMitC3ZgniEgTmHuISBcUpVzFoq0einiVgPMPY/I8eyUnaTIB5x/G4GF0gkrHJSLSF8zPpM+K0qwIdWKeICJNYO4hIl1Q1HIVi7Z6KCBUAgOxSC1jG4hF2Hq5aM5qISJSN+Zn0mdFaVaEOjFPEJEmMPcQkS4oarmKRVs9FHw/WuVXHeTSZAKCH0SrZWwioqKO+Zn0VVGbFaFOzBNEpAnMPUSkC4parmLRVs+8T5FCoubbBCWxSUhMkap1H0RERQ3zM+mzojYrQl2YJ4hIE5h7iEgXFMVcxaKtnomMTYR6rjn8RwDwJDZRzXshIipamJ9JnxW1WRHqwjxBRJrA3ENEuqAo5ioWbfVMqlRWpPZDRFRUMD+TviqKsyLUhXmCiDSBuYeIdEFRzFUs2uoZI8PC+ScvrP0QERUVzM+kr4rirAh1YZ4gIk1g7iEiXVAUcxWzop6xszaHejrG/Uf0//0QEVHuMT+TviqKsyLUhXmCiDSBuYeIdEFRzFUs2uoZc2ND2FqZqXUfttZmMDc2VOs+iIiKGuZn0ldFcVaEujBPEJEmMPcQkS4oirlK949eKc/cHGzUukKzW00btYxNRFTUMT+TPiqKsyLUiXmCiDSBuYeIdEFRy1Us2uqh/s62al2heUBTW7WMTURU1DE/kz4qirMi1Il5gog0gbmHiHRBUctVhVq0TUyRIvx5PG5I3iD8eXyRWMVXF9mXtUCrGqVVfvXBQCxCqxqlUcPGQqXjEhHpC+Zn0ldFbVaEOjFPEJEmMPcQkS4oarlKJAiCWhfsjXiVgIBQCYLvR0MSl6S0OrAIgK2VGdwcbNDf2Rb2ZZmoC0tUXBLaLT+LFBUuymFsKMbJia6orObZMkRERRnzM+mjiFcJaL/inNrGPznRpUgVBJgniEgTmHuISBcUpVyltpm2UXFJGLghFO1XnIN/aCQiPyvYAoAAIDIuCf6hkWi/4hwGbghFVFySukKidCpbmWGel6NKx5zv5chftkREBcT8TPqoqM2KULfKVmbw7cI8QUSFi8coRKQLilKuUkvRdkeYBO2Wn0XIo1gAyLGfhPzxkEexaLf8LHaESdQRFn3Gx8kWkzvUVMlYUzo4oI8T+xAREamCj5MtJrSp/ukvBbwhhvmZdMXC7nVhqOKiraFYhIXd66p0TG1x//BGvDn7p0rGYp4gotziOSQR6YKikqtUXrT1C47AtL23kCKV5bn5b5pMQIpUhml7b8EvOELVoVEmxrjZY3GPujA2FOf5w2AgFsHYUIwlPepitFsNtcRHRKSvIvatxNvjq2FkIMrz7EPmZ9JF6pgV4etZu0jO4Fq7di3mzZuH6V4NFMdxzBNEVFjSn0My9xCRtioKuUqlPW13hEkwbe8tVQ2HJT3q8spbIYmKS4LHvK14Z1YRBmJRtgV3+eOtapTGwu51i+TJEBGRJm3btg39+/fH2rVr0bnXQMwIvIXzD2OYn0kv+AVHYNnxBwUeJ/6cP7ram2Djxo0Qiwt17V21CgwMhLe3N0aPHo1ff/0VIpEIUXFJzBNEVOjS5x6xCMhuzhZzDxFpSlRcEoatPYGIBEOdy1UqK9oWpUa/+uj169eoWLEipi1aAbG9K4IfREMSm8nCcdZmcKtpgwFNbYtcfzgiIm1w584dODk5oWfPntiyZQtEok9XhRULezI/kx7YESbB3KBwSGVCnu7cMhCLYCgWYb6XI9IiLmDAgAEYOXIk/Pz8FN8lXXbhwgW0a9cOXbt2xfbt2zMUo5kniEgTIl4lYPLv+3Ht+QcUK1WeuYeItE6nTp0QKzVC57ELdeo4SWVF24EbQhHyKDbPLRGyYyAWoXk1a/gPd1bZmJS5X375BdOnT8fz589hbW0NAEhMkeJJbCJSpTIYGYphZ20Oc2NDDUdKRFR0vX//Hk5OTjAwMEBoaCjMzc0zfR7zM+kDVcweXb9+Pb7++mtMmTIFS5Ys0enCbXh4OFq2bIn69evj6NGjMDY2zvb5zBNEVJh69uyJ2NhYHDp2krmHiLTK06dPYWtri3Xr1uGrr74CoDvHSSqJKOJVAs4/jFHFUErSZALOP4zBw+gEralyF0WCIGDDhg3o3r27omALAObGhnCsYKnByIiI9IcgCPj666/x9OlTXL16NcuCLcD8TPqhspUZ/Ic7F2j26FdffYXExERMmDABxYsXx5w5cwr1PahKVFQUOnbsCFtbW+zbty/Hgi3APEFEhUcQBFy8eBHDhw9n7iEirfPnn3/CxMQEvXv3VmzTlVylkqJtQKgkxxkQ+WUgFmHrZQl8VbwwBf0nNDQUd+7cwYoVKzQdChGR3lq9ejV27NiBnTt3wsHBQdPhEGkN+7IW8PVyhC8c8zUrYvz48UhMTMTMmTNhbm6OSZMmFVLkqhEXFwd3d3cYGBjgyJEjsLTU/hMMItIv//77L169eoUWLVpoOhQiIiWCIGDTpk3o1asXSpQooelw8kwlRdvg+9FqKdgCn2bbBj+Ihi9YtFWX9evXo0qVKmjbtq2mQyEi0ktXrlzBxIkTMW7cOKUrwESkLL+zImbMmIHExERMnjwZZmZmGDlypBqiU73k5GR4eXkhOjoaFy9eRIUKFTQdEhFRBhcvXgQANGvWTMOREBEpu3DhAh4+fIg//vhD06HkS4GLtu9TpJDEJakilixJYpOQmCItlP4SutLXQlXev3+PnTt3YvLkyUVqZWUiIl0RGxuLXr16oVGjRli6dKmmwyEqshYsWID3799j1KhRMDc3x6BBgzQdUrakUin69u2LGzdu4PTp05yBT0Ra6+LFi6hTpw5KlSql6VCIiJRs2rQJVatWhYuLi6ZDyZcCVyMjYxOhnjm2/xEAPIlNVFu/CUWvtPvRkMRl0ivNygxuDjbo72wL+7JFq7furl27kJiYiKFDh2o6FCIivSOTyTBgwAAkJiZi165dMDIy0nRIREWWSCTCihUrkJSUhKFDh8LMzAze3t6aDitTgiBg9OjROHjwIPbv3w9nZy7KS0Ta6+LFi2jVqpWmwyAiUvL+/Xvs2rULU6dO1dlJigUu2qZKZaqIQyP7yc2qxAKAyLgk+IdGYvOlJxlWJdZ1GzZsQIcOHWBra6vpUIiI9M7ChQtx7NgxHDlyBJUrV9Z0OERFnkgkwtq1a5GUlIS+ffvC1NQUHh4emg4rg/nz52PdunXYuHGjVsZHRCQXFxeHO3fuYNq0aZoOhYhIyV9//YWkpCQMHjxY06HkW4FLzUaGhVOtVvV+doRJ0G75WYQ8igWAHHvyyh8PeRSLdsvPYkeYRKXxaMLdu3cREhKC4cOHazoUIiK9c/LkScyZMwdz5syBu7u7psMh0hsGBgbYvHkzPD090bNnT5w6dUrTISlZt24dfH19sXDhQt4JRURa79KlSwDARciISOts2rQJbdu21elJigWuhNpZm0OkikiyJSD+6UMIgmoaMfgFR2Da3ltIkcryvIBamkxAilSGaXtvwS84QiXxaMqGDRtgbW0NLy8vTYdCRKRXnj59ir59+6J9+/aYPXu2psMh0jvFihXDjh070Lp1a3h5eSkW0dG0ffv2YeTIkRgzZgxnrRGRTrhw4QLKlSuHqlWrajoUIiKFiIgInD9/HsOGDdN0KAVS4KKtubEhbNXcKiDt7Uu0cG6MihUrYvjw4dizZw/i4+PzNdaOMAmWHX+gkriWHX+AnTo64zY1NRV//vknBg4cCGNjY02HQ0SkNz5+/Ig+ffrAxMQEW7duhYGBgaZDItJLxsbG2Lt3Lxo3bozOnTvj2rVrGo3nwoUL6Nu3L3r06IEVK1ZAJFL/tAgiooK6ePEiWrRowZxFRFpl8+bNsLS0RLdu3TQdSoGopOeAm4MNDMTqSdIGYhGGuDfB6dOn0b9/f1y+fBne3t4oXbo0WrdujZ9++gm3b9/O1SzcqLgkzA0KV2l8c4LCERWXpNIxC8PBgwfx+vVrtkYgIipkU6dOxZUrV7Br1y6UKVNG0+EQ6TUzMzMcPHgQtWrVQocOHXD79m2NxBEeHo4uXbrA2dkZ/v7+vJhDRDohNTUVYWFhaNmypaZDISJSSEtLw5YtWxTrF+gykaCCngMRrxLQfsU5VcSTqZMTXVDDxkLx9ydPnuDIkSM4fPgwTp06heTkZFSuXBmdO3dG586d0aZNGxQvXjzDOAM3hCLkUWyeWyJkx0AsQvNq1vAfrlur+np4eCA2NhaXL1/WdChERHpjz5498Pb2xq+//opx48ZpOhwi+r83b97Azc0NL1++xPnz52Fvb19o+46KikLz5s1hZWWFc+fOwdLSstD2TURUEJcvX0azZs1w5coVODk5aTocIiIAwLFjx9CxY8cikZtUMtPWvqwFWtUorfLZtgZiEVrVKK1UsAUAOzs7jBw5EgcOHEBcXByOHj2K7t274+TJk+jatSusra3RoUMHrFixAg8ePIAgCIh4lYDzD2NUWrAFPvW4Pf8wBg+jE1Q6rjo9ffoUR48e5SxbIqJC9ODBAwwdOhS9evXC2LFjNR0OEaVTqlQpHD9+HKVKlULbtm0RGRlZKPt98+YNOnbsCLFYjCNHjrBgS0Q65eLFizAzM0P9+vU1HQoRkcLGjRvh6OiIxo0bazqUAlPJTFvgU+uBdsvPIkUqU8VwAABjQzFOTnRF5Tz0zI2IiMDhw4dx+PBhnD17FikpKahevToqd/sOkcWqQHXR/cdALMJA5yrw9XJUw+iqt2DBAixatAgvXrxAiRIlNB0OEVGRl5SUhKZNmyIlJQVhYWHMvURa6tmzZ3BxcYFIJMK5c+dQoUIFte0rOTkZHTp0wJ07d3Dx4kV88cUXatsXEZE69OjRA2/evEFwcLCmQyEiAgDExcWhfPnyWLhwISZNmqTpcApMJTNtAaCylRnmqbhoOd/LMU8FWwCwt7fH+PHjcezYMcTGxuLAgQPo0KEDHn0wU0vBFvg02zb4QbSaRlctmUyGjRs3onfv3iwaEBEVAkEQMGrUKDx8+BB79uxh7iXSYhUrVsSpU6eQkpKC9u3b4/Xr12rZT1paGvr164dr167h0KFDLNgSkc4RBAEXLlxAixYtNB0KEZHC9u3bIZPJMGDAAE2HohIqK9oCgI+TLSZ3qKmSsaZ0cEAfJ9sCjWFubg5PT0/8tPw3iC3Uu9iLJDYJiSlSte5DFc6cOYPHjx/jq6++0nQoRER6YcOGDdiyZQt+//131KlTR9PhEFEO7OzscPLkScTExMDd3R1v377N8TWJKVKEP4/HDckbhD+Pz/aYUBAEjB49GgcOHMCuXbvQtGlTFUZPRFQ4Hj58iNevX7NoS0RaZePGjfDw8EDZsmU1HYpKGKp6wDFu9ihd3Bhzg8IhlQl56iFrIBbBUCzCfC/HAhds04uMTYRqO9lmJAB4EpsIxwra3Ytsw4YNcHBwQPPmzTUdChFRkXfjxg2MGTMG33zzDQYOHKjpcIgolxwcHHDy5Em4urqiU6dOOHHiRIZFbiNeJSAgVILg+9GQxCUpHWuKANhamcHNwQb9nW1hX/a/9Rl++OEH/P7779iwYQM8PT0L5w0REanYxYsXIRKJ0KxZM02HQkQEALh58yauX7+OOXPmaDoUlVFZT9vPRcUlYUbgLZx/GAMDsSjb4q388VY1SmNh97p5bomQkxuSN+i+JkSlY2YmcGRzNLAtpfb95NebN29Qvnx5/PDDD5gyZYqmwyEiKtLevn2LRo0aoWTJkrh48SJMTEw0HRIR5VFYWBjatm2Lxo0b49ChQzA1NS3QMe7RPQEYMWIEFixYgJkzZxbiOyEiUq2vv/4aoaGhuHnzpqZDISICAEycOBHbtm3D06dPUaxYMU2HoxIqn2krV9nKDP7Dnf+bhfAgGpLYTGYhWJvBraYNBjS1RQ0bi6yGKxAjQ5V2gdD4fvIrICAAaWlpGDRokKZDISIq0gRBwODBgxEXF4eTJ0+yYEuko5ycnHD48GG4u7ujZ8+e6DfLDz8cuQ/p/wu1Od1RJn885FEs2iw7jVdH9mD06NGYMWOG2mMnIlKnixcvwtXVVdNhEBEBAFJTU7F161YMHjy4yBRsATUWbeXsy1rA18sRvnBEYooUT2ITkSqVwchQDDtrc5gbqz0E2FmbQwSotUWC6P/70QZZ/Zw3bNiALl26FJneHkRE2mrZsmUICgpCUFAQqlatqulwiKgAWrZsif3796PPvI24c+BuvsZIkwmQCoBVx7Go2c4eIpFIxVESERWe2NhY3L17lxegiEhrHDx4EDExMRg6dKimQ1Ep9VdM0zE3NtRIz1dzY0PYWpkhMi5JffsQknH9yiU0a9YMhoaF+mMFkHNftbLmBogs3QSj+3oVemxERPrk3LlzmD59OqZNm4YuXbpoOhwiUoEYy5qwaNm/QGPIC7W/nIxA2RImKl2/gYioMIWEfGo9yEXIiEhbbNq0CU2aNIGjo6OmQ1Ep7b6fX4XcHGxgIFbPrAaRIMPbuxfh4uICGxsb9O/fH9u3b0dcXJxa9pdeVFwSBm4IRfsV5+AfGonIzwq2wKcZxi8T02DR0AMLrgMDN4QiSo0FbCIiffXy5Uv06dMHrVq1wg8//KDpcIhIBaLikjA3KFylY84JCuexGBHprIsXL6J8+fKws7PTdChERHjx4gUOHz5c5GbZAnpUtO3vbJtj37H8EkRinP59Hq5cuYKxY8fi3r176NevH2xsbODq6oqlS5fi7t27UPWabzvCJGi3/CxCHsUCyLmvmkhsAOBTX7V2y89iR5hEpfEQERVFiSlShD+Pxw3JG4Q/j0diijTT50mlUvTt2xcAsH37do3cdUFEqjcj8Jaih62qSGUCZgTeUumYRESF5eLFi2jZsiVbvRCRVvD394eRkRF8fHw0HYrKiQRVVxK12MANoQh5FKvS4q2BWITm1azhP9xZafuzZ89w+PBhHDx4ECdOnEBycjKqVasGT09PeHp6wsXFBcbGxvner19wBJYdf1DQ8DG5Q02McbMv8DhEREVJTi1nbK3M4OZgg/7OtrAv+2kRzenTp2Pp0qU4ffo0XFxcNBI3EalWxKsEtF9xTm3jn5zooraFeImI1CElJQWWlpZYsmQJxo8fr+lwiEjPCYKA2rVro2HDhggICNB0OCqnV0XbqLgktFt+FilSmcrGNDYU4+REV1S2MsvyOcnJyThz5gwOHjyIAwcOICoqCsWLF0eHDh3g6emJzp0752lxsB1hEkzbq7rZGUt61GVfNSIifPo9MSPwFs4/jIGBWJTtRT75461qlEabEtEY1tsLS5Yswffff1+IEROROvkGhcM/NFItd2sZiEUY6FwFvl5Fq/caERVtly5dQvPmzREWFobGjRtrOhwi0nPynHTixAm0a9dO0+GonF4VbQHNFzwFQcDt27dx8OBBHDx4EJcuXYIgCGjSpIliFm79+vWzvNVEU4VnIqKibkeYBHODwiGVCXkq0BiIAOnHVNi+voxzmxbzVkGiIsR1abBaF7KtYm2Gs5Pd1DY+EZGqLVu2DHPnzsXbt29RrFgxTYdDRHpuxIgROHbsGB4/fgyxuOh1gC167ygHPk62mNyhpkrGmtLBIc8zVEUiEerWrYvp06fj4sWLePXqFf7880/Y2dlh2bJlaNiwISpVqoRvvvkGBw4cQFKS8okC+6oREameX3AEpu29hRSpLM8z6tIEQGRQDFHlXbDqzEM1RUhEhe19ihQSNS8WJolNyrJPNhGRNrp48SKcnZ1ZsCUijUtMTMSOHTswZMiQIlmwBQC9XCVljJs9Shc3zt+MKrEIhmIR5ns5qqSlQJkyZTBw4EAMHDgQHz9+xIULFxRtFNatWwcTExO0adMGnp6ecGzWFucfxhR4n59Lkwk4/zAGD6MT2FeNiPTOjjBJwXuE/3927bLjD1CmuDFbzhAVAZGxiVD37WgCgCexiXCsYKnmPRER5U9iihRPYhORKpXByECMi6FX8c3wIZoOi4j0jFIuMhTDztoce/fuRUJCAoYMGaLp8NRG79ojpJff3oULu9ctlFYCDx48ULRROH/+PCxaD0OJhh6A2EDl+2JfNSLSR2w5Q0RZuSF5g+5rQtS+n8CRzdHAtpTa90NElFvZLcgqCAJsTEXwaGCntCArEZGq5bQ4tMGHNzB78xh7lkwosrlIr4u2cooPwoNoSGIzWSXc2gxuNW0woKmtxmaivn37Fm2Xn0dsqvqmfLOvGhHpm4EbQhHyKFaliwwZiEVoXs0a/sOdVTYmERW+8Ofx8Fh5Qe37OTS2JWfaEpFW0PZJTUSkH/KSi8QQIIOoyOYiFm0/k9mUa3NjzXeReJ8iRV3fY2q9TU8E4Lavu1a8XyIidYt4lYD2K86pbfyTE13YcoZIhyWmSFFHzcdegICNnUqhZVMnGBkZqXVPqqKtx8pEVDD5XpD1/+0D53k5woftoYiogJiLlPEI6zPmxoZaOduBfdWIiFQrIFSS45Xb/DIQi7D1soQtZ4h0mLmxIWytzBCpxsXI0t6+RFvXLjA1NUWLFi3g6uqK1q1bw8nJCcbGxmrbb17ldHuirZUZ3BxseKs0kY7yC47Id3//tP8XVqbtvYWY9ykY42av4uiISF8wF2VUNJdXK4JSVdhvURv2Q0SkacH3o9VSsAU+HTQEP4hWy9hEVHjcHGxgIBapZWwDsQhDOzrjypUrmDdvHoyNjbF06VK0atUKJUuWRNu2bTF//nycO3cOKSkpaokhJ1FxSRi4IRTtV5yDf2gkIj8r2AKfLvpHxiXBPzQS7Vecw8ANoYhSY6GbiFRLJQuy/t+y4w+wM0yikrGISL8wF2WO7RF0BPuqERGpDlvOEFFuFHYblbS0NPz99984c+YMzp49i3PnziE+Ph7GxsZo1qyZYiZu06ZNYWJiora4AN6eSKQPuCArEWkD5qKscaatjrCzNod65nn8R/T//RARFXWF2XKGiHSXfVkLtKpRWuWzbQ3EnxbM+LzvtYGBARo1aoRJkyYhKCgIsbGxuH79OhYtWgRLS0v89ttvcHNzQ8mSJeHq6oq5c+fi9OnTSE5OVml8fsERmLb3FlKksjzfkZAmE5AilWHa3lvwC45QaVxEpFozAm9BquK7jqQyATMCb6l0TCIq2piLssairY6Q91VTJ1trM84IIyK9wJYzRJRbC7vXhaGKi7aGYhEWdq+b4/MMDAzQoEEDTJw4Efv27UNMTAz+/vtvLFmyBNbW1vDz80Pbtm1RsmRJtGrVCrNnz8bJkyeRmJj/C0a8PZFIP0S8SsD5hzEqbxWVJhNw/mEMHkYnqHRcIiqamIuyx6KtDlF3XzW3mjZqGZuISNsYGRbOr7/C2g8RqU9lKzPMU/GigvO9HPN1u55YLEa9evUwfvx47N27F69fv8bNmzexbNkylC1bFmvXrkX79u1RqlQptGjRAjNnzsTx48fx/v37XI0fFZeEuUHheY4rO3OCwtnjlkgLyRdkVQf5gqxERDlhLsoezyZ1SH9nW7UumjOgKfuOEZF+YMsZIsqLEq9vI/78VpWMNaWDA/qoqNerWCxG3bp1MXbsWOzevRuvXr3C7du3sXz5clSsWBF//PEH3N3dUapUKTRr1gzTp0/H0aNHkZCQ+awT3p5IpD+4ICsRaQPmouyxaKtDCruvGhFRUcWWM0SUW2FhYejVqxdcSydjYTdHGBuK83wsZiAWwdhQjCU96mK0Ww01RfqpiOvo6IjRo0dj165dePXqFe7cuYPffvsNVapUwaZNm9CpUyeUKlUKzs7OmDp1Kg4fPox3797x9kQiPfI+RQqJmmfAS2KTkJgiVes+iEi3MRfljEVbHaPJvmpEREUJW84QUU7+/fdfeHh4oF69eti2bRv6Odvh5ERXNK9mDQA55hD5482rWePkRFeVzbDNLZFIhFq1amHkyJHYsWMHXrx4gXv37mHVqlWoXr06/P394eHhgVKlSsHru6UQCerpw10Ubk8kKkq4ICsRaQPmopxxCpCOkfdVm7ZXdbeZ5bevGpEuSUyR4klsIlKlMhgZimFnbc5ZkHquv7MtNl96opax2XKGSPe9fv0aHTt2RKlSpXDgwAGYmX06VqpsZQb/4c6IeJWAgFAJgh9EQxKbpHTSIcKn2fZuNW0woKmt1tzNJBKJ4ODgAAcHB3zzzTcQBAEPHz7EmTNn8PN9Cwgi9cznkN+e6AvV9gYmovzR1QVZeTxPVLToai4qTMxwOsjHyRYx71MKtrKvIAAikUr7qhFpG8UJ9f1oSOIyOaG2MoObgw36O9vCvqx2nFBT4ZG3nAl5FKvS24ENxCI0r2atNUUaIsq7xMREeHp6IiEhAZcuXYK1tXWG59iXtYCvlyN84aizhQSRSAR7e3uUt62Khb7H1Lov+e2JuvBzISrqdGlBVh7PExVdupSLNIVHTTpqjJs9Shc3xtygcEhlQp4KDiIISJOmYtiX5mrtq0akKVFxSZgReAvnH8bAQCzK9PshAIiMS4J/aCQ2X3qCVjVKY2H3upx1rmcWdq+LdsvPqrRoy5YzRLpNKpXCx8cHd+7cwdmzZ1G1atUcX2NubAjHCpaFEJ16FObtibr8cyIqKuQLsqrze1/QBVl5PE9U9OlCLtI03S03E3ycbPPVV61FjdKwC/fH+hnD8ebNG7XHSVSYdoRJ0G75WYQ8igWAHItx8sdDHsWi3fKz2BHGnnv6RN5yRpXYcoZIdwmCgFGjRuHo0aPYvXs3GjZsqOmQCgVvTyTSL9q+ICuP54n0g7bnIm3Aoq2Ok/dVOzHBBQOdq6CKtRk+L92KAFSxNsNA5yo4OdEFW4c3xbY/ViIxMREjR46EIKh7bgVR4fALjsC0vbeQIpXleeZkmkxAilSGaXtvwS84Qk0RkjbycbLF5A41AaDA+ZAtZ4h024IFC/DHH39g/fr1cHd313Q4hYa3JxLpH3UuyCrI0vD86gn8/PPPePnyZZ5ey+N5Iv3CxaGzJxJYsStycttXbefOnfDx8cHWrVvRv39/DURKpDo7wiQqXaBvSY+6LL7pkb///huth8+EVftvIDIwzNNJgoFYBEOxCPO9HPmZIdJhmzZtwrBhw7BgwQLMnDlT0+EUqsQUKer4HlP77Ym3fd11erYLUVES8SoB7VecU9v4X0btx7G/tiAtLQ2dO3fGkCFD4OnpCSMjoyxfw+N5Iv2j7lx0cqKLTq81wsvdRZC8r1oD21JwrGCZ5cFxnz590L9/f4waNQqRkZGFHCWR6kTFJWFuULhKx5wTFI6ouCSVjknaKSEhAb1790Y14QWOj2+V55YzzatZ4+REV54UEOmwI0eO4Ouvv8a3336LGTNmaDqcQsfbE4n0j3xBVlXPcDMQi9CqRmkEbV2HFy9eYOXKlXj16hV69uyJChUqYPz48bhx40aGu5t4PE+kn9Sdi3S5YAuwaKv3/Pz8YGlpicGDByMtLU3T4RDly4zAW5CqcCEpAJDKBMwIVN2VftJOgiDgm2++wcuXL7Fr1y7UKF8qzy1n/Ic7s4ctkQ67evUqevXqBQ8PD/j5+UEkUs8tetqOtycS6Z/p7atCJv2o0nZ56RdktbKywsiRIxEaGorbt29j2LBh2LVrFxo2bIj69etj+fLliI6OBsDjeSJ9trB7XRiq+BikqCwOzfYIhDNnzqBNmzZYvHgxvv/+e02HQ5QnvJ2CCmL9+vX4+uuvsX37dvj4+GT6nNy2nCEi3fPo0SM0a9YM1apVw6lTp2Bmpr8XYPj7lEi/REVFwcPDAy/MqsLcbYTKxs2pJYFUKsWxY8ewefNmBAUFQSaToW33/rhXvZfKYvgc8w+R9mN7lMxxpi2hdevWmDx5MmbNmoW///5b0+EQ5UlAqEStM4O2Xubqs0XVrVu3MHbsWIwYMSLLgi2Q+5YzRKRbXr9+jY4dO8LS0hIHDhzQ64ItwNsTifTJtWvX4OzsjHfv3uHMhoWKBVkLKjcLshoaGsLDwwN//fUXnj9/juXLl+NJscoQZOq565PH80S6If3i0AVVlBaH5kxbAgCkpKTA2dkZHz9+xNWrV2FqaqrpkIhyxXVpMCLV2KuqirUZzk52U9v4pBmJiYlo3LgxihUrhtDQUOY8Ij2TlJSENm3a4PHjx7h06RKqVaum6ZC0QlRcEtotP4sUqUxlYxobinFyoivbyBBpif3796Nfv36oU6cOgoKCULZsWQCfZrnNDQqHVCYU+oKsPJ4nIjlN5iJtxJm2BAAwNjZGQEAA/v33X0yfPl3T4RDlyvsUKSRqXlxAEpuExBSpWvdBhW/06NGIiorCrl27WLAl0jNSqRQ+Pj64ffs2Dh8+zIJtOpWtzDDPy1GlY873cmTBlkgLCIKAX375Bd27d0enTp0QHBysKNgCn2a5nZzoWugLsvJ4nojSk+eiRhWLAwDEyL5wW9QXh+Y9nqTg6OiIJUuWYMKECfDw8ED79u01HRJRtiJjE3NI4QUnAHgSmwjHCpZq3hMVli1btmDLli34888/8cUXX2g6HCIqRIIgYPTo0Th8+DAOHjyIRo0aaTokrePjZIuY9ylYdvxBgccqSrcnEukyqVSKsWPHYu3atZg6dSoWLlwIsTjj/K3KVmbwH+6MiFcJCAiVIPhBNCSxSUrH2yIAttZmcKtpgwFNbQvc+oTH80T0ucpWZmieeh0HN/2G8X67EfIkXu25SFuxaEtKxo4di4MHD2LIkCG4efMmrK2tNR0SUZZSVXj7pjbsh9Tv7t27GDVqFIYMGYKBAwdqOhwiKmQ//vgj1q1bh02bNqFjx46aDkdrjXGzR+nixrw9kagIePfuHXr37o1Tp07hjz/+wFdffZXja+zLWsDXyxG+cFT7gqw8nieizOzYsQPtnetikXdDAPq7OHTRf4eUJ2KxGJs3b0bdunXx7bffYteuXRCJ1LPIE1FBGRkWToeXwtoPqVdycjJ69+6NKlWqwM/PT9PhEFEh27x5M2bPno0ffvgBQ4YM0XQ4Ws/HyRYtqpfGjMBbOP8wBgZiUbbFW/njzatZY2H3umyJQKQFJBIJPD09IZFIcOTIEbRr1y7PY8gXZFUXHs8T0ecePXqE0NBQbN++XbFN3blIW7FoSxlUrFgRv//+O3r37g1/f38MGjRI0yERZcrO2hwiQK23VIn+vx/SfePHj8e///6LK1euwNyc/6ZE+uTo0aP46quvMGLECMycOVPT4egMTdwqTUSqcfXqVXTp0gUmJiYICQlB7dq1NR1Spng8T0Sf27lzJ8zMzNClSxdNh6JxIkEQ1N1ChnTU4MGDERgYiJs3b8LOzk7T4RBliqvNUm5s374d/fr1w/r16zF8+HBNh0NEhejatWtwdXVFmzZtsHfvXhgacs5CQejr7YlEuiQwMBD9+/fHl19+if379ystOKaNeDxPROnVq1cPtWvXVpppq694jwBl6bfffoOVlRUGDhyItLQ0TYdDlCk3B5scV7fNLwOxCG41bdQyNhWeiIgIjBgxAv369cOwYcM0HQ4RFaJHjx6hc+fOqFOnDnbs2MGCrQrIb09sYFsKjhUsWbAl0iKCIGDZsmXo2bMnPD09ERwcrPUFW4DH80T0nzt37uDmzZvw8fHRdChagUVbypKlpSX8/f1x8eJFLF26VNPhEGWqv7NtnhZIyYs0mYABTbmIii778OEDevfujfLly2Pt2rXs0U2kR2JiYtCxY0dYWlriwIEDMDNjj1UiKro+fvyIkSNHYsqUKZg2bRp27NgBU1NTTYeVKzyeJyK5HTt2wNLSkgvG/h+LtpStVq1aYerUqZg9ezauX7+u6XCIMrAva4FWNUqr/Oq8gViEVjVKszefjps8eTLu3r2LXbt2wcKC/5ZE+iIpKQldunRBfHw8jh49ijJlymg6JCIitYmPj4enpyc2bNiADRs2YOHChRCLdedUn8fzRAR8ultg+/bt6NGjB4yNjTUdjlbQnUxOGjNv3jzUrVsX/fv3R1KS+noNEeXXwu51YajigzxDsQgLu9dV6ZhUuPbs2YNVq1Zh+fLlqF+/vqbDIaJCIpVK0bdvX9y6dQuHDh1CtWrVNB0SEZHaREZGokWLFggNDcXRo0d1thUUj+eJ6Pr163j48CH69u2r6VC0Bou2lCMjIyMEBATgyZMnmDp1qqbDIcqgspUZ5nk5qnTM+V6OqGzFW2l11aNHjzB8+HB4e3vj22+/1XQ4RFRIBEHAmDFjcOjQIfz1119o3LixpkMiIspWYooU4c/jcUPyBuHP45GYIs31a69cuQJnZ2ckJSXh0qVLaNu2rRojVS8ezxPR9u3bYWNjAzc3LhwoJxIEQT3NY6jI8fPzw9ixY3HkyBH2FyGt5BccgWXHH+T79YIgQCQSYUoHB4x2q6HCyKgwpaamomXLloiJicGNGzdgaWmp6ZCIqJD8+OOPmDVrFjZu3IihQ4dqOhwiokxFvEpAQKgEwfejIYlLQvoTchEAWyszuDnYoL+zLezLZn5r/549ezBw4EDUr18f+/btg41N0Vhsq6DH83I8nifSLTKZDFWqVEHXrl3h5+en6XC0Bou2lGuCIKBTp074559/cOvWLZQuXVrTIRFlsCNMgrlB4ZDKhDwtaGAgFgGyNEQf8cPOhePRuXNnNUZJ6vTdd9/Bz88PISEhnGVHpEc2b96MoUOHYv78+Zg9e7amwyEiyiAqLgkzAm/h/MMYGIhF2R6ryh9vVaM0Fnavq5gxKggCli1bhu+//x59+vTBpk2bdGbBsdwqyPG8oViE+V6O6OPExceIdMn58+fh4uKC8+fPo2XLlpoOR2uwaEt58vz5c3z55ZdwdXXF7t27uRI7aaX8HhAv6OqIMUP74sKFC7h27Rr7IOqgoKAgdO3aFStWrMD48eM1HQ4RFZJjx47B09MTQ4cOxe+//87jEyLSOgUtRM7zckTP+uUxZswYrFu3DjNmzMAPP/ygUwuO5YUqCtxEpDtGjx6NoKAgREZGFtm8lh8s2lKe7d27Fz179uSth6T1FLeePYiGJDaTW8+szeBW0wYDmtoqVpV9+/YtGjdujOLFiyMkJARmZjzo0xUSiQT169eHi4sLAgMDWbQh0hPXr1+Hi4sL3NzcEBgYCENDQ02HRESkRFW3/Fs/C8HN7T9h3bp1enMelp/jeSLSLVKpFBUqVMDgwYOxdOlSTYejVVi0pXwZNmwY/vrrL/zzzz+cjUg6ITFFiiexiUiVymBkKIadtTnMjTM/sb958yaaNm0Kb29vbNmyhcU/HfDx40e4urri2bNnuHHjBqysrDQdEhHlQ15yNQA8fvwYzZo1Q5UqVXD69GmYm5sXYrRERDnbESbBtL23VDbe0NpGmDuwvcrG0yV5/R1BRLrh+PHjcHd3x9WrV9GoUSNNh6NVWLSlfElISEC9evVQvnx5nD17NtNZLfylSrps27Zt6N+/P/z8/DB69GhNh0M5mDZtGn7++WecO3cOzZo103Q4RJQH+V2QJyYmBi1atIBMJkNISAjKlClT6LETEWUnKi4J7ZafRYpUprIxjQ3FODnRlS0AiKjIGDp0KC5evIj79+9zwtRnWLSlfLt48SJcXFwwf/58zJw5E4BqVkIl0hbjx4/H6tWrcfbsWTRv3lzT4VAWjhw5gs6dO+Onn37ClClTNB0OEeVSQfoVWpsA7dq1w7///ouQkBBUr169ECMnIsqdgRtCEfIoNk89bHNiIBaheTVr+A93VtmYRESakpKSAhsbG4wfPx7z58/XdDhah0VbKpBZs2ZhyZIl2HfyAnY8BBvFU5GSmpqKNm3a4PHjx7h27RrKlSun6ZDoM8+ePUP9+vXh5OSEgwcPsmk9kY4o6II8NlFncP0vP5w5cwZOTk7qC5SIKJ8iXiWg/Ypzahv/5EQX9nAlIp23b98+dO/eHeHh4ahdu7amw9E6LNpSgXz8+BENvMfgvUMnGBgWQ1oePk3pV0L1cbJVX5BEBfDixQs0bNgQNWvWxMmTJ1GsWDFNh0T/J5VK0bZtWzx8+BD//PMPSpcuremQiCgXCrwgjyAAIhG6VAFWfuuhusCIiFTINygc/qGRKp1lK2cgFmGgcxX4ejmqfGwiosLk4+ODu3fv4p9//tF0KFqJU5KoQH6/8ATva3tBEBvmqWALAGkyASlSGabtvQW/4Aj1BEhUQOXLl8dff/2FkJAQTJ06VdPhUDrz58/HhQsXsH37dhZsiXTEjjBJwVdQ/3+vswORwM4wiQqiIiJSveD70Wop2AKfzqOCH0SrZWwiosLy/v17BAUFoW/fvpoORWuxaEv5lv7Eq6DNopcdf8ATL9JaLVu2xM8//4zly5dj586dmg6HAJw8eRILFizA/Pnz4eLioulwiCgXouKSMDcoXKVjzgkKR1RckkrHJCIqqPcpUkjUnJsksUlITJGqdR9EROp04MABJCcno0+fPpoORWuxaEv5whMv0jdjx45Fv379MHz4cISHq/azT3nz8uVLDBgwAO3atcP06dM1HQ4R5dKMwFuQqnjWmVQmYEbgLZWOSURUUJGxiVB3D0IBwJPYRDXvhYhIfbZv346mTZuiatWqmg5Fa7FoS/nCEy/SNyKRCOvWrUO1atXQvXt3xMfHazokvZSWloYBAwZAJBLB39+fC48R6YiIVwk4/zBG5bcKp8kEnH8Yg4fRCSodl4ioIFKlsiK1HyIiVXvz5g2OHj0KHx8fTYei1Xi2S3nGEy/SV+bm5ti7dy+io6MxePBgyGQ8UC5sixYtwunTpxEQEICyZctqOhwiyqWAUAkMxAVrpZQVA7EIWy+zxRIRaQ8jw8I5zS6s/RARqVpgYCCkUil69+6t6VC0GrM85RlPvEif1ahRA1u3bsX+/fuxZMkSTYejV86ePYu5c+di9uzZaNOmjabDIaI84II8RKRP7KzNoZ6zpf+I/r8fIiJdtH37drRu3Rrly5fXdChajUVbyjOeeJG+8/T0xOzZszFr1iycOHFC0+HohdevX6Nfv35wcXHBnDlzNB0OEeUBF+QhIn1jbmwIWyszte7D1toM5saGat0HEZE6vHr1CqdPn0bfvn01HYrWY9GW8oQnXkSfzJ07Fx06dEDfvn0RGRmp6XCKNJlMhkGDBuHjx48ICAiAgYGBpkMiojzggjxEpI/cHGzUeneiW00btYxNRKRuf/31F8RiMXr06KHpULQei7aUJzzxIvrEwMAAAQEBsLCwQM+ePfHhwwdNh1RkLV26FEePHoW/vz8qVKig6XCIKI+4IA8R6aP+zrZqvTtxQFNbtYxNRKRuO3bsgLu7O6ytrTUditZj0ZbyhCdeRP+xsrLC3r17ER4ejjFjxmg6nCIpJCQEM2fOxLRp0+Du7q7pcIgoH7ggDxHpI/uyFmhVo7TKZ9saiEVoVaM0athYqHRcIqLCIJFIcPHiRfj4+Gg6FJ3Ao1vKE554ESlr0KAB1q5diw0bNuCPP/7QdDhFSmxsLHx8fNC0aVP88MMPmg6HiPKJC/IQkb76sVsdQJYGQVDdjFtDsQgLu9dV2XhERIVp586dMDExQdeuXTUdik5gZYzyhCdeRBkNHjwYI0eOxJgxY3DlyhVNh6P1ElOkCH8ejxuSNwh/Hp9pD2tBEDB06FAkJiZi+/btMDTkQhtEuooL8hCRPvr48SPmTx2P6MN+EIlUdwY138sRldWcU4mI1GX79u3w9PSEhQXvFsgNHt1SnshPvCLVuBgZT7xIFy1fvhzXr1+Ht7c3rl27hjJlymg6JK0S8SoBAaESBN+PhiQuSak3tgiArZUZ3Bxs0N/ZFvZlLbBixQocOHAAQUFBqFy5sqbCJiIVcXOwgX9opFr6O3JBHiLSNm/fvkXPnj1x/vx5rF+/Hu8q18Sy4w8KPO6UDg7o48RetkSkm+7fv48bN25g1qxZmg5FZ7AyRnnGEy+ijIyNjbF79240atQIPj4+OHbsGGeHAoiKS8KMwFs4/zAGBmJRpnlDABAZlwT/0EhsvvQEX5YphhMLf8Z3332HLl26FH7QRKRy/Z1tsfnSE7WMzQV5iEibPH78GB4eHnj58iVOnDgBV1dXAEDp4saYGxQOqUzI03mUgVgEQ7EI870cWbAlIp22Y8cOWFhYoFOnTpoORWewPQLlGVdCJcpcpUqVsHPnTpw9e5ZXDwHsCJOg3fKzCHkUCwA55g354zdffUDZYX5o4D1a7TESUeHggjxEpA8uX74MZ2dnpKam4tKlS4qCLQD4ONni5ERXNK/2abX0nPKh/PHm1axxcqIrC7ZEpNMEQcCOHTvQvXt3mJqaajocncGiLeUZT7yIsta6dWssWbIES5Yswd69ezUdjsb4BUdg2t5bSJHK8n6RR2wAGBTDrKC78AuOUE+ARFToFnavC0MVHztwQR4i0hZ//fUX3Nzc4ODggMuXL8PBwSHDcypbmcF/uDNOTHDBQOcqqGJtlmG9EBGAKtZmGOhcBScnusB/uDN72BKRzvvnn39w7949+Pj4aDoUnSISVLmUJemNqLgktFt+FilSmcrGNDYU4+REVx6UkM4TBAF9+vTBkSNHEBYWhi+++ELTIRWqHWESTNt7S2XjLelRl7NLiIoI5gciKmoEQcDixYsxY8YM9OvXDxs3boSxsXGuX5+YIsWT2ESkSmUwMhTDztqc63sQUZEzbdo0rF+/Hi9evECxYsU0HY7OYNGW8o0nXkRZS0hIgLOzMwAgNDQ0y9Uxi9qBOi/oEFF2BEFAp8m/4p6RPT51tM7/zNspHRww2q2GymIjIsqr1NRUjBw5Ehs3bsTcuXMxd+5ciESqvaOAiEjXCYKAqlWromPHjli7dq2mw9EpulsZII3zcbJFzPsUroRKlAkLCwsEBgbCyckJw4YNw65duxQH8RGvEhAQKkHw/WhI4pKQ/sqZCICtlRncHGzQ39kW9mV1q13IjMBbkKq457VUJmBG4C34D3dW6bhEVPiWLVuGY798j2+X/onT8aW5IA8R6ay3b9+iZ8+eOH/+PP78808MHDhQ0yEREWmly5cvIzIykq0R8oEzbanAdoRJ8rUSqiBLgxgCFnk3gA9PvKiICgwMRI8ePbB06VL0GTYKMwJv4fzDGBiIRdl+X+SPt6pRGgu719WJWaYRrxLQfsU5tY1/cqILe14T6bBt27ahf//+mDVrFn744QdExSUV6ZxIREXX48eP4eHhgZcvXyIwMFBpwTEiIlI2btw47NmzBxKJBAYGBpoOR6ewaEsqkZ8Try9KAqcWDceqn+ZjxIgRhRcsUSGbPn06Vh29gbKdx0KG7L8fn5PPKpvn5aj1Fzd8g8LhHxqZ94XHcsFALMJA5yrw9XJU+dhEpH6nT59Gx44d0b9/f2zcuFHp9mHF3QcPoiGJzeTuA2szuNW0wYCmtrxwQ0Qad+nSJXTt2hWWlpY4dOgQatasqemQiIi0VlpaGipWrIi+ffti+fLlmg5H57BoSyqV1xOvESNGYNu2bbh58yaqVaumqbCJ1Oq3Uw/wy8kICIJQoD5nkzvUxBg3exVGplquS4MRGZektvGrWJvh7GQ3tY1PROpx8+ZNtGrVCs2aNcOBAweyXXyiqPX5JqKiZdeuXRg0aBCcnJwQGBiI0qVLazokIiKtdvr0abRt2xahoaFo0qSJpsPROSzaktrk5sQrISEBX375JSpXrowzZ85ALBZrKFoi9dCXBfvep0hR1/cY1PkLRQTgtq87CzhEOkQikaBZs2YoV64czpw5k+WijERE2kwQBCxatAgzZ85E//79sWHDBhgbG2s6LCIirZJZDWjCmJE4ffo0Hj58yIUa84FnvqQ25saGcKxgme1zLCwssHnzZri5uWHFihX47rvvCik6IvWLikvC3KBwlY45JygczauX1rp+jpGxiWot2AKf1pl/EpuYY14hIu3w5s0bdOrUCUZGRjh06BALtkSkUfmdyZ+amopvv/0WmzZtwty5czF37lwWHoiI/i+nRbalxZqhdm9nPIx+r3OLbGsDzrQlrfDdd99h9erVuH79OmrXrq3pcIhUYuCGUIQ8ilVpj1cDsQjNq1nDf7izysYsiOTkZLx48QLn70Rh7sX3at9f4MjmaGBbSu37IaKC+fDhA9zd3XH79m2EhITAwcFB0yERkR7KqZhga2UGNwcb9He2zbSY8ObNG3h7e+P8+fPYsGEDBg4cWGixExFps7ysayQWATIBXFA2H1i0Ja2QnJyMhg0bwtzcHJcuXcq23x2RLoh4lYD2K86pbfyTE13UuiCPvBj74sULPH/+XPHf9P//4sULvHnzBgBQzKYqKgxbqbZ45A6NbcmZtkRaTiaTwcfHBwcOHMDp06fRrFkzTYdERHomP4skf15MePToETw8PBAdHY3AwEC4uLgUVvhERFptR5gEc4PCIZUJRXaRbW3B9gikFUxNTfHnn3+iWbNmWLRoEebMmaPpkIgKJCBUkuNJQn4ZiEXYelkCXy/HPL/2w4cPmRZgsyrGypmYmKBChQooX748KlSoAEdHR8X/V6hQASVLl4X3jki197S1szZX4x6ISBUmTZqEPXv2YM+ePSzYElGhS19MAJDjsZj88ZBHsWi3/CzmeTmiivQZunbtCktLS1y+fBn29tq7ECwRUWHyC47AsuMP8vXatP8XeaftvYWY9ylavci2tmDRlrSGk5MTZsyYgR9++AEeHh5o1KiRpkMiyrfg+9FqKdgCn37ZBT+Ihi/+K9pmVYz9fFtOxdjatWsrCrHpi7KWlpY59m+zPf4akXFJannPAGBrbcZFyIi03C+//IIVK1Zg1apV6Natm6bDISI9o6piQsLFbajt4IB9+/bB2tpaxVESEemmHWGSfOfYzy07/gBlihtr5SLb2oTtEUirpKamomnTpkhJScG1a9dgYmKi6ZCI8ux9ihR1fY+pd2EuQUD1G6vx6pkEz58/z7QYm77oKv//z7eVLFlSZYtp+AaFwz80Um2ziwc6V8nX7GIiKhw7d+6Ej48Ppk2bhkWLFmk6HCLSMzvCJJi295bKxlvgVQsDmlVT2XhERLosKi4J7ZafRYpUprIxjQ3FODnRlT1us8GiLWmd27dvo1GjRhg3bhyWLl2q6XCI8iz8eTw8Vl5Q+34cJUGoUTrz4qwqi7G5pet9fIko/86cOQN3d3f06dMHW7Zs4crqRFSoWEwgIlIvfVhkWxvxPlPSOnXq1MEPP/yAadOmoWvXrmjZsqWmQyLKk1QVnjBkZ8GixWhgW6pQ9pUb9mUt0KpGabX9MmfBlkg73bp1C926dYOLiwvWr1/Pgi0RFboZgbcUPWxVRSoTMCPwFosJRKT3Il4l4PzDGJWPmyYTcP5hDB5GJ/BcLwtiTQdAlJlJkyahefPmGDx4MN6/f6/pcIjyxMiwcFJrYe0nLxZ2rwtDsWoLNoZiERZ2r6vSMYlINZ4+fYpOnTqhatWq2LNnD4yMjDQdEhHpGXkxQdXtmdIXE4iI0ktMkSL8eTxuSN4g/Hk8ElOkmg5JreSLbKuDfJFtyhxn2pJWMjAwwObNm1GvXj1MmTIFa9as0XRIRLlmZ20OEaDWnrai/+9H21S2MsM8L0eV9pSb7+XIWxOJtNDbt2/RqVMnGBgY4NChQyhRooSmQyIiPSQvJqirp/7WyxL21CciRLxKQECoBMH3oyGJS1I61xMBsLUyg5uDDfo728K+bNGaNVrYi2zTf7RvmhbR/9WoUQNLly7F2rVrcezYMU2HQ5Rr5saGsFVzkdHW2gzmxtp53c3HyRaTO9RUyVhTOjhwRVEiLZSSkoLu3bvj2bNnOHr0KCpUqKDpkIhITxVGMYGI9FdUXBIGbghF+xXn4B8aicjPCrbAp8k6kXFJ8A+NRPsV5zBwQyii4pI0Ea7KvU+RQqLm9yKJTSrys5Xzi0Vb0mojR45E+/btMXz4cLx580bT4RDlmpuDjVpvIXGraaOWsVVljJs9FveoC2NDcZ5/DgZiEYwNxVjSoy5Gu9VQU4RElF8ymQxDhgzBpUuXEBQUhFq1amk6JCLSUywmEJE67QiToN3yswh5FAsAOV4gkj8e8igW7ZafxY4w3b/tPzI2Ua13kAKfit5PYhPVvBfdxKItaTWRSISNGzfi/fv3GDdunKbDIcq1/s62ap31MaCp9s8+9XGyxcmJrmhezRoAcizeiv9/ONC8mjVOTnTlDFsiLfX9999j586d2LZtGxcLJSKNYjGBiNTFLzgC0/beQopUlufzujSZgBSpDNP23oJfcISaIiwchbXIdmHtR9do5721ROlUqlQJK1euxKBBg9C9e3f06NFD0yER5ci+rAVa1SiNkEexKi3eGohFaF7NWmdW16xsZQb/4c7/9YB6EA1JbMYeUAbJcSiR+BS7fhyjM++NSB/9+uuv+Pnnn/Hbb7/x9zERaRyLCUSkDjvCJFh2/IFKxlp2/AHKFDfW+gkpKSkpiIyMxKNHj/D48WPFfx+8TgKaj1b7/rVxkW1tIBIEQd0XJ4kKTBAE9OjRAxcuXEB4eDhsbLT71nAi4FP/o3bLzyJFhQf6xoZinJzoqtMLcyWmSPEkNhGpUhmMDMWwszbHlg3rMG7cOLx48QJlypTRdIhElIndu3ejd+/emDx5Mn766SdNh0NEhPDn8fBYeUHt+zk0tiUcK1iqfT9EpHlF9RxOEAS8fPlSUZBNX5x99OgRnj17Bnl50NDQEFWqVEG1atVQuWoNnLLqjE9TbdRDBOC2r7vWrtmiSSzaks6Ijo5GnTp10KJFC+zduxcikfqSBpGq7AiTYNreWyobb0mPulp/lTY/Xr9+jfLly2PlypUYOXKkpsMhos+cP38e7du3R48ePbB161aIxZwNQUSal5giRR3fY2ptkcBiApF+GbghVG13S/oPd1bZmJl5//59lkXZJ0+eIDk5WfFcGxsbVKtWDVWrVlX6b7Vq1VCxYkUYGv6X81yXBiNSjf3DjVPjMdXxAzp37oxy5cqpbT+6iEVb0imBgYHo0aMHtmzZgkGDBmk6HKJc8QuO+P/tNQIKcoVySgeHIr0wV+fOnfH+/XucO3dO06EQUTp37txBixYt0LBhQxw+fBjGxsaaDomISEHdxYQq1mY4O9lNbeMTkfaIeJWA9ivUdy5ycqJLgVrBSaVSPH36NENRVv7f169fK55ramqqVIhN//92dnYoXrx4rvfrGxQO/9BItazZIoIAs6dXcTdgPgRBQKNGjeDh4QEPDw80btxYayYKZHa3aGFczGPRlnTOoEGDsH//fty+fRuVK1fWdDhEubLy8DUsC46E2KAYBFHuf/EYiEUwFIsw38uxSM6wTc/f3x+DBg1CZGQkbG2L9nsl0hXPnj1Ds2bNUKpUKZw7dw6Wlrw9mIi0izqLCQZiEQY6V4Gvl6PKxyYq6jRV5CoITecTQRAQFxeXYZas/P8lEgmkUimAT4u2V6pUKUNBVv7/ZcuWVdndyYVRzC4pTsHRo0dx6NAhHD16FG/fvoWNjQ06d+4MDw8PtG/fvtCPQxXrstyPhiQu47ostlZmcHOwQX9nW9iXVc+6LCzaks55+/Yt6tSpg1q1auH48eNsk0BaTxAEdO7cGeFPXqLFxFW49PgNDMSibA8G5I+3qlEaC7vX1eketrmVkJAAGxsbzJ8/H1OmTNF0OER6Lz4+Hi4uLnjz5g0uXbqEihUrajokIqIMtH1mHJE+0YYiV0EUxsz9Y2Oa4cmTJ1m2MUhISFA839LSEtWrV8+0KGtra1uodz8VZtsIqVSKkJAQHDp0CIcOHUJ4eDgMDQ3RqlUrxSxcBwcHtdWCouKSMCPwFs4/jNH4eTuLtqSTjh8/Dnd3d6xatQqjRo3K8LguXtWjomvXrl3o06cP9u/fDy8vr/8OZh5EQxKbycGMtRncatpgQFNbvTtJ6NOnDx48eIAbN25oOhSiIicvvxtTU1PRuXNnXLt2DRcuXICjI2eZEZH2GrghFCH/xiBNhWe2hdWDkqgo0KYiV369T5Girpp7ZAuCgKhfekH4+AEAUKxYMcWCX5nNmC1VqpQao8kbTS7Q9uTJExw+fBiHDh3C6dOn8eHDB1SrVg2enp7w8PCAq6urygrYO8IkmBsUDqlMyFOBWn6H7DwvR/io8A5ZFm1JZ40aNQpbtmzBP//8gxo1auj8VT0qmuLj4/HFF1+gadOmCAwMzPA4LzAo27dvH7p37447d+6gVq1amg6HSOfl53ejTCbDoEGD8Ndff+H48eNwdXXVSOxERLl19uptDN4ZAcGgmMpmXmnDau9EukDbilz5Ff48Hh4rL6h9P9/YxqFZrcqKBb8MDAzUvk9V0YZFtpOSknD69GnFLNyoqCiYm5ujXbt28PDwQOfOnfN9d9h/a9EUzOQONTHGzb7A4wAs2pIOe//+PerVq4fSdg6w7zsHF/6N1dmrelR0jRkzBlu2bMGdO3fYgzkXUlJSULZsWYwbNw7z58/XdDhEOqsgM15W/TQfS5Yswc6dO9G7d+9CjJqIKO+Cg4PRo0cPlGnaDan1vVU2bn6KCUT6RhuLXPl1Q/IG3deEqH0/gSObo4Gt9sygzStV/ZurYpFtQRBw+/ZtRQE3JCQEMpkMDRo0ULRRcHJyylVhXBsK0plh0ZZ02g/bTuOP629hUMwIAnJ/VV3brupR0XTlyhU0bdoUv/zyCyZMmKDpcHTG8OHDcfbsWURERLBnNVE+FGTGi0iWhleH/fDD0E7MW0Sk9TZt2oQRI0bAzc0Nu3btwtYbr7WmmEBU1GlrkSs/4uLicOTyLcw8/17t+zo0tiUcK+j2wq4FnV2trkW24+LicOzYMRw8eBBHjx5FXFwcSpcujU6dOsHDwwPu7u4oWbJkhtdpsvVDTli0JZ2luMIjCEABCjvacFWPih6pVAonJyeIRCJcuXIFhob62/Igr06dOoV27dohNDQUTZo00XQ4RDqloLMfBEGASCTi70Yi0moymQyzZs3CokWLMGLECPj5+aFYsWIAtLeYQFSUaHORKyuCIODFixe4c+cO7t69q/Tn1atXEBUzQeXv/lLrpBERgNu+7kWiHZ629zFOS0vD5cuXFbNwb968CQMDA7Rs2VIxC7dWrVoQiUSFushaXrFoSzqpKF3Vo6Jp+fLlmDRpEkJDQ+Hk5KTpcHRKWloaKlWqBB8fHyxfvlzT4RDpDP5uJCJ9kJycjMGDB2P37t1YunQpvvvuuwxFFm0vJhDpOm0ucqWlpeHJkyeKgmz6Iu27d+8AAEZGRqhZsyZq1aqFWrVqoXbt2qhVqxbGHI2B5E2yKt5OpqpYm+HsZDe1ja8JurLItkQiUSxmdurUKSQnJ8POzg4uXfrgrHkrte335ESXAr1vFm1J5+jiVT3SLxKJBLVr18bQoUOxcuVKTYejkyZMmICdO3fi6dOnOtWcn0hT+LuRiPTBq1ev0LVrV9y8eRMBAQHo3r17ts/XlWICkS6JeJWA9ivOqW383Ba5UlJSEBERoTRj9s6dO3jw4AE+fPgAAChevLiiMJu+QFu1atVM74T0DQqHf2ikSovRcgZiEQY6V4Gvl6PKx9YWiSlShEui0crVDYt/XIAR/bpr5azi5ORknDlzBocOHcLBF6YQarSCSKz6c05V/JuzaEs6R5uv6hEBQLdu3XDlyhXcvXsXlpa63a9IU65cuQJnZ2ecOnUKbdq00XQ4RFqPvxuJqKgLDw+Hh4cHUlJScODAATRu3DhPr09MkeJJbCJSpTIYGYphZ22ulcUEIm1X2IXNhIQE3Lt3L0NLg3///RdpaWkAgNKlSytmy6b/U6lSpTy1O9CWgrQuEwQBxYoVw8qVKzFy5EhNh5Mjl6XBkMQlqW38gs6u5m8p0ikRrxJw/mGMysdNkwk4/zAGD6MTinwSJfXav38/9u/fj127drFgWwBOTk6oXr06tm3bxqItUQ74u5GIiroTJ07A29sbdnZ2OHjwICpXrpznMcyNDXV+8R8ibRB8P1otBVvg07HH/rAIvD39h6I4GxUVpXi8cuXKqFWrFjp37qxUnC1durRK9m9f1gKtapRW24VwfTieEolEsLKyQlxcnKZDydH7FCmi1FiwBQBJbBISU6T5vkgoVnE8RGoVECqBgVg9jcENxCJsvSxRy9ikH96/f4+xY8eiU6dO8Pb21nQ4Ok0kEqFv377YvXs3UlJSNB0OkVbj70YiKsrWrVuHTp06oUWLFrhw4UK+CrZEpBrvU6RqnZUIAHEfDXH4+CmYmpqif//++PPPPxEWFoaEhARIJBIcO3YMy5cvx4gRI9CqVSuVFWzlFnavC0MVH1cZikVY2L2uSsfUZtbW1jpRtI2MTYS6Ww8IAJ7EJub79Szakk5R91W94AfRahmb9MPcuXMRExODVatWqXXVUX3Rr18/xMfH4+jRo5oOhUir8XcjERVFaWlpmDx5Mr755huMHDkSQUFBsLAo+rPUiLRZYRS5RCIRAk9ewL59+7Bo0SIMHDgQjRs3RvHixdW8508qW5lhnor7zs73ctSrNQJ0ZaZtqgrXglDXfli0JZ1RGFf15FPXifLq77//xq+//oo5c+agatWqmg6nSKhVqxbq1auHbdu2aToUIq3F341EVBQlJibC29sby5cvx2+//YaVK1dmumgQERUuXShyqYKPky0md6ipkrGmdHBAHydblYylK6ysrBAbG6vpMHJkZFg4JdGC7IdFW9IZujB1nfRTWloavvnmG9SqVQuTJk3SdDhFSr9+/RAUFISEhARNh0Kklfi7kYi0UWKKFOHP43FD8gbhz+PzdOHn+fPncHV1xcmTJxEUFISxY8eqMVIiygtdKHKpyhg3eyzuURfGhuI8t6EyEItgbCjGkh51Mdqthpoi1F66MtPWztoc6r4/VvT//eQXL1eSztCXq3qke37//XdcuXIFFy9eRLFixTQdTpHi4+ODqVOnYv/+/RgwYICmwyHSOvzdSETaIuJVAgJCJQi+Hw1JXJLSBSURAFsrM7g52KC/sy3sy2be5uCff/6Bp6cnAODChQuoV6+e+gMnolyTF7nUecG4oEUuVfJxskWL6qUxI/AWzj+MgYFYlG1LKvnjzatZY2H3unrVEiE9KysrhIWFaTqMHJkbG8LWygyRarxrzdbaLN+LkAEs2pIO0aereqQ7nj9/junTp+Prr79G8+bNNR1OkWNra4uWLVti27ZtLNoSZYK/G4lI06LiknIsaAgAIuOS4B8aic2XnqBVjdIZChqHDh2Cj48PatasiQMHDqBChQqF+C6IKDd0ocilapWtzOA/3Pm/C1MPoiGJzeTClLUZ3GraYEBTW9Sw0e/+27oy0xYA3Bxs4B8aqZb1IQzEIrjVtCnQGNrzTSDKgb5d1SPdMHHiRBgbG2Px4sWaDqXI6tevH8aOHYvXr1+jTJkymg6HSKvwdyMRadKOMAnmBoVD+v+T3ZxOeuWPhzyKRbvlZzHPyxE+TrZYuXIlJkyYAE9PT2zbtg3m5sw5RNpK24tc6mJf1gK+Xo7whSMSU6R4EpuIVKkMRoZi2Fmba1WhWdPkRVtBELR+ge7+zrbYfOmJWsZOkwkY0LRg/Yw5bYJ0hvyqnjpp21U90m5Hjx7Frl278Msvv8DKykrT4RRZ3t7eAIDdu3drOBIi7cPfjUSkKX7BEZi29xZSpLI8F2/SZAJSpDJM23sLHSetwLhx4zBhwgTs3buXBVsiLdff2VYtBVtANUWuwmBubAjHCpZoYFsKjhUseZz0GSsrK3z8+BGJidq/JoJ9WQu0qlE6z32Lc2IgFqFVjdIFnnXNoi3pFDcHG5V/meS0+aoeaZ+kpCSMGjUKbdu2Rf/+/TUdTpFWpkwZdOjQAdu2bdN0KERaib8biaiw7QiTYNnxByoZ656RPb796U/8/PPPMDAwUMmYRKQ+2l7kIs2TT2jSlRYJC7vXhaGKP8+GYhEWdq9b4HFYtCWdwqt6pC0WLFiAZ8+eYfXq1Vp/y0dR0K9fP1y4cAESiUTToRBpHf5uJKLCFPW/9u4/OKry3uP45+wuxMnKoNkYRykbhBDqxKhgIXEr0DhKUabhRov8kDBaZAan4phqKkYqwcskInj1jkzxtjKpIMJtR3IvKtjKj2BAzWClyk2vJOiFTUdLSuJoTEzSze79o6YghoTAOfucTd6vf0/me74Dz/N8c77nyXOa27R8W62tMXd/kaoGB8/IBGAvNze5YF6iNW1HpiRrRX6WrTEfz8+y5UN0NG2RUHirBzeora3V6tWrVVJSoszMTNPpDAozZ87UBRdcoC1btphOBXAdaiOAeCqpPPTPM2ztEonGVFJ5yNaYAJzj5iYXzEu0pq0kzZkY1EPT7Hm2L542TrMn2rPpgaYtEg5v9WBSNBrV4sWLNXr0aC1dutR0OoPGsGHDlJ+fzxEJwBlQGwHEQ/3xFlUfOWH77v6uaEzVR07oSGOLrXEBOMetTS6Yl4hNW0m6L2+snrgtW0k+T783Q3g9lpJ8Hq26LVs/zcuwLSeatkg4vNWDSRUVFdq3b5+ee+45JSUlmU5nUJk3b57ef/99/fnPfzadCuA61EYA8bCpJuzoGdovvsMxSEAicWOTC+YNHz5cHo8n4Zq20j9eRuwsmqrQ6IAk9Tmuu6+HRge0s2iq7S8faNoiIfFWDyY0NjaquLhYCxYsUF5enul0Bp3p06froosu0ubNm02nArhS5/9W6bO9G2yJRW0E0JM9hxsdPUN7T12jI7EBOMdtTS6Y5/F4dPHFFydk01b6x2aIjQtz9MYDU1SYk670QLJOH9WWpPRAsgpz0rWzaIo2LsxxZLODFYvFnKm6QBxsORDW8m21ikRj/fsFMtolxaL613/JVmFojHMJYkBZsGCBXnvtNX344Ye65JJLTKczKN1zzz2qqqpSfX09H4ADTvH666/rRz/6ke666y7duOgXKn2l/7XR67Hk81h6PD+LBygA3/JlR0TZpb+Xkw+PlqT/Kf2h/Ek+B+8CwCn1x1u0qSasPXWNCje1fWO9sCQFA8nKy0zT/NwgZ+YPcJmZmZo5c6ZWr15tOhVbtHZEdLSpVZ2RqIb6PBoV8MelVlENkdDmTAzq+2NSVVJ5SNVHTsjrsXp9QO2+PmHEhfr9yrv0QWy2FBoYiwictXv3bm3cuFHPP/88DVuD5s6dq/Xr1+vAgQOaNGmS6XQAV3j33Xf14x//WNOnT9e6devk8/l0Q0b/a2NodEBlBdkciQCgR8eaWh1t2EpSTNLRplZlXT7c4TsBcMLYS4epND9Lpcoy1uSCO6SkpCTsTtue+JN8RmoTO20xYPT3rd6qVatUUlKi/fv3Kzc311TaSAAdHR26+uqrlZaWpr1798rj4WQZU7q6uvSd73xHs2fP1jPPPGM6HcC4jz76SKFQSFdccYV27dolv9//jevseAFgl4Phz1Sw7i3H71N5b0jjgxc7fh8AgHNuvfVWJSUlqbKy0nQqCY2mLQaks3mrF4lEFAqF1NLSooMHD+qCCy4wlC3cbsWKFVq5cqX+9Kc/KSvL3g/9oP+Kioq0ZcsW/eUvf5HX6zWdDmBMY2OjQqGQvF6v9u/fr9TU1F5/nh0vAM5H7Sefa8az+xy/z2tLbmCnLQAkuPnz56uhoUF79+41nUpCY7sYBqTurevjgxcr6/LhPT6U+nw+VVRU6OOPP9aKFSsMZIlEUFdXp7KyMhUXF9OwdYm5c+fqr3/9q6qqqkynAhjz5ZdfasaMGWptbdXrr7/eZ8NWOrvaCABnMirg/9aHWOxmfX0fAEBiG2jHI5hC0xaDWlZWlpYvX64nn3xSBw4cMJ0OXCYWi+nee+/ViBEjtGzZMtPp4GsTJ07UmDFj9NJLL5lOBTDi73//u2bNmqXDhw9r+/btuuKKK0ynBGAQ8Cf5FHT4zOtgIJkXSgAwANC0tQdNWwx6xcXFuvbaa3X33Xero6PDdDpwkU2bNmn37t365S9/qeRkPszjFpZlad68eXr55ZeZsxh0YrGYFi1apF27dmnr1q0aP3686ZQADCJ549Lk9Tiz39brsZSXmeZIbABAfNG0tQdNWwx6Q4YMUUVFherq6rRy5UrT6cAlmpub9bOf/UyzZ8/W9OnTTaeD08ydO1eff/65duzYYToVIK4effRRvfDCC/rNb36jm266yXQ6AAaZO3OC6oo680mUrmhM83ODjsQGAMRXSkqK2tvb9dVXX5lOJaHRtAUkXX311Vq2bJnKy8v13nvvmU4HLrB06VJ1dHTo6aefNp0KenDllVfq2muv1ebNm02nAsTN2rVrVV5erjVr1mjevHmm0wEwCI29dJgmZ6TavtvW67E0OSNVGWnDbI0LADAjJSVFkthte55o2gJfe+SRR3TVVVfp7rvvVmdnp+l04KDWjohqP/lcB8OfqfaTz9XaEfnG9X379unXv/61ysvLddlllxnKEn2ZO3eutm3bppaWFtOpAI57+eWXdf/996uoqEgPPvig6XQADGJlBdmyYlHFYvbtuPV5LJUVZNsWDwBgFk1be1gxO6stkOAOHjyoSZMmadmyZVq+fLnpdGCj+uMt2lQT1p7DjQo3t+nUhc+SFExJVt64NN0x4TLdMX2K/H6/3nrrLXm9XlMpow/hcFjp6enasGGDCgsLTacDOKa6ulo333yzCgoKtGnTJnk8vHMHYEYsFlNpaan+7b/eUeDW+22Lu+q2bM2eyNEIADBQ1NXVady4caqqqtLUqVNNp5OwaNoCp3nsscdUXl6ud999V9dcc43pdHCeGprbVFJ5SNVHTsjrsXo9h637evv/HdQL903XD2/4XhwzxbmYMmWKLrzwQm3fvt10KoAjamtrdcMNN2j8+PHasWOHkpKSTKcEYJBqb2/XT37yE23evFnl5eXyT7pdT71Rd95xi6eN00/zMmzIEADgFk1NTUpNTdXWrVtVUFBgOp2ERdMWOE1nZ6euu+46DRkyRDU1NRoyZIjplHCOthwIa/m2WkWisX59NMOKRTV0iE8r8rM0h10frrZu3TotWbJEn376qS655BLT6QC2amhoUCgUUkpKit58800NHz7cdEoABqm//e1vKigo0B//+Edt2LBBs2bNknTuv2t5PZZ8HkuP52exwxYABqCuri75fD49//zzWrhwoel0EhZ/XwecZujQoaqoqNAHH3ygJ5980nQ6OEdr99Rr6dZD6ohE+/2V45jlUUckqqVbD2ntnnqHMoQdZs2aJcuy9Lvf/c50KoCtPvvsM91yyy3yer3asWMHDVsAxnz44YfKzc1VfX29qqqq/tmwlaQ5E4PaWTRVodEBSerzA2Xd10OjA9pZNJWGLQAMUF6vVxdddJGamppMp5LQ2GkLnEFJSYnWrFmj9957T1dddZXpdNAPWw6EtXTrIdvicc6au82YMUNffPGFqqurTacC2KK9vV3Tpk1TbW2t9u/fr+9+97umUwIwSO3evVu33367RowYoVdffVWjRo0648/+8/sBdY0KN/Xw/YBAsvIy0zQ/N6iMtGFOpw4AMGzMmDGaNWuWnnjiCdOpJCyatsAZtLe3a8KECfL7/Xr77bfl8/m+cb21I6KjTa3qjEQ11OfRqIBf/iTfGaIhXhqa23TT03vVEYnaFjPJ59HOoqkamZJsW0zY58UXX1RhYaGOHTumYDDI3ITr9GdMdnV16Y477tCOHTu0a9cuXX/99XHOFsBAcj41cf369Vq8eLFuvPFG/fa3v+3Xjn9qMQDge7nfV8b46/XwI49SC84RTVugFzU1NQqFQiorK9PDDz98cgfB4UaFm3vYQZCSrLxxabozJ6ixl7KDwITC9TV66+Omfh+J0Buvx1JodEAbF+bYFhP2aWlp0eVXXqcf3LNMX/hHMjfhCudSL2KxmJYsWaJ169apsrJS+fn5RnIHkNjO9/fVaDSqkpISrVq1SosXL9azzz77rc0LAAD05NQadKypVbJOHpvDc1n/0bQF+vDzn/9cays264e/eEEHP/1KXo/Va0Ow+/rkjFSVFWSzOzOO6o+36OZn3nQs/s6iKfw5n8s0NLeppPKQqo+ckKJdksd7xp9lbiIeTh2T/a0XL/7Hv6ukpES/+tWvtGjRojhmDWAgOJ/1p7smtrW1acGCBdq6daueeuopPfDAA7Ks3s+pBQDAjhqEb6NpC/Rhw/6P9Nh/H5IsT68NodN1fxV3RX6W5nAealyUbqvVxppjtu6y7eb1WCrMSVdpfpbtsXFuzveL1cxN2O18xqQVi+r4a8/qwYLrVVpa6lySAAYkO2riD0YOVX5+vmpra/XSSy9p5syZDmYMABgoeC5zDk1boBdr99RrzR/qzjvOQ9MydV/eWBsyQm+mrt6jY81tjsVPDyRr70N5jsXH2WNuwm3Od0zGYjFZlqUHb87UkhsZkwDOnl01Mfb+NkXef1WvvPKKJkyYYENmAICBjucyZ3lMJwC41ZYDYVsWH0la84c6/eeBsC2x0LMvOyIKO9iwlaRwU5taOyKO3gN9Y27CbewYk91/fvzUG4xJAGfPzppoXZOvX1Rsp2ELADgrPJc5j6Yt0IOG5jYt31Zra8zHttWqweGm4mB2rKlVTv/ZQEzS0aZWh++C3jA34TaMSQCmOLH+PF39CesPAKBP/A4cHzRtgR6UVB5SxOZzUSPRmEoqD9kaEyd1RqID6j7oGXMTbsOYBGAK6w8AwBRqUHzQtAVOU3+8RdVHTtj+MauuaEzVR07oSGOLrXHxD0N98VnO4nUffBtzE27DmARgCusPAMAUalD80H0ATrOpJiyvx3Ikttdj6cV3OKfFCaMCfjnzv3aS9fV9YAZzE27DmARgCusPAMAUalD80LQFTrPncKPtb4y6dUVj2lPX6Ejswc6f5FMwJdnRewQDyfIn+Ry9B86MuQm3YUwCMIX1BwBgCjUofmjaAqf4siOisMMHX4eb2tTaEXH0HoNV3rg0R9/45WWmORIbfWNuwm0YkwBMYf0BAJhCDYovmrbAKY41tcqZ90UnxSQdbWp1+C6D0505QUff+M3PDToSG31jbsJtGJMATGH9AQCYQg2KL5q2wCk6I9EBdZ/BZuylwzQ5I9X23bZej6XJGanKSBtma1ycPeYm3IYxCcAU1h8AgCnUoPiiaQucYqgvPlMiXvcZjMoKsuWzuWnr81gqK8i2NSb6h7kJt2FMAjCF9QcAYAo1KL74VwBOMSrglzMnop5kfX0fOGNkSrJW5GfZGvPx/CyNdPgjZ+gdcxNuw5gEYArrDwDAFGpQfNG0BU7hT/Ip6HBzLhhIlj/J5+g9Brs5E4N6aFqmLbGKp43T7ImcZWsacxNuw5gEYArrDwDAFGpQfNG0BU6TNy7N9jNRu3k9lvIy0xyJjW+6L2+snrgtW0k+T7//P70eS0k+j1bdlq2f5mU4lCH6i7kJt2FMAjCF9QcAYAo1KH5o2gKnuTMnqK6oM99D7IrGND+XXZvxMmdiUDuLpio0OiBJfRaW7uuh0QHtLJrKDluXYW7CbRiTAExh/QEAmEINih+atsBpxl46TJMzUm1/c+T1WJqckaqMtGG2xkXvRqYka+PCHL3xwBQV5qQrPZD8rTN4LEnpgWQV5qRrZ9EUbVyYwxm2LsTchNswJgGYwvoDADCFGhQ/ViwWc6Y9DiSwhuY23fT0XnVEorbFTPJ5tLNoKs1AF2jtiOhoU6s6I1EN9Xk0KuDnzJwEwdyE2zAmAZjC+gMAMIUaFB/stAV6MDIlWSvys2yN+Xh+FouPS/iTfMq6fLjGBy9W1uXDadgmEOYm3IYxCcAU1h8AgCnUoPigaQucwZyJQT00LdOWWMXTxnE+KmAT5ibchjEJwBTWHwCAKdQg53E8AtCHLQfCWr6tVpForF+HbXs9lnweS4/nZ7H4AA5gbsJtGJMATGH9AQCYQg1yDk1b4Cw0NLeppPKQqo+ckNdj9boQdV+fnJGqsoJstvcDDmJuwm0YkwBMYf0BAJhCDXIGTVugH+qPt2hTTVh76hoVbmrTqZPHkhQMJCsvM03zc4N88RCII+Ym3IYxCcAU1h8AgCnUIHvRtAXOUWtHREebWtUZiWqoz6NRAT8ftAJcgLkJt2FMAjCF9QcAYAo16PzRtAUAAAAAAAAAF/GYTgAAAAAAAAAAcBJNWwAAAAAAAABwEZq2AAAAAAAAAOAiNG0BAAAAAAAAwEVo2gIAAAAAAACAi9C0BQAAAAAAAAAXoWkLAAAAAAAAAC5C0xYAAAAAAAAAXISmLQAAAAAAAAC4CE1bAAAAAAAAAHARmrYAAAAAAAAA4CI0bQEAAAAAAADARWjaAgAAAAAAAICL0LQFAAAAAAAAABf5f3/PMBNEMP7lAAAAAElFTkSuQmCC",
      "text/plain": [
       "\u001b[1m<\u001b[0m\u001b[1;95mFigure\u001b[0m\u001b[39m size 140\u001b[0m\u001b[1;36m0x300\u001b[0m\u001b[39m with \u001b[0m\u001b[1;36m5\u001b[0m\u001b[39m Axes\u001b[0m\u001b[1m>\u001b[0m"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# maybe create a reporting library? or just copy paste for now\n",
    "\n",
    "def draw_topn_patterns(df, n):\n",
    "    analysisdf = df\n",
    "\n",
    "    subcol = [c for c in df.columns if '#Sub' in c][0]\n",
    "    fig = plt.figure(figsize=(3*n-1,3))\n",
    "    fig.suptitle(f'Top {n} patterns')\n",
    "    gs = fig.add_gridspec(nrows=1, ncols=n)\n",
    "\n",
    "    for i in range(n):\n",
    "        ax1 = fig.add_subplot(gs[0, i])\n",
    "        most = analysisdf.sort_values(subcol, ascending=False).iloc[i]\n",
    "        nx.draw(patterns[most.pattern].graph)\n",
    "        plt.title(f'Rank {i+1}: {most[subcol]:.3E} times')\n",
    "    plt.tight_layout(w_pad=5)\n",
    "    \n",
    "draw_topn_patterns(result_df, 5)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.10.13"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}