{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Handling the log files : solution to exercise on N2 molecule \n",
    "\n",
    "Compare the values of the HOMO and HOMO-1 eigenvalues for the LDA and the HF run.\n",
    "Change the values of the *hgrid* and *crmult* to find the converged values. Note that a thorough description of those two (essential) parameters is provided later on when studying a [CH4 molecule](./CH4.ipynb).\n",
    "\n",
    "Note that, both in the LDA and in the HF calculation, a norm-conserving PSP is used.\n",
    "\n",
    "The results can be compared to all-electron calculations, done with different basis sets, from references (units are eV)\n",
    "(1) S.&nbsp;Hamel <i>et&nbsp;al.</i> J. Electron Spectrospcopy and Related Phenomena 123 (2002) 345-363 and (2) P.&nbsp;Politzer, F.&nbsp;Abu-Awwad, Theor. Chem. Acc. (1998), 99, 83-87:\n",
    "\n",
    "eigenvalues          | LDA(1) | HF(1) | HF(2) | (Exp.)\n",
    "---------------------|--------|-------|-------|-------\n",
    "3&sigma;<sub>g</sub> | 10.36  | 17.25 | 17.31 | (15.60)\n",
    "1&pi;<sub>u</sub>    | 11.84  | 16.71 | 17.02 | (16.98)\n",
    "2&sigma;<sub>u</sub> | 13.41  | 21.25 | 21.08 | (18.78)\n",
    " \n",
    "The results depends, of course, on the precision chosen for the calculation, and of the presence of the pseudopotential.\n",
    "As it is well-known, the pseudopotential appoximation is however much less severe than the approximation induced by typical XC functionals. We might see that, even in the HF case, the presence of a LDA-based pseudopotential (of rather good quality) does not alter so much the results. \n",
    "\n",
    "Here you can find the values from BigDFT calculation using a very good precision (*hgrid=0.3*, *crmult=7.0*). \n",
    "Note that 1 Ha=27.21138386 eV.\n",
    "\n",
    "eigenvalues          | LDA    | HF\n",
    "---------------------|--------|------\n",
    "3&sigma;<sub>g</sub> | 10.41  | 16.82\n",
    "1&pi;<sub>u</sub>    | 11.88  | 17.37\n",
    "2&sigma;<sub>u</sub> | 13.42  | 21.18\n",
    "\n",
    "How much do these values differ from the calculation with default parameters? Do they converge to a given value?\n",
    "What is the *correlation* for the N2 molecule in (PSP) LDA?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "from BigDFT.Systems import System\n",
    "from BigDFT.Fragments import Fragment\n",
    "from BigDFT.IO import XYZReader\n",
    "\n",
    "N2 = System()\n",
    "with XYZReader(\"N2\") as ifile:\n",
    "    N2[\"N:0\"] = Fragment(xyzfile=ifile)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "from BigDFT.Calculators import SystemCalculator\n",
    "from BigDFT.Logfiles import Logfile\n",
    "from BigDFT.Inputfiles import Inputfile\n",
    "\n",
    "HtoeV = 27.21138386 #Conversion Hartree to meV\n",
    "\n",
    "inp = Inputfile()\n",
    "inp['psppar.N']={'Pseudopotential XC': 1}\n",
    "study = SystemCalculator(skip=True, verbose=False) #Create a calculator\n",
    "\n",
    "#Run the code with the name scheme LDA\n",
    "inp.set_xc('LDA')\n",
    "LDA = study.run(input=inp, name=\"LDA\",posinp=N2.get_posinp(), run_dir=\"work\") \n",
    "\n",
    "#Run the code with the name scheme HF\n",
    "inp.set_xc('HF')\n",
    "HF = study.run(input=inp, name=\"HF\",posinp=N2.get_posinp(), run_dir=\"work\")\n",
    "\n",
    "#Run the code with the name scheme HF \n",
    "inp.set_xc('PBE0')\n",
    "PBE0 = study.run(input=inp, name=\"PBE0\",posinp=N2.get_posinp(), run_dir=\"work\") "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The variables *first*, *LDA*, *HF* and *PBE0* are instances of the class *BigDFT.Logfiles.Logfile* which contain all information as the total energy.\n",
    "We should also use directly this call loading the corresponding output file as"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "PBE0 = Logfile(\"work/log-PBE0.yaml\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We compare the values for LDA and HF."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "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>eigenvalues</th>\n",
       "      <th>LDA</th>\n",
       "      <th>HF</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>3$\\sigma_g$</td>\n",
       "      <td>-10.39</td>\n",
       "      <td>-16.80</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>1$\\pi_u$</td>\n",
       "      <td>-11.86</td>\n",
       "      <td>-17.36</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2$\\sigma_u$</td>\n",
       "      <td>-13.41</td>\n",
       "      <td>-21.17</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   eigenvalues    LDA     HF\n",
       "0  3$\\sigma_g$ -10.39 -16.80\n",
       "1     1$\\pi_u$ -11.86 -17.36\n",
       "2  2$\\sigma_u$ -13.41 -21.17"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from pandas import DataFrame, options\n",
    "\n",
    "df = DataFrame(columns=[\"eigenvalues\", \"LDA\", \"HF\"])\n",
    "options.display.float_format = '{:,.2f}'.format\n",
    "\n",
    "lda_evals = LDA.evals[0][0]\n",
    "hf_evals = HF.evals[0][0]\n",
    "\n",
    "df.loc[0] = [\"3$\\sigma_g$\", lda_evals[-1]*HtoeV, hf_evals[-1]*HtoeV]\n",
    "df.loc[1] = [\"1$\\pi_u$\", lda_evals[-3]*HtoeV, hf_evals[-3]*HtoeV]\n",
    "df.loc[2] = [\"2$\\sigma_u$\", lda_evals[-4]*HtoeV, hf_evals[-4]*HtoeV]\n",
    "\n",
    "display(df)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Modifications of the calculation parameters\n",
    "Then we do a convergence curve varying *hgrid* which controls the grid step of the Daubechies basis set and *crmult* the extension. The default values are:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "hgrids [0.45, 0.45, 0.45]\n",
      "rmult [5.0, 8.0]\n"
     ]
    }
   ],
   "source": [
    "print('hgrids',LDA.log['dft']['hgrids'])\n",
    "print('rmult',LDA.log['dft']['rmult'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*hgrids* is an array of 3 values for the x, y, and z direction. A simple scalar can be indicated for the input.\n",
    "*rmult* is composed into two multiplied factors, one for the coarse grid, and the second one for the fine grid.\n",
    "We build our script for LDA and run it (on one core, it takes 10 minutes roughly for all calculations)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "Hgrids = [0.45, 0.40, 0.35, 0.30, 0.25]\n",
    "Crmult = [3.0, 5.0, 7.0, 9.0]\n",
    "log_LDA = {}\n",
    "\n",
    "inp.set_xc('LDA')\n",
    "\n",
    "emin_LDA = 0.0\n",
    "for crmult in Crmult:\n",
    "    log_LDA[crmult] = []\n",
    "    for hgrid in Hgrids:\n",
    "        inp.set_hgrid(hgrid)\n",
    "        inp['dft']['rmult'] = [ crmult, 8.0]\n",
    "        name = \"LDA-%4.2f-%04.1f\" % (hgrid,crmult)\n",
    "        log = study.run(input=inp, name=name, posinp=N2.get_posinp(), run_dir=\"work\")\n",
    "        log_LDA[crmult].append( log )\n",
    "        emin_LDA =min(emin_LDA,log.energy)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We do the same loops to run the Hartree-Fock calculations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "log_HF = {}\n",
    "inp.set_xc(\"HF\")\n",
    "\n",
    "emin_HF = 0.0\n",
    "for crmult in Crmult:\n",
    "    log_HF[crmult] = []\n",
    "    for hgrid in Hgrids:\n",
    "        inp.set_hgrid(hgrid)\n",
    "        inp['dft']['rmult'] = [ crmult, 8.0]\n",
    "        name = \"HF-%4.2f-%04.1f\" % (hgrid,crmult)\n",
    "        log = study.run(input=inp, name=name, posinp=N2.get_posinp(), run_dir=\"work\")\n",
    "        log_HF[crmult].append(log)\n",
    "        emin_HF = min(emin_HF,log.energy)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA4UAAAG5CAYAAAA595FfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nOzdd3hc53Xv++9CJ0ASJNhAAgTA3imKgkQVqlJdoopFNUu2JRfFzomTe2+unbjk2FasY+fkHsd24hLFPVavlGxZnRLVKBIsEnsHCLB3gujArPvHHggDcNBIAANgfp/nmUeYXdfes0Fhzfu+6zV3R0REREREROJTQqwDEBERERERkdhRUigiIiIiIhLHlBSKiIiIiIjEMSWFIiIiIiIicUxJoYiIiIiISBxTUigiIiIiIhLHlBSKSFwzs1+a2T/F4Lx/MbPPnea+eWZ20swSuzqu/sjMLjKzreF7dksHti8wMzezpJ6I70yZ2e/M7Pvhny82s80xjGWKma02s3Iz+9tuOP5lZlYW8X69mV0W/tnM7LdmdtTMloeXfcXM9oc/+2FdHU9v1teeYxGJLSWFItJvmVmxmVWF/0A9Zmbvm9mXzeyTf/vc/cvu/s89HZu7X+fuv+/ItuHruDJi313uPtDdG7ovwn7lQeA/wvfs+ZYrW97frhROYtzMftZi+btmdl/45xvC74+Z2T4z+y8zG3Q653P3d9x9SheEfrq+Drzl7oPc/afdfTJ3n+Hub4XfzgeuAnLd/TwzSwZ+BFwd/uwPd3c8kbrzuTodZvaWmX0x1nGISO+kpFBE+ruF7j4IyAd+CPwD8OvYhiQ93MqZD6zvwfO1VAF81swKWlmfCXwfGANMA3KBf+2RyE5TuFUu2t8Qp32vu6BFKx8odveK8PtRQNoZxNOjLfFq0RORWFJSKCJxwd2Pu/sLwJ3A58xsJpzS9W64mf0p3GJzxMzeafzD18z+wcx2h1sdN5vZgvDyVDP7sZntCb9+bGapjec1s5vNbI2ZnTCz7WZ2bXj5J9/am9kEM3vTzA6b2SEze8TMhoTX/TeQB7wY7gL39ZbdwsxsjJm9EI55m5l9KeL83zWzJ83sD+HY15tZYWv3ycymmtlr4WNtNrM7Itb9zsx+ZmZ/Dh/rQzOb0Il9f2FmL5lZBXC5mc21pq6GT5nZExGfxTozWxixf3L43sxpJe4vha/9SPhejAkv3w6Mj7h/qS32O+X+Rqy+x8x2hc/7rYh9EszsH8Of5+Hw/c1q7Z4Cx4DfAd+JttLdH3X3l9290t2PAv8FXNTawczsbDNbFb5vTxAkPo3rWnavLDazr5nZx2ZWYWa/NrNRFnRfLjez181saMT251vQon7MzD6ycNfM8Lq3zOwhM3sPqAzf18i43gQuB/4jfC8nm1lm+Nk7aGYlZvbtiN+p+8zsPTP7NzM7Anw3yrUOCD87R81sA3Bui/XFZnalmX0B+BVwQfjcjwGN3WiPhWM7nWd0jJk9E45/p0V0iW3rd6ud56rZZ2XBvy37gN+Gj/mUmf0xfMy14fv4DTM7YGalZnZ1y+tvEdMfo5zrIeDiiM/mP1puIyJxzt310ksvvfrlCygGroyyfBfwlfDPvwO+H/75B8AvgeTw62LAgClAKTAmvF0BMCH884PAMmAkMAJ4H/jn8LrzgOMEXdoSgBxganjdW8AXwz9PDG+TGj7GUuDHrV1H+PwOJIXfvw38nCA5mAMcBBaE130XqAauBxLD17islfuVEb7O+4EkYC5wCJgRca+OhK8rCXgEeLwT+x4nSHYSgMFACfB34Xv9KaA24rP4OvBERGw3A2tbifuK8Lnmhu/hvwNL23sOOnB//wsYAJwF1ADTwuv/r/Bnnhs+338Cj7Vy7MuAMiAbOAFMCS9/F7ivlX1+3Hhfo6xLCd+3/zt83xYBdRH37TKgrMW1LSNoNcsBDgCrgLPDsb8JfCe8bQ5wOPysJBA8k4eBERHP7C5gRvgzTo4S31uEn+vw+z8Ai4FB4fu6BfhCeN19QD3w1fDxBkQ53g+Bd4AsYCywLsr1XRlxvHfb+D3p7DOaDqwE/mf4vo8HdgDXdOR3i/afu8vC1/8v4c9iQMQxrwnH+AdgJ/Ct8Of9JWBnG8/ud4E/tnL9zT4bvfTSS6/Il1oKRSQe7SH4I7OlOmA0kO/udR6Mz3KggeCPtulmluzuxe6+PbzPPcCD7n7A3Q8C3wM+E173BeA37v6au4fcfbe7b2p5UnffFt6mJnyMHwGXduRCzGwswViqf3D3andfQ9Bi8pmIzd5195c8GIP43wRJTjQ3EnS/+62717v7KuAZgsSj0bPuvtzd6wmSwjmd2Hexu7/n7qHwfknAT8P3+llgecS2fwSuN7PB4fefCccezT0E93mVu9cA3yBoMSpoZfuO+p67V7n7R8BHNN23vwK+5e5l4fN9F1hkbXT/c/d9BF84PNjWCc3sKuBzBIlINOcTJAc/Dt+3p4EV7VzHv7v7fnffTZBgfejuq8OxP0eQIALcC7wUflZC7v4aUESQ9DT6nbuvD3/Gde1cSyJBy/w33L3c3YuB/0PzZ3OPu/97+HhVUQ5zB/CQux9x91LgTMYpdvYZnUWQED/o7rXuvoPgi4K7Irbv6O9Wa0IESXlNxPW/4+6vhH/HniL4ouiH4fv9OFBg4Z4EIiJdRUmhiMSjHIIWr5b+FdgGvGpmO8zsHyFI2ghah74LHDCzxxu7JxKMAyuJOEZJeBkELRvbaYeZjQwfc7eZnSBIiIZ38FrGAEfcvbxFDDkR7/dF/FwJpLWSwOQD88JdB4+Z2TGChCu7jWMN7MS+pS3i3h1Ouk9Z7+57gPeA28J/AF9HkIRG0+wzcPeTBC1cOa1s31FtXetzEde5keCLg1HtHO9fgGvMLGriYGbnA48Ci9x9SyvHiHbfSlrZttH+iJ+roryPvK7bW3yG8wm+KGkU+Rm2ZzhNLZuRsUZ+Lu0db0yLbdq71rZ09hnNB8a02P6bNP+cO/q71ZqD7l7dYlnLz+eQNxWVakwcByIi0oU0qFlE4oqZnUvwR+m7LdeFE6u/B/7ezGYAS8xshbu/4e6PAo+GW67+k+AP/M8QtDpGFtfICy+D4A/MCbTvBwTdvGa7+2ELpk2IHPPj0XeD8LmyzGxQRGKYB+zuwHlbKgXedverumnfyOvYC+SYmUUkOC2T6N8DXyT4f9UH4ZauaBo/AwDMLAMYRsfvQVv3N5pS4PPu/l5ndgp/tj8GTql2a2ZnAy+Ej/tGG4eJdt/y6MCXDx1QCvy3u3+pjW06c68OEbS+5wMbwstaPpvtHW8vwXMR+ft1ujr7jJYSdNWcdJrn68i96uyz11IFQTfXRtmtbdgF5xKRfkwthSISF8xssJndSND96o/uvjbKNjea2UQzM4LxXw1AgwVzr11hQZGSaoJv6xu/uX8M+LaZjTCz4QTd/hoLPfwauN/MFlhQnCTHzKZGCW8QcJKgIEYO8LUW6/fToqhHo3CXuveBH5hZmpnNJui22lqrWlv+BEw2s89YUNgl2czONbNp3bDvBwT38G/MLMnMbiYYqxjpeYJxX39HMLaqNY8S3Oc54c/ofxF0kSzuQNzQxv1txS+Bh8wsHyD82d/cwX1/BFxIUGWU8P4zgZeBr7r7i+3s/wHBOLS/Dd+3T3HqfTtdfwQWmtk1ZpYYfp4uM7Pc0zlYuHXrSYJ7NSh8v/4fmn4/OuJJ4BtmNjQcx1dPJ5awzj6jy4ET4UIwA8L3ZGb4i6WO6OxzdTrWAHeFr6WQ5l1hYxGPiPRRSgpFpL970czKCb71/xbBH+X3t7LtJOB1ggTtA+DnHsyBlkpQ8OIQQXexkQTdyCCYSqAI+BhYS1DE4/sA7r48fK5/Iyhg8TYRLVoRvkeQ/BwH/gw822L9DwgSz2Nm9v9G2f9ugqISewjGiH0nPB6sU8ItjVcTjJnaE77WxiIYXbqvu9cSFJf5AkF1znsJ/mividimimDM1zhOvSeRx3oD+KfwtnsJWmfvam37KNq7vy39hKBV79Xws7UMmNeRE7n7CeB/03xM698TjBv7dbgy5EkzizqNQsR9uw84SjBmr9V70xnhLxhuJni2DxL8znyNM/tb4asErVk7CFrnHwV+04n9v0fQZXQn8Cqtjytt12k8ow3AQoLxrzsJfv9/RTCFSEd09rk6Hf9E8LwfJbhXj7ax7U8Ixr4eNbNun0NSRPoWaz4sQUREJDbM7EPgl+7+24hl/xOY7O73xi4yERGR/k0thSIiEhNmdqmZZYe7QX4OmE3QjbJxfRZBS+LDsYpRREQkHigpFBGRWJlCMNXDcYIulIvcfS8Ek9ETdF/8i7svjV2IIiIi/Z+6j4qIiIiIiMQxtRSKiIiIiIjEsbiYp3D48OFeUFAQ6zBERERERERiYuXKlYfcfUS0dXGRFBYUFFBUVBTrMERERERERGLCzEpaW6fuoyIiIiIiInFMSaGIiIiIiEgcU1IoIiIiIiISx5QUioiIiIiIxLF+nRSa2UIze/j48eOxDkVERERERKRX6tdJobu/6O4PZGZmxjoUERERERGRXqlfJ4UiIiIiIiLSNiWFIiIiIiIicUxJoYiIiIiISBxTUigiIiIiIhLHlBSKiIiIiIjEMSWFIiIiIiIicUxJoYiIiIiISBxTUigiIiIiInKGHnnkEQoKCkhISKCgoIBHHnkk1iF1mJLCGOjLD4yIiIiISLwJhULU1NRQUVHB8ePHOXz4MPv372f37t0UFxfzox/9iC9+8YuUlJTg7pSUlPDAAw/0mb/zzd1jHUO3Kyws9KKioliHAQQJ4QMPPEBlZeUny9LT03n44Ye55557YhiZSPd65JFH+Na3vsWuXbvIy8vjoYce0jMv/Z6eexHpb9yd+vr6T151dXXN3re3vK/uc7o5U35+PsXFxV37IZwmM1vp7oVR1ykp7FkFBQWUlJScsnzIkCF885vfxMw+Wdb4c7Rlrf0c6/W9KZa+FGtviqU7Yn3ppZf4/ve/T3V19SfbpKWl8Z3vfIcbb7yRhIQEzKzHXqd7PpHO0JeAEs/i4QsRd6ehoaFXJjDduU8oFIrpfU9OTiYpKemUV2vL21rXlft8/vOfjxqvmcX8njVSUtiLksKEhITT/qZBRGKvJ5PQnk56db6uPd/3vvc9jh49esozlJWVxfe//33MWv8i5UzXxXr/3hybrq37j/3ss8/yta99jaqqqk+Wp6Wl8eCDD3Ldddf1+qSno/vU19cTS4mJiT2aDPVUctXWusTExJje87YUFERv+FFLYS/Sm5LC1h6YsWPHsmHDhk/eN34ukZ9Pez/39vWKJfbnilWst956a9QvQ8yMJ554AncnFArh7j320vl0vtZeItK3JCQk9KnEpiuOl5iY2CwZl9jrC71D2koKk3o6mHj30EMPRX1gfvCDHzBw4MAYRibSffLy8qJ+GZKXl8ftt98eg4hE2tYVSehZZ51FWVnZKcfOzc1lxYoVzc4V7b+nuy7W+/fm2HRtPXPsL3/5y0RjZjz55JNdmlwlJSWRkKC6iRJ7jYlfX+02rZbCGIiHfvYikfrCt2ciXU3PvcSrvtCNTiQetdVSqK9WYuCee+6huLiYUChEcXGx/jiQfu+ee+7h4YcfJj8/HzMjPz9ffxhLv6fnXuLVQw89RHp6erNl6enpPPTQQzGKSETao5ZCEREREelS6hUl0vuo0IySQhERERERiWPqPioiIiIiIiJR9bnqo2aWAfwcqAXecvdHYhySiIiIiIhIn9UrWgrN7DdmdsDM1rVYfq2ZbTazbWb2j+HFnwKedvcvATf1eLAiIiIiIiL9SK9ICoHfAddGLjCzROBnwHXAdOBuM5sO5AKl4c0aejBGERERERGRfqdXJIXuvhQ40mLxecA2d9/h7rXA48DNQBlBYghtxG9mD5hZkZkVHTx4sDvCFhERERER6fN6RVLYihyaWgQhSAZzgGeB28zsF8CLre3s7g+7e6G7F44YMaJ7IxUREREREemjenOhGYuyzN29Ari/p4MRERERERHpj3pzS2EZMDbifS6wpzMHMLOFZvbw8ePHuzQwERERERGR/qI3J4UrgElmNs7MUoC7gBc6cwB3f9HdH8jMzOyWAEVERERERPq6XpEUmtljwAfAFDMrM7MvuHs98DfAK8BG4El3Xx/LOEVERERERPqbXjGm0N3vbmX5S8BLPRyOiIiIiIhI3OgVLYXdRWMKRURERERE2tavk0KNKRQREREREWlbv04KRUREREREpG1KCkVEREREROKYkkIREREREZE41q+TQhWaERERERERaVu/TgpVaEZERERERKRt/TopFBERERERkbYpKRQREREREYljSgpFRERERETiWL9OClVoRkREREREpG39OilUoRkREREREZG29eukUERERERERNqmpFBERERERCSOKSkUERERERGJY0oKRURERERE4piSQhERERERkTjWr5NCTUkhIiIiIiLStn6dFGpKChERERERkbb166RQRERERERE2qakUEREREREJI4pKRQREREREYljSgpFRERERETimJJCERERERGRONavk0JNSSEiIiIiItK2fp0UakoKERERERGRtvXrpFBERERERETapqRQREREREQkjikpFBERERERiWNKCkVEREREROKYkkIREREREZE4pqRQREREREQkjikpFBERERERiWNKCkVEREREROJYv04KzWyhmT18/PjxWIciIiIiIiLSK/XrpNDdX3T3BzIzM2MdioiIiIiISK/Ur5NCERERERERaZuSQhERERERkTimpFBERERERCSOKSkUERERERGJY0oKRURERERE4piSQhERERERkTimpFBERERERCSOKSkUERERERGJY0oKRURERERE4piSQhERERERkTimpFBERERERCSOKSkUERERERGJY/06KTSzhWb28PHjx2MdioiIiIiISK/Ur5NCd3/R3R/IzMyMdSgiIiIiIiK9Ur9OCkVERERERKRtSgpFRERERETimJJCERERERGROKakUEREREREJI4pKRQREREREYljSgpFRERERETimJJCERERERGROKakUEREREREJI4pKRQREREREYljSgpFRERERETimJJCERERERGROKakUERERERE5ExkZ4PZqa/s7FhH1iFKCkVERERERM7E/v0ALJ48jLR/HcbiycOaLe/tkmIdgIiIiIiISJ9y7Bhs3tz0Ag4MSGbRX9dRn1bOoq8MZvc3kxlZVRfjQDtGSaGIiIiIiEhL9fWwc2eQ9G3a1DwJPHCgabvERAAK/24M9Sm7wZz61ArO/bscSn5YEqPgO0dJoYiIiIiIxK9Dh5onfI2v7duhLqKlb/hwmDIFbrwx+O+UKfjkyawbmMgX/td1lOaXQEJDsG1iPbvy9vD5mwr4TWyuqlOUFIqIiIiISP9WWxskedFa/Y4cadouORkmToSpU+GWWz5J/pgyhZ3U8JetH/Du7lWsO/I+paW/53jZLjz5JJwd5ZyJdfz+yoNKCkVERERERHqEe1DYJVqr386d0NDQtG12dpDsLVrULPGjoID91Sf4y9ZlvFO6ko8OvUTJqh9zdFUJDcnHmvZvSGMQY5mScjEzhsxg1zsvsGLaTkiMaFlsSOb+d8fCV3vuFpwuc/dYx9DtCgsLvaioKNZhiIiIiIjImaqqgq1boyd/J040bZeWBpMnN0/6pkwJlmVmcqL6JK9uX87bJStZdfAjdlRs5nComLqUQ03HaEgmvT6X0cnjmZo5nfNGz+Gq8fOYlzuNhITmEznk/eJiSlkGifXQkEQeF1DylaU9dFPaZ2Yr3b0w2ro+11JoZuOBbwGZ7r4o1vGIiIiIiEgXc4fdu6MnfiUlwfpGublBsnfvvc2Tv7w8SEigpr6WJcUreWPHclZu/iXbVm7kQP1OapL3g4WPE0oktWE0oxOnM2nAVM7NnsOCcfO4OH8WqcnJHQq56N4Xyfn9BOoTjpDUkMmKz73QDTeme/RoUmhmvwFuBA64+8yI5dcCPwESgV+5+w9bO4a77wC+YGZPd3e8IiIiIiLSjU6ehC1bTk38tmyBioqm7TIygkTvggvgvvuat/plZABQ31DPst3reG3Hhyxf9jybX93IvrrtVCXubSoA40Zy3SiGJ4xjQuq1zB05m8vyz+PKCXMZlDbgjC5l5KAhPL1gMXe+fhdPXPk4IwcNOaPj9aSebin8HfAfwB8aF5hZIvAz4CqgDFhhZi8QJIg/aLH/5939ACIiIiIi0jeEQrBrV/RWv7Kypu3MID8/SPYuvrh5q19OTrAeCHmIj/dv4dXtH/LBX37PxqPr2V27jZOJuyGhaUxfYs1wsshn5oBLOGvYLC7JO4drJp7HyEGZ3XapN0+bT/W0svY37GV6NCl096VmVtBi8XnAtnALIGb2OHCzu/+AoFXxtJjZA8ADAHl5ead7GBERERER6Yjjx6Mnflu3QnV103aZmUGid/nlzRO/iRNhQFNrnbuz/eguXt7+Lu+tWMW6I+sprdnKCduFJ9Z8sl1CzRAyPZ+zUgqZmTWTi3PP4eqJ5zFu2MievPo+rTeMKcwBSiPelwHzWtvYzIYBDwFnm9k3wsnjKdz9YeBhCArNdF24IiIiIiJxqr4eioujJ3/79jVtl5gI48YFyd5VVzUlflOnwsiRn7T6Ndpzcj+vbH6Vd8qK+PjQOoqrNnOUEkKJlZ9sY7UDGdiQx9S0a5iROYMLxszl6onzmJGdi7U4nnROb0gKo32CrSZx7n4Y+HL3hSMiIiIiEucOH46e+G3b1nxC92HDgmTvuuuat/pNmAApKacc9mj1MV7f9Cpvlaxg9cG17KjYxOFQMfVJEVVD6waQXp/L+ORLmZIxjXmjz+aq8fM4N28CiS0qfkrX6A1JYRkwNuJ9LrAnRrGIiIiIiMSHurqmCd1bvg5FTMuQnBwkeVOmwMKFTS1+U6YESWEUJ2tP8s72Zby5czlFB9awtXwzBxt2UJsUMVF8QwqptTmMSTyXSalNFT/nF0wlLaVjFT+la/SGpHAFMMnMxgG7gbuAT3fFgc1sIbBw4sSJXXE4EREREZG+xR0OHIie+O3Y0XxC91GjgkTv1lubt/qNGwdJ0dOGmoYaPtyzhtd3LGf5vjVsPhFU/KxOPBgx3UMSydWjGZ44gwlp05g7chaX5Z/HFRNmk5me1gM3QdrTo5PXm9ljwGXAcGA/8B13/7WZXQ/8mKDi6G/c/aGuPK8mrxcRERGRfq26OujaGS35O3asabvUVJg0qXlrX+PUDkNan0KhPlTPRwc28NqO5Szbu4oNx9azp3Y7FQl7wULBRqEEEqtHkWX5FKRPYc7wWVySV8iVE+aSPWRQN98AaU+vmbze3e9uZflLwEs9GYuIiIiISJ/iDnv2RE/8ioubT+iekxMke3fffeqE7omJrZ4i5CE2H93Gazs+5P3dq1h3ZB2l1Vspt914Qn04DsOqh5Pp+cwecAmzh85kfu45XDWxkHHDs1T0pQ/qDd1HRURERESkUWVl9AndN28OJntvlJ4etPDNmwef/WzzVr+BA9s8hbuzq7yU13cu452ylXx8aB0lVVs4xi5CCbVNG1YPZVBDHlNTC5meOYMLx8zlqonnMn10toq+9CP9OinUmEIRERER6ZVCISgtjZ74lUbM1mYWtO5NmQL333/qhO7tJGbuzr7KfSzdVcRbu4pYffAjdlZs4bAX05BQ1bRhzWDS63IZl3INUzOnMy97DldOmMfcsWNJTe7XKYPQw2MKY0VjCkVEREQkJsrLT036Nm0KJnSvikjKBg9unvA1viZNajahe1uOVB/h/T0rWbJzBUUHPmJb+SYOhHZSn1DetFFdOik1OYxKGs/kQdMpHDWbBePP44KCiQxMO3UKCek/es2YQhERERGRfqehofUJ3ffubdouIaFpQvcFC5onf9nZp0zo3pry2nJW7F/Dm8XLWb43qPi5v34HNQlHmzaqTyOpagwjEs9lwsCpzB1xFpcXFHLRuKmMGJzetdcvfZ6SQhERERGRjjh6NHqr37ZtUBsxDi8rK0j0rrnm1AndU1M7fLqq+io+PrSO13d+yId717Dx2Hp2126nKuFg00YNySRUZZPFTGalT+bsEbO5ZOw5XDp+FrnDBqnoi3RIv04KNaZQRERERDqlri6Yvy9aq9/BiGQsKalpQvcbbmg+xcPw4Z06ZW1DLRuPbmJJyfJwxc/1lFZv5aTti5jrLxGrHEmmT2BS2rXMGjaL+blzuWLCHCaMGqqiL3JGNKZQREREROKLOxw6FL3Vb8cOqK9v2nbkyOhj/caNg+TkTp22PlTPtuPbeKdsJe/sKuLjI+soqdzCcXbjFp5E3hOgcjgDG8aSmzqJGUOnc0HOOVwx7mym54xQ0Rc5bRpTKCIiIiLxp6am9Qndj0aMv0tJCQq6zJwJt93WPPkbOrTTpw15iJLyEt7fs5Klu4pYfXAtOys2c9h34VYXbOQG1Vmk1eYwLrmQqUOmM2/0HK4Yfw5zxo5R0RfpUUoKRURERKTvcod9+05t8Wuc0D0Uatp2zJgg0bvzzuaJX35+mxO6t35qZ0/lHlbsW81bu4pYuf8jtp3cxMGGYhqsumnD6iEkV48hJ/kaJg+axrmj5rBg3LnMzc9l2KCOVRYV6U5KCkVERESk96usDKZxiNbqVx4x5cKAAcHk7YWFcO+9zSd0HzTotE9/sOogqw99xJLiFazY11Txs84qmjaqHURi5WiGJ17CxIHTOGfEbC4rKOS8gnGMGTpQRV+k11JSKCIiIiK9QygEu3ef2uK3eTPs2tV828YJ3T/3ueatfrm57U7o3pZjNcf4+PBa3t5VxLI9q9h4fCN7ardRY8ebNqobgFWOJot5jEufzJzhs7gkv5ALCiYxbmSmir5In9Ovk0JVHxURERHpYdnZsH//qctHjQq6eULQsrdly6ktflu2BC2CjQYODBK9+fObKns2TuiefmZz7Z2sO8n6I+t5t2wV7+9eyfqjQcXPSjvctFF9KlSOYnBoFpPTJjNr2Ezm557DRQVTmZozjJSkznc5FemNVH1URERERLpOW10kL788SP727GlalpAABQXRK5JrCzEAACAASURBVHyOHt3hCd1bU11fzaZjm1i2dzXvlhbx0eG1lFRtpZx9TRs1JEFlNul1OUHFz6wZXJgzl0sKZjE9d5iKvki/oOqjIiIiItI9amqCsX4bNwavtlRVwVVXNU/8Jk7s1ITurakL1bH1+FaK9q8JKn4e+pidFVs4GiqLmOsvASpHklIzhoKUS5g2ZAbzRp/FpQVzmDV2pIq+SNxSUigiIiIi7Tt5Mhjjt3EjbNjQlARu3w4N4Tn22mvV++CDMw6jIdTAzvKdrD74EW+XFrFq/8dsO7mJQw27cAvPL+gGVcNJqhpNTtI5TBo8nfNGncWlBWczJ380o4dkqOiLSAQlhSIiIiLS5PDhUxO/jRubF3pJSgrG9c2aBXfcAdOmwfTpQYXPjIwuCcPd2XVyFx8fWcs7u1ayYt9qtpRvYn/dThqstmnDqiwSqrIZYdcwYeBUzhl5Fpfmz+Xs/BwKRgxW0ReRDlBSKCIiIhJv3INxfS0Tvw0b4ODBpu0GDAgSvosvbkr8pk2DCRMgObnNUyyePIw7vwRP/BfcvOVwq9u5O/uq9rH28Fre372KZXtXs+nYBvbUbqfOqpo2rMnEKkYxhEsZlz6FOcNnc0neXM7Jz2fy6KEq+iJyBvp1UqjqoyIiIhLXGhpg585TE79Nm+DEiabthgwJEr6bbmqe/OXldX56h1GjOHDiCIv+uo76tHIWfWUwu7+ZzMjBWRyuPsy6I+v4cN9q3i9bxbqj6yir2UYNEfMM1mZgldlk1M9j8oDJzMqayfyxczkvfwLTc4aRkdZ2MioinafqoyIiIiJ9XctiL40tgJs3B+sajR4dJHuRid+0acF0EV04xm7sLy6izD6EhAYIGUmhQaQkJVLJ0aaN6tOgIpu02qDi58yhM7kgdy4X5k9ieu4wsgaq6ItIV1L1UREREZH+oKPFXgoKgmTv6qubEr9p04IWwS7m7pRWlPLOnvf40/YlvFjyAhUJ+6Exx0xw6q2chuMF5Pl1TBsynXmj53JR3lRm5o1Q0ReRXkBJoYiIiEhvc/hw9PF+paVN27RV7OUMJ3ZvS3V9NUUHi3h559u8vmsp606sooJDwcqGJLBQU0LYyBxLP0TJ/3ik2+ISkdOnpFBEREQkFqIVe2n8OVqxl0su6XSxlzMPMagAuqTsHf60fQnLD35IWe3mpqkfqrJIOjmOiak3cuGoC7hxynyeKX6cJw7/CyRGVAhtSOH+3K93a6wicvo0plBERESkO7Us9hKZBJZHFFhpLPbScszf6RR7OU1V9VUUHSzixW1v8mbZu2wsX00l4cqhDclQPpashinMHlLI1fmXcPW0WcweO5zkFpU/835xMaUsg8R6aEgijwso+crSHrkGEYlOYwpFREREultksZfIxK+1Yi+f/Wy3Fntpj7tTXF7MG6VL+dP2JRQd+pA9dVtxC49NrBpGSsU4pqXdwvzRF3Lj5Iu4cHIewwe1XwCm6N4Xyfn9BOoTjpDUkMmKz73QzVcjImeiXyeFmpJCREREulxjsZeWY/5iWOylIyrrK/lw/3Je2Pomb+1+l80Va6hqrAbakIKVj2WEX8ucoedy3fhLuWrqTKaOyTqtyd9HDhrC0wsWc+frd/HElY8zclBsrllEOqbT3UfNLAOodveG7gmp66n7qIiIiHRaZ4q9tOz22c3FXtrj7uws38nLO9/iLzvfYuXh5eyr39bUClg5nLSq8UxKP4tLx8znpqnzmTchh8HpqTGLWUS61xl1HzWzBOAu4B7gXKAGSDWzg8BLwMPuvrUL4xURERHpGe6we3f08X69pNhLR1TUVfDBvg95fssbLN3zPlsr11Btx4KVDSkklOeRbTdwzrBzuX78ZVw1dQYTRg3RVBAiAnSs++gS4HXgG8A6dw8BmFkWcDnwQzN7zt3/2H1hioiIiJyBzhZ7uemmmBV7aY+7s6N8B3/a/iZ/2fEWa44UsT+0LZgKAqByOBnVU5k78Gwuy7mIhdPmc964MaSnxj55FZHeqSNJ4ZXuXtdyobsfAZ4BnjEz/SsjIiIisdfHir10xMm6k7y7ZxnPbXmDd/e+x/aqj6mx48HK+lQSTuYxNvEmzh0+jxsnBq2AucMGxTZoEelTOpIU/puZPeru77e2QbSkUURERKTbRCv2smED7NjRq4u9tMfd2Xp8Ky9sfZOXi9/io2NFHArtjGgFHMHg2hmcNWguV+TO5+Zp85lbkE1KiykhREQ6oyNJ4Vbg/5jZaOAJ4DF3X9O9YYmIiIgAhw6dWuhl48boxV5mz4Y77+w1xV464mTdSZaUvsvzW97kvX3vs7N6LbUJJ4KV9akkVRQwLvkW5o04n4WTLufKKdMZmdm7r0lE+p52k0J3/wnwEzPLJyg481szSwMeAx539y3dHKOIiIj0Z/2k2Et73J3Nx7bw7ObXeLV4KWtPrOSI7wQLV4KvGMnQhrM4b/BcFoy9mFumX8SssSNPa0oIEZHO6PSUFABmdjbwG2C2u/f6/gqakkJERKQXiCz20nKqh2jFXiKneOhlxV46ory2nNdK3uH5La/zwf5llNSupS7hZLCyPo3kigIKUmdxwagLuHnSFVwxZQpDMtJiG7SI9FtnNCVFxEGSgWsJWgsXAG8D3+uSCEVERKT/aCz20jLx68PFXtrj7qw/spGnN73G6yXvsL58Jcco+aQV0CpHMSw0l1mZhVyVfym3Tr+IKaOzNCWEiPQKHZmn8CrgbuAGYDnwOPCAu1d0c2xnzMwWAgsnTpwY61BERET6n/LyoNhLy26ffbzYS0ecqD3BS9vfYvHWN1h+8EN21a2jPiH8p1F9GqlV45mWeicXZV/ArVMWcOmkyWSk9f4uriISn9rtPmpmS4BHgWfC01D0Oeo+KiIicgYii71Etv5FK/bSsttnHyj20p6Qh/j44Hqe3PgaS8reZWP5Ko7brqAV0A2rHMVIpnHW0HO4tuAybp1xEfnDM9UKKCK9yhl1H3X3y8MHMTO7Fxjv7g+aWR6Q7e7LuzZcERER6XGRxV5advvsR8VeOuJYzTEWb32TF7e+SdGh5ZQ1rKchoTJYWTeA9JoJzE77NBeNvoBF06/kovETSE3u8IgcEZFepzP/gv0cCAFXAA8C5QST15/bDXGJiIhId2gs9tIy8Wut2MtNN/XpYi/tCXmIlfvW8uTGV3ir7F02V6yhPKHsk1bAhKpssu0Czh56LteNu5RbZ17EmKGaGF5E+pfOJIXz3H2uma0GcPejZpbSTXGJiIjImaipgS1bTp3jb8uWflvspSOOVB3lmc2v86dtS1h1ZDl7GjYSSmxqBRxYO5G56Zdxyej53DZtAeePG0dSYv9JgkVEoulMUlhnZomAA5jZCIKWQxEREYmV1oq9bN8OofD/piOLvVxzTb8p9tKekIf4YPcqntrwOkv3vMeWqjVUJJYFK91IrBpNbuJFnDP4PG6YeDm3zLiAYQP79vhHEZHT0Zmk8KfAc8BIM3sIWAR8u1uiEhERkeY6U+xl9my4885+VeylIw5WHuaJ9a/yl51vseboCvaFNhFKrApW1qUzuG4S52VcyWU587lj5pWcnZtHQkL/aw0VEemsDieF7v6Ima0kmKPQgFvcfWO3RSYiItKXZWfD/v2nLh81Cvbti75Pa8VeNmwIksJGcVDspT0NoQbe3lXE0xtf5d29H7C9+iMqk/YEK91IqsohP/lSzs06j4UTr+Dm6eczaEBqbIMWEemlOjJPoXl43gp33wRsamsbERERIXpC2Lg8WrGXDRuCbqDRir3cfHO/LvbSEXtPHuCxda/wSvHbfHxsJft9E55YHaysy2Bo/WTOTruOK8ZezB0zFjBjdI6mhBAR6aCOtBQuMbNngMXuvqtxYbjIzHzgc8AS4HfdEqGIiEh/k5ER18Ve2lPfUM9rOz/k2U2v8/7+D9hR/THVyXuDlW6k1OQyIflyzss6n1smLeDGqecxIDU+WkhFRLpDR5LCa4HPA4+Z2TjgGDAASABeBf7N3dd0X4giIiK9kHvQ6ldcDCUlp/63LV/9atwUe+mI4mN7eHzdq7xWspS1x4s4ZFvwxCBptroMskJTmZdyI1fmXcKds65k0ojsGEcsItK/dGTy+mqCOQp/bmbJwHCgyt2PdXdwIiIiMRMKwd690ZO+4mLYtQuqq5vvM3RoUOVz8uSgO2hr/vVfuy3s3q62vo4/b32f57e8wbIDyyip/Zia5HBXW08grX4sU1Kv4vwRF3DrlCu5dtJcUjQxvIhIt+rUv7LuXgfs7aZYREREek59PezZc2qy1/hzaSnU1jbfZ/jwIOmbPRsWLgx+zs9v+u/gwU3bxmG3z2i2Hi7lsXWv8Maud1h/YhWHE7ZAYnBfrW4QI3wqswbcytX5l3LXrCvJGzo8xhGLiMQfffUmIiL9U10dlJWdmuxFJn0NDc33yc4OkrvCQrjttiDZa0z48vODsYAdNWpU69VH+6nquloWb3qXxVvfYPnBZZTWraM25UCw0hMYEMpjRvJ1XJh9PrdNvYorJ5xNoiaGFxGJuY5UH53h7ut7IpiuZmYLgYUTJ06MdSgiItLVamqCxK61pG/37qbJ2yFouRszJkjuLrzw1Fa+vLxgqoeu0tq0E/3Iuv07eWztqywpe4eNJ1dzLHHbJ62ACfWDGWXTOGvg7VxbcBl3zbqSUYPie+ykiEhvZe3NJGFmq9x9bvjnL7r7ryLWpbt7ZTfHeMYKCwu9qKgo1mGIiEhnVFUFCV608XwlJcF4v8j/hyUkQG7uqcleY2vf2LGQkhKDC+kfKmqreXr92/xp2xKKDn1IWcN66lMOBitDCWTU5zNxwBwuHn0ht0+7mvn5M0mIs2kzRER6MzNb6e6F0dZ1pPto5KCIvwZ+FfH+HeCcM4hNRETi1cmTzZO+lq19LbteJiUFiV1BAVx9dfOunQUFkJMTNxO3dzd3Z9Xu7Ty+/hXe3v0emytWcyJpOyTWAZAYymR04nTmDL6b68dfzp0zryArfXA7RxURkd6qI0lhZFNiy1Hz+gpQRESiO3Gi9a6dxcVw+HDz7VNSgi6cBQVBEZeWrX1jxkBiYs9eQ5w4VlnJU+uW8OcdS1h1ZAV7GjbQkHooWBlKZJAVcE7aTVwy5iLumHY188ZO18TwIiL9SEeSwmwzuw/4iFOTwrb7noqISP/kDseOtZ30HWsxc1FaWlOSV1jYvGtnfn5Q5EXdDbudu/N+8Wae3PAq7+x9j62VH3EyeccnrYBJPoSc5BmcM/Qz3DDxcm6fdjmD0wbGOGoREelOHUkKvwsUAvcDuWa2HtgUfqlutIhIf+QOhw61Pp6vuBjKy5vvk5HRlOBdeOGpSd/IkZqmIQYOlpfz5Lol/GXHW6w+uoJ9vpFQariVNpRIZuJ45qXfyuW587lj+tXMyZ6sVkARkTjTkcnrH458b2a5wGxgFvB2N8UlIiLdyT0Ys9faeL7iYqhsUUds8OCmJO+yy04t6DJsmJK+GGsIhXh72wae3vQq7+57n+1VH1OZshMS6wFIThhKftJMzh0+j4UTr+DWqZeSkZIe46hFRCTWOj1PobuXAWVmdgLI6/qQRETkjIVCQXXOaElfcTHs2gXV1c33GTo0SPCmTGkq5BLZ2jdE0wn0NruPHuPxtW/waslSPjpWxAE24qlHg5WhJIYmT+DsgYu4YuzF3DXjGqaPmBDbgEVEpFfqVFJoZnOATwN3APuBqcD/6Ia4RESkLfX1sGdP6+P5Skuhtrb5PiNGBEne7NlBIZeWE7MPVvXI3qy2voHXt3zMs5tf5/39H7Cz+mOq00ogIWgFTEkcxoSU2cwbMY+bJy1g4eRLSEtKi3HUIiLSF3Rk8vrJwF0EyeBJ4CngMncvNrOd3RyfiEh8qquDsrLWx/OVlQWJYaTs7KYiLosWNe/amZ8fjPmTPmP7gUM8ue5NXtu1lLXHizicsBlPbSzek8Sw1InMG3wnC/Iu5u7p1zAxqyCW4YqISB/WkZbCTcAKYJG7r2uxTtVHRUROR01N0JrX2ni+3buDLqCNzIIpGQoKgiIuLbt25uUF1T2l11m88V3ufP0unrjycW6eNj/qNlW1dfxl0xoWb3mDZfuXUVz7MbUDdkFCAwCpycOZnHo2F4y6gFunLODqcReqFVBERLpMR5LC2whaCl8zs9eBJ4GX3b2uWyMTEenLqqqCcXvRxvOVlATj/Tzie7WEBMjNbV7EJTLxGzs2mMdP+pQD5cdY9MbN1KccZdEbt7A7dxsjBmayce9+nlq/hNd3LWV9+SqOJm6G1OMAWGIywwdMYvaQT3N1/qXcMeNqCgaPjfGViIhIf9aR6qPPAc+ZWQZwC/BXwK/M7CVAA1BEJD6dPBkkd61N2bB/f/Ptk5KCxK6gAK655tSJ2XNyIDm5p69CulnhHxdSn3gCzKlPOsqY304hoSGdugGln7QCpqWOYHpaIRdlX8htU67i8vx5pCTqCwAREek5HS404+4VwCPAI2aWBdwOFHRTXCIisXXiROvj+UpKgjn8IqWkNI3dW7jw1Dn6xoyBxMSevgqJkbIjx7n+yS9Rau9DQrgbcEKIhtQDJNdnc3nWZ7h23KXcPu0qcgfmxDZYERGJe52ekgLA3Y8A/xl+iYj0Le5w7Fjr4/mKi4P1kdLSmhK8wsJT5+jLzg66gEpcCoWcd7Zv5xcrnuD1fS9zOHUVJFeeuqFBLRW8csdvez5IERGRVpxWUigi0qu5w+HDrY/nKy6G8vLm+2RkNCV4LQu55OfDyJGamF2aqaiu47FV7/O79U9SdGIJNRlbIKGBpPQMzk6/hCRPY0XlnyExYgh+Qwr35349dkGLiIhEoaRQRPoe92DMXmvj+YqLobJFK83gwUGCN25cUyGXyC6eWVlK+qRdxQeP8fMPX+S5Hc+zveEDPGMvAIMycrhi+Gf58py7uH78FSQlBP97zfvFxZQ2LIPEemhIIo95/Ormb8fyEkRERE7R4aTQzP4GeMTdj3ZjPCLSH2Vnn1p4BWDUKNi379TloVBQnbO18XwlJVBd3XyfrKwgyZsyBa6++tTqnUOGdPVVSRxoCIVYumUHvyh6ijf2v8yR1FWQchJSE8hNnMUNYz/P3xTew8xh06LuX3Tvi+T8fgL1CUdIashkxede6OErEBERaV9nWgqzgRVmtgr4DfCKu2ueQhFpX7SEsHH5f//3qa19paVQW9t82xEjggRv1qzmhVwai7sMVjFk6RonKmt4YtWH/Hb9U6w6sYSaQZshoZ6kjHTmZFzEPVMXcf/M2xg2YFi7xxo5aAhPL1j8yTyFIwfpywkREel9rDN5nZkZcDVwP1BIMGfhr919e/eE1zUKCwu9qKgo1mGI9E3uwUTrFRVNr8rK5u/belVWwtNPt3+e7OxTx/FFTsyekdG91ylxbdv+I/znspd4budidoQ+wAfuBmCgj+bi4VfxV3Pu5vpxC0hO0LQhIiLSN5nZSncvjLauU2MK3d3NbB+wD6gHhgJPm9lr7q6R8yKx0tAQPVHrbPLW2rpQqHPxDBgA6elBItdeMrd5c5D0paWd/vWLdFJ9Q4i3N23nFyufYcn+lzmSthpST8AAIydxJjfmfZa/nvtpZg2bgWmsqYiI9HOdGVP4t8DngEPAr4CvuXudmSUAW4EeSwrN7BbgBmAk8DN3f7Wnzi1yWtyD7pCnk5B1ZNuW4+vak5jYlLBFJm8ZGUGVzcj3ka+W27a2TcupGdr6o3ry5M7fT5HTcLSimseLlvP79c+w+uSb1A7aDIl1JA5OY0560C30vlmLGJ42PNahioiI9KjOtBQOBz7l7iWRC909ZGY3dvQgZvYb4EbggLvPjFh+LfATIBH4lbv/sLVjuPvzwPNmNhT4/wAlhXLmQqH2k7MzSd4aGjoXT1pa9IRrxIigS2VHkrPW1qWkqNKm9HvuzuY9R/iv5S/z3M4XKPYP8EGlkAQDM0eyYMSdfGn2nVw/7ipSE1NjHa6IiEjMdCYpPA7c1qIbzXFgpbuv6cRxfgf8B/CHxgVmlgj8DLgKKCMoaPMCQYL4gxb7f97dD4R//nZ4P4kXbbW2nWnyVlXVuVgSElpPukaO7FirWmsJXHp60JrXX4wa1Xr1UZEuVFvfwJsbdvDwqudYsv9ljg1YDWnHIMPISZrGjXlf56/m3M2c4WepW6iIiEhYZ5LCcwiKy7wYfn8DsAL4spk95e7/uyMHcfelZlbQYvF5wDZ33wFgZo8DN7v7DwhaFZsJF7z5IfAXd18V7Txm9gDwAEBeXl5HQpOu0Nja1lVj2VpuV1/fuXhSU6MnYS1b204ngUtNVWtbR0WbdkKkixwqr+KJoiL+sOEZVpcvoS5zEyTWkjgklbMyLuDTU2/js9MXkZ2eHetQRUREeqXOJIXDgLnufhLAzL4DPA1cAqwEOpQUtiIHKI14XwbMa2P7rwJXAplmNtHdf9lyA3d/GHgYguqjZxBb1+rsfG3dIbK1rauTt9NtbYuWlA0ffvrj2vpja5uIAEG30HWlh/jNitd5fudiSnwZPngXJDsZWcNZMHwRX5p9J9cVXMWApAGxDldERKTX60xSmAdEThxWB+S7e5WZ1ZxhHNGaW1pN5Nz9p8BPz/CcsdHWfG2NQqEgueqOKpIVFWfe2taYlA0b1jRVwOkWJlFrm4h0QHVtPa+v386vVi/mrf2vcjx9NQw4AoNhTOJUbsz/e740+y7OGTFX3UJFREQ6qTNJ4aPAMjNbHH6/EHjMzDKADWcYRxkwNuJ9LrDnDI/Z94wc2ZTUdYZZ60lXVtbpj2trXJbUqZlLRES6xN6jJ3lq1Wr+uP5Z1lS8TV3mRkiqJiErhdkZ5/Ppqbdx77TbyMnIiXWoIiIifVqH/9p39382s5eA+QQte19298YZ4e85wzhWAJPMbBywG7gL+PQZHhMzWwgsnDhx4pkeqmfcdtvpJXBpaWptE5E+z91ZXXyA369YwvPFi9nFcsjcCWlOetpQFoy4mc/PuoPr868hIzkj1uGKiIj0G+be/nC7cGGXXHcvbXfj9o/1GHAZwRQX+4HvuPuvzex64McEFUd/4+4Pnem5GhUWFnpRUVH7G/aEtpK3DnwWIiL9SWVNHa+t28GvV7/IW/tfoTxjDaQfAmB00iRuzF/IF2bewbkjzyXBEto5moiIiLTGzFa6e2G0dR1qKXR3N7PnCSqQnhF3v7uV5S8BL53p8UVEpHcrPVzO00VreGTj83xU8Rb1QzZCchUJw5OYNXAen57yTe6e8inyB+XHOlQREZG40JnBYsvM7Fx3X9Ft0cQDzdcmInEmFHJW7NjHH4uW8nzxC+y25XjmDkgPMSA9kwUjbuDzs+7kurxrGJQyKNbhioiIxJ3OJIWXE8xJWAxUEIwrdHef3R2BdYVeOaZQ87WJSBwor6rllbU7+N2aP/P2/lc5OfAjyNgPWTAqaTwL87/KfTNu5/yR55OYoKljREREYqlDYwoBzCxqPx53L+nSiLpBrxpTKCLST+08cJynV37MoxsWs7ZqKQ1D10NyJeaJTB94LvdMuY07J3+K8YPHxzpUERGRuHPGYwrDdhFUGR3v7g+aWR6QDfT6pFBERLpefUOIZdv28ujK91i88wX2Jq7Ah2yHQQ2kDRrEFSOu4b4Zd3BDwXVkpmTGOlwRERFpRWeSwp8DIeAK4EGgHHgGOLcb4hIRkV7oWEU1L3+0kz+seZm3D75C5aCPYeBeGAEjk/JZmP8VPjv9di4cdSFJCZrjVEREpC/ozP+x57n7XDNbDeDuR80spZvi6hK9ckyhiEgfs2XvUZ5duZ7HNixmXfVSQlkbIOUklp3A9IGF3D3l77hj4q1MHjI51qGKiIjIaehMUlhnZomAA5jZCIKWw17L3V8EXiwsLPxSrGMREekr6uobeG/zbh5b+SEvlLzA/qQifOg2GFJPKgOZP/IK7pt+O9fnXU9WWlaswxUREZEz1Jmk8KfAc8AoM3sIWAR8u1uiEhGRHnW4vIqX1mznj2teY+mhV6kevBYG7YZRMDwpl5vyH+Ceqbdx8eiLSU5IjnW4IiIi0oU6nBS6+yNmthJYEF50i7tv7J6wRESkO7k7G3Yf5vmVG3hs45/YWPNO0C009QSWk8C0jLO5a8pfs2jCLUwbMg0zi3XIIiIi0k06nBSaWSowF8gM73e7meHuD3ZXcCIi0nVq6upZurGMJ1at4IXiP3EopQgfuhWy6khhAPNHXsG9UxdxY/4NjBgwItbhioiISA/pTPfRxcBxYCVQ0z3hiIhIV9p/vII/r9rBYx+/wTuHXqMmcy0MLoUxkJU0mpvy7+fuKbdx6ehLSU1MjXW4IiIiEgOdSQpz3f3aboukG6j6qIjEG3fno5KDLF61icc3/pnNte/iwzbAgGMw1piScRZ3Tf4inxp/C7OyZqlbqIiIiHQqKXzfzGa5+9pui6aLqfqoiMSDqto63lxXylOrVvKnkj9zJG0VnrUZRtSSTBoXjbyMe6bcxo35N5Kdnh3rcEVERKSX6UxSOB+438x2EHQfNcDdfXa3RCYiIq0qO1zOn1fv4ImP3+K9w69TO2QtDN4FY52hSaNYmP8Z7pz0KS4fczkDkgbEOlwRERHpxTqTFF7XbVGIiEibQiGnaOc+Xly1mSc3vczWuvfxYeth4BEYCJMy/v/27jy8qurs+/j3TgiEQAwzAoEEiKJCEREV55FBIBGclRYVFKG2YAeltfVR21qH9ulDfdXGgKi0TnVABsMsg4qgyCQIGASRBGQeEwIZ1vtHDhggkMA5yd7n5Pe5rlyes/bea99nsbM9d9baa/2E2854lH5tbuC8hudpWKiIiIhU2Mkkhd8D/YE2zrk/mVkr4HRgfaVEJiJSze3LP8iMr9bz3uIlfLg+k121F5cMC22aTw1qcWnjK7i93Y2kJqXSok4Lr8MVERGRMHUySeGLQDFwDfAnYC/wHnBBJcQlIlItfbd1Nx8uXst/YI48ygAAIABJREFUl37MZztnUlBvOSSsgyRHQo2GpCbdzi1t+3Fti2upE1PH63BFREQkApxMUniRc66zmS0GcM7tNLOalRRXSGj2URHxu6LiYuZnbWLSom94Z9VU1hbPLxkWmrANEqBtnXO4NeX39G19A10adyHKorwOWURERCLMySSFBWYWDTgAM2tMSc+hb2n2URHxo915B5i27DvGLV5G5vrJ7KmzBNdwFbTYTzQxXNL4cm4780b6tOpDUnyS1+GKiIhIhDuZpPA5YBzQxMyeBG4G/lgpUYmIRJisH3by4eK1vLt0Hgt2zaSo/gpcvbXQppj46Pr0SbqJm9r0pXtid+JrxnsdroiIiFQjFU4KnXOvm9mXwLWULEfR1zm3stIiExEJYwWFRcz7ZiMTF2Xx3qrprGdBySLyDTZDA2gd145bzniIG5LSuKjJRURHRXsdsoiIiFRTJ9NTiHNuFbCqkmIREQlr2/fuZ+qy7/hg0XImb5hKbt3AsNBWuURRg4sbX8qtKb8nNSmVNqe18TpcEREREeAkk0IREfmRc46VOTvIXLyWd5fNZ+HuWRQ3WIGr9y2kFFE3OoHerdK4sU1ferTsQULNBK9DFhERETmGkkIRkZNwoKCQj1flMGnRGt5bNZOcqM9LhoU23gSNoVVcW25p+yBpyWlc0vQSakTpNisiIiL+pm8rIiLl2LI7j8lL1jF+8QqmZk9nf/wyaPg1rvU+oojmokYXc3PKb0htlcqZ9c70OlwRERGRk1JuUmhmewksQ3H0JsA5504LeVQhonUKReRUOOdY9v02Mhev5b2ln7N43xxcgxW4+mvgzELqRMdzfcvr6dv6Bq5veT0NYht4HbKIiIjIKTPnysr3IkuXLl3cwoULvQ5DRHxs/8ECZn+9gYlfruGD1bP5ocbCkmGh8TkAJNZO4qa2fUlLSuPyZpcTExXjccQiIiIiFWdmXzrnupS17aSGj5pZfeAMIPZQmXNubnDhiYh4Y+POfWQuXsuExSuZkT2T/ISvoNHXuJQ9GFFc0OgibmrzS1KTUjm73tmYmdchi4iIiIRchZNCM7sXGA4kAkuArsBnwDWVE5qISGgVFzsWfbeZDxetZdyyRXyVOxfXcAU0WIM7+yC1o+rQs2UPbmidRq+WvWhcu7HXIYuIiIhUupPpKRwOXADMd85dbWZnAU9UTlgiIqGRm1/AjOXrmbToWyas/phtNb/ENVyBS9wAQLPYRG5scx9pyWlc2exKakXX8jhiERERkap1MklhvnMu38wws1rOuVVm1q7SIhMROUXrt+7hw8Vrmbh4NbM2zeJgwldYo5UUt9sJwPkNL6Bf6/tJTUrlJw1+omGhIiIiUq2dTFKYbWb1gA+A6Wa2E9hYOWGJiFRcUXExC9b8QObitXywdAkr8z8pmSSm4Te4+geIjapNt8Ru3JCcRu9WvTk97nSvQxYRERHxjQonhc65foGXj5vZLCABmFwpUYmIlGNP3gGmfbWeSYu/ZdKqeeyMXYRr9DUuaT2Yo2lsM/q2vpu0pDSubn41tWvU9jpkEREREV86mYlmnnHOjQBwzs05VAaMqKTYRESOsOaHnXy4eB2TFq9mzqa5FNZfgTX+muJztgPQqUFn+ra+m9SkVM5reJ6GhYqIiIhUwMkMH+3GsQng9WWUiYiERGFRMfO+yeHDxWuZsPQrvjk4D9fwa6zhNxQ32k/NqFpc1+Ja0pLS6JPUhxZ1WngdsoiIiEjYKTcpNLOhwM+BNma2rNSmeODTygosFMwsFUhNSUnxOhQRAVZkb+PO5z7kjWG9aZ/YqMx9duzbz9Sl35G5eC2TVi9gT9xiaPQ1xW2+AyumUa0m3JB8J6lJqVzX4jrqxNSp2g8hIiIiEmHMOXfiHcwSgPrAU8DvSm3a65zbUYmxhUyXLl3cwoULvQ5DpFrLzS/grEf/Snbzf9Fy01BW/ukR6sTG4Jxj1cYdZC5ey8TFWczb/CmF9VcQ1XglRbFbAOhQ/yf0Tb6B1KRUujTuQpRFefxpRERERMKLmX3pnOtS1rZyewqdc7uB3cAdZnYucHlg08dAWCSFIuK9uzPGk9Pyn1BrF9mJ/6T331vTOaklE5YsZ13R57iGK4hqtJqipnnEWE2ubnE1aUmp9GnVh6T4JK/DFxEREYlYJzPRzDBgMPB+oOg/ZpbhnPt/lRKZiESMV2cv54N9f8XV3wvmcDF7mVP793y8KQ535jqcFdGgZkNSk28htVUq3RO7E18z3uuwRURERKqFk5lo5l7gIudcLhyeefQzQEmhiJzQ8OnPUJi4AqILSwqiCyF+I+7gaYw47yFSW6VyUZOLiI6K9jZQERERkWroZJJCA4pKvS8KlImIHNfn325iX/PxEH3wmG3xtWvw1IVPeRCViIiIiBxyMknhK8ACMxsXeN8XeDn0IYlIJDhYWMRf3p/PMxM+p36bfuxs/jbFUQcOb492tRh52d89jFBEREREAMqdws/MagA45/4B3EPJ5DI7gXuccyMrNzwRCUdL12+l66Ov89T4BfzsinP49ncZ9G2dCq7klmPFMfRtncY97e7xOFIRERERqciSFIucc52rKJ5KoSUpRKpGYVExf5v4BX9+/zMaxtcmfVA3enduA0BuQS5NxjYlrzCX5nEt+eb2lVpjUERERKSKBLUkBXpuUEQqYNXGHQxMn8IX3/7ArV3b8dzd19Awvvbh7XVi6pCWlMo7695hWu/JSghFREREfKIiSWFjM/v18TYGhpWKSDVVXOx4buoiHn37E+JqxfD6L3pz68Xtyty3QWwD6tWsR/sG7as4ShERERE5nookhdFAXdRjKCJHWbdlN4NemsLHq3Lo07kN/xrUjdPrHb8H0OEw061ERERExE8qkhRucs79qdIjEZGw4Zxj9KyveOg/c4iOMkYP7sGAK85RwiciIiIShvRMoYiclJwdexk8ahrTlq3n2g6tyLivO60anVahY51zmG4pIiIiIr5SkaTw2kqPopKYWSqQmpKS4nUoImHPOccbn67kwddmcbCoiOfuvob7rz2XqKiKJ3kaPioiIiLiP+Umhc65HVURSGVwzk0EJnbp0uU+r2MRCWdbdufxwJgZfLBwDZec2ZyX7+9Byun1vQ5LREREREKgIj2FIlKNjfsii5+/PIO9+Qd55s4rGH59Z6Kjok6pLg0fFREREfEfJYUiUqadufkMf+0j3vx0FZ1bN+WVIT05J7Fh0PUqKRQRERHxFyWFInKMKUvWcf/oaWzZs5/HbrqYEWkXElMjOuh6HS4E0YmIiIhIKCkpFJHD9u4/yG//M5sxs5fTPrEh437Tl86tm4asfk00IyIiIuI/SgpFBIA5X2/g3oypfL9tLw+lXsBjN11MrRjdIkREREQinb7xiVRzeQcK+MPbn/D81MWkNK3HrP+5jUvObF4p59JEMyIiIiL+o6RQpBqbn7WRQS9N5ZtNO3mgeyeevO1y6sTGVNr5NHxURERExH+UFIpUQwcKCvnz+/P528QvSGxQl2mP3MzV7Vt5HZaIiIiIeEBJoUg1s/i7LQxMn8LyDdu456oO/L3/lZwWV6tKzq3hoyIiIiL+o6RQpJooLCrmmQmf85dx82kcX5vxv+1Lr/PaVHkcSgpFRERE/EVJoUg1sDJnOwPTp7Bw7WZuv7gd/7z7GhrUrV3lcWidQhERERH/UVIoEsGKiot5bvIiHn3nU+rWiuHNYX24+aIzPYtHE82IiIiI+I+SQpEI9e3mXdz70lQ+WZ1D2vlteXHQdTRNqON1WCIiIiLiM0oKRSKMc46MmcsY8cZcakRHMWZIT3562dm+6KHTRDMiIiIi/qOkUCSCbNi+l8EZ05ixfD3XdUgiY3B3WjaM9zqswzR8VERERMR/lBSKRADnHP/++Gt+NXYWRcWO5++5lsHXdlQCJiIiIiLlUlIoEuY2787l5y/PYMKX33JZuxaMvr8HbZvW8zqsMmn4qIiIiIj/KCkUCWPvLfiGB8bMYN+BAp7tfyXDep5HdFSU12GJiIiISBgJu6TQzM4GhgONgJnOuX95HJJIlduxbz/DX/2Itz5bzfltmvLKkJ6c3aKh12GVy6GeQhERERG/qdKk0MzGAH2ALc65DqXKewL/BKKB0c65p49Xh3NuJTDEzKKAUZUcsojvZC5ey/2jp7Nt736euPkSHk67kBrR4dE7qIlmRERERPynqnsKXwWeB8YeKjCzaOAFoBuQDXxhZhMoSRCfOur4gc65LWaWBvwuUJdItbAn7wC/+c9sXp2zgg4tGzHxoX50Sm7idVgiIiIiEuaqNCl0zs01s+Sjii8E1jjn1gKY2VvADc65pyjpVSyrngnABDP7EHijrH3MbDAwGKBVq1YhiV/EK7NWfM+9L00le8c+RqRdyKM3dqVWTNiN/tZEMyIiIiI+5IdvlS2ADaXeZwMXHW9nM7sKuBGoBWQebz/nXAaQAdClSxcXikBFqlpufgGPvP0xL05bwpnN6jPnsdvoekZzr8MKioaPioiIiPiLH5LCsr4hHjeJc87NBmZXVjAifjHvm40MSp/Cms27+GXP8/jLrZcRVyvG67CC4o7/qy0iIiIiHvFDUpgNtCz1PhHY6FEsIp47UFDIE+99xv9OWkjLhvHM+MMtXHlOy/IPDAMaPioiIiLiP35ICr8AzjCz1kAOcDtwZygqNrNUIDUlJSUU1YlUukXrNjMwfQorsrcz6Oqf8Lf+VxJfu6bXYYmIiIhIBKvSeezN7E3gM6CdmWWb2SDnXCHwC2AqsBL4r3NuRSjO55yb6JwbnJCQEIrqRCpNQWERf37vMy597E127Mtn4kP9SL+3W8QlhFqnUERERMR/qnr20TuOU57JCSaNEYlkK7K3MTB9KovWbeaOS89i5ICraVC3ttdhVQqtUygiIiLiP34YPipSLRUVFzMy80v+5515nFa7Jm8P78ONF57pdVgiIiIiUs1EdFKoZwrFr9b8sJNBL01l3jcbuaFLW14c2I0mCXFeh1XpNNGMiIiIiP9EdFLonJsITOzSpct9XsciAlBc7EifsZTfvzWXmOhoXh3akzsvPbtaDamsTp9VREREJBxEdFIo4iffb9vDfRnT+GjF93TvmMRL93YnsWG812FVKa1TKCIiIuI/SgpFKplzjtfmruA3/55NUbHjxUHXce/VP6mWPWYaPioiIiLiP0oKRSrRD7tyGfrydCYtWsvlZ7Vg9P09aNOkntdhiYiIiIgcFtFJoSaaES/997PV/PLVmeQdKODvP72KX/Y4j6io6t1LpnUKRURERPwnopNCTTQjXti2dz/DXpnJOwu+4YK2pzNmSE/Oat7A67B8ozoOmxURERHxs4hOCkWq2qRF3zJk9HR27Mvnz7deym/7XECN6Civw/INTTQjIiIi4j9KCkVCYHfeAX7979mMnbuCn7RqxIcjbuLcpMZeh+U7mmhGRERExH+UFIoEaeby9dyXMY2cHfv4/Q0X8ccbu1KzRrTXYYmIiIiIVEhEJ4WaaEYqU25+Ab97cy7pM5bSrll95j5+OxelNPM6LF9zOD1TKCIiIuIzEZ0UaqIZqSyfrs5h0EtTWLtlN8Ov78yfb72U2jVjvA7L9zR8VERERMR/IjopFAm1/IOFPP7uPP6RuZDkRgnM+MMtXHF2S6/DEhERERE5ZUoKRSpo4dofGJQ+la9ztnPfNR155s4riK9d0+uwworWKRQRERHxHyWFIuUoKCzirx8s4KnxCzg9oQ4fjriR7h2TvQ4rbOmZQhERERF/UVIocgLLN2zjnvQpLPluC/0vO5v/G3A19evEeh1W2NI6hSIiIiL+E9FJoWYflVNVVFzMPz5cyOPvfkZCXE3eeTCVvhec4XVYYU8TzYiIiIj4T0QnhZp9VE5F1g87GZg+hflZm+h3QQovDLyOxqfFeR2WiIiIiEiliOikUORkFBc7Xpy+hEfe+phaNaIZ+/Pruf2Ss/QMXAhpnUIRERER/1FSKAKs37qHezOmMvvrDfQ8N5mX7utO8/p1vQ4rImn4qIiIiIi/KCmUas05xytzlvPbf8/B4Ui/txsDr+qg3qxKoolmRERERPxHSaFUWxt37mPI6OlMXrKOK89OZPT9PUhunOB1WBFNE82IiIiI+I+SQql2nHO8/dlqhr06k/0HC/nHz67ige7nERWlZKUqqBdWRERExF+UFEq1snVPHr98ZSbvfZ7FRSnNePn+HrRr3sDrsKoNDR8VERER8Z+ITgq1TqGUNmHhGoa+PIOdufk8edtl/Lp3F2pER3kdVrWi4aMiIiIi/hPRSaHWKRSAXbn5/Orfs/nPx19zblJjpvz+Jn7SqrHXYYmIiIiI+EJEJ4Ui05d9x+BR09i0K5c/9LuIR/p2pWaNaK/Dqra0TqGIiIiI/ygplIi0L/8gI96YS8bMZZzdvAH/fTyNC9qe7nVYgtYpFBEREfEbJYUScT5Zlc2gl6aybutuftXrfJ645RJq14zxOixBE82IiIiI+JGSQokY+w8W8Ng78xg5+UtaN07goz/eymVnJXodlpSiiWZERERE/EdJoUSEL779gYHpU1i1cQf3X3cuT99xOXVja3odloiIiIiI7ykplLB2sLCIJ8fN55kJn9OsXh0yR9xIt47JXoclx6GJZkRERET8R0mhhK2vvt/KPelTWLp+Kz+7/Bz+8bOrqFcn1uuwpBwaPioiIiLiL0oKJewUFhXzv5MW8sR786hfJ5b3fpVGWpcUr8OSCnBOE82IiIiI+E1EJ4VmlgqkpqQoYYgUqzfuYGD6FD7/9gduuvAMnh94HY3ia3sdllSQho+KiIiI+E9EJ4XOuYnAxC5dutzndSwSnOJix/PTFvPHtz8hNiaa//yiF7d2bacEIwxp+KiIiIiIv0R0UiiRYd2W3dyXMZU5K7Pp1ak16fd2o1n9ul6HJafAoSUpRERERPxGSaH4lnOOl2d9xUOvz8EwMu7rzt1XtlfvYBhzzhFlUV6HISIiIiKlKCkUX8rZsZcho6czZel3XN2+JaPu60FS49O8DktEREREIlBBQQHZ2dnk5+d7HUrQYmNjSUxMJCYmpsLHKCkUX3HO8ea8VQx/9SMOFBYx8q6rGXpdJ6Ki1DsYCTTRjIiIiPhRdnY28fHxJCcnh/V3Fecc27dvJzs7m9atW1f4OCWF4htb9+TxwJgZjPtiDV3PaMaYIT054/T6XoclIaZnCkVERMRv8vPzwz4hBDAzGjZsyNatW0/qOCWF4gsffJHFz8fMYHfeQZ66/XJ+1ft8oqP07Fmk0TqFIiIi4lfhnhAeciqfQ0mheGpnbj4PvjaLNz5dSafkJkx7pCcdWjbyOiypJBo+KiIiIuI/6ooRz0xduo5OI17j7c9W8eiNXZn3xB1KCKsBDR8VERGRcPf666+TnJxMVFQUycnJvP7660HXmZ+fz4UXXsi5555L+/bteeyxx47ZxznHsGHDSElJoWPHjixatCjo84J6CsUDe/cf5OE35jD6o684p0VD3v9NX85v3dTrsKQKODR8VERERMLb66+/zuDBg8nLywNg/fr1DB48GID+/fufcr21atXio48+om7duhQUFHDZZZdx/fXX07Vr18P7TJ48maysLLKysliwYAFDhw5lwYIFwX0glBRKFZu7cgODXprK+m17+E3vLjx+8yXE1tRlWF04p8XrRURExN8efPBBlixZctzt8+fP58CBA0eU5eXlMWjQIEaNGlXmMZ06dWLkyJEnPK+ZUbduXaBkiYyCgoJjHrsZP348AwYMwMzo2rUru3btYtOmTTRr1qwiH+24NHxUqsT+gwX85t+zue7Jd4iOMmY9ehtP33mFEkIRERERCStHJ4TllZ+MoqIiOnXqRJMmTejWrRsXXXTREdtzcnJo2bLl4feJiYnk5OQEfV59I5dKt2DNJgalT2H1pp0M7XYuT91+BXViK76YpkQOTTQjIiIifldej15ycjLr168/pjwpKYnZs2cHde7o6GiWLFnCrl276NevH8uXL6dDhw6Ht5c1k3sovltFdE+hmaWaWcbu3bu9DqVaOlhYxKP//YQrHn+LvIOFTPn9TTx397VKCKs5DR8VERGRcPbkk08SFxd3RFlcXBxPPvlkyM5Rr149rrrqKqZMmXJEeWJiIhs2bDj8Pjs7m+bNmwd9vohOCp1zE51zgxMSErwOpdpZun4rXR99nafHf87PrjiHxU8P4NoOSV6HJR7TOoUiIiIS7vr3709GRgZJSUmYGUlJSWRkZAQ1yQzA1q1b2bVrFwD79+9nxowZnHXWWUfsk5aWxtixY3HOMX/+fBISEoJ+nhA0fFRCrLComL9N/II/v/8ZDerGMu43N9Cnc1uvwxKf0PBRERERiQT9+/cPOgk82qZNm7jrrrsoKiqiuLiYW2+9lT59+pCeng7AkCFD6NWrF5mZmaSkpBAXF8crr7wSknMrKZSQWbVxB/f8azIL127m1q7teO7ua2gYX9vrsEREREREfK9jx44sXrz4mPIhQ4Ycfm1mvPDCCyE/t5JCCVpxseO5qYt49O1PiKsVw+u/6M2tF7fzOizxIYeWpBARERHxGyWFEpS1W3Zx70tT+XhVDr3Pa0P6vd04vV4dr8MSH1NSKCIiIuIvSgrllDjnGPXRVzz8+hyio4zRg3sw4Ipz9LyYnJAmmhERERHxHyWFctKyt+/l/tHTmLZsPde0b8Wowd1p1eg0r8OSMKCJZkRERET8R0mhVJhzjtc/WcmDY2dRUFTEc3dfw/3XnktUlL7kS8Vp+KiIiIiIvygplArZsjuPn4+ZzviF33LJmc15+f4epJxe3+uwJMxo+KiIiIiI/ygplHK9//k3PDBmJnv2H+SZO69g+PWdiY6K8josCUOafVRERETC3umnw+bNx5Y3bQo//BBU1cnJycTHxxMdHU2NGjVYuHDhEdudcwwfPpzMzEzi4uJ49dVX6dy5c1DnBCWFcgI79u1n+GuzeGveKjq3bsqYIT1on9jI67BERERERLxTVkJ4ovKTNGvWLBo1Kvs79+TJk8nKyiIrK4sFCxYwdOhQFixYEPQ5lRRKmSYvWcv9o6azde9+HrvpYkakXUhMjWivw5Iwp4lmRERExPcefBCWLDm1Y6+6quzyTp1g5MhTDumQ8ePHM2DAAMyMrl27smvXLjZt2kSzZs2CqldjAOUIe/IOcP+oaaT97QMa1I3l0yfu4I83XqyEUEJGw0dFREREymZmdO/enfPPP5+MjIxjtufk5NCyZcvD7xMTE8nJyQn6vOoplMNmr/ieezOmsWH7Xh5KvYDHbrqYWjG6RCR0NNGMiIiI+F55PXonGvU0e3ZQp/70009p3rw5W7ZsoVu3bpx11llcccUVh7eX9V0qFKOw1FMo5B0o4FdjZ9Htr+8SEx3FrP+5jb/efrkSQgk5DR8VEREROb7mzZsD0KRJE/r168fnn39+xPbExEQ2bNhw+H12dvbhY4KhpLCa+yxrI10e+TfPT13MA907sfCvP+OSM4O/sESOR8NHRUREJKw1bXpy5RWUm5vL3r17D7+eNm0aHTp0OGKftLQ0xo4di3OO+fPnk5CQEPTzhKDho9XWgYJC/vTeZ/x90kISG9Rl2iM3c3X7Vl6HJRFOw0dFREQk7AW57MTxbN68mX79+gFQWFjInXfeSc+ePUlPTwdgyJAh9OrVi8zMTFJSUoiLi+OVV14JybmVFFZDi7/bwsD0KSzfsI17rurA3/tfyWlxtbwOS6oJ9RSKiIiIHKtNmzYsXbr0mPIhQ4Ycfm1mvPDCCyE/t5LCaqSgsIhnJ37BX8bNp1F8bcb/ti+9zmvjdVhSjTjUUygiIiLiN0oKq4mvs7cz8KUpfLl2M7df3I6Rd11Dw/jaXocl1YwmmhERERHxn7BMCs2sDjAXeMw5N8nrePysqLiY5yYv4tF3PqVurRjeHNaHmy860+uwpBrT8FERERERf6nSpNDMxgB9gC3OuQ6lynsC/wSigdHOuafLqWoE8N9KCzRCfLt5F4Nemsqnq3NI7dyWf917HU0T6ngdllRjmmhGRERExH+quqfwVeB5YOyhAjOLBl4AugHZwBdmNoGSBPGpo44fCHQEvgZiqyDesOSc46WZyxjxxhxioqMZM6QnP73sbA3bE89p+KiIiIiI/1RpUuicm2tmyUcVXwiscc6tBTCzt4AbnHNPUdKreAQzuxqoA5wD7DezTOdccRn7DQYGA7RqVX2WWtiwfS+DM6YxY/l6ruuQRMbg7rRsGO91WCIiIiIiEWVF9jbufO5D3hjWm/aJjbwOJyh+WLy+BbCh1PvsQFmZnHN/cM49CLwBjCorIQzsl+Gc6+Kc69K4ceOQBuxHzjnGzl1BpxGv8VnWRp6/51oyf3ejEkLxFeecnikUERGRsJebX0Das+NYuXE7N/ztA3LzC4Kuc/Xq1XTq1Onwz2mnncbIkSOP2Mc5x7Bhw0hJSaFjx44sWrQo6POCPyaaKesbYrkPHjnnXg19KOHph125/PzlGUxc9C2XtWvB6Pt70LZpPa/DEimTkkIREREJd/dmTGXLnjycg827c7lv1FTe+OUxgxxPSrt27ViyZAkARUVFtGjR4vBi9odMnjyZrKwssrKyWLBgAUOHDmXBggVBnRf8kRRmAy1LvU8ENoaiYjNLBVJTUlJCUZ0vvbvgG34xZgb7DhTwbP8rGdbzPKKj/NABLHIsrVMoIiIifvfrf89i6fqtx93+w859rNm8m+LABHr5BUW8tyCL9t+N4fT6dcs85tykxvzjZ1dXOIaZM2fStm1bkpKSjigfP348AwYMwMzo2rUru3btYtOmTTRr1qzCdZfFD9nDF8AZZtbazGoCtwMTQlGxc26ic25wQkJCKKrzle179/PT5z/kjucmkdwkgS+e/Cm/6nW+EkLxNU00IyIiIuFu3dY9hxPCQ4qdY93WPSE7x1tvvcUdd9xxTHlOTg4tW/58HYltAAARDklEQVTYn5aYmEhOTk7Q56vqJSneBK4CGplZNiXrDL5sZr8AplIy4+gY59yKqowr3GQuXsv9o6ezbe9+nrj5Eh5Ou5Aa0UoGJTxo+KiIiIj4WXk9eq/OXs7wsR+Rd6DwcFlcrRo8d9c13HVlhxMcWTEHDx5kwoQJPPXU0QsxlL28Vyj+4F7Vs48em+6WlGcCmVUZSzjak3eA3/xnNq/OWUGHlo2Y+FA/OiU38ToskQrTOoUiIiIS7u6+qgNTl33HpEXfkl9QRGxMNL3PaxOShBBKnhvs3LkzTZs2PWZbYmIiGzb8OEdndnY2zZs3D/qcEd29ZGapZpaxe/dur0MJ2kfLv+e8341l7NyvGZF2IfP/fKcSQglLGj4qIiIi4W704B40OS0OA5om1GHUfT1CVvebb75Z5tBRgLS0NMaOHYtzjvnz55OQkBD084QQ4UlhJDxTmJtfwPDXPqLHU+8SW7MGcx67jb/cdhm1YvwwR5DIydFEMyIiIhIJ6sTGMOHhfpzdoiHjH+pLndiYkNSbl5fH9OnTufHGGw+Xpaenk56eDkCvXr1o06YNKSkp3Hfffbz44oshOa8yCx+b981GBqVPYc3mXfyy53n85dbLiKsVmgtOxAtap1BEREQiRfvERix99q6Q1hkXF8f27duPKBsyZMjh12bGCy+8ENJzgpJCX8o/WMgT783jHx9+ScuG8Ux/5Gauat/K67BEQkJJoYiIiIi/KCn0mUXrNjMwfQorsrcz6Oqf8Lf+VxJfu6bXYYmEhIaPioiIiPhPRCeF4bR4fUFhEU+P/5y/jl9A4/jaTHyoHz07tfY6LJGQ0jqFIiIiIv6jiWY8siJ7G+c+/BorsrexInsblz72Jn96/zNu6XomS54ZoIRQIpaGj4qIiIj4S0T3FPpVbn4Bac+OY8OOvVz1xNvkHiggIa4Wbw/vw40Xnul1eCKVRusUioiIiPhPRPcU+tW9GVPZvDsP52BX3gEan1abpc/cpYRQqgUNHxUREZFIsGLHCjq804EVO1Z4HUrQlBRWsVdnLydzyVoOFBYdLtuVe4DMxWs9jEqkamiiGREREYkEuQW59JrSi693fk3vKb3JLcgNSb3//Oc/6dChA+3bt2fkyJHHbHfOMWzYMFJSUujYsSOLFi0KyXkjOik0s1Qzy9i9e7fXoRz2h7c/Ie9A4RFleQcL+cPbn3gUkUjV0TqFIiIiEgkGzhnIlv1bcDg279/MoDmDgq5z+fLljBo1is8//5ylS5cyadIksrKyjthn8uTJZGVlkZWVRUZGBkOHDg36vBDhzxQ65yYCE7t06XKf17Ec8uRtlzF87EdHJIZxtWrw19sv8zAqkaqjpFBERET87MF5D7Jk+5Ljbt+Uu4k1e9dQ7IoByC/K552177D4rcU0q9OszGM6NezEyEuO7fkrbeXKlXTt2pW4uDgArrzySsaNG8fDDz98eJ/x48czYMAAzIyuXbuya9cuNm3aRLNmZZ+3oiK6p9CP7r6qA706tSE2JhqA2Jhoep/Xhruu7OBxZCKVT8NHRUREJNyt27fucEJ4SDHFrNu3Lqh6O3TowNy5c9m+fTt5eXlkZmayYcOGI/bJycmhZcuWh98nJiaSk5MT1HkhwnsK/Wr04B50fPhVNmzfS9OEOoy6r4fXIYlUGU00IyIiIn5WXo/emNVjGPbpMHILf3yOMK5GHM9f+jz3tLvnlM979tlnM2LECLp160bdunU599xzqVHjyHStrJncQ/HdSj2FHqgTG8OEh/txdouGjH+oL3ViY7wOSaRKqKdQREREwt3AdgPp3ao3sdGxAMRGx5LaKjWohPCQQYMGsWjRIubOnUuDBg0444wzjtiemJh4RO9hdnY2zZs3D/q8Sgo90j6xEUufvYv2iY28DkWkymiiGREREYkEY64cQ5PaTTCMprWb8vKVL4ek3i1btgDw/fff8/7773PHHXccsT0tLY2xY8finGP+/PkkJCQE/TwhaPioiFQxDR8VERGRcFcnpg6ZPTO5beZtvH3t29SJqROSem+66Sa2b99OTEwML7zwAvXr1yc9PR2AIUOG0KtXLzIzM0lJSSEuLo5XXnklJOeN6KTQzFKB1JSUFK9DERE0fFREREQiR/sG7Vl+y/KQ1vnxxx8fUzZkyJDDr82MF154IaTnhAgfPuqcm+icG5yQkOB1KCKCho+KiIiI+FFEJ4Ui4j9KCkVERET8RUmhiFQZDR8VERERvypruYdwdCqfQ0mhiFQpTTQjIiIifhMbG8v27dvDPjF0zrF9+3ZiY2NP6riInmhGRPxFPYUiIiLiR4mJiWRnZ7N161avQwlabGwsiYmJJ3WMkkIRqTKaaEZERET8KCYmhtatW3sdhmcieviomaWaWcbu3bu9DkVEAjR8VERERMRfIjop1JIUIv6i4aMiIiIi/hPRSaGI+I+Gj4qIiIj4i4X7DDsVYWZbgfVex1GGRsA2r4OoptT23lHbe0dt7x21vbfU/t5R23tHbe8dv7Z9knOucVkbqkVS6FdmttA518XrOKojtb131PbeUdt7R23vLbW/d9T23lHbeycc217DR0VERERERKoxJYUiIiIiIiLVmJJCb2V4HUA1prb3jtreO2p776jtvaX2947a3jtqe++EXdvrmUIREREREZFqTD2FIiIiIiIi1ZiSQhERERERkWpMSWGImFlPM1ttZmvM7HdlbO9vZssCP/PM7NxS274zs6/MbImZLSxV3sDMpptZVuC/9avq84STU217M2sXaPNDP3vM7MHAtsfNLKfUtl5V/bnCQQXa/oZAuy8xs4Vmdll5x+q6r7hTbX8za2lms8xspZmtMLPhpY7RtV8BQV77uucHIYjrXvf8IJXX9qX2u8DMiszs5vKO1XVfMafa9rrfBy/I6z587vfOOf0E+QNEA98CbYCawFLgnKP2uQSoH3h9PbCg1LbvgEZl1Pss8LvA698Bz3j9Wf32E2zbH1XPD5Qs6gnwOPBbrz+fn38q2PZ1+fHZ5Y7AqvKO1XVfJe3fDOgceB0PfFOq/XXtV2LbB97rnu9R2x9Vj+75IW77Uvt9BGQCN5d3rK77Sm973e89avtAedjc79VTGBoXAmucc2udcweBt4AbSu/gnJvnnNsZeDsfSKxAvTcArwVevwb0DVG8kSRUbX8t8K1zbn2lRhtZKtL2+1zgjgfUAVwFjtV1XzGn3P7OuU3OuUWB13uBlUCLKos8/AVz7Z+Irv3yhartdc8/eeW2fcAvgfeALRU8Vtd9+U657XW/D1ow1/2J+O66V1IYGi2ADaXeZ3PiX7hBwORS7x0wzcy+NLPBpcqbOuc2QckvNdAkRPFGkmDb/pDbgTePKvtFYAjSGF906/tPhdrezPqZ2SrgQ2BgBY7VdV8xwbR/6e3JwHnAglLFuvZPLNi21z3/1IXkukf3/FNRbtubWQugH5B+Esfqui9fMG1fep9kdL8/WcG2fdjc75UUhoaVUVbmX4XN7GpKEpMRpYovdc51pmRo4wNmdkXoQ4xYwbY9ZlYTSAPeKVX8L6At0AnYBPxvKIKNMBVqe+fcOOfcWZT8FezPJ3OsnFAw7V9SgVldSv6y+aBzbk+gWNd++YJte93zT10ornvd809NRdp+JDDCOVd0CsfK8QXT9iUV6H5/qoJt+7C53yspDI1soGWp94nAxqN3MrOOwGjgBufc9kPlzrmNgf9uAcZR0lUNsNnMmgWObUbFu6Srk6DaPuB6YJFzbvOhAufcZudckXOuGBjFj/8m8qMKtf0hzrm5QFsza1TOsbruKyaY9sfMYij5gvC6c+79Uvvp2i9fUG2ve35Qgmr7AN3zT01F2r4L8JaZfQfcDLxoZn3LOVbXffmCaXvd74MTVNuH0/1eSWFofAGcYWatA3+BvB2YUHoHM2sFvA/8zDn3TanyOmYWf+g10B1YHtg8Abgr8PouYHylforwdMptX8odHDWM6NAvakA/fvw3kR9VpO1TzMwCrztT8pD29nKO1XVfMafc/oGyl4GVzrl/HHWMrv3yBdP2uucHJ5j7ziG655+actveOdfaOZfsnEsG3gV+7pz7oJxjdd2X75TbXvf7oAXT9mF1v6/hdQCRwDlXaGa/AKZSMvvQGOfcCjMbEtieDvwP0JCSvx4AFDrnugBNgXGBshrAG865KYGqnwb+a2aDgO+BW6rwY4WFINseM4sDugH3H1X1s2bWiZIhAt+Vsb3aq2Db3wQMMLMCYD9wW2ACiDKPDVSt674Cgml/K5mi/2fAV2a2JFDlI865THTtlyvIttc9PwhB3nd0zw9CBdv+pI4NbNZ1X45g2h64FN3vT1mQbR9W9/tDUzaLiIiIiIhINaThoyIiIiIiItWYkkIREREREZFqTEmhiIiIiIhINaakUEREREREpBpTUigiIiIiIlKNKSkUEZGIYGZNzewNM1trZl+a2Wdm1u84+zY3s3ePs222mXWp4Dn7mtk5wcR9VH3NzGxS4PVVZrbbzJaY2TIzm2FmTco5/nEz+20FzzXDzOqHIm4REQlvSgpFRCTsBRZo/gCY65xr45w7n5JFhhPL2LeGc26jc+7mEJy6LxCypBD4NTCq1PuPnXOdnHMdKVlE+YFgT2AlooB/Az8Ptj4REQl/SgpFRCQSXAMcLL2QsHNuvXPu/wGY2d1m9o6ZTQSmmVmymS0PbKttZm8FeuPeBmqXdQIze9rMvg7s93czuwRIA/4W6M1rG/iZEuip/NjMzgoc+6qZpQfKvjGzPsf5HDcBU44uDCS98cDOwPsGZvZBIJb5Ztax1O7nBHo715rZsMD+yWa20sxeBBYBLYEJwB0VbWAREYlcNbwOQEREJATaU5LsnMjFQEfn3A4zSy5VPhTIc851DCRXx9RjZg2AfsBZzjlnZvWcc7vMbAIwyTn3bmC/mcAQ51yWmV0EvEhJwgqQDFwJtAVmmVmKcy6/1DlaAzudcwdKnfpyM1sCNARygUcC5U8Ai51zfc3sGmAs0Cmw7SzgakqSyNVm9q9AeTvgHufc4d5BM6tlZg2dc9vLaTsREYlg6ikUEZGIY2YvmNlSM/uiVPF059yOMna/AvgPgHNuGbCsjH32APnAaDO7Ecgr45x1gUuAdwKJ3EtAs1K7/Nc5V+ycywLWUpK8ldYM2HpU2aHhoy2BV4BnA+WXUTL8E+fcR0BDM0sIbPvQOXfAObcN2AI0DZSvd87NP6r+LUDzMj6viIhUI+opFBGRSLCCkqGXADjnHjCzRsDCUvvknuB4d6LKnXOFZnYhcC0lzyr+gh97AA+JAnY55zodffxxznH0+/1A7AnCmAC8F3htJ6i/dE9jET/+v76szx8bOK+IiFRj6ikUEZFI8BEQa2ZDS5XFVfDYuUB/ADPrAHQ8eodAL2CCcy4TeJAfh2rupWSYJs65PcA6M7slcIyZ2bmlqrnFzKLMrC3QBlh91Gm+oWSI6fFcBnxbRsxXAdsC56+wwHOKpwPfncxxIiISedRTKCIiYS/wnF9f4P/M7GFKhmHmAiMqcPi/gFfMbBmwBPi8jH3igfFmFktJL92vAuVvAaMCE7rcTEmi9i8z+yMQE9i+NLDvamAOJcM5h5R+njDwGXLN7NvAs4ZrAsWHnik0YDdwb6D88VIx5wF3VeBzHu18YL5zrvAUjhURkQhizp1wxIyIiIgEycxepdSENCfYrx9wvnPuj1UQ0z+BCc65mZV9LhER8Tf1FIqIiPiEc26cmTWsotMtV0IoIiKgnkIREREREZFqTRPNiIiIiIiIVGNKCkVERERERKoxJYUiIiIiIiLVmJJCERERERGRakxJoYiIiIiISDX2/wEJkviNxR7V+gAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 1080x504 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "%matplotlib inline\n",
    "# Lists of markers and colors (for matplotlib)\n",
    "colors = ['#74a9cf', '#2b8cbe', '#045a8d', '#009900', '#FF8000']\n",
    "colors = ['#000000', '#ff0000', '#045a8d', '#009900', '#FF8000']\n",
    "markers = ['o','s','d','d','d']\n",
    "\n",
    "plt.figure(figsize=(15,7))\n",
    "# Plot with matplotlib\n",
    "for i,crmult in enumerate(Crmult):\n",
    "    im = i%len(colors)\n",
    "    ener = [ HtoeV*(l.energy-emin_LDA) for l in log_LDA[crmult] ]\n",
    "    plt.plot(Hgrids, ener, marker=markers[im], \n",
    "             ls='-', label=str(crmult),color=colors[im])  \n",
    "\n",
    "plt.yscale('log')\n",
    "plt.xlabel('Grid step (Bohr)')\n",
    "plt.ylabel('Total energy $\\Delta E$ (eV)')\n",
    "plt.title('Dissociation energy of the N2 dimer for different rmult')\n",
    "plt.legend(loc=4)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA4UAAAG5CAYAAAA595FfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nOzdd3hc133n//d3Bh1gARs6CIIEe5UgqpCS1YtFSrIiW7JkR8W2NtmNs/ltNskm9m4Sb7R2slnH9qatYstykdwdq1iS1a1GSmIXKfYCEmAnAbAAIMp8f3/cgTADDkCARJ/P63nmIWbuufeee+cCnM+cc88xd0dERERERESSU2iwKyAiIiIiIiKDR6FQREREREQkiSkUioiIiIiIJDGFQhERERERkSSmUCgiIiIiIpLEFApFRERERESSmEKhiCQdM/tXM/vvg7Df583s/vNct9TMTplZuK/rNRKZ2RIz2x49Z3f0oHyZmbmZpQxE/S6UmT1uZn8T/flKM9s6iHWZYWZrzeykmf1hP2z/ajOrjnm+ycyujv5sZvZdM6s1s/eir/2+mR2Kvvfj+7o+w1XsNSMi0plCoYiMKGa2x8waox9Q68zsHTP7PTP76O+du/+eu//Pga6bu9/i7t/rSdnocVwfs+5ed89x97b+q+GI8hXgH6Pn7FedF3Y+v30pGmLczP6p0+tvmdkD0Z9vjT6vM7ODZvZvZjbqfPbn7m+6+4w+qPr5+lPgdXcf5e7f6u+dufscd389+nQpcANQ7O6LzSwV+DpwY/S9P9bf9YnVn9dVX+octEVEFApFZCRa7u6jgMnA14A/A74zuFWSAW7lnAxsGsD9dXYa+F0zK+ti+Rjgb4BCYBZQDPzvAanZeYq2yiX63HDe57oPWmYnA3vc/XT0eR6QcQH1GRIt8d2caxGRfqE/OCIyYrl7vbs/DdwN3G9mc+GsrncTzOzZaIvNcTN7s/3DmJn9mZnVRFsdt5rZddHX083sG2a2P/r4hpmlt+/XzG43s3VmdsLMdprZzdHXXzezz0d/nmpmr5rZMTM7amZPmNnY6LIfAKXAM9EucH/auXujmRWa2dPROu8wsy/E7P+vzOynZvb9aN03mVllV+fJzGaa2UvRbW01s0/FLHvczP7JzH4d3da7Zja1F+v+i5k9Z2angWvM7CLr6Gr4MzP7Scx7sdHMlsesnxo9Nwu7qPcXosd+PHouCqOv7wTKY85feqf1zjq/MYvvM7O90f1+KWadkJn9t+j7eSx6fsd1dU6BOuBx4C8TLXT3J939BXdvcPda4N+AJV1tzMwWmdma6Hn7CUHwaV/WuXvlHjP7EzPbYGanzew7ZpZnQfflk2b2spnlxpS/zIIW9TozW2/RrpnRZa+b2SNm9jbQED2vsfV6FbgG+MfouZxuZmOi194RM6sysy/H/E49YGZvm9k/mNlx4K8SHGtm9NqpNbMPgUs6Ld9jZteb2eeAbwOXR/f9I6C9G21dtG7nc40WmtkvovXfbTFdYrv73TrHdRVb/+7+PsSdawt+5/+jBV2hT5rZ/7Tgb8eK6Po/NbO0mHP7Vqd9uZlN6/RaNvA8UBit56n23x0RSWLuroceeugxYh7AHuD6BK/vBX4/+vPjwN9Ef/4q8K9AavRxJWDADGAfUBgtVwZMjf78FWAlMAmYCLwD/M/ossVAPUGXthBQBMyMLnsd+Hz052nRMunRbbwBfKOr44ju34GU6PPfAv9MEA4WAkeA66LL/gpoAj4OhKPHuLKL85UdPc4HgRTgIuAoMCfmXB2PHlcK8ATw416sW08QdkLAaKAK+M/Rc30n0BzzXvwp8JOYut0OfNBFva+N7uui6Dn8v8Ab57oOenB+/w3IBBYAZ4BZ0eV/FH3Pi6P7+3/Aj7rY9tVANZAPnABmRF9/C3igi3W+0X5eEyxLi563/y963u4CWmLO29VAdadjW0nQalYEHAbWAIuidX8V+Mto2SLgWPRaCRFck8eAiTHX7F5gTvQ9Tk1Qv9eJXtfR598HngJGRc/rNuBz0WUPAK3AF6Pby0ywva8BbwLjgBJgY4Ljuz5me29183vS22s0C1gN/I/oeS8HdgE39eR3i3Nfd+f6+xB3rqPH8jTB784cgmvylWi9xgAfAvcnOhfR1xyYluDv3tWx51QPPfTQQy2FIpIs9hN8yOysBSgAJrt7iwf3ZznQRvABeraZpbr7HnffGV3nPuAr7n7Y3Y8Afw18Nrrsc8Bj7v6Su0fcvcbdt3TeqbvviJY5E93G14GP9eRAzKyE4F6qP3P3JndfR9Bi8tmYYm+5+3Me3IP4A4KQk8gygu5333X3VndfA/yCIHi0+6W7v+furQShcGEv1n3K3d9290h0vRTgW9Fz/UvgvZiyPwQ+bmajo88/G617IvcRnOc17n4G+HOCFqOyLsr31F+7e6O7rwfW03He/gPwJXevju7vr4C7rJvuj+5+kOALh690t0MzuwG4nyCIJHIZQUD4RvS8/Rx4/xzH8X/d/ZC71xAErHfdfW207v9OEBABPgM8F71WIu7+ErCKIPS0e9zdN0Xf45ZzHEuYoGX+z939pLvvAf4P8dfmfnf/v9HtNSbYzKeAR9z9uLvvAy7kPsXeXqPzCALxV9y92d13EXxRcE9M+Z7+biVyrr8Pic7137r7CXffRBCQX3T3Xe5eT9DitwgRkQukUCgiyaKIoMWrs/8N7ABeNLNdZvbfIAhtBK1DfwUcNrMfx3SxKiRouWlXFX0NgpaNnZyDmU2KbrPGzE4QBKIJPTyWQuC4u5/sVIeimOcHY35uADK6CDCTgUujXQfrzKyOIHDld7OtnF6su69TvWuiofus5e6+H3gb+B0LutLeQhBCE4l7D9z9FEELV1EX5Xuqu2P995jj3EzwxUHeObb3t8BNZpYwOJjZZcCTwF3uvq2LbSQ6b1VdlG13KObnxgTPY4/rk53ew6UEX5S0i30Pz2UCHS2bsXWNfV/Otb3CTmXOdazd6e01OpmgW2Vs+b8g/n3u6e9WIuf6+5Do3PT0vRQROW/DYuhtEZELYWaXEHwofavzsmiw+mPgj81sDvCamb3v7q+4+5PAk9GWq/9H8AH/swStjrGDa5RGX4PgQ91Uzu2rBF275rv7MQumTfjH2Kp1s+5+YJyZjYoJhqVATQ/229k+4LfufkM/rRt7HAeAIjOzmIDT+UPy94DPE/z/tCLa0pVI+3sAfHSf1Hh6fg66O7+J7AMecve3e7NS9L39BnDWaLdmtoiga+BD7v5KN5tJdN5K6cGXDz2wD/iBu3+hmzK9OVdHCVrfJxN0bYSzr81zbe8AwXUR+/t1vnp7je4Ddrt7xXnu71zHdq6/D729LmOdJuj+CoCZ5XdT9kL2IyIjkFoKRWTEMrPRZrYM+DHwQ3f/IEGZZWY2zcyM4P6vNqDNgrnXrrVgkJImgm/k26eD+BHwZTObaGYTCLr9/TC67DvAg2Z2nQWDkxSZ2cwE1RsFnCIYEKMI+JNOyw/RaVCPdtEude8AXzWzDDObT9AtratWte48C0w3s89aMLBLqpldYmaz+mHdFQTn8A/MLMXMbie4xyrWrwju+/rPBPemdeVJgvO8MPoe/S+CLpJ7elBv6Ob8duFfgUfMbDJA9L2/vYfrfh24gmCUUaLrzwVeAL7o7s+cY/0VBPfh/WH0vN3J2eftfP0QWG5mN5lZOHo9XW1mxeezsWiXyp8SnKtR0fP1X+j4/eiJnwJ/bma50Xp88XzqEtXba/Q94IQFg0xlRs/J3OgXSz1xruuqp38fzsd6YE70dyKDBIP4dKrneDMb00f7FpFhTqFQREaiZ8zsJMG38l8i+FD+YBdlK4CXCQLaCuCfPZgDLZ1gwIujBN3FJhF0I4NgKoFVwAbgA4JBPP4GwN3fi+7rHwgGlPgtMS1aMf6aIPzUA78Gftlp+VcJgmedmf3XBOt/mmBQjf0E94j9ZfR+sF6JtjTeSHDP1P7osf4twfH36bru3kwwuMznCEbn/AzBh/YzMWUaCe75msLZ5yR2W68A/z1a9gBB68s9XZVP4Fznt7NvErTqvRi9tlYCl/ZkR+5+Avg74u9p/WOCAYa+EzMCZMJpFGLO2wNALcE9e12em96IfsFwO8G1fYTgd+ZPuLDPB18kaLXaRdA6/yTwWC/W/2uCLqO7gRfp+r7SczqPa7QNWE5w/+tugt//bxMM6tIT3V5Xvfj70GvR7sdfIfh7tp0EPSNiym4h+HJrV7SuGn1UJMlZ/C0KIiIiA8fM3gX+1d2/G/Pa/wCmu/tnBq9mIiIiyUMthSIiMmDM7GNmlh/tBnk/MJ+gG2X78nEELYmPDlYdRUREko1CoYiIDKQZBPc+1RN0obzL3Q9AMBk9QffF5939jcGrooiISHJR91EREREREZEkppZCERERERGRJJYU8xROmDDBy8rKBrsaIiIiIiIig2L16tVH3X1iomVJEQrLyspYtWrVYFdDRERERERkUJhZVVfL1H1UREREREQkiSkUioiIiIiIJLERHQrNbLmZPVpfXz/YVRERERERERmSRnQodPdn3P3hMWPGDHZVREREREREhqQRHQpFRERERESkewqFIiIiIiIiSUyhUEREREREJIkpFIqIiIiIiCQxhUIREREREZEkplAoIiIiIiKSxBQKRUREREREkphCoYiIiIiISBJTKBQREREREUliCoUiIiIiIiIX6IknnqCsrIxQKERZWRlPPPHEYFepxxQKB8FwvmBEzpeuexERERmO3J2Wlhaampo4deoUdXV1HD16lIMHD1JdXc2ePXv4+te/zuc//3mqqqpwd6qqqnj44YeHzecdc/fBrkO/q6ys9FWrVg12NYDgg/HDDz9MQ0PDR69lZWXx6KOPct999w1izUT6j657SVZPPPEEX/rSl9i7dy+lpaU88sgjuuYlKejaHxncndbW1o8ebW1tcc87P0bq8kgkct7ncPLkyezZs6fv3pQLYGar3b0y4TKFwoFVVlZGVVXVWa/n5ORw//33f/TczOKWxz7vblkylR2KdVLZxMv+7M/+jGPHjtHZ+PHj+frXv46Z9fmjvR4D/Ris/V7ovqXv6csQSVYj4drvHIYGO5QM1vILCUN9KRwOk5KS0uVjKCx/6KGHEtbdzIbMeTSFwqETCkOhEF2d8/HjxwOctTz2eXfLRnpZEelfyRaE+3u/L7zwAo2NjWed56ysLJYvXx63neH071Cow0j9dyjUoS/+/dznPseRI0forP2LwKEUerpaPlQ+xA+l0DNYy0OhUNzvyFBVVpa44We4tBSmDHRlkl1paemQv2CGi6EQVFW2Z2WvuOIKampq6KywsJA33ngDd+/TR/t+B/oxWPvVMffdIxKJ9Nl+EwVCgIaGBtauXRv3+zJc/h0qdZDh6dixY3G9orpzocEiLS2NrKysHq8/FEPRcAlDEnjkkUcStpA/8sgjg1irnlMoHGDD/YIZShJ9sypD09/+7d8mvO7/7u/+jqlTpw5izUT6T3ffGm/dunUQajSyDHY4HUlBu6//XbZsGQcPHqSzwsJC3nzzzXMGI4UhGY7au0YP13tpFQoH2HC/YETOh657SUb6ErB/Jep6KUPD3//933f5RWB5efkg1kykf913333D9rPNsLun0MyygX8GmoHX3f2c47wOpXsKRUQkeWgERklWuvZFhp7u7ikcEqHQzB4DlgGH3X1uzOs3A98EwsC33f1rZvZZoM7dnzGzn7j73efavkKhiIiIiIgks+5C4VCZvP5x4ObYF8wsDPwTcAswG/i0mc0GioF90WJtA1hHERERERGREWdIhEJ3fwM43unlxcAOd9/l7s3Aj4HbgWqCYAjd1N/MHjazVWa2KtGwyCIiIiIiIjJEQmEXiuhoEYQgDBYBvwR+x8z+BXimq5Xd/VF3r3T3yokTJ/ZvTUVERERERIapoTz6aKLhxNzdTwMPDnRlRERERERERqKh3FJYDZTEPC8G9vdmA2a23Mwera+v79OKiYiIiIiIjBRDORS+D1SY2RQzSwPuAZ7uzQbc/Rl3f3jMmDH9UkEREREREZHhbkiEQjP7EbACmGFm1Wb2OXdvBf4A+A2wGfipu28azHqKiIiIiIiMNEPinkJ3/3QXrz8HPDfA1REREREREUkaQ6KlsL/onkIREREREZHujehQqHsKRUREREREujeiQ6GIiIiIiIh0T6FQREREREQkiY3oUKh7CkVERERERLo3okOh7ikUERERERHp3ogOhSIiIiIiItI9hUIREREREZEkplAoIiIiIiKSxEZ0KNRAMyIiIiIiIt0b0aFQA82IiIiIiIh0b0SHQhEREREREemeQqGIiIiIiEgSUygUERERERFJYgqFIiIiIiIiSWxEh0KNPioiIiIiItK9ER0KNfqoiIiIiIhI90Z0KBQREREREZHuKRSKiIiIiIgkMYVCERERERGRJKZQKCIiIiIiksQUCkVERERERJLYiA6FmpJCRERERESkeyM6FGpKChERERERke6N6FAoIiIiIiIi3VMoFBERERERSWIKhSIiIiIiIklMoVBERERERCSJKRSKiIiIiIgkMYVCERERERGRJKZQKCIiIiIiksRGdCjU5PUiIiIiIiLdG9GhUJPXi4iIiIiIdG9Eh0IRERERERHpnkKhiIiIiIhIElMoFBERERERSWIKhSIiIiIiIklMoVBERERERCSJKRSKiIiIiIgkMYVCERERERGRJKZQKCIiIiIiksQUCkVERERERJKYQqGIiIiIiEgSG9Gh0MyWm9mj9fX1g10VERERERGRIWlEh0J3f8bdHx4zZsxgV0VERERERGRIGtGhUERERERERLqnUCgiIiIiIpLEFApFRERERESSmEKhiIiIiIhIElMoFBERERERSWIKhSIiIiIiIklMoVBERERERCSJKRSKiIiIiIgkMYVCERERERGRJKZQKCIiIiIiksQUCkVERERERJKYQqGIiIiIiEgSUygUERERERG5EPn5YHb2Iz9/sGvWIwqFIiIiIiIiF+LQIQCemj6ejP89nqemj497fahTKBQRERERETkfJ0/Chg0AHM5M5a7/2MKZUce56/dbOZyZOsiV67mUwa6AiIiIiIjIkNTWBjU1sGtX4seRIx8VvfiPCmhNqwFzWtNPc8l/LqLqa1WDWPmeUygUEREREZHkVV8Pu3cnDn179kBLy0dFT2Sl8eZF01g5ZwIbrp3LzlGNHEyv41haFaTuBYsWDLeyt3Q/D91WxmODclC9M+xCoZmVA18Cxrj7XYNdHxERERERGcJaW6G6uuvWvmPH4oqfzBvPmxdNY8WNU/hgwhR2Zp/mQOox6lMO05paC/bhR2WtNYvsSAGEWzoCYbtwC9+7/ohCYWdm9hiwDDjs7nNjXr8Z+CYQBr7t7l/rahvuvgv4nJn9vL/rKyIiIiIiw0BdXdehr6oqCIbtUlI4OW0yb86bzLsfW8KG3DZ2ZJzgQPgo9eGDtKYeB3v3o+LWmkVWWwHF4flMzipnVu4MKgvm8LHJC5k2oQiAh+6ZzXev2gEpHa2KtKXy4Fsl8MWBOgnnz9x94HZmdhVwCvh+eyg0szCwDbgBqAbeBz5NEBC/2mkTD7n74eh6P+9pS2FlZaWvWrWqbw5CREREREQGVksL7NvXdfCrrY0vP2ECJyvKeGt2ASsLM9kwupkdaXUcsMPUhw4ELX4xrCWbrEgBE1OKmZw9lVm506ksmMuVk+dTMaEIs87NgGcr/Zcr2cdKCLdCWwqlXE7V77/Rl2fhgpjZanevTLRsQFsK3f0NMyvr9PJiYEe0BRAz+zFwu7t/laBV8byY2cPAwwClpaXnuxkREREREelv7kGw6yr07d0bDPrSLi0Nyso4VTGZt66cwcpJqWzIbmRH+CgHOEQd+2lNjW8UspZsstoKKLKFTE4vZ9a4IPgtLZ3PjIk9C37dWfWZZyj63lRaQ8dJaRvD+/c/fUHbG0hD4Z7CImBfzPNq4NKuCpvZeOARYJGZ/Xk0PJ7F3R8FHoWgpbDvqisiIiIiIr3W3ByEu66CX319fPlJk6C8nJNLLuGdz17PynEh1qefYocd5kDbfuqooTVlW9wq1pxDViSfopQg+M3MnU5l/hyWTl7AzEkXHvy6M2nUWH5+3VPc/fI9/OT6HzNp1Nh+21dfGwqhMNE702WIc/djwO/1X3VERERERKTX3INBW7oKffv2QSTSUT49HaZMCYLf0ktYUZrDitERNqTWsqPtEPtb9lHnW2hNWRm3G2seRWZbPoUpi5icdnbwC4X6L/idy+2zltI0q3rQ9n++hkIorAZKYp4XA/v7YsNmthxYPm3atL7YnIiIiIhIcjtzJhi4pavgd/JkfPn8fCgvhyuv5ER5Ee8WZbAy6wzr7Bg7zlSz/0wVdf42rSkxrYQtQPNostryKEy5iNLUcmbmVlCZP5elk+czK29wg99INBRC4ftAhZlNAWqAe4B7+2LD7v4M8ExlZeUX+mJ7IiIiIiIjmnswIXtXoa+6OijTLiMjCH3l5fCxj1E/pYD381J4J62B9ZFD7Di9m/1Ne6jzZ2gNR4NfQ3TdM6PJasunMKWS0pQpzBw3nYvz57C0dD6z8gsJh0IDfvjJaqCnpPgRcDUwwcyqgb909++Y2R8AvyEYcfQxd980kPUSEREREUkaTU3BpOxdBb/Tp+PLFxYGoe+aa6C8nPqyfFZNgHdCJ1jfuI8dJ3ZQ07STOn+D1vAJiG0sPDOazLZ8CsKVlKZPYWbudC7Om8OSyfOZnV9ISljBbygY0CkpBoumpBARERGRpOEOhw51HfpqauLLZ2V1tPZFH/Vleawa28qKtmOsq90WBL/GPdT5/iD4xTozhszWPCaklFCSNYUZ0Xv8lpTOZ3Z+Aakp4YE7dunSkJmSYqDpnkIRERERGZEaGrpv7Wts7ChrBkVFQeC74YaPgl/t5AmsHdXCilN7WXdkM9vrd1DT9AJ1kRpaD5yEAzH7OzOGjNZ88sOXUJJa/lGL3xWlc5lTWEiagt+wppZCEREREZGhJhKBgwe7Dn0HDsSXz8k5q7WP8nKOl4xjfWYD7xzawrrDH7KjfgfVTbuD4Bc+Fb+NprFktOYxPhy0+M3Knc5F+dHgV1BAeuqIbk8a8ZK2pVBEREREZMg6fRp2704c+nbvDu79a2cGJSVB2Lvllrjgd6xoLBsih1lR8wHrDm9me/0mqht/TV1VDa37YoKfG5wZS0bLJPJSLqMkawozxwbB7/KSucwpzCczLXXgz4MMOoVCEREREZH+EInA/v1dt/YdOhRffvToIOjNmgW33vpR6PMpUziWn8PGE7t4p3oD6w59yPb6N6k5+X1q19bQuj5mYJho8EtvyWNSuCP4LcqfzRUl85hdkE92hoKfxBvRoVD3FIqIiIhIvzp5svvWvubmjrKhEJSWBmFv+fK41j6fMoWjWRE2Hd/KiuoPWHf4Q7bX/4rqqj3U7q6hNdQp+DXlkt46iYmhy4PglzudRXlzuKJkLrML88jJSBv4cyHDlu4pFBERERHpSltbMDdfe8jrHPyOHIkvP2YMTJ2a8P4+LynhSFsdHx7bEnT1PLSZ7fXbqW7cTW2khtZQQ8d2osEvrWUS48MlFGdOYebYChblzeHy0jnMLshjdFb6wJ4LGdZ0T6GIiIiISFfq688OfDt3Bv9WVUFLS0fZcBgmTw6C3ic+cXbwGzuWw42H2Vy7lRXVG1h/eAPbDv6C6t17qI1U0xqKGRXUDZrGkdY8iQnhJUHwy40Gv5I5zCqYxNjsjIE/H5J0FApFREREZGRrbYV9+7ru5nnsWHz5ceOCkHfRRXDXXfHBr6QED4c51HiILce3snL/B6w79BpbNz9KzZo9QYufxQa/EDTlknomj/HhpcE8fmMrWDhpNpeVzGF24STG5WQO7PkQ6UShUERERESGv7q6xIFv507YuzcIhu1SUqCsLAh5n/xkfOibMgXGjsXdOdh4kK2121hZ8wHrDv+CbW9vp7pxD8cj1bRZzMigHoLGcaQ2T2J8+EqKM8uYHu3qeVnxbGYVTmJ8TgZmNuCnRaQnRnQo1EAzIiIiIiNES0vQ2tfVSJ61tfHlJ0wIQt7ixXDPPfHBr7gYwmHcnQMNB9hat4339n/AuppX2LppWzT41cQHv0gImsaRcmYS48NXUZxZxoyx01k4aTaXFs9iVuEkJo7OVPCTYUkDzYiIiIhI38nPP3uqBYC8vGAy9q64B8Guq9C3d28w6Eu7tLSO1r7OjylTgukdgIhHONBwgG1124Pgd3gzW2u3xbT4nenYZiQETeNJOTOR3FAxxZkdXT0vLZnFzIKJ5I3JUvCTYUkDzYiIiIjIwEgUCNtfb24Owl1Xwa++Pn6dSZOCkHf55XDfffHBr7AwGPSFIPjtP72fbXXbWXXwDda++a9srYsGv7Zq2ixmWohIGJrGEW6aRG7oY5RkTmH62AoW5s3mkqKZzCyYQGFujoKfJBWFQhEREREZGJmZwYTu7dLTg1a98nJYsuTs1r6cnI+KRjxCzekattdvZ9XBf2ftBx+ytTaYzuF4W83Zwa9xPOEzE8kNXUNxRhkVYytYOGkWi4tnMaNgAkW5OYRCCn4ioFAoIiIiIr3lHszdt3UrbNkS/2/UU9PHc/cX4Cf/Brdvi47u+eUvxwe/goJgQveoiEeoPl3NtrptrNn6GmsPfRi0+DXs6SL4TSDUNJFcu5qi9ha/SbO5pHgGM/InUDw+h3DM9kUkMd1TKCIiIiKJNTTAtm1nh79t2+D06Y5yo0bBzJkwYwb88Icczkyl6KuZtGacJKVxNDV/0cCkxhZwpy3SFgS/+m2sO7yZtQc/ZEvdNvY1BC1+EYuZEzCSAo3jCTVNZKwVBsFvTBD8KoumM6NgAqUTRin4ifRAd/cUjuhQGDP66Be2b98+2NURERERGXrcoaYmcavf3r0d5cyCSdvbw1/sv/n5wfJoudL/Npl9pTUQboVIiKxTo5hxCPbNyOV4WzURi5keoi0l2uI3gTFWFNPVMwh+FfnjKJs4hpSwgp/IhUjaUNhOLYUiIiKS9BoaYPv2xK1+p051lMvJ6Qh8seGvoiK4J7CTlkgLm2s3szTevIQAACAASURBVOLg+7xatZLnN/2Uk6NPQOzteg6cGY2dmsxYCinKLKNiTAULJs7m4qIKpheMZ8rE0aSmhPv9NIgkK40+KiIiIpIM3GH//iDwdQ5/e/cGyyFo1SstDQLf0qXxrX4FBR2tfp2caD7BuqPr+G31u7xZ/R6b6jZwsGVXR8tfWyqMaosPhAAGFnYa/2QN6an6+Cky1Oi3UkRERGS4aWxM3Oq3dWt8q192dhD2liyBhx6Kb/XLyupy8+7O/ob9vH9oFa9UreTdg6vZfmoTdb6/o1BzDnaqkLGR65k+ag6X5VdyTflFPLnjB/z0+N9BOGZQmLY0Hir+UwVCkSFK3UdFREREhiL3YLL3zqFvyxaoqupo9YOOVr/Y7p4zZkBRUZetfu3aIm1srd/KO/uD7p9rjq6lqmkLTcTMGdgwgXBDEQUp05kzdj5XFi/mqvLZLCidyOis9LO2WfovV7KPlcE9hW0plHI5Vb//Rl+dGRE5D+o+KiIiIjJUNTWd3erX/vPJkx3lsrKCoHf55fDAA/GtftnZPdrV6ZbTbDi2gVf3ruDN6vfZVLeBAy07OqZ6iIThdAGZZ+YwM3MmF01YyLWTL+ey8jJmFIzr8WAvqz7zDEXfm0pr6DgpbWN4//6ne3lSRGQgjehQGDP66GBXRURERJKZOxw6lLjVb8+e+Fa/kpIg7N1/f3zLX1FR3Jx+53Ko4RDvH17Ny7tXsPLgarad2khtpBosuq+WTOx0Ebl+NRXZc7gs/2KuK6/k4ilFFIzNxs7RwtidSaPG8vPrnuLul+/hJ9f/mEmjxp73tkSk/6n7qIiIiEhfOXMGduxIHP5OnOgol5UF06efPcrn9Ok9bvVrF/EIO0/s5K2a93h5zwrWHl3HnsbNNNrxjkJNuaQ0FFMQrmBO7gKWFl3CNVPnM790IjkZaX108CIylKn7qIiIiEhfcYfDh8/u6rl1K+zeDZFIR9ni4iDwffaz8a1+xcW9avVr19TaxAfHP+DlPSt4s/o9NtZt4EDLdlqtKSgQCUFDHlnN05mVMZOLJizimrJLWVJeQUX+WE3yLiIJKRSKiIiIJHLmDOzcmbjVrz5mEJaMjCDoXXwx3HdffKtfTs557/5403HeO7Sa3+x+h5UHVrH91CaORarAoqGzNR07Xcg4X0pFzmwuza/khvLFVE4pJm9M71obRSS5KRSKiIhI8nKHI0cST+2wa1d8q19RURD47r03fqTPkpLzavXrqIJTdaqK1/et5NU9K1hzdB27GzfTYEc6Cp0ZQ0pDMSUpt300+ud1Uxcxr2QiWempF3ACREQUCkVERCQZNDcHrX6Jwl9tbUe5jIyghW/RIvj0p+Nb/UaNuvBqtDWz6fiHvLj7Hd6sfo8Patezv3UbrdYQFHCDhklkt5QxO/0WFk1YxLWTL+WqaTMpnzSWUOj8B38REemKQqGIiIiMDO5w9GjXrX5tbR1lCwqCsHf33fGtfqWlF9TqF6u+uZ73Dq7h+Z1v8e7B1Ww7uYljvge31qBAW2rQ/ZNLqciew6X5F3Nj+eUsLi9lwqjMPqmDiEhPKBSKiIjI8NLS0nWr3/GYETfT04MWvgULgvDXPsrnjBkwenSfVcfdqTldw2t7V/LynhWsObqGPY1bOGUHOwo155DSWExJyi3MGRuM/nnDtErmFU8iI00fx0RkcOmvkIiIiAxNnVv92n/etQtaWzvK5ecHrXyf/OTZrX7hcJ9WqTXSyubaLbyw823eqH6PjbXr2d+yneZQzHQTDRPIbpnM7IzrWBSd/P3qqXOYMmnMBc39JyLSX0Z0KNTk9SIiIkNcS0sQ8hK1+h071lEuLQ0qKmDePLjrrvj5/caM6ZeqnW45zcqDq3lhx9usPLiarSc3ctR34dYSFIiEsdMFjGMRFTlB98+byy/n0vIp5GZn9EudRET6gyavFxERkf537NjZc/pt2RJ0A41t9cvLO3tC9xkzoKysz1v9Yh08fZBXqlbw8p6VrD6yhj2Nmzlp+8Gin5NaMkltLKEgZTpzxs5jadFibqpYzLziPNJS+q9eIiJ9RZPXi4iISP9rbe261e/o0Y5yaWkwbRrMmQN33hkfAseO7dcqRjzCluPbeG7nW7yxL5j8vaZlK82huo5CjePIaZ3M7PSrPur+eX3FAkrGj1L3TxEZkRQKRUREpHdqaxNP6L5zZ9AdtN2kSUHQ+8Qnzm71S+n/jyBNrU2sPLCG53e8zYoD77P11CaORnYSCZ0JCkRChBrzGcdcpmXPYXHexdw89QqWlE9jdFZ6v9dPRGSoUCgUERGRs7W2wu7dibt8HomZVD01NWj1mzUL7rgjPvzl5g5YdY82HuWl3St4afc7rDm6lt2NWzhh+8Cik8+3ppPaWEJxyrXMHjuPK4sWc0vFZcwrLiAl3DdTUIiIDFcKhSIiIsmstjZx8NuxI77Vb+LEIOjdfnt88JsyZUBa/dq5O9vqdvLrbW/xRvW7fFC7npqWbZwJxwxKc2YMOS2TmZNxKQvHB5O/3zj9Yopy1f1TRCSRXv8VN7NsoMnd285ZWERERAZfWxvs2ZO4y+fhwx3lUlKCVr8ZM+C22+LD37hxA17t5rZmVu5fy3Pb32JF++ifkZ20hRuCAm6EGicxjhksyp7D4ryLuGnqEq6aOoOcjLQBr6+IyHB1zlBoZiHgHuA+4BLgDJBuZkeA54BH3X17v9ZSREREzq2+PvEgL9u3Q3NzR7kJE4Kgt3z52a1+qamDUvXapjp+s+ttXt4djP65u3Ez9bYXQtHvoNvSSGsqoihlKXOy57O0aDG3TL+c+UWFhEPq/ikiciF60lL4GvAy8OfARnePAJjZOOAa4Gtm9u/u/sP+q6aIiIgAQatfVdXZwW/rVjh4sKNcSgpMnRqEvVtvjQ9/48cPWvXdnZ11VTyz7c1o988N1DRvpSklpsWyOSfa/fMTLBy7kGsmX8Yt0y+hMHf0oNVbRGQk60kovN7dWzq/6O7HgV8AvzCzwflaUUREZKQ6caLrVr8zZzrKjRsXhL1bbukIfTNnQnn5oLX6tWuNtLKiZh3PbX+blQdXseXERo74TtrCJz8qE2qcyDjKWZi2LDr65xKuKp9Ftrp/iogMmJ6Ewn8wsyfd/Z2uCiQKjSIiIkktPx8OHTr79by8jha9tjbYuzdx+DtwoGOdcLij1e+WW+Jb/SZMGJjjOYcTZ07y/I53eHnPiqD7Z8Nm6kNVeCj6ESGSQlpTIcUplzIrZx5Liy7h1ulLmF9YQiikwV9ERAZTT0LhduD/mFkB8BPgR+6+rn+rJSIiMswlCoTtr3/yk0Hw27YtvtUvNzcIezfddHarX9rQaTnbXVfNU1vf4I19Hd0/G1MOgnlQoCWTnNYpzEpfxsLchVxbejkfn3kpBWPU/VNEZCg6Zyh0928C3zSzyQQDznzXzDKAHwE/dvdt/VxHERGRoe/kyWAah+3bg3+7s25dEPZuvDE+/E2YAENoyoSIR1hRvSEY/fPAKrac3MiRyA5aU+o/KhM6M45xPpUFmTexOO8ibpxyBddOnU9muu4sEREZLszde7+S2SLgMWC+u4f7vFZ9xMyWA8unTZv2he3bNUCqiIhcoBMn4oNf7M9dtQwmch7/9/a3U2caeH7HO7y4+x3WHFnLrobN1If24OFoS2YkRNqZQvJTKpg9Zh5LihazrGIpC4pKNfefiMgwYGar3b0y4bKehsLoYDI3E7QWXgf8lqAr6a/6qqL9pbKy0letWjXY1RARkeGgvv7swNf+b+ycfgCFhVBREczt1/5v+yMnp+t9DHIo3FN3gKe3vsFv973LB8c3UN2ylcaU/WCRoEBrOjmtZZSmz2BBdPL3ZTOuIH/0mEGtt4iInL/uQmFP5im8Afg0cCvwHvBj4GF3P92ntRQRERko9fVnB772f48ciS9bVBSEvNtuiw9/U6dCdvbg1L+HIpEIK6s38eto98+tJzdyOLKD1tTjH5UJNY9lnJezIPUaLsm7iJumXsF1UxeRMcgjl4qIyMDpyUAzfwE8CfzX6DQUIiIiQ19d3dmBr/3no0fjyxYXB0HvjjvODn5ZWee3/7y8rkcf7QenzzTy3PYVvLR7BauPrGV344fU2R48pTEo4EZaawGF4bnMyp7H0qLFLJt+JQsKytT9U0QkyfVkoJlrACzwGaDc3b9iZqVAvru/19+VFBERSai2NvH9fdu3w7Fj8WVLSoKgd+edHV08KyqCkT3PN/h1J3Yi+T62t+4wT20Jun9uqF1PTfM2GlKqIdQWFGhLIydSyqz065g/fgHXlF7GbTOWkD86t9/qJCIiw1dPWgrb/TMQAa4FvgKcJJi8/pJ+qJeIiEjg+PHE3Ty3bw+WxSopCYLe7/xO/L1+5eWQmTko1d9UfZR7v/VrnvzDW5lT3Ls5BSORCO9Vb+XZbW+w4sBqtpz8IOj+mdbR0hlqGcU4pjIvfUnQ/bN8CTdMvZh0df8UEZEe6k0ovNTdLzKztQDuXmtmQ2fSJBERGb6OHes6+NXWdpQz6wh+n/xkfPCbMmXQgl9XTje1cPM3/4nqwn/hlm/tYfNX/oLsjMRhraH5DM9ve5cXdweTv+9q+JA6242ndtzCnxbJoyA8g9k5n2RJ4SUsm34lCwumqvuniIhckN6EwhYzCwMOYGYTCVoORUREuufeEfwShb/Owa+0NAh6d98dP6pneTlkZAzecfTSA48+RU3JNyG9jurib/Lgv83hp1+8i6rjR3hm61u8vu9dPqhdT3Xz1qD7Z7glWDGSQjalzEy7igXjFnB16WXcPvNK8keNG9wDEhGREak3ofBbwL8Dk8zsEeAu4Mv9UisRERl+2oNfV6N61tV1lDWDyZODoHfPPfGDu5SXQ3r64B1HH3n89Y386tT/wnNPgjmeVs/Pmr9Ayrf+E23pR8CCaSlCbVnkhqYyN/MTLJ5UyY3lV3DDtEoyUob/ORARkeGhx6HQ3Z8ws9UEcxQacIe7b+63momIyNDjHozc2dWonvX1HWVDoY7gd++98cFvypQREfwS2Vt7lG+veZr/tervaRu/5aPwRygC6XVEThdz44QvsKSwkmUzrmJR/nR1/xQRkUHVk3kKzaMz3Lv7FmBLd2VERGSYcw/m6utqVM8TJzrKhkJQVhYEvc98Jj74lZWN2OAX6+CJOr6z5lme2fUiG0+t5HTarmAU0LEEX6HGMhg1poHf3Pv/BqOqIiIiCfWkpfA1M/sF8JS7721/MTrIzFLgfuA14PF+qaGIiPQ9dzh8OHE3zx074oNfONwR/C67LH5wl7IySEuuMceOnTrF42uf51c7X2DDiZWcSNsGoVZwY3R4Kldl38sd024mNaeOP3rnv9BmZz5aN+zpfGPp3w9i7UVERM7Wk1B4M/AQ8CMzmwLUAZlACHgR+Ad3X9d/VRQRkfPiHkye3lXwO3myo2x78KuogCuuiA9+kycnXfCLdaKxkR+se4lfbn+eNXXvUJe6BcLN4EZOeDKXZ93F8vIbeWjhbeTljI9b97X9r/Hvu5/CQy1YJJU7ym/jwRkPDtKRiIiIJGa96fVpZqnABKDR3evOVX6oqKys9FWrVg12NURE+p57MEl6ovv7duyAU6c6yobDwb18sYEvNvhpXjsAGpubeXL9q/xs63OsOv42x1I+hJQmADJbipiTfSkfL7uRzy+6nZIx+d1u63TLaSp+NJMDjdUUZpWw7Z7NZKdmD8RhiIiIxDGz1e5emWhZb0Yfxd1bgAN9UisREemZ9uDX1aiepzvmsSMlpSP4XXVVfPhT8EuouaWVn218i59sfpZ3j73J4dAmiM4NmB7OZ0HGDdw8+To+v+gTTBtX2qttZ6dm89KyF7j7lbv5yXU/USAUEZEhqVehUERE+ok7HDjQdfBraOgom5ISTNswbRpcffXZwS9Ff9q709LaxtOb3+XJD59hxZE3OcgGPC3oSpsaHs/s9KXcVHI9n1t0B3MmTrvg/c0ZN4eNn9x4wdsRERHpLz0ZfXSOu28aiMqIiIxokUh88Ovc5TM2+KWmdgS/a6+ND36lpQp+vRCJOC9sXcsPNj7N24d+S42vI5Ie3AGREhpLRdolXFd0LQ8t/AQX58/S9BAiIpJ0evKp4gfARQBm9nl3/3b7AjPLcveGLtcUEUk2kQjs39/14C6NjR1l09I6gt9118UHv5ISBb/z5O68vmMT3/vgad448DpVbWuJZBwFIBzKoSz1Iq4uuJoH59/BkpKFCoEiIpL0evKJI/Z/y/8IfDvm+ZvAxX1aIxGRoS4SgZqaxMFv586zg9/UqUHQu/76+MFdSkqCwV/kgrg771bt4LF1T/Ha/lfZ3bqGtoxDAIRCmRSnLuSqvM/zu/Nu57opiwlZaJBrLCIiMrT0JBTGDk/a+evUQfmf1czuAG4FJgH/5O4vDkY9RGQEi0SgujrxqJ47d0JTU0fZ9uBXUQE33hgf/IqLFfz6mLuzvmYvj617ipf2vcrO5lW0ZNYAYJZOYdo8rphwH5+Zexu3TltKOKTzLyIi0p2ehMJ8M3sAWM/ZobDn81lEmdljwDLgsLvPjXn9ZuCbQBj4trt/rattuPuvgF+ZWS7w9wTzJYqI9E578Es0uMvOnXCmY9Jx0tM7gt/NN8cHv6IiBb9+tvngfr695mle3PcK25repzljL5iDpZKXPpvLxt3JvbNv5Y4Z15AWTt45FUVERM5HT0LhXwGVwINAsZltArZEHxPOY5+PA/8IfL/9BTMLA/8E3ABUA++b2dMEAfGrndZ/yN0PR3/+cnQ9EZHE2tq6Dn67dsUHv4yMjuD38Y/H3+NXXAwhdTscKLuOHuE7q5/h+aqX2dzwHk0ZuyEUAQ8zMW0ml4x9mE/N/Dh3z7mRjJSMwa6uiIjIsHbOUOjuj8Y+N7NiYD4wD/htb3fo7m+YWVmnlxcDO9x9V3QfPwZud/evErQqxrFgVICvAc+7+5pE+zGzh4GHAUpLezevlIgMM21tsG9f4lE9d+6E5uaOshkZQcibMQOWLYsPfkVFCn6DpKa2jsfWPMezu19k46mVNGTsgFAbeIjc9AqWjrmfu2Z8nPvm3kJOmub6ExER6Uu9HtrO3auBajM7AfRV2ioC9sU8rwYu7ab8F4HrgTFmNs3d/zVBPR8FHgWorKzsdTdXERli2tpg797Eg7vs2hUf/DIzg5A3cyYsXx783B7+CgsV/IaAIydP8fiaF/jVjhfYcHIlp9K3QbgF3BiTNoWPjbqHT1Tcwu/Ov5XcjLGDXV0REZERrVeh0MwWAvcCnwIOATOB/9QH9Ug0HniXQc7dvwV8qw/2KyIDIT8fDh06+/W8PDh4sON5a2v3wa+lpaNsVlYQ9GbPhttui7/Hr6BAwW+IqTvdyPfXvswvtz/PuvoV1KdthpSg625OWglX5HyC26bexAMLlpOXPXGQaysiIpJcejJ5/XTgHoIweAr4GXC1u+8xs919VI9qoCTmeTGw/0I3ambLgeXTpk270E2JyIVIFAjbX/+jP+oIf7t3Jw5+c+fCHXecHfw0v9yQdaqxmSfXv8bPtz3Pqtq3qU3ZBKnBVB1ZaQVckvNxlk25gQcX3E7J6MJBrq2IiEhy60lL4RbgfeAud9/YaVlfdct8H6gwsylADR0h9IK4+zPAM5WVlV+40G2JSCfu0NAAx4+f/Th2LP55d7797SDozZ8Pd94ZH/zy8xX8homGM838/IN3+OmWX/Pesbc4Et4IaacASE+dyMKsa7ll8vU8tPAOpuWWDW5lRUREJE5PQuHvEIS0l8zsZeCnwAvu3tL9aomZ2Y+Aq4EJZlYN/KW7f8fM/gD4DcGIo4+5+6bz2b6I9JI7nD6dOMx1FfLaH7Ejd3aWng7jx8O4cd3v/+RJBb9hqLm1jac2vseTHz7LyqNvcsg24On1AKSl5jIn8wpuKrmehxbezpwJ0we5tiIiItIdc+9ZY5+ZZQN3AJ8GLgGeA5a7+/lMSzGgKisrfdWqVYNdDZH+5R4ErJ6Euc6vt3TzHU9mZke46/zo7vXMzI5tdBf6evg3SAZXa1uE5zav5cmNz/LW4d+yn/V4RtAKnNI2iqlplVxffC0PLriDi/LmYAr6IiIiQ4qZrXb3ykTLejzQjLufBp4AnjCzccAngbI+qWE/0T2FMiy5Q31971rs2l9va+t6u9nZ8aFtzpxzh7zc3PhwJ0mjLRLh1e0f8v0NT/Pmwd+yL7KOSGYwRWw4NZvy1EVcU3g1D86/g8uLLlIIFBERGcZ63FI4nKmlUAZFJNIR7nrTLbO2tvtwN2pU71rs2sNdevrAHXtnPR19VAZNJOK8s3sH3133NK8feI2q1jW0ZR0AINSWQUnKAj6WfzW/O+82ri69lHAoPMg1FhERkd7ok5ZCkaTV1gZ1db3vlllb233XyNGj48Pc5MnnDnnjxkFq6sAde19R8Bty3J1VVVV8d90zvFLzKruaV9OaVQ3mWEoahWlzuXLSPdw3dzk3l19JSkj/XYiIiIxU+l9ekkdraxDUetsts66u++2OHRsf5qZOPXdL3tixwzPcybDl7nxQs5/vrnuWl/a9wramVbRkVYFFIJxCftYsLp94G/fOXsbyimtIDw9iy7KIiIgMqB6HwujooE+4e20/1qdP6Z7CEaqlpSPc9aZbZn1919s0C7pYxoa4iopzd8scOxbC6kYnQ9PWg4f57tpf83zVy2xteI8zWbsh1AahEBMzZ7B43Oe4e9at/M6MG8hKzRrs6oqIiMgg6U1LYT7wvpmtAR4DfuND/IZEzVM4xDU3934glePHgxE2uxIKBeGuPbTl5cGsWefuljl2bLCuyDC2+8hxHl/zAr+ueokPT71LY+Z2CLeCGeOypnJV7mf41IxbuXv2zYxKGzXY1RUREZEhojejj37ZzP47cCPwIPCPZvZT4DvuvrO/KijDwJkzvZ/f7tixYG68roTD8aGtsBDmzTt3t8zRoxXuJGnsrz3B99a+xNM7X+CDk+9yOmMbpARzR47JLOPSsZ/izopb+MzcW8nNyB3k2oqIiMhQ1at7Ct3dzewgcBBoBXKBn5vZS+7+p/1RQRlAjY29n9/u+HFoaOh6mykp8cGtpAQWLOh+IJXx44MRNjXEvUicIydO84O1r/CrHb9hff0KTqRvgdRGAHIyi1gy+jZun3oTn523jPzsvEGurYiIiAwXvbmn8A+B+4GjwLeBP3H3FjMLAdsBhcKe6O+h+d2DkHY+3TKbmrrebmpqR5AbPx6mTIGLLz53t8ycHIU7kfNUe6qRH69/k19se57VdW9Tl/ohpAUt7FmZeSzOuZFl5TfwwPw7KBlVNMi1FRERkeGqNy2FE4A73b0q9kV3j5jZsr6tVt8YkgPNJAqEiV53D7pX9rZb5vHjQXfOrqSnxwe4adPOPZjKuHGQlaVwJ9LPTjSc4RcbVvLTrc/x/vG3OJayEdJPAJCeMZ5FWVfx8bIbeWDB7UwbO2WQaysiIiIjRY8nrzezPwY6F64HVrv7ur6uWF8aUpPXdxesliyJD3ctLV2XzczsWZjr/HpmZt8fk4icl4YzLTy9cTVPfvgs7x59k8PhDyAjGOA5tW0MMzIv4ebS63lg/u3MHj8D0xczIiIicp76avL6i4FK4Jno81uB94HfM7OfufvfXVg1hfR0mDPn3CEvN1fhTmQYampu5fkPN/DEpmd558gbHLT1eOZRAFIyc5iZfjE3FF/HA/NvZ9GkeQqBIiIiMiB6EwrHAxe5+ykAM/tL4OfAVcBqQKHwQr3yymDXQET6UHNrGy9v2cwTG5/lzUOvU+PriWQF9w6HMzKZmrqIa4v+A/fPu4PLCi8iZBo5V0RERAZeb0JhKdAc87wFmOzujWbWzU1sIiLJobUtwhs7tvODDb/m9QOvsbdtLZGs/WBOKD2NkpQFXF3wWX533m1cVXwZKaFeDQAtIiIi0i9684nkSWClmT0Vfb4c+JGZZQMf9nnN+sCQHGgmL6/r0UdFZFiJRJyVu/bwvfXP8mrNa+xuXUNb9l4wx9JSKArP5cq8u7hv7m3cMHkpaeG0wa6yiIiIyFl6PNAMgJldDCwFDHjL3YfI6C3dG1IDzYjIsOXurKmq4fF1v+bl6lfZcWYVrdlVEGoDD5NvM7li0lXcO3sZt5ZfS0ZKxmBXWURERATog4FmLBjtoNjdVxPcPygiMuK5Ox/UHOIH617gN1WvsLXpPZqzd0GoFVKNiWnTuXT8A9wzaxl3VNxAdmr2YFdZREREpNd6FArd3c3sVwQjkIqIjEjuzraDx/jBupd4bs+LbD79Hk3ZOyDcDKkwLqWcynGf5lMzbuWumTczJm3MYFdZRERE5IL15p7ClWZ2ibu/32+1EREZYHuO1PHDta/wzK4X2XhqJQ1Z2yClCcIwJqeUy3Pv5M6KW7hn1seZkDlhsKsrIiIi0ud6EwqvIZiTcA9wmuC+Qnf3+f1RMRGR/rD/+EmeWPs6v9r5GzacWMmpzM2Q2gAhyMkpYMnoW/lExU3cO3sZBdkFg11dERERkX7Xm1B4S7/VQkSknxyuP81P1r/DL7c9z5q6tzmRsRnSTgKQmTORxaOu47byG/nMnNuYPLp0kGsrIiIiMvB6Ewr3AvcB5e7+FTMrBfKBqn6pWR8YklNSiEi/On6qkV+sf4+fbX2OVbVvU5u6CTLqAEjPyWVR9hXcOuUmfnfucqaNmUowjpaIiIhI8urxlBRm9i9ABLjW3WeZ/f/t3Xd8VFX6x/HvkxCJAYwCUoO0KCgIiCioKFhQASmhSFvByoIFWXXV1bWv664/USwoYkFZEawQkCAISFEBC0VBxCCKJCBNmkAg5fz+yIABAkmYSe6dyef9euW1M+eW88zxcneeOeeeYydJmu6cO6c4AwwFlqQAIteO3XuV/O1ivfPDFC3c8pk2l/lOOn6LJCkm5wQ1PL6Frqx9ma5p3FmNK51BEggAAEqloJekCGjpnGtuZoslAADXIQAAIABJREFUyTm31cxYiRlAidqVkamU5d9p3PdT9MWmudoQ9a0Ut1GSVKZcnBrGnq3La12qaxp3UvMqzRRlUR5HDAAA4G9FSQozzSxakpMkMztZuT2HAFBs9uzL1IwVP+itZSn6bMNsrbdv5cqtkyRFxZVV4nFn6ZKaN6j/mV3UqloLRUdFexwxAABAeClKUvicpAmSqprZ45J6SPpnsUQFICIlr/hMvWb01juXjVeX01vnu8/ezCzNXvmTxn6Xojm/faq0nCXKKZ8mmZOVi1HtMk3UtkYfXdO4iy6s2UoxUTEl/CkAAAAiS6GTQufcWDP7RtKlgaKuzrkVxRMWgEizcec29ZjZRVnHbVWPmV2VnrBKVSqcqMysbH2+6heN/e5jzUqfpV+yFiun/BopKkcWF62a0Y10UbUk9Wt0lS495SKVjS7r9UcBAACIKIVOCs2srKTmkuIDx/U0MznnHi2u4ABEjhZvdVJW9A7JnLKit6v+a+erhjXRT5nfKLv8z1JUtnR8lKpFNdAFVQaqzxkd1b7uJYorE+d16AAAABGtKMNHkyVtl/SNpL3FEw6ASHR98mNaqy+l6Kzcgugs/RG3Qj+6FTr5+FPV6uQB6nN6R3Wsd5lOOO4Eb4MFAAAoZYqSFCY4564stkgARKw30/5Pitl3WHlUVnltvOVHDyICAADAfkWZq/0LMzuz2CIpBmbWycxGbd++3etQgFJr7ZadOum3JCn7kBVsso/TdQn3eBMUAAAADihKUtha0iIzW2lm35rZd2b2bXEFFgrOucnOuYHx8fFehwKUSvNT1+m8B8YqZ905OnlfcyknsHB8dhmdopZ6tQsTGAMAAHitKMNH2xdbFAAizlvzvtdfX/1ECRXLa/p9PVU5/i+q9lY1ObdXZbLj9dWASV6HCAAAABWtp/BXSRdKGuCcW6PcReyrFktUAMJWdk6O7hs/T9eN/FjnnVpdXzzaV2ckVFKVCifqyprtJZnev3SiqlQ40etQAQAAoKL1FL4oKUfSJZIelbRT0geSzimGuACEoZ179qn/iyn6aNFq3XRJEz074GLFlIk+sL1uxRqqtLXiEReuBwAAQMkrSlLY0jnX3MwWS5JzbquZHVfQQQBKh583blfSsIn6Yd3venbAJRrcrqnM7KB9nNxhZQAAAPBWUZLCTDOLVu6wUZnZycrtOQRQys1bkaaewycpO8dpyj3ddGnj2vnu55yTiaQQAADAT4qSFD4naYKkKmb2uKQekpg6ECjlXvv0O906eqbqV4nXh3d21WnVTzrivi73NyUAAAD4SKGTQufcWDP7RtKlkkxSV+fcimKLDICvZWXn6O635+j5jxer3Zm19fZtHXViudgCj2P4KAAAgL8UpadQzrkfJP1QTLEACBPbdmWo7/NT9Ml3a3TblWfpyb5tVCa64MmMGT4KAADgP0VKCgHgx/Vb1W3YRP20cbteuuEy3XhJk0If60RSCAAA4DckhQAKbeayNer97EeKjjJN+0d3XXR6La9DAgAAQJBICgEUyDmnlz5Zqjv+96ka1qioCXd2Vd0q8UU/D0tSAAAA+E6BSaGZ7ZTynTLQJDnn3AkhjypEzKyTpE6JiYlehwKErcysbA0d86lGzfxWHc+qpzE3t9cJcWWP6Vw8UwgAAOA/BSaFzrkKJRFIcXDOTZY0uUWLFjd5HQsQjrbs3KNez07WnBVp+nunc/TY1RcoOqrgCWUAAAAQPoo0fNTMTpJ0qqQD88475+aGOigA3vs+bYuShk1U2u9/aPSgK/WXC88I+pwMHwUAAPCfQieFZnajpNslJUhaIqmVpPmSLime0AB4ZeqS1er3fIriypbRjH/21Hmn1gjJeRk+CgAA4D9FGQd2u6RzJK1xzl0s6SxJm4olKgCecM7pmSlfq8tTE1Wvary+eLRvyBJCAAAA+FNRho9mOOcyzExmVtY594OZNSi2yACUqL2ZWbr59ZkaM3e5up1zql4fdKXKxcaEtA7WKQQAAPCfoiSFaWZ2oqSJkj4xs62S1hVPWABK0sbtu9XjmUman7pO/0xqpQe6naeoqNAnbzxTCAAA4D+FTgqdc0mBlw+b2aeS4iVNLZaoAJSYpWs2qduwidq4Y7fG3tpRV5/HAAAAAIDSpNDPFJrZf/e/ds7Ncc5NkvSvYokKQIlI/nqV2jwyXlk5Ofr0wV7FnhAy0QwAAID/FGWimXb5lLUPVSAASo5zTv9JXqgez0zS6TUrav5j/dSiXrXir5fhowAAAL5T4PBRMxss6WZJ9czs2zybKkj6vLgCA1A89uzL1MBR0zV+/kr1Pr+hRt3UTscfF9oJZQAAABA+CvNM4dvKfXbwCUn35inf6Zz7vViiAlAs1m39Q92fTtbXqzfoX1e31t2dzynRnjuGjwIAAPhPgUmhc267pO2S+phZU0kXBjbNk0RSCISJr1f/pu5PJ2v77n16/2+d1aVFYonHwJIUAAAA/lOUiWaGSBorqUrg7y0zu624AgMQOu/OX6mLH31HZaKiNOeh3p4khAAAAPCnoqxTeKOkls65XdKB2UjnS3q+OAIDELycHKdHP/xCj09YqPNPq6H3hnZWlfg4z+JhohkAAAD/KUpSaJKy87zPDpQB8KFdGZm6buRUTfhqlQa0aaQR112qsjFF+ScfejxTCAAA4D9F+YY4WtJCM5sQeN9V0muhDwlAsH7dvEPdnk7Wd79u1pP92mho++a+6KFzcl6HAAAAgEMUZkmKMs65LOfc02Y2W1Jr5fYQXuecW1zcAQIomvmp69TzmUnasy9LE+/qovbN6nkd0kH8kJwCAADgT4XpKfxSUnNJcs4tkrSoWCMqgJmdLul2SZUlzXTOveRlPICf/G/e9xr06ieqVam8pt/XU2ckVPI6pIMwfBQAAMB/CjP7aMi+wZnZ62a20cyWHVJ+pZmtNLNVZnbvkY6XJOfcCufcIElXS2oRqtiAcJadk6N/jJur60d+rPNPq6HPH+nru4RQYkkKAAAAPypMT+HJZnbHkTY6554uQn1vSHpB0pj9BWYWLWmEpHaS0iR9ZWaTJEVLeuKQ4693zm00s86S7g2cCyjVduzeq/4vTtWUxas18NImGt7/YsWUifY6LAAAAISJwiSF0ZLKKwQ9hs65uWZW55DicyWtcs6tliQzGy+pi3PuCUlXHeE8kyRNMrMpkt4ONi4gXK3euE3dhiXrh3W/67lrL9Hgds28DumoWJICAADAfwqTFK53zj1ajDHUlLQ2z/s0SS2PtLOZtZXUTVJZSSlH2W+gpIGSdMopp4QiTsBX5q5Yq6uHT1aOc0q5p7suaez/65xnCgEAAPynMElhcX+Dy+/8R5y33jk3W9Lsgk7qnBslaZQktWjRgnnwEVFenfWtbntjlupXideEu7rq1GoneR0SAAAAwlRhksJLizmGNEm18rxPkLSumOsEwlJWdo7+PnaOXpi2WJc3qa2xt3bUieVivQ6r0Bg+CgAA4D8FJoXOud+LOYavJJ1qZnUlpUvqLalvKE5sZp0kdUpMTAzF6QBPbd2Vob7PTdGMZWt0e/vm+k+fi1QmujATCPsHw0cBAAD8p0S/UZrZOEnzJTUwszQzu8E5lyXpVknTJK2Q9K5zbnko6nPOTXbODYyPjw/F6QDP/Lh+q1o/NE5zVqzVyze101N/aRt2CSEAAAD8qTDDR0PGOdfnCOUpOsqkMUBpNuO7Nerz3EcqEx2l6ff1UOuGCV6HdMxYpxAAAMB/SjQpLGkMH0U4c87pxU+W6M7/zVbDGhU18a6uqnNyePd680whAACA/0T0+DOGjyJcZWZl65bXZ2rom5+qQ7N6mvdwn7BPCAEAAOBPEd1TCISjLTv3qNezkzVnRZru7nSOHru6taKiIqN3jYlmAAAA/IekEPCR5Wmb1W1YstK3/qE3Bl+pfq3P8DqkkGL4KAAAgP+QFAI+kbJ4tf7yQoriypbRjPt7qtWpNbwOCQAAAKVARCeFTDSDcOCc0zMp3+jecXPVtHYVfXhHF9WqVMHrsIoFw0cBAAD8h4lmAA/tzczSDS9P0z1vz1W3c07V7Ad6RWxCKLEkBQAAgB9FdE8h4Gcbtu9Sz2cma37qOj3Y7Tzdn9QqYiaUAQAAQPggKQQ8sHTNJnUbNlGbdu7RuCFXqUfL07wOqUQw0QwAAID/RHRSyDOF8KOJX6VqwEtTdVJcrGY/2EvN61b1OqQSwzOFAAAA/sMzhUAJcc7piYkL1XP4ZDVOqKz5j/UtVQnhfvQUAgAA+EtE9xQCfrFnX6ZuGjVd78xfqT4XNNSoGy9X7HGl75+fk/M6BAAAAByi9H0rBUrYuq1/qNvTyVr08wY93qu1/t7pnFLbW8bsowAAAP5DUggUo69++k09nknW9t379P7QzurconQ/38ozhQAAAP5DUggUk3fm/6AbX56mqvHlNPfh3mpyyslehwQAAAAcJqKTQmYfhRdycpwe+eAL/XviQl3QoKbeG9pJJ58Q53VYvsCSFAAAAP7D7KNACO3KyFTv5ybr3xMX6to2jTT9vh4khHkwfBQAAMB/IrqnEChJv27eoaRhyVq2drOe+ktbDbnyLHrFAAAA4HskhUAIfPHjOvV8ZpIyMrM06e9ddUXTul6H5EsMHwUAAPAfkkIgSGPmLtfg12bolEoVNOOfPXV6zUpeh+RbLEkBAADgPySFwDHKzsnR/eM/07ApX+viRrU0fshVqlj+eK/DAgAAAIqEpBA4Bjt279U1I1KUsuRnDbqsqZ6+pq1iykR7HZbvMdEMAACA/0R0UsiSFCgOqzduU9JTyVq5/nc9d+0lGtyumdchhQ2eKQQAAPAflqQAimDO92t1/gNva/22PzT13u4khAAAAAh7Ed1TCITSK7O+1ZA3Zimx6omacGcXJVY7yeuQwg7DRwEAAPyHpBAoQFZ2ju56a7ZGTF+iK5rU0djbOio+rqzXYYUlho8CAAD4D0khcBRbd2Wo73NTNGPZGg1tf7b+0/dCRUdF9KhrAAAAlDIkhcAR/Lh+q5KGTdTPG7dr1E2X67q2jb0OKeyxTiEAAID/kBQC+Zjx3Rr1ee4jlYmO0vT7e6p1g5pehxQRnHOKMnpaAQAA/ISkEMjDOacR05forrdm6/SalTThzi6qczKz1wIAACBykRQCAfuysnX7m7P06qzv1Kl5fb15c3tVOP44r8OKKEw0AwAA4D8RnRSyeD0Ka/POPer17GTNXZGmezqfq0d7XqCoKJKXUGNJCgAAAP+J6Id7WLwehbE8bbPOf+BtLVy1Xm/e3F7/6tWahLAY0VMIAADgLxHdUwgUZMqi1bpmRIrKxcZo5j+vVsvE6l6HFNGcnNchAAAA4BAkhSiVnHN6esrX+sf4eWpWu4o+vKOLEipV8DqsiMeSFAAAAP5DUohSZ29mlga/NkP/m/e9erQ8Ta/99QrFlY3xOqxSg6QQAADAX0gKUaps2L5LPZ6ZpAWp6/VQ9/N0f1IrnnErQc4xfBQAAMBvSApRaiz5ZaO6PZ2szTv3aPyQq9S95Wleh1TqsCQFAACA/5AUolSY8FWqrn1pqiqWi9XsB3uped2qXodUKrEkBQAAgP+QFCKiOef0xMSFeuj9L3Ru/Wp6/2+dVf2k8l6HBQAAAPgGSSEi1p59mbpp1HS9M3+l+l5wul6+sZ1ij+OS9xLDRwEAAPyHb8iISOm/71T3ZyZp0c8b9O/erXXXVeeQjPgAS1IAAAD4D0khIs5XP/2m7k8na2fGPn3wty7qdHZ9r0MCAAAAfIukEBFl/Bc/6KZR01TtxHKack9vnXnKyV6HhDyYaAYAAMB/IjopNLNOkjolJiZ6HQqKWU6O0yMffKF/T1yo1g1q6t2hnXTyCXFeh4VD8EwhAACA/0R5HUBxcs5Nds4NjI+P9zoUFKM/Mvap17OT9e+JC3V928aadl8PEkIAAACgkCK6pxCRb82mHUp6eqKWr92iYde01W1XnEVPlI8xfBQAAMB/SAoRtj5fma6ewydpX1aOJt+dpMub1PE6JBSA2UcBAAD8h6QQYenNOcs1+LVPVLvyCZpwV1c1rFHR65AAAACAsERSiLCSnZOjf4ybp2dSvtEljU7RuCEdVbH88V6HhUJiohkAAAD/ISlE2Ni+e6+uGZGiqUt+1s2XN9NT/doopky012GhCHimEAAAwH9IChEWftqwTUlPTVTqhm164bpL9dfLmnodEo4RPYUAAAD+QlII35u9/Ff1eu4jSdLUe7qpbaNTPI4Ix8rJeR0CAAAADkFSCF8bNfNb3f7mLJ1a9URNuKur6lc90euQEASGjwIAAPgPSSF8KSs7R3e+NVsvTl+i9s3q6n+3dFB8XFmvw0IIkBQCAAD4C0khfGfrrgz1ee4jzVz2q/7W4Ww90edCRUdFeR0WQoDhowAAAP5DUghf+WHd7+o2bKJ+2bRDrwy8XNe2aex1SAghlqQAAADwH5JC+Mb0b39R3+en6LgyUfrk/p66oEFNr0NCMWD4KAAAgL+QFMJzzjk9P22x/v7WHDWqVUkT7uiq2ief4HVYKAbOMXwUAADAb0gK4al9Wdm6bfRMvT57mbq0qK83BrdX+djjvA4LxYThowAAAP5DUgjPbN65R1cPn6R5P6TrH11a6uEe5ysqioQhkrEkBQAAgP+QFMITy9ZuVtKwiVq/bZfG3NxefS443euQAAAAgFIpLOf5N7NyZvaNmV3ldSwouo8W/aQLHx6nvZnZmvXA1SSEpYgTPYUAAAB+U6JJoZm9bmYbzWzZIeVXmtlKM1tlZvcW4lT3SHq3eKJEcXHO6anJX6nb08k6rXpFzX+sr86tX93rsFCCeKYQAADAf0p6+Ogbkl6QNGZ/gZlFSxohqZ2kNElfmdkkSdGSnjjk+OslNZH0vaTYEogXIZKxL0uDXvtEYz9boatbNdArAy9XXNkYr8MCAAAASr0STQqdc3PNrM4hxedKWuWcWy1JZjZeUhfn3BOSDhseamYXSyon6QxJe8wsxTmXk89+AyUNlKRTTjkllB8DRfTbtl3q8cwkLVy1Xg/3OF/3dW1Jb1EpxUQzAADAjzIzM5WWlqaMjAyvQwlabGysEhISFBNT+A4YP0w0U1PS2jzv0yS1PNLOzrn7JcnMrpW0Ob+EMLDfKEmjJKlFixYsjuaRxb9sVLdhE/X7rgy9c/tV6nbuaV6HBA8xfBQAAPhRWlqaKlSooDp16oT1dxXnnLZs2aK0tDTVrVu30Mf5ISnMr9ULTOKcc2+EPhSE0odf/qjrRn6sSuWP1+wHe+usOlW8DgkAAAA4TEZGRtgnhJJkZqpUqZI2bdpUpOP8kBSmSaqV532CpHUexYIQcM7p8QkL9MgH89Uysbre/1tnVTuxnNdhwQcYPgoAAPwq3BPC/Y7lc/ghKfxK0qlmVldSuqTekvqG4sRm1klSp8TExFCcDoWwe2+mbnx5mt5b+KP+cuEZeun6yxR7nB8uM/gBS1IAAAD4T0kvSTFO0nxJDcwszcxucM5lSbpV0jRJKyS965xbHor6nHOTnXMD4+PjQ3E6FCD99526+LF39f6XP+o/fS7U63+9goQQAAAAEWfs2LGqU6eOoqKiVKdOHY0dOzboc2ZkZOjcc89V06ZN1ahRIz300EOH7eOc05AhQ5SYmKgmTZpo0aJFQdcrlfzso32OUJ4iKaUkY0FoffnTevV4epJ2ZuzThDu6qmPzel6HBB9iohkAABDuxo4dq4EDB2r37t2SpDVr1mjgwIGSpH79+h3zecuWLatZs2apfPnyyszMVOvWrdW+fXu1atXqwD5Tp05VamqqUlNTtXDhQg0ePFgLFy4M7gPJH8NHiw3DR0vGuM9X6KZXpqvGieWVcm93Na5V2euQ4FM8UwgAAPxu6NChWrJkyRG3L1iwQHv37j2obPfu3brhhhv0yiuv5HtMs2bNNHz48KPWa2YqX768pNwlMjIzMw/7MT05OVn9+/eXmalVq1batm2b1q9fr+rVqxfmox1RiQ4fLWkMHy1eOTlOD7z7mfq/OFXn1q+mLx7rS0KIAtFTCAAAwtmhCWFB5UWRnZ2tZs2aqUqVKmrXrp1atjx4pb709HTVqvXnHJ0JCQlKT08Put6I7ilE8fkjY58GvDhVk775Sde3baznr7tUx5WJ9jos+JwreLUZAAAATxXUo1enTh2tWbPmsPLatWtr9uzZQdUdHR2tJUuWaNu2bUpKStKyZcvUuHHjA9udO/y7VCh+cI/onkIUjzWbduiiR8bro0Wr9fQ1bTXyxnYkhCgUho8CAIBw9/jjjysuLu6gsri4OD3++OMhq+PEE09U27Zt9fHHHx9UnpCQoLVr1x54n5aWpho1agRdX0QnhWbWycxGbd++3etQIsZnK9N13oNj9evmnfro7iTddmVzhgOiSEgKAQBAOOvXr59GjRql2rVry8xUu3ZtjRo1KqhJZiRp06ZN2rZtmyRpz549mjFjhho2bHjQPp07d9aYMWPknNOCBQsUHx8f9POEUoQPH3XOTZY0uUWLFjd5HUskeGPOMt382gzVOfkETbizqxrUqOh1SAgzDB8FAACRoF+/fkEngYdav369BgwYoOzsbOXk5Ojqq6/WVVddpZEjR0qSBg0apA4dOiglJUWJiYmKi4vT6NGjQ1J3RCeFCI3snBzd+/Y8DZ/6jS5tfIrGDblKJ5WL9ToshCGWpAAAAMhfkyZNtHjx4sPKBw0adOC1mWnEiBEhr5ukEEe1ffde/eWFKfp46S+65fJmeuovbVUmOqJHHaOYMXwUAADAX0gKcUSrftuqbsOSlbphm0Zcf5kGXtrE65AQ5vKbMQsAAADeiuikkMXrj93s5b+q13MfSZI+vre72pxRq4AjgIIxfBQAAMB/InocIIvXH5uXZyxV+/9+qGrx5fTFo31JCBFSDB8FAADwl4juKUTRZGZl6863ZuulT5aqfbO6euuWDjohrqzXYSGCMHwUAADAf0gKIUn6/Y896vPcFM1a/qvu6Hi2/t37QkVHRXRHMjzgxOL1AAAAfsO3fuiHdb/rggfH6bOV6Xp14BX6b982JIQoFjxTCAAAwl61apLZ4X/VqgV96jp16ujMM89Us2bN1KJFi8O2O+c0ZMgQJSYmqkmTJlq0aFHQdUr0FJZ605b+rH4vpKhsmWhNv6+HLmhQ0+uQAAAAAP/asKFo5UX06aefqnLlyvlumzp1qlJTU5WamqqFCxdq8ODBWrhwYdB1RnRSyOyjR+ac03MfL9bdY+eoca3K+vCOLqp98gleh4UI5xzDRwEAgM8NHSotWXJsx7Ztm395s2bS8OHHHNJ+ycnJ6t+/v8xMrVq10rZt27R+/XpVr149qPNG9BhBZh/N376sbA169RPd9dZsdT67vuY81IuEECWC4aMAAABHZma6/PLLdfbZZ2vUqFGHbU9PT1etWn+uDJCQkKD09PSg643onkIcbtOO3bp6+GR9tjJd93VtqYe6n6+oKL6kAwAAAJIK7tE72g/cs2cHVfXnn3+uGjVqaOPGjWrXrp0aNmyoiy666MD2/GZyD8UP7hHdU4iDfffrJp3/4Nv6evVv+t8tHfRIzwtICFGiGD4KAABwZDVq1JAkValSRUlJSfryyy8P2p6QkKC1a9ceeJ+WlnbgmGCQFJYSk7/5SRc9Ml57M7M164Fe6n1+Q69DQinEkhQAACDsVa1atPJC2rVrl3bu3Hng9fTp09W4ceOD9uncubPGjBkj55wWLFig+Pj4oJ8nlBg+GvGcc3rqo690/zufqXndqvrgb51Vs2IFr8NCKcYzhQAAIKz99luxnHbDhg1KSkqSJGVlZalv37668sorNXLkSEnSoEGD1KFDB6WkpCgxMVFxcXEaPXp0SOomKYxgGfuy9NdXP9Hbn6/Q1a0a6JWBlyuubIzXYaEUczp8HDwAAACkevXqaenSpYeVDxo06MBrM9OIESNCXndEJ4WleUmK37btUvenk/XlT7/pkR7n6x9dW9JDA8/xTCEAAID/RPQzhaV1SYpFP2/QeQ+M1bK0zXp3aCfdl9SKhBC+wbUIAADgLxHdU1gafbDwR1038mNVrnC85jzYW83qVPE6JOAAho8CAAD4D0lhhHDO6fEJC/TIB/N13qk19N7fOqlqfDmvwwIOwvBRAAAA/yEpjAC792bqxpen6b2FP+qaC8/QSzdcprIx/KeFP5EUAgAA+EtEP1NYGqRt2am2j76j97/8Uf/te5Fe++sVJITwLYaPAgCASLE8bbOa3v2mlqdt9jqUoJEUhrGFq9brvAff1qrftmninV11R8cWTOIBX3NyXKMAACDs7crIVOcnJ2jFui3q8n8TtSsjM+hzrly5Us2aNTvwd8IJJ2j48OEH7eOc05AhQ5SYmKgmTZpo0aJFQdcrkRSGrbc/X6FL//Wu4o4ro3mP9FaHs+p5HRJQKAwfBQAA4e7GUdO0ccduOSdt2L5LN70yLehzNmjQQEuWLNGSJUv0zTffKC4u7sBi9vtNnTpVqampSk1N1ahRozR48OCg65V4pjDs5OQ4Pfje5/rvpC910ekJeuf2Tqpc4XivwwIKxTmGjwIAAH+743+faumaTUfc/tvWP7Rqw3blBL7XZGRm64OFqWr0y+uqdlL5fI9pWvtkPX3NxYWOYebMmapfv75q1659UHlycrL69+8vM1OrVq20bds2rV+/XtWrVy/0ufMT0T2FZtbJzEZt377d61BCYueefeo5fJL+O+lL3XjJmZp6b3cSQoQVho8CAIBw9/OmHQcSwv1ynNPPm3aErI7x48erT58+h5Wnp6erVq1aB94nJCQoPT096PoiuqfQOTdZ0uQWLVrc5HUswfpl03YlDUvWivQtGj7gYt3crhlfrhGWGD4KAAD8rKAevTdmL9PtY2Zp996sA2VxZcvouQGXaECbxkHXv2/fPk2aNElPPPHEYdvyG3UVipwgonsKI8VnK9N13gNlFjfWAAAQNUlEQVRvK23LTn10dzfdcvlZJIQISwwfBQAA4e7ato3VoVk9xcZES5JiY6LV8ax6IUkIpdznBps3b66qVaseti0hIUFr16498D4tLU01atQIuk6SQp8bPXuZLn/8PVUsH6vPH+2ry86sXfBBgE85sXg9AAAIf68OvEJVToiTSaoaX06v3HRFyM49bty4fIeOSlLnzp01ZswYOee0YMECxcfHB/08oRThw0fDWXZOju59e56GT/1G7c6srbG3ddRJ5WK9DgsIGr3cAAAg3JWLjdGku5PU97kpentIR5WLjQnJeXfv3q1PPvlEL7/88oGykSNHSpIGDRqkDh06KCUlRYmJiYqLi9Po0aNDUi9JoQ9t371X/Z6fomnf/qJbrzhL/9evjcpE06mL8Mfi9QAAIFI0SqispU8OCOk54+LitGXLloPKBg0adOC1mWnEiBEhrVMiKfSdVb9tVdKwZK3asE0v3nCZbrqkidchASHjHMNHAQAA/Iak0Ec+Xf6rej07WVFmmvaP7rro9FoFHwSEEZakAAAA8B+SQp8YOWOphr45Sw1rVNSHd3ZRvSoneh0SAAAAgFKAB9U8sjxts5re/aaW/rJRt42eqdtGz9SVTetq7kO9SQgRsRg+CgAA4D/0FHpgV0amOj85QWt/36nWD49TRma27rqqhf7Vq7Wio8jTEblYkgIAAMB/yEA8cOOoadqwfbeckzIys3Vu/Wp6os9FJIQoFXimEAAARILlvy9X4/caa/nvy70OJWhkISXsjdnLlLJktfZmZR8oW5a2WW/MXuZhVEDJYEkKAAAQCXZl7lKHjzvo+63fq+PHHbUrc1dIzvvss8+qcePGatSokYYPH37YduechgwZosTERDVp0kSLFi0KSb0RnRSaWSczG7V9+3avQzng/nc+0+69WQeV7d6bpfvf+cyjiICSwzOFAAAgElw/53pt3LNRTk4b9mzQDXNuCPqcy5Yt0yuvvKIvv/xSS5cu1UcffaTU1NSD9pk6dapSU1OVmpqqUaNGafDgwUHXK0X4M4XOucmSJrdo0eImr2PZ7/FerXX7mFkHJYZxZcvo371bexgVUHIYPgoAAPxs6BdDtWTLkiNuX79rvVbtXKUclyNJysjO0Hur39Pi8YtVvVz1fI9pVqmZhp9/eM9fXitWrFCrVq0UFxcnSWrTpo0mTJigu++++8A+ycnJ6t+/v8xMrVq10rZt27R+/XpVr55/vYUV0T2FfnRt28bq0KyeYmOiJUmxMdHqeFY9DWjT2OPIgOLH8FEAABDufv7j5wMJ4X45ytHPf/wc1HkbN26suXPnasuWLdq9e7dSUlK0du3ag/ZJT09XrVp/rmWekJCg9PT0oOqVIryn0K9eHXiFmtz9htZu2amq8eX0yk1XeB0SUGIYPgoAAPysoB6911e+riGfD9GurD+fI4wrE6cXLnhB1zW47pjrPf3003XPPfeoXbt2Kl++vJo2baoyZQ5O15w7/Af2UIzCoqfQA+ViYzTp7iSdXrOSkv/eVeViY7wOCSgxJIUAACCcXd/genU8paNio2MlSbHRsep0SqegEsL9brjhBi1atEhz585VxYoVdeqppx60PSEh4aDew7S0NNWoUSPoekkKPdIoobKWPjlAjRIqex0KUCLy+2ULAAAgHL3e5nVVOb6KTKaqx1fVa21eC8l5N27cKEn69ddf9eGHH6pPnz4Hbe/cubPGjBkj55wWLFig+Pj4oJ8nlBg+CqCE7H+ekIlmAABAuCsXU04pV6ao18xeeufSd1QuplxIztu9e3dt2bJFMTExGjFihE466SSNHDlSkjRo0CB16NBBKSkpSkxMVFxcnEaPHh2SekkKAZQoho8CAIBI0KhiIy3rGdq1xufNm3dY2aBBgw68NjONGDEipHVKDB8FUEIYPgoAAOBPJIUASgTDRwEAAPyJpBBAiWL4KAAA8KNIGdV0LJ+DpBBAiWDhegAA4FexsbHasmVL2CeGzjlt2bJFsbGxRTqOiWYAlIj9N1l6CgEAgN8kJCQoLS1NmzZt8jqUoMXGxiohIaFIx5AUAihRPFMIAAD8JiYmRnXr1vU6DM8wfBRAiWD4KAAAgD+RFAIoEQwfBQAA8CeSQgAliuGjAAAA/mLhPsNOYZjZJklrvI4jH5UlbfY6iFKKtvcObe8d2t47tL23aH/v0Pbeoe2949e2r+2cOzm/DaUiKfQrM/vaOdfC6zhKI9reO7S9d2h779D23qL9vUPbe4e29044tj3DRwEAAACgFCMpBAAAAIBSjKTQW6O8DqAUo+29Q9t7h7b3Dm3vLdrfO7S9d2h774Rd2/NMIQAAAACUYvQUAgAAAEApRlIIAAAAAKUYSWGImNmVZrbSzFaZ2b35bO9nZt8G/r4ws6Z5tv1iZt+Z2RIz+zpPeUUz+8TMUgP/e1JJfZ5wcqxtb2YNAm2+/2+HmQ0NbHvYzNLzbOtQ0p8rHBSi7bsE2n2JmX1tZq0LOpbrvvCOtf3NrJaZfWpmK8xsuZndnucYrv1CCPLa554fhCCue+75QSqo7fPsd46ZZZtZj4KO5bovnGNte+73wQvyug+f+71zjr8g/yRFS/pJUj1Jx0laKumMQ/Y5X9JJgdftJS3Ms+0XSZXzOe+Tku4NvL5X0n+9/qx++wu27Q85z2/KXdRTkh6WdJfXn8/Pf4Vs+/L689nlJpJ+KOhYrvsSaf/qkpoHXleQ9GOe9ufaL8a2D7znnu9R2x9yHu75IW77PPvNkpQiqUdBx3LdF3vbc7/3qO0D5WFzv6enMDTOlbTKObfaObdP0nhJXfLu4Jz7wjm3NfB2gaSEQpy3i6Q3A6/flNQ1RPFGklC1/aWSfnLOrSnWaCNLYdr+Dxe440kqJ8kV4liu+8I55vZ3zq13zi0KvN4paYWkmiUWefgL5to/Gq79goWq7bnnF12BbR9wm6QPJG0s5LFc9wU75rbnfh+0YK77o/HddU9SGBo1Ja3N8z5NR/8Hd4OkqXneO0nTzewbMxuYp7yqc269lPuPWlKVEMUbSYJt+/16Sxp3SNmtgSFIr/uiW99/CtX2ZpZkZj9ImiLp+kIcy3VfOMG0f97tdSSdJWlhnmKu/aMLtu255x+7kFz34p5/LApsezOrKSlJ0sgiHMt1X7Bg2j7vPnXE/b6ogm37sLnfkxSGhuVTlu+vwmZ2sXITk3vyFF/gnGuu3KGNt5jZRaEPMWIF2/Yys+MkdZb0Xp7ilyTVl9RM0npJw0IRbIQpVNs75yY45xoq91ewx4pyLI4qmPbPPYFZeeX+sjnUObcjUMy1X7Bg2557/rELxXXPPf/YFKbth0u6xzmXfQzH4siCafvcE3C/P1bBtn3Y3O9JCkMjTVKtPO8TJK07dCczayLpVUldnHNb9pc759YF/nejpAnK7aqWpA1mVj1wbHUVvku6NAmq7QPaS1rknNuwv8A5t8E5l+2cy5H0iv78b4I/Fart93POzZVU38wqF3As133hBNP+MrMY5X5BGOuc+zDPflz7BQuq7bnnByWotg/gnn9sCtP2LSSNN7NfJPWQ9KKZdS3gWK77ggXT9tzvgxNU24fT/Z6kMDS+knSqmdUN/ALZW9KkvDuY2SmSPpR0jXPuxzzl5cyswv7Xki6XtCyweZKkAYHXAyQlF+unCE/H3PZ59NEhw4j2/0MNSNKf/03wp8K0faKZWeB1c+U+pL2lgGO57gvnmNs/UPaapBXOuacPOYZrv2DBtD33/OAEc9/Zj3v+sSmw7Z1zdZ1zdZxzdSS9L+lm59zEAo7lui/YMbc99/ugBdP2YXW/L+N1AJHAOZdlZrdKmqbc2Yded84tN7NBge0jJT0oqZJyfz2QpCznXAtJVSVNCJSVkfS2c+7jwKn/I+ldM7tB0q+SepbgxwoLQba9zCxOUjtJfz3k1E+aWTPlDhH4JZ/tpV4h2767pP5mlilpj6RegQkg8j02cGqu+0IIpv0td4r+ayR9Z2ZLAqe8zzmXIq79AgXZ9tzzgxDkfYd7fhAK2fZFOjawmeu+AMG0vaQLxP3+mAXZ9mF1v98/ZTMAAAAAoBRi+CgAAAAAlGIkhQAAAABQipEUAgAAAEApRlIIAAAAAKUYSSEAAAAAlGIkhQCAiGBmVc3sbTNbbWbfmNl8M0s6wr41zOz9I2ybbWYtCllnVzM7I5i4DzlfdTP7KPC6rZltN7MlZvatmc0wsyoFHP+wmd1VyLpmmNlJoYgbABDeSAoBAGEvsEDzRElznXP1nHNnK3eR4YR89i3jnFvnnOsRgqq7SgpZUijpDkmv5Hk/zznXzDnXRLmLKN8SbAWWK0rS/yTdHOz5AADhj6QQABAJLpG0L+9Cws65Nc655yXJzK41s/fMbLKk6WZWx8yWBbYdb2bjA71x70g6Pr8KzOw/ZvZ9YL+nzOx8SZ0l/V+gN69+4O/jQE/lPDNrGDj2DTMbGSj70cyuOsLn6C7p40MLA0lvBUlbA+8rmtnEQCwLzKxJnt3PCPR2rjazIYH965jZCjN7UdIiSbUkTZLUp7ANDACIXGW8DgAAgBBopNxk52jOk9TEOfe7mdXJUz5Y0m7nXJNAcnXYecysoqQkSQ2dc87MTnTObTOzSZI+cs69H9hvpqRBzrlUM2sp6UXlJqySVEdSG0n1JX1qZonOuYw8ddSVtNU5tzdP1Rea2RJJlSTtknRfoPwRSYudc13N7BJJYyQ1C2xrKOli5SaRK83spUB5A0nXOecO9A6aWVkzq+Sc21JA2wEAIhg9hQCAiGNmI8xsqZl9laf4E+fc7/nsfpGktyTJOfetpG/z2WeHpAxJr5pZN0m786mzvKTzJb0XSORellQ9zy7vOudynHOpklYrN3nLq7qkTYeU7R8+WkvSaElPBspbK3f4p5xzsyRVMrP4wLYpzrm9zrnNkjZKqhooX+OcW3DI+TdKqpHP5wUAlCL0FAIAIsFy5Q69lCQ5524xs8qSvs6zz66jHO+OdnLnXJaZnSvpUuU+q3ir/uwB3C9K0jbnXLNDjz9CHYe+3yMp9ihhTJL0QeC1HeX8eXsas/Xn/9fn9/ljA/UCAEoxegoBAJFglqRYMxucpyyukMfOldRPksyssaQmh+4Q6AWMd86lSBqqP4dq7lTuME0553ZI+tnMegaOMTNrmuc0Pc0syszqS6onaeUh1fyo3CGmR9Ja0k/5xNxW0uZA/YUWeE6xmqRfinIcACDy0FMIAAh7gef8ukp6xszuVu4wzF2S7inE4S9JGm1m30paIunLfPapICnZzGKV20v3t0D5eEmvBCZ06aHcRO0lM/unpJjA9qWBfVdKmqPc4ZyD8j5PGPgMu8zsp8CzhqsCxfufKTRJ2yXdGCh/OE/MuyUNKMTnPNTZkhY457KO4VgAQAQx5446YgYAAATJzN5QnglpjrJfkqSznXP/LIGYnpU0yTk3s7jrAgD4Gz2FAAD4hHNugplVKqHqlpEQAgAkegoBAAAAoFRjohkAAAAAKMVICgEAAACgFCMpBAAAAIBSjKQQAAAAAEoxkkIAAAAAKMX+HxZGhTOqDlTmAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 1080x504 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(15,7))\n",
    "# Plot with matplotlib\n",
    "for i,crmult in enumerate(Crmult):\n",
    "    im = i%len(colors)\n",
    "    ener = [ HtoeV*(l.energy-emin_HF) for l in log_HF[crmult] ]\n",
    "    plt.plot(Hgrids, ener, marker=markers[im], \n",
    "             ls='-', label=str(crmult),color=colors[im])  \n",
    "\n",
    "plt.yscale('log')\n",
    "plt.xlabel('Grid step (Bohr)')\n",
    "plt.ylabel('Total energy $\\Delta E$ (eV)')\n",
    "plt.title('Dissociation energy of the N2 dimer for different crmult')\n",
    "plt.legend(loc=4)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In order to converge the result, you need to decrease the grid step and also increase the extension of the mesh. For a given *crmult*, the curve are almost flat. For a *hgrid* value of of 0.35, there is few difference between the values with *crmult*=5.0 and 7.0 but for a *hgrid* value of 0.20 it is important. Now we give the HOMO-1 and HOMO eigenvalues both for LDA and HF functionals"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "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>LDA</th>\n",
       "      <th>HF</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>10.41</td>\n",
       "      <td>16.82</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>11.88</td>\n",
       "      <td>16.82</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>11.88</td>\n",
       "      <td>17.37</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>13.42</td>\n",
       "      <td>21.18</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>28.35</td>\n",
       "      <td>40.21</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "    LDA    HF\n",
       "0 10.41 16.82\n",
       "1 11.88 16.82\n",
       "2 11.88 17.37\n",
       "3 13.42 21.18\n",
       "4 28.35 40.21"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "ih = Hgrids.index(0.3)\n",
    "lda_evals = log_LDA[7.0][ih].evals[0][0]\n",
    "hf_evals = log_HF[7.0][ih].evals[0][0]\n",
    "    \n",
    "df = DataFrame(columns=[\"LDA\", \"HF\"])\n",
    "for idx, (l,h) in enumerate(zip(reversed(lda_evals),reversed(hf_evals))):\n",
    "    df.loc[idx] = [-HtoeV*l, -HtoeV*h]\n",
    "    \n",
    "display(df)"
   ]
  }
 ],
 "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.7.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}