{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "a096f8df",
   "metadata": {},
   "source": [
    "# How to handle errors during optimization\n",
    "\n",
    "## Try to avoid errors\n",
    "\n",
    "Often, optimizers try quite extreme parameter vectors, which then can raise errors in your criterion function or derivative. Often, there are simple tricks to make your code more robust. Avoiding errors is always better than dealing with errors after they occur.  \n",
    "\n",
    "- Avoid to take ``np.exp`` without further safeguards. With 64 bit floating point numbers, the exponential function is only well defined roughly between -700 and 700. Below it is 0, above it is inf. Sometimes you can use ``scipy.special.logsumexp`` to avoid unsafe evaluations of the exponential. Read [this](https://en.wikipedia.org/wiki/LogSumExp) for background information on the logsumexp trick.\n",
    "- Set bounds for your parameters that prevent extreme parameter constellations.\n",
    "- Use the ``bounds_distance`` option with a not too small value for ``covariance`` and ``sdcorr`` constraints.\n",
    "- Use `estimagic.utilities.robust_cholesky` instead of normal\n",
    "  cholesky decompositions or try to avoid cholesky decompositions.\n",
    "- Use a less aggressive optimizer. Trust region optimizers like `fides` usually choose less extreme steps in the beginnig than line search optimizers like `scipy_bfgs` and `scip_lbfgsb`. \n",
    "\n",
    "## Do not use clipping\n",
    "\n",
    "A commonly chosen solution to numerical problems is clipping of extreme values. Naive clipping leads to flat areas in your criterion function and can cause spurious convergence. Only use clipping if you know that your optimizer can deal with flat parts. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4c551530",
   "metadata": {},
   "source": [
    "## Let estimagic do its magic\n",
    "\n",
    "Instead of avoiding errors in your criterion function, you can raise them and let estimagic deal with them. If you are using numerical derivatives, errors will automatically be raised if any entry in the derivative is not finite. \n",
    "\n",
    "### An example\n",
    "\n",
    "Let's look at a simple example from the Moré-Wild benchmark set that has a numerical instability. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "5ec31d93",
   "metadata": {},
   "outputs": [],
   "source": [
    "import warnings\n",
    "\n",
    "import estimagic as em\n",
    "import numpy as np\n",
    "from scipy.optimize import minimize as scipy_minimize\n",
    "\n",
    "warnings.simplefilter(\"ignore\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "fec56a0b",
   "metadata": {},
   "outputs": [],
   "source": [
    "def jennrich_sampson(x):\n",
    "    dim_out = 10\n",
    "    fvec = (\n",
    "        2 * (1.0 + np.arange(1, dim_out + 1))\n",
    "        - np.exp(np.arange(1, dim_out + 1) * x[0])\n",
    "        - np.exp(np.arange(1, dim_out + 1) * x[1])\n",
    "    )\n",
    "    return fvec @ fvec\n",
    "\n",
    "\n",
    "correct_params = np.array([0.2578252135686162, 0.2578252135686162])\n",
    "correct_criterion = 124.3621823556148\n",
    "\n",
    "start_x = np.array([0.3, 0.4])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "13c144d7",
   "metadata": {},
   "source": [
    "### What would scipy do?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "21383146",
   "metadata": {},
   "outputs": [],
   "source": [
    "scipy_res = scipy_minimize(jennrich_sampson, x0=start_x, method=\"L-BFGS-B\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "36d8e926",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "scipy_res.success"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "40511eb9",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([0.2578, 0.2578]), array([0.3384, 0.008 ]))"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "correct_params.round(4), scipy_res.x.round(4)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ca245e3b",
   "metadata": {},
   "source": [
    "So, scipy thinks it solved the problem successfully but the result is far off. (Note that scipy would have given us a warning, but we disabled warnings in order to not clutter the output).\n",
    "\n",
    "### Estimagic's error handling magic"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "617108b1",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([0.25782521, 0.25782521]), array([0.25782521, 0.25782522]))"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "res = em.minimize(\n",
    "    criterion=jennrich_sampson,\n",
    "    params=start_x,\n",
    "    algorithm=\"scipy_lbfgsb\",\n",
    "    error_handling=\"continue\",\n",
    ")\n",
    "\n",
    "correct_params, res.params"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7fba61dd",
   "metadata": {},
   "source": [
    "### How does the magic work\n",
    "\n",
    "When an error occurs and `error_handling` is set to `\"continue\"`, estimagic replaces your criterion with a dummy function (and adjusts the derivative accordingly). \n",
    "\n",
    "The dummy function has two important properties:\n",
    "\n",
    "1. Its value is always higher than criterion at start params. \n",
    "2. Its slope guides the optimizer back towards the start parameters. I.e., if you are minimizing, the direction of strongest decrease is towards the start parameters; if you are maximizing, the direction of strongest increase is towards the start parameters. \n",
    "\n",
    "Therefore, when hitting an undefined area, an optimizer can take a few steps back until it is in better territory and then continue its work. \n",
    "\n",
    "Importantly, the optimizer will not simply go back to a previously evaluated point (which would just lead to cyclical behavior). It will just go back in the direction it originally came from.\n",
    "\n",
    "In the concrete example, the dummy function would look similar to the following:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "dbf49b7b",
   "metadata": {},
   "outputs": [],
   "source": [
    "def dummy(params):\n",
    "    start_params = np.array([0.3, 0.4])\n",
    "    # this is close to the actual value used by estimagic\n",
    "    constant = 8000\n",
    "    # the actual slope used by estimagic would be even smaller\n",
    "    slope = 10_000\n",
    "    diff = params - start_params\n",
    "    return constant + slope * np.linalg.norm(diff)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5958751d",
   "metadata": {},
   "source": [
    "Now, let's plot the two functions. For better illustration, we assume that the jennrich_sampson function is only defined until it reaches a value of 100_000 and the dummy function takes over from there.  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "061ba6c5",
   "metadata": {},
   "outputs": [],
   "source": [
    "from plotly import graph_objects as go\n",
    "\n",
    "grid = np.linspace(0, 1)\n",
    "params = [np.full(2, val) for val in grid]\n",
    "values = np.array([jennrich_sampson(p) for p in params])\n",
    "values = np.where(values <= 1e5, values, np.nan)\n",
    "dummy_values = np.array([dummy(p) for p in params])\n",
    "dummy_values = np.where(np.isfinite(values), np.nan, dummy_values)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "2556c2fb",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAArwAAAH0CAYAAADfWf7fAAAgAElEQVR4XuydB5hU1f2Gv3tndwFpYotRY4u9JWLUv0aNsUSMDUuMGgshGrG3YAODBtTYW1QExR4UDImJvSVRE1s0sUWjib0D0payZWb+z53dRcDClO/euTv7zvPwADPnfOfc93fQl8OZe4N8Pp8XLwhAAAIQgAAEIAABCNQogQDhrdHKclkQgAAEIAABCEAAAgUCCC8LAQIQgAAEIAABCECgpgkgvDVdXi4OAhCAAAQgAAEIQADhZQ1AAAIQgAAEIAABCNQ0AYS3psvLxUEAAhCAAAQgAAEIILysAQhAAAIQgAAEIACBmiaA8NZ0ebk4CEAAAhCAAAQgAAGElzUAAQhAAAIQgAAEIFDTBBDemi4vFwcBCEAAAhCAAAQggPCyBiAAAQhAAAIQgAAEapoAwlvT5eXiIAABCEAAAhCAAAQQXtYABCAAAQhAAAIQgEBNE0B4a7q8XBwEIAABCEAAAhCAAMLLGoAABCAAAQhAAAIQqGkCCG9Nl5eLgwAEIAABCEAAAhBAeFkDEIAABCAAAQhAAAI1TQDhrenycnEQgAAEIAABCEAAAggvawACEIAABCAAAQhAoKYJILw1XV4uDgIQgAAEIAABCEAA4WUNQAACEIAABCAAAQjUNAGEt6bLy8VBAAIQgAAEIAABCCC8rAEIQAACEIAABCAAgZomgPDWdHm5OAhAAAIQgAAEIAABhJc1AAEIQAACEIAABCBQ0wQQ3pouLxcHAQhAAAIQgAAEIIDwsgYgAAEIQAACEIAABGqaAMJb0+Xl4iAAAQhAAAIQgAAEEF7WAAQgAAEIQAACEIBATRNAeGu6vFwcBCAAAQhAAAIQgADCyxqAAAQgAAEIQAACEKhpAghvTZeXi4MABCAAAQhAAAIQQHhZAxCAAAQgAAEIQAACNU0A4a3p8nJxEIAABCAAAQhAAAIIL2sAAhCAAAQgAAEIQKCmCSC8NV1eLg4CEIAABCAAAQhAAOFlDUAAAhCAAAQgAAEI1DQBhLemy8vFQQACEIAABCAAAQggvKwBCEAAAhCAAAQgAIGaJoDw1nR5uTgIQAACEIAABCAAAYSXNQABCEAAAhCAAAQgUNMEEN6aLi8XBwEIQAACEIAABCCA8LIGIAABCEAAAhCAAARqmgDCW9Pl5eIgAAEIQAACEIAABBBe1gAEIAABCEAAAhCAQE0TQHhrurxcHAQgAAEIQAACEIAAwssagAAEIAABCEAAAhCoaQIIb02Xl4uDAAQgAAEIQAACEEB4WQMQgAAEIAABCEAAAjVNAOGt6fJycRCAAAQgAAEIQAACCC9rAAIQgAAEIAABCECgpgkgvDVdXi4OAhCAAAQgAAEIQADhZQ1AAAIQgAAEIAABCNQ0AYS3psvLxUEAAhCAAAQgAAEIILysAQhAAAIQgAAEIACBmiaA8NZ0ebk4CEAAAhCAAAQgAAGElzUAAQhAAAIQgAAEIFDTBBDemi4vFwcBCEAAAhCAAAQggPCyBiAAAQhAAAIQgAAEapoAwlvT5eXiIAABCEAAAhCAAAQQXtYABCAAAQhAAAIQgEBNE0B4a7q8XBwEIAABCEAAAhCAAMLLGoAABCAAAQhAAAIQqGkCCG9Nl5eLgwAEIAABCEAAAhBAeFkDEIAABCAAAQhAAAI1TQDhrenycnEQgAAEIAABCEAAAggvawACEIAABCAAAQhAoKYJILw1XV4uDgIQgAAEIAABCEAA4WUNQAACEIAABCAAAQjUNAGEt6bLy8VBAAIQgAAEIAABCCC8rAEIQAACEIAABCAAgZomgPDWdHm5OAhAAAIQgAAEIAABhJc1AAEIQAACEIAABCBQ0wQQ3pouLxcHAQhAAAIQgAAEIIDwsgYgAAEIQAACEIAABGqaAMJb0+Xl4iAAAQhAAAIQgAAEEF7WAAQgAAEIQAACEIBATRNAeGu6vFwcBCAAAQhAAAIQgADCyxqAAAQgAAEIQAACEKhpAghvTZeXi4MABCAAAQhAAAIQQHhZAxCAAAQgAAEIQAACNU0A4a3p8nJxEIAABCAAAQhAAAIIL2sAAhCAAAQgAAEIQKCmCSC8NV1eLg4CEIAABCAAAQhAAOFlDUAAAhCAAAQgAAEI1DQBhLemy8vFQQACEIAABCAAAQggvKwBCEAAAhCAAAQgAIGaJoDw1nR5uTgIQAACEIAABCAAAYSXNQABCEAAAhCAAAQgUNMEEN6aLi8XBwEIQAACEIAABCCA8LIGIAABCEAAAhCAAARqmgDCW9Pl5eIgAAEIQAACEIAABBBe1gAEIAABCEAAAhCAQE0TQHhrurxcHAQgAAEIQAACEIAAwssagAAEIAABCEAAAhCoaQIIb02Xl4uDAAQgAAEIQAACEEB4WQMQgAAEIAABCEAAAjVNAOGt6fJycRCAAAQgAAEIQAACCC9rAAIQgAAEIAABCECgpgkgvDVdXi4OAhCAAAQgAAEIQADhNayBD6bONaQQ4SSwVO8GzZnXqnktOWcsWQYCyy/VQ59Mm6tc3hBGhI1AGAZarm83fTRtni2TIA+B7g0ZLdEto09nNXsCSbESWGHpHtY8wuIhgPAauCK8BojmCITXDNQYh/AaYRqjEF4jTHMUwmsGao5DeM1AY4pDeA1gEV4DRHMEwmsGaoxDeI0wjVEIrxGmOQrhNQM1xyG8ZqAxxSG8BrAIrwGiOQLhNQM1xiG8RpjGKITXCNMchfCagZrjEF4z0JjiEF4DWITXANEcgfCagRrjEF4jTGMUwmuEaY5CeM1AzXEIrxloTHEIrwEswmuAaI5AeM1AjXEIrxGmMQrhNcI0RyG8ZqDmOITXDDSmOITXABbhNUA0RyC8ZqDGOITXCNMYhfAaYZqjEF4zUHMcwmsGGlMcwmsAi/AaIJojEF4zUGMcwmuEaYxCeI0wzVEIrxmoOQ7hNQONKQ7hNYBFeA0QzREIrxmoMQ7hNcI0RiG8RpjmKITXDNQch/CagcYUh/AawCK8BojmCITXDNQYh/AaYRqjEF4jTHMUwmsGao5DeM1AY4pDeA1gEV4DRHMEwmsGaoxDeI0wjVEIrxGmOQrhNQM1x9WS8L75zofK5XL65qormikVF/f7ex/TFeMmqbm5VTtsvYnOOOFgZTJhcZ0X0wrhNWBEeA0QzREIrxmoMQ7hNcI0RiG8RpjmKITXDNQcV0vCe+Ho27XhOqtrp203NVNafNxb736kwSeep1uuGKZll+mnU0Zdo43WW12D9h2w+M5FtEB4i4C0uCYI7+IIJf85wps882JHRHiLJZVsO4Q3Wd6ljIbwlkIr+ba1Irx3P/ykzr7sZvXuuYT69O6pC84Yogcf/Yd69eyhvz/zkp765yu6ffQI9evbW6eeM0bRbnBra1b7DdxOh/1k1wL46L1fXXKj3nj7QzU01OvaC4dqlZW+pkeffF4XXzNR85qatdnG62j4cQcVPl/wdd34ezSrcY6OP2yfwtuv/vcdnXH+OE0cc6alqAivASPCa4BojkB4zUCNcQivEaYxCuE1wjRHIbxmoOa4SoT3iWdy+mRK3jyjxcdtsWlGyy3z+XbRrup2W/Wfv8M79ta7dMvvHtSoUw7V1ptvWOiQy+X1r5dfV/8N19KUT2do14NP0z23nKe+vXtq4E+H6ejBe2qnbTfT9BmNBXH+eMo0HXjUKN3ym2FafrmldPq51+qbq66gQw/YZaEJ/PKCcYXMgQO2Krzf1NyiLXc7Ss/eP2bxF1REC4S3CEiLa4LwLo5Q8p8jvMkzL3ZEhLdYUsm2Q3iT5V3KaAhvKbSSb1uJ8F5ydatefjV54T1+SJ02WDf4HKwvEt7X33hP558x5EvB/uSoURp23IEKgkCnnztWvx83aqG2N028Xx98PFWnHn1A4f2n//lq4ZzuzVecvlC7k0eOLsj2gO9vNv/99bcdpJf+fH0hu9IXwlspQUkIrwGiOQLhNQM1xiG8RpjGKITXCNMchfCagZrjKhHezrDD29TUUti17Xi9+8EnuubmP+njydMUeeiLr7yhsRcN1dRPZ2r8Hx7S6PNOWojwBVffpjvv+5v69ulZeD/6Uly0G3zb6BELtRtx4fWF88P77Pq9wvuNs+dqmz2P1XMPjLVUDOE1YER4DRDNEQivGagxrjMK72uvB7plfEZrrZnTgfvnjDTSE4XwpqcWi84E4U1vbaKZVSK8abuyU86+Rtt9d+EjDc3NLTrqp58J776Hn6mD99lJu+64RWH6P/r5mRpx0iGFX//y/HGadN3IhS7rxon3a8rUGTppyL5febk3TLhPk6dM19Aj9yu0e/HVNzXigs/nlcsM4S2X3AL9EF4DRHMEwmsGaozrjML77HOB7rwro+/0z2n3XRFe43IgqggCCG8RkKrYpJaE99wrbi18KW3IwbsXiEZneBcV3h32PVGXnHW0Nlx3dT37wmv62UnnF87nrrPGytr1oNM09Ij9tP3W/Qs7tNFfpGfMnK2Djj1Ho399otZYbcXC2dzJU6drpa8vu1DV3v9oig459hzddPnp7XdpGK01V19JRxy8h6W6CK8BI8JrgGiOQHjNQI1xnVF4//LXUI/8NdT3v5cr/KjFFzu86a0qwpve2kQzqyXhje6ycOzwy5XN5XTxmUfpsade+JzwRndzuGTMRPVcors2Xn9NNbe06IC9dtAGa6+m//zvXY285Ca98/7H6tG9m645/ySt+o3l9cQ/XlZ0y7NIdBvq63TkoIHa64fbfK6w9zz8lC4afbvmNjVp68030sihgz93N4dyVwPCWy45dngN5OKLQHjjY1tpcmcU3j/dFeqZ50LtvktW39kk+S+YVMq8mP4IbzGUqtMG4a0O92JHrSXhLfaaO2M7hNdQNXZ4DRDNEQivGagxrjMK7y3jQ732eqgD989qrTURXuNyIKoIAghvEZCq2AThrSL8EoZGeEuA9WVNEV4DRHMEwmsGaozrjMJ79dg6ffihdMTPW/X15Y0wUhTFDm+KirHIVBDe9NYmmhnCm+76dMwO4TXUCeE1QDRHILxmoMa4zii8519cp8ZG6eQTW9WrlxFGiqIQ3hQVA+FNbzG+YGYIb+coF8JrqBPCa4BojkB4zUCNcZ1NePN5acTIOoWhNGJYa+G+k7X4QnjTW1V2eNNbG3Z4012bBWeH8BpqhfAaIJojEF4zUGNcZxPeGTMDXXRpRn365PWL47NGEumKQnjTVY8FZ4Pwprc2CG+6a4PwmuuD8JqBGuIQXgPEmCI6m/C+936gMddltNKKef38ZwhvTMuC2K8ggPCme3lwpCHd9emYHTu8hjohvAaI5giE1wzUGNfZhPeVV0ONnxBqvXVy2m/f2rwHb1RedniNi9wchfCagZrjEF4z0JjiEF4DWITXANEcgfCagRrjOpvwPvVMqLvvDbX5pjntsjPCa1wKRBVJAOEtElSVmiG8VQJf4rAIb4nAvqg5wmuAaI5AeM1AjXGdTXgfeiTUo4+H2mG7nLbZCuE1LgWiiiSA8BYJqkrNEN4qgS9xWIS3RGAIrwFYAhEIbwKQyxyiswnvpDsz+tfzgfbcI6uNv1WbD52ISsmRhjIXdALdEN4EIFcwRC0Jb/Ro4Vwup2+uumIFRMrvms/nNfrmP+qmCffribuuKj/oC3oivAac7PAaIJojEF4zUGNcZxPeG2/J6H9vBDrkwKy+uTrCa1wKRBVJAOEtElSVmtWS8F44+nZtuM7q2mnbTROn2dKa1dBfXa3llllSdz30hP7+xyutc0B4DTgRXgNEcwTCawZqjOtswvubq+v0yWTp6CNatdyyRhApi2KHN2UFWWA6CG96axPNrFaE9+6Hn9TZl92s3j2XUJ/ePXXBGUP04KP/UK+ePfT3Z17SU/98RbePHqF+fXvr1HPGKNoNbm3Nar+B2+mwn+xaKFL03q8uuVFvvP2hGhrqde2FQ7XKSl/To08+r4uvmah5Tc3abON1NPy4gwqfL/p67KkX9d1NN9AWux2pp+6+2lp4hNeAE+E1QDRHILxmoMa4zia8555fp7nzpNNOblWP7kYQKYtCeFNWEIQ3vQVZZGaVCG/zo/cr99F7iV9rwzYDFC7/+WMLp4y6Rttt1X/+Du/YW+/SLb97UKNOOVRbb75hYZ65XF7/evl19d9wLU35dIZ2Pfg03XPLeerbu6cG/nSYjh68p3badjNNn9FYEOePp0zTgUeN0i2/Gabll1tKp597rb656go69IBdvvC6W7NZfXf3oxHexFdFEQMivEVASrgJwpsw8BKG60zCm81JZ41qe8ramcNbS7jKztcU4U1vzdjhTW9toplVIryNZ5+o1uefTvwCe51+keq+vfnnxv0i4X39jfd0/hlDvnSOPzlqlIYdd6CCINDp547V78eNWqjtTRPv1wcfT9WpRx9QeP/pf76qK8ZN0s1XnI7wJl75CgdEeCsEGEN3hDcGqKbIziS8Uz8NdNlvMlp6qbyOO7p2HzoRlRbhNS3wGGIQ3higGiMrEd7OsMPb1NRS2LXteL37wSe65uY/6ePJ0wqPWn/xlTc09qKhmvrpTI3/w0Mafd5JC9G94OrbdOd9f1PfPj0L70dfiot2g28bPQLhNa7DRKIQ3kQwlzQIwlsSrkQbdybhfevtQONuzGjVVfIafAjCm+hCYbD5BBDedC+GSoQ3bVd2ytnXaLvvLnykobm5RUf99DPh3ffwM3XwPjtp1x23KEz/Rz8/UyNOOqTw61+eP06Trhu50GXdOPF+TZk6QycN2beoy+VIQ1GYqtMI4a0O968aFeFNX006ZtSZhPellwJNmJTRBhvkte9eCG96V1VtzwzhTXd9a0l4z73i1sKX0oYcvHsBenSGd1Hh3WHfE3XJWUdrw3VX17MvvKafnXR+4XzuOmusrF0POk1Dj9hP22/dX42z5xb+5WjGzNk66NhzNPrXJ2qN1VZUU3OLJk+drpW+/sXfAkZ4U7zeEd70FQfhTV9NOqPw/v2JUPc9GGqLLXLaecfafehEVBuONKT3zwzCm97aRDOrJeGN7rJw7PDLlc3ldPGZR+mxp174nPBGd3O4ZMxE9VyiuzZef001t7TogL120AZrr6b//O9djbzkJr3z/sfq0b2brjn/JK36jeX1xD9eVnTLs0h0G+rrdOSggdrrh9t8YWER3hSvd4Q3fcVBeNNXk84ovPc9EOrvT4YasGNOW26B8KZ3VdX2zBDedNe3loQ33aQrmx23JauMX6E3wmuAaI5AeM1AjXGd6UjDxN9l9OLLgX60V1YbblC7D51gh9e4wGOIQnhjgGqMRHiNMGOMQngNcBFeA0RzBMJrBmqM60zCG31hLfriWvSFteiLa7X84khDequL8Ka3NtHMEN5016djdgivoU4IrwGiOQLhNQM1xnUm4Y1uSRbdmuy4o7JaemmE17gMiCqBAMJbAqwqNEV4qwC9jCER3jKgLdoF4TVANEcgvGagxrjOJLxnjqpTLieNGN6qTGiEkMIodnhTWJT2KSG86a0NO7zprs2Cs0N4DbVCeA0QzREIrxmoMa6zCO+8edI559epe3fp9JNr+ylrUXkRXuMiN0chvGag5jh2eM1AY4pDeA1gEV4DRHMEwmsGaozrLML7yWTpN1fXadll8jrmyNq+By/Ca1zgMUQhvDFANUYivEaYMUYhvEXAff7f/9PF10zQjZed9oWtEd4iICbcBOFNGHgJw3UW4f3fG4FuvCWj1VfLa9BBCG8JJaapmQDCawZqjkN4zUBjiutywvvXJ54v3DC5qblZvXstoeHHHaSN1vtmAe/v731MV4ybpObmVu2w9SY644SDlcmEQnhjWn0xxiK8McKtMLqzCO8/nw/0+zsz+tZGee09EOGtsOx0r4AAwlsBvAS6IrwJQDYM0aWEt6U1q+/vfbx+e9Vwrbzi1/Tkc//WOZfdoj/eeI7eevcjDT7xPN1yxTAtu0w/nTLqGm203uoatO8AhNew0JKOQHiTJl78eJ1FeB99PNRDj4TaZqucdtiuth86EVWPM7zFr+GkWyK8SRMvbTyEtzRe1WrdpYR39px5+uGBp+iRiZcUdm4/nT5L+xz2y8Lvrxt/j2Y1ztHxh+1TqMWr/31HZ5w/ThPHnLmQ8M6Z26RBx5+rIQftru226l9oy5GGai3fLx8X4U1fTTpm1FmE9+57Qz31TKhdBuS0+WYIb3pXVO3PDOFNd40R3nTXp2N2XUp4o4s+/8rx+uDjqTpgzx009ta7tOfOW+uH22+uX14wTv03XEsDB2xVYNPU3KItdztKz94/Zr7wXn/JqTpm+GXaYpP1deDeO86vMMKbvsWO8KavJp1NeG+bEOrfr4ba70c5rbcuwpveFVX7M0N4011jhDfd9emywvvmOx/q8JMvKuzwfn25pXXhiCO11JK9dfLI0YUd2wHf32x+5dbfdpBe+vP1euGVN3Tp2Ilad41VlMvnderRByxU3aaW2v+fYedYzp/Nsj4TKJvLK1fbzwrobGUpzLehPlRLS05pL82Fl+f05jt5/eKYUKutEnRK1qVMOgik+kyo5lb+e1YKtyTahqGUCQK1ZNP+pyYJGukbo1t9jd+kO33Iy5pRl9rhnTptpg48epQuH3Ws1lxtJf3hvsc17rZ79buxZ2nUZTdrw3VW1z67fq8AsnH2XG2z57F67oGxhR3eI069WJkw1P4Dt9eRgwYuBHvqzOay4NMpPgK9l6hTU3NWza38DyI+yuUlL9WnQdNmNSuf8tKcc2GgGTOkU0/Kq9+S5V1rZ+oVhFK/ng36dBb/PUtb3aK/JEZSNWtO7d8POm3si5nP0n0aimlGmyoT6FLCe/fDT+pvT7+oc047bD72A44cqV+eeEjhC2yTp0zX0CP3K3z24qtvasQF4zTpupEF4T3roht07UUna9+fj9ClI4/RBmuvNj+DIw1VXsVfMDxHGtJXk44ZdYYzvJGMn3V213nKWlQbvrSW3j8zHGlIb22imXGkId316ZhdlxLel//zlk4eNVq/vfIM9e3TU++8/4kOOe4c/fGGczSzcY4OOfYc3XT56e13aRitNVdfSUccvMdCX1p77KkXdN6V43XH2LPUvVvb3+oQ3vQtdoQ3fTXpTMLb2Cidf3GdevaUTjmpa+yqIbzp/TOD8Ka3Nghvumuz4Oy6lPBGF37rpId02x8eLpwf7NZQX7grw9abb1Rgcs/DT+mi0bdrblNT4b2RQweroaH+c7cli3Z76+oyGnbcQQhvStc6wpvSwkjqDDu8H34kXT2mTssvLx35c4Q3vaupa8wM4U13ndnhTXd9uuQOb1wlYYc3LrLl5yK85bOLu2dnEN7/vBbo1tsyWmvNvA7cv/YfOsGRhrhXfWX5CG9l/OLujfDGTdiT3+V2eD3YFk5BeOOgWlkmwlsZvzh7dwbh/cezgf54d0ab9M9rj10R3jjXA9mLJ4DwLp5RNVsgvNWkX/zYCG/xrL60JcJrgGiOQHjNQI1xnUF4//zXUNGPbb+X03bf6xq36eIMr3GRm6MQXjNQcxzCawYaUxzCawCL8BogmiMQXjNQY1xnEN4778ro2ecC7bZLVptukvL7p5lqg/CaQMYQg/DGANUYifAaYcYYhfAa4CK8BojmCITXDNQY1xmE9+bxGb3+eqCf7JfV2mshvMbyE1UGAYS3DGgJdkF4E4RdwVAIbwXwOroivAaI5giE1wzUGNcZhDe6Q0N0p4Yhh7ZqhRWMF5/iKHZ401schDe9tYlmhvCmuz4ds0N4DXVCeA0QzREIrxmoMa4zCG90D97oXrwnn9iqXr2MF5/iKIQ3vcVBeNNbG4Q33bVZcHYIr6FWCK8BojkC4TUDNcalXXizOemsUXWFKz7rjFYFgfHiUxyF8Ka3OAhvemuD8Ka7NgivuT4IrxmoIQ7hNUCMKSLtwjtjRqCLLsuoT++8fnFC17glWVRqhDemBW+IRXgNEGOM4EhDjHCN0ezwGmAivAaI5giE1wzUGJd24X33/UBjr8topRXy+vmhCK+x9ESVSQDhLRNcQt0Q3oRAVzgMwlshwKg7wmuAaI5AeM1AjXFpF95/vxLqtomh1l0np/337Rr34GWH17jAY4hCeGOAaoxEeI0wY4xCeA1wEV4DRHMEwmsGaoxLu/A+9XSou+8LtdmmOe26M8JrLD1RZRJAeMsEl1A3hDch0BUOg/BWCJAdXgPAGCIQ3higmiLTLrwPPRLq0cdDbf/9nL63NcJrKjsxFRBAeCuAl0BXhDcByIYhEF4DRHZ4DRDNEQivGagxLu3CO+nOjP71fKA9d89q4293jYdOROXlS2vGRW6OQnjNQM1xCK8ZaExxCK8BLMJrgGiOQHjNQI1xaRfeG27O6I03Ax38k6zW+CbCayw9UWUSQHjLBJdQN4Q3IdAVDoPwVggw6o7wGiCaIxBeM1BjXNqF94qrMpo8JdDRQ1q13HLGC095FDu86S0Qwpve2kQzQ3jTXZ+O2SG8hjohvAaI5giE1wzUGJd24T3n/DrNmyedOrRVS/QwXnjKoxDe9BYI4U1vbRDedNdmwdkhvIZaIbwGiOYIhNcM1BiXZuHteMpaGEpnDm81XnX6oxDe9NYI4U1vbRDedNcG4TXXB+E1AzXEIbwGiDFFpFl4p04NdNmVGfXrl9cJx3Sdh05EpUZ4Y1rwhliE1wAxxgiONMQI1xjNDq8BJsJrgGiOQHjNQI1xaRbet94ONO7GjFZZOa+fDUJ4jWUnqgICCG8F8BLoivAmANkwBMJrgIjwGiCaIxBeM1BjXJqF98WXAk2clNEG6+e1794Ir7HsRFVAAOGtAF4CXRHeBCAbhkB4DRARXgNEcwTCawZqjEuz8P7t76HufyjUFv+X084/6DoPnYjKy5EG4yI3RyG8ZqDmOITXDDSmOITXABbhNUA0RyC8ZqDGuDQL770PhnriiVA77ZjTd7dAeI1lJ6oCAghvBfAS6IrwJgDZMATCa4CI8BogmiMQXn4KhFYAACAASURBVDNQY1yahXfCpIxeeinQPntltdEGXeehE+zwGhd4DFEIbwxQjZEIrxFmjFEIrwEuwmuAaI5AeM1AjXFpFt7rbsjo7XcC/fTgrFZbFeE1lp2oCgggvBXAS6ArwpsAZMMQCK8BIsJrgGiOQHjNQI1xaRbeS6/I6NNpgY49KqtllkZ4jWUnqgICCG8F8BLoivAmANkwBMJrgIjwGiCaIxBeM1BjXJqF98xRdcrlpDNOa1V9vfGiO0EUX1pLb5EQ3vTWJpoZwpvu+nTMDuE11AnhNUA0RyC8ZqDGuLQK79x50rnn16lbgzTs1K71lLWovAivcZGboxBeM1BzHMJrBhpTHMJrAIvwGiCaIxBeM1BjXFqF95NPpN+MrtOyy+R1zJFd6x68CK9xgccQhfDGANUYifAaYcYYhfAa4CK8BojmCITXDNQYl1bh/d8bgW68JaPVV8tr0EEIr7HkRFVIAOGtEGDM3RHemAGb4hFeA0iE1wDRHIHwmoEa49IqvP96PtCkOzP69kZ57TUQ4TWWnKgKCSC8FQKMuTvCGzNgUzzCawCJ8BogmiMQXjNQY1xahfevj4V6+M+httoypx/s0LUeOhGVlzO8xkVujkJ4zUDNcQivGWhMcQivASzCa4BojkB4zUCNcWkV3rvuDfX0M6F23imnLTZHeI0lJ6pCAghvhQBj7o7wxgzYFI/wGkAivAaI5giE1wzUGJdW4R0/IdQrr4b68Y9yWn9dhNdYcqIqJIDwVggw5u4Ib8yATfEIrwEkwmuAaI5AeM1AjXFpFd4x12X03vuBDvtpVt/4Rtd66ERUXo40GBe5OQrhNQM1xyG8ZqAxxSG8BrAIrwGiOQLhNQM1xqVVeC+8JKOZswKdeGxWSy6J8BpLTlSFBBDeCgHG3B3hjRmwKR7hNYBEeA0QzREIrxmoMS6NwpvPSyNG1hWucsTwVmVC4wV3kih2eNNbKIQ3vbWJZobwprs+HbNDeA11QngNEM0RCK8ZqDEujcLb2Cidf3GdllhCOvUXXe8pa1F5EV7jIjdHIbxmoOY4hNcMNKY4hNcAFuE1QDRHILxmoMa4NArvBx9Io6+t09e+Jh11OMJrLDdRBgIIrwFijBEIb4xwjdEIrwEmwmuAaI5AeM1AjXFpFN7/vBbo1tsyWnONvA46oOs9dIIdXuMCjyEK4Y0BqjES4TXCjDEK4TXARXgNEM0RCK8ZqDEujcL7zLOB/nR3Rv03zmvgbgivsdxEGQggvAaIMUYgvDHCNUYjvAaYCK8BojkC4TUDNcalUXgf+Wuov/w11Pe2yWn7bbvePXjZ4TUu8BiiEN4YoBojEV4jzBijEF4DXITXANEcgfCagRrj0ii8d96V0bPPBdptl5w23QThNZabKAMBhNcAMcYIhDdGuMZohNcAE+E1QDRHILxmoMa4NArvLeMzeu31QD/ZL6e110J4jeUmykAA4TVAjDEC4Y0RrjEa4TXARHgNEM0RCK8ZqDEujcJ71Zg6ffSRNOTQVq2wgvFiO1EUtyVLb7EQ3vTWJpoZwpvu+nTMDuE11AnhNUA0RyC8ZqDGuDQK768vrNOcOdLQE1rVu7fxYjtRFMKb3mIhvOmtDcKb7tosODuE11ArhNcA0RyB8JqBGuPSJrzZnHTWqLanrJ11RquCwHixnSgK4U1vsRDe9NYG4U13bRBec30QXjNQQxzCa4AYU0TahHf6jEAXX5ZRn955/eKErnlLsqjUCG9MC94Qi/AaIMYYwZGGGOEao9nhNcBEeA0QzREIrxmoMS5twvvue4HGjstoxRXyOvxQhPejafOM1SbKQQDhdVCMLwPhjY+tMxnhNdBEeA0QzREIrxmoMS5twvvyK6FunxhqnbVzOuDHXfMODezwGhd4DFEIbwxQjZEIrxFmjFEIrwEuwmuAaI5AeM1AjXFpE94nnw51z31h4f670X14u+qLIw3prTzCm97aRDNDeNNdn47ZIbyGOiG8BojmCITXDNQYlzbhffDhUI/9LdR2389p260RXo40GBe7KQrhNYGMKQbhjQmsORbhNQBFeA0QzREIrxmoMS5twjvpDxn964VAA3fPqv+388Yr7VxR7PCmt14Ib3prww5vumuz4OwQ3iJq9fy//6eLr5mgGy877QtbI7xFQEy4CcKbMPAShkub8N5wc0ZvvBnooAOyWnMNhJcd3hIWc0JNEd6EQJc5DDu8ZYJLuFuXFN5bfvegxv/hYTU1t2jA9zfTL4b8uID99/c+pivGTVJzc6t22HoTnXHCwcpkQiG8Ca9Kw3AIrwFiTBFpE94rrspo8pRARx7equW/FtNFd4JYdnjTWySEN721YYc33bXp0ju8N0y4T0//8xWdfeqh6tf3s0cqvfXuRxp84nm65YphWnaZfjpl1DXaaL3VNWjfAQhv51nP82eK8Ka3aGkT3nPOq9O8JunUoa1aokd6ucU9M4Q3bsLl5yO85bNLoic7vElQrnyMLrXDm83mtNMBQ3XHmLO0ZN9eC9G7bvw9mtU4R8cftk/h/Vf/+47OOH+cJo45cyHhnTO3SYOOP1dDDtpd223Vv9CWIw2VL0R3AsLrJurLS5PwtrRII8+tUxhKZw5v9V1kJ0xCeNNbNIQ3vbVhhzfdtemyO7yvvfGehv36Wm3ef1099tQL6t1zCZ14+I/Uf8O19MsLxhV+HjhgqwKf6LjDlrsdpWfvHzNfeK+/5FQdM/wybbHJ+jpw7x3nc0R407fgEd701aRjRmkS3ilTA11+ZUb9lszrhGO77kMnotogvOn9M4Pwprc2CG+6a9Nlhfexp17UcWdcrnNPP0w7brOpXnn9bR097FLdffN5OvPC6ws7ttGZ3o7X+tsO0kt/vl4vvPKGLh07UeuusYpy+bxOPfqAhSo8rbG581S8i8y0V/c6NbXk1JLtureZSmupl+zVoBmNzUrD18P+96Y0+tpAq64sHXV4GmZUvaoFQaC+S9Rr+mz+e1a9KnzxyPV1obrVhWqc17X/FSJtdemYT79eDWmdGvNagECXOtLwxD9e1lU33qmbrzh9PoLBJ5yn4w7bR5PueVQbrrO69tn1e4XPGmfP1TZ7HqvnHhhb2OE94tSLlQlD7T9wex05aOBCi2huU9feGUrjn6iG+lDRERZ8N33ViXarmlqyyqfAL//xz7xuHJ/XxhtJgw8M0wcrwRkFgdStPtS8Zv6SmCD2oobKhIEymUDNLdSmKGAJN+rRLZPwiAxXDoEuJbzvvP+Jhpxyke655bz5rA48+mwNO+5APfXPVzR5ynQNPXK/wmcvvvqmRlwwTpOuG1kQ3rMuukHXXnSy9v35CF068hhtsPZq8zM40lDO0ou3D0ca4uVbSXqajjT87YlQ9z8YaovNc9p5p64tExxpqGRVx9uXIw3x8q00nS+tVUowmf5dSngjpAcdc7Z23+m7+tGu2+qfL72uk0ddo7tvOleTP52hQ449Rzddfnr7XRpGa83VV9IRB++x0JfWorO/5105XneMPUvdu7X9MwbCm8xiLWUUhLcUWsm2TZPw3vtAqCeeDPWDHXLaakuEd7m+3cR9eJP981DMaAhvMZSq1wbhrR77UkbucsIb7fIO+/VYRT8vt0w//fKEg7XhuqsXmN3z8FO6aPTtmtvUpK0330gjhw5WQ0P9525LFu321tVlNOy4gxDeUlZbgm0R3gRhlzhUmoR3wu8yeunlQHvvmdW3NkzBGYsSWTqbs8PrpOnNQni9PN1pCK+baDx5XU5448DIDm8cVCvLRHgr4xdn7zQJ7x2/D/XCi6EGH5LVqqsgvOzwxrnyy89GeMtnl0RPhDcJypWPgfBWzpAjDQaG7giE103Ul5cm4T3vojrNni2ddEJWfXsjvAivb507kxBeJ01/FsLrZxpHIsJroMoOrwGiOQLhNQM1xqVFeJubpVG/rlNDvTT8NG73xJEG4yI3RyG8ZqDmOITXDDSmOITXABbhNUA0RyC8ZqDGuLQI77vvBhp7fUYrfD2vIYdxa0GE17jIzVEIrxmoOQ7hNQONKQ7hNYBFeA0QzREIrxmoMS4twvvcvwL94Y+ZwpfVoi+tdfUXwpveFYDwprc20cwQ3nTXp2N2CK+hTgivAaI5AuE1AzXGpUV4H3go1ON/D7Xd93PaduuufUuyqLwIr3GRm6MQXjNQcxzCawYaUxzCawCL8BogmiMQXjNQY1xahPe3t4d69T+h9t0npw3WQ3gRXuMiN0chvGag5jiE1ww0pjiE1wAW4TVANEcgvGagxri0CO9lV2Y0dWqgo4a06mvLGS+wk0YhvOktHMKb3tpEM0N4012fjtkhvIY6IbwGiOYIhNcM1BiXBuHN5aSzzq5TPi+NGNaqTMZ4gZ00CuFNb+EQ3vTWBuFNd20WnB3Ca6gVwmuAaI5AeM1AjXFpEN7JUwJdcVVGS/XL6/hj+MJaVF6E17jIzVEIrxmoOY4dXjPQmOIQXgNYhNcA0RyB8JqBGuPSILyvvBpq/IRQa66Z10H7I7wIr3GBxxCF8MYA1RiJ8BphxhiF8BrgIrwGiOYIhNcM1BiXBuF97G+hHnw41Jb/l9OAH/CFNYTXuMBjiEJ4Y4BqjER4jTBjjEJ4DXARXgNEcwTCawZqjEuD8E66M6N/PR9o912y+s4mXfuRwh2l5UiDcZGboxBeM1BzHMJrBhpTHMJrAIvwGiCaIxBeM1BjXBqEd8x1Gb33fqDBh2S16ioILzu8xgUeQxTCGwNUYyTCa4QZYxTCa4CL8BogmiMQXjNQY1wahHfUuXVqbpFOOalVPXsaL64TR7HDm97iIbzprU00M4Q33fXpmB3Ca6gTwmuAaI5AeM1AjXHVFt7GRun8i+vUUC8NP63VeGWdOwrhTW/9EN701gbhTXdtFpwdwmuoFcJrgGiOQHjNQI1x1RbeN98KdP1NGX3jG3kd9lPu0NBRWoTXuMjNUQivGag5jh1eM9CY4hBeA1iE1wDRHIHwmoEa46otvM88G+hPd2e08bfz2nN3hBfhNS7umKIQ3pjAmmIRXhPImGMQXgNghNcA0RyB8JqBGuOqLbz3PhDqiSdD/WCHnLbakluSIbzGxR1TFMIbE1hTLMJrAhlzDMJrAIzwGiCaIxBeM1BjXLWF9+bfZvT6fwMd8OOc1lkb4UV4jYs7piiENyawpliE1wQy5hiE1wAY4TVANEcgvGagxrhqC+8ll2c0bXqg447MaulluCUZwmtc3DFFIbwxgTXFIrwmkDHHILwGwAivAaI5AuE1AzXGVVN4s1nprLPrFATSiGGtCkPjhXXyKL60lt4CIrzprU00M4Q33fXpmB3Ca6gTwmuAaI5AeM1AjXHVFN6PPpauuqZOyy6T1zFH8oW1BcuK8BoXuTkK4TUDNcchvGagMcUhvAawCK8BojkC4TUDNcZVU3hfejnQhN9ltO7aOe3/Y87vIrzGhR1jFMIbI1xDNMJrgJhABMJrgIzwGiCaIxBeM1BjXDWF9y+PhnrkL6G2/m5OO26P8CK8xoUdYxTCGyNcQzTCa4CYQATCa4CM8BogmiMQXjNQY1w1hfeOSRm98FKgvfbI6tvf4gtrCK9xYccYhfDGCNcQjfAaICYQgfAaICO8BojmCITXDNQYV03hvXpsnT78UDrsZ1l9Y0WEF+E1LuwYoxDeGOEaohFeA8QEIhBeA2SE1wDRHIHwmoEa46opvGeNqlM2Jw0/tVUNDcaLqoEovrSW3iIivOmtTTQzhDfd9emYHcJrqBPCa4BojkB4zUCNcdUS3hkzA110aUa9ekknn9hqvKLaiEJ401tHhDe9tUF4012bBWeH8BpqhfAaIJojEF4zUGNctYT3f28EuvGWjFZdJa/Bh3BLskVLivAaF7k5CuE1AzXHscNrBhpTHMJrAIvwGiCaIxBeM1BjXLWE96mnQ919X6hN++e0267coQHhNS7qmKMQ3pgBVxiP8FYIMKHuCK8BNMJrgGiOQHjNQI1x1RLeu+4N9fQzoQb8IKct/w/hRXiNizrmKIQ3ZsAVxiO8FQJMqDvCawCN8BogmiMQXjNQY1y1hPeGmzN6481ABx2Q1ZprcIcGhNe4qGOOQnhjBlxhPMJbIcCEuiO8BtAIrwGiOQLhNQM1xlVLeC+8JKOZswKdcGxW/ZZEeBFe46KOOQrhjRlwhfEIb4UAE+qO8BpAI7wGiOYIhNcM1BhXDeHNZqWzzq5TJpRGDOcODV9UTr60Zlzk5iiE1wzUHIfwmoHGFIfwGsAivAaI5giE1wzUGFcN4f3gA2n0tXVa/mvSkYcjvAivcUEnEIXwJgC5giEQ3grgJdgV4TXARngNEM0RCK8ZqDGuGsL7/IuBfvf7jNZfP68f780tyRBe44JOIArhTQByBUMgvBXAS7ArwmuAjfAaIJojEF4zUGNcNYT34T+H+utjobbdJqfttuUODQivcUEnEIXwJgC5giEQ3grgJdgV4TXARngNEM0RCK8ZqDGuGsJ7+8RQL78Sap89s9poQ76whvAaF3QCUQhvApArGALhrQBegl0RXgNshNcA0RyB8JqBGuOqIbxXjq7Tx59IQw5t1QorGC+mhqL40lp6i4nwprc20cwQ3nTXp2N2CK+hTgivAaI5AuE1AzXGVUN4zxpVp2xOGjGsVZmM8WJqKArhTW8xEd701gbhTXdtFpwdwmuoFcJrgGiOQHjNQI1xSQvvp9MCXXpFRn365PWL4/nC2peVEuE1LnJzFMJrBmqOY4fXDDSmOITXABbhNUA0RyC8ZqDGuKSF97XXA90yPqNvrp7XIQcivAivcTEnFIXwJgS6zGEQ3jLBJdwN4TUAR3gNEM0RCK8ZqDEuaeH9+xOh7nsw1Oab5bTLAO7QgPAaF3NCUQhvQqDLHAbhLRNcwt0QXgNwhNcA0RyB8JqBGuOSFt4778ro2ecC7bpzTpttivAivMbFnFAUwpsQ6DKHQXjLBJdwN4TXABzhNUA0RyC8ZqDGuKSFd9yNGb31dqBBB2W1+mrckgzhNS7mhKIQ3oRAlzkMwlsmuIS7IbwG4AivAaI5AuE1AzXGJS28511Up9mzpZNOyKpvb4QX4TUu5oSiEN6EQJc5DMJbJriEuyG8BuAIrwGiOQLhNQM1xiUpvM3N0qhf16mhXhp+WqvxKmovirs0pLemCG96axPNDOFNd306ZofwGuqE8BogmiMQXjNQY1ySwvvuu4HGXp/RCl/Pa8hh3KHhq8qI8BoXuTkK4TUDNcchvGagMcUhvAawCK8BojkC4TUDNcYlKbzP/SvQH/6Y0bc2zGvvPRFehNe4kBOMQngThF3GUAhvGdCq0AXhNUBHeA0QzREIrxmoMS5J4X3goVCP/z3Udt/PadutuUMDwmtcyAlGIbwJwi5jKIS3DGhV6ILwGqAjvAaI5giE1wzUGJek8P729lCv/ifUvvvktMF6CC/Ca1zICUYhvAnCLmMohLcMaFXogvAaoCO8BojmCITXDNQYl6TwXnZlRlOnBjpqSKu+tpzxImowijO86S0qwpve2kQzQ3jTXZ+O2SG8hjohvAaI5giE1wzUGJeU8OZy0lln1ymfl0YMa1UmY7yIGoxCeNNbVIQ3vbVBeNNdmwVnh/AWUavn//0/XXzNBN142Wlf2BrhLQJiwk0Q3oSBlzBcUsI7eUqgK67KaKl+eR1/DF9YW1yJEN7FEare5whv9dgXMzI7vMVQqn6bLiu8Uz6dod0POV3Djz9YP9x+80Ilfn/vY7pi3CQ1N7dqh6030RknHKxMJhTCW/2FWuoMEN5SiSXXPinhfeXVUOMnhFpzzbwO2h/hXVyFEd7FEare5whv9dgXMzLCWwyl6rfpssJ7zLDLNHvOPO2z67YF4X3r3Y80+MTzdMsVw7TsMv10yqhrtNF6q2vQvgMQ3uqv05JngPCWjCyxDkkJ72N/C/Xgw6G22CKnnXfkC2uLKzDCuzhC1fsc4a0e+2JGRniLoVT9Nl1SeP/4wN/03Auvq1u3en1rvTUKwnvd+Hs0q3GOjj9sn0JVXv3vOzrj/HGaOObMhYR3ztwmDTr+XA05aHdtt1X/QluONFR/IS86A4Q3fTXpmFFSwjvpzoz+9Xyg3XfJ6jub8Ejhxa0IhHdxhKr3OcJbPfbFjIzwFkOp+m26nPBOnjpdQ065WDddfrouGTNR/TdcqyC8v7xgXOHXAwdsVahKU3OLttztKD17/5j5wnv9JafqmOGXaYtN1teBe+84v3oIb/UXMsKbvhp82YySEt4x12X03vuBBh+S1aqrILyLWyEI7+IIVe9zhLd67IsZGeEthlL123Q54T3q9Et18I920uYbr6tRl948X3hPHjm6sGM74Pubza/K+tsO0kt/vl4vvPKGLh07UeuusYpy+bxOPfqAhSqXy/E/0+ov5YVnEARSoSqUJm2lURAEyke3Toj5dcwprWpqli4aWafevWIerEbiC39u4i9NjdBK8DICKYj+c0ZtEoRe/FDRXxZ5pZ9AlxLe6Etp/37tLQ077qBCZRYU3hEXXq8N11ld++z6vcJnjbPnaps9j9VzD4wt7PAecerFyoSh9h+4vY4cNHChyn746dz0V7qLzbBfrwbNndeqea2c3Uxb6b/Wr4cmT5+rOP+e2NgonXdRnRrqpTNOb00bglTOJ/qf9rJ9uunj6fNSOb+uPKnu9Rn16JbRtMbmrowhtdf+9aV6pHZuTOwzAl1KeA/7xYV6/t//LewwRa+mpmZlMhntMWArrbzicpo8ZbqGHrlf4bMXX31TIy4Yp0nXjSwI71kX3aBrLzpZ+/58hC4deYw2WHu1+RQ50pC+P1Kc4U1fTTpmlMSRhrfeDjTuxoy+8Y28Dvspd2goZjVwpKEYStVpw5GG6nAvdlSONBRLqrrtupTwLop6wR3e9z+aokOOPadwtrftLg2jtebqK+mIg/dY6Etrjz31gs67crzuGHuWundrKEQivNVdxF80OsKbvpokKbzPPBvoT3dntPG389pzd4S3mNWA8BZDqTptEN7qcC92VIS3WFLVbYfwtn9pLSrDPQ8/pYtG3665TU3aevONNHLoYDU01H/utmTRbm9dXWb+0QiEt7qLGOFNH/+vmlESO7z3PhDqiSdD/WCHnLbakmMtxawQhLcYStVpg/BWh3uxoyK8xZKqbrsuLbwu9Aivi6Qvhx1eH0t3UhLCe/NvM3r9v4EO+HFO66yN8BZTQ4S3GErVaYPwVod7saMivMWSqm47hNfAH+E1QDRHILxmoMa4JIT3ksszmjY90HFHZrX0Mny1vZjyIbzFUKpOG4S3OtyLHRXhLZZUddshvAb+CK8BojkC4TUDNcbFLbzZrHTW2XWKvps6YlirwtA4+RqOQnjTW1yEN721iWaG8Ka7Ph2zQ3gNdUJ4DRDNEQivGagxLm7hnfxJoCtGZ7TsMjkdcyTHGYotHcJbLKnk2yG8yTMvZUSEtxRa1WuL8BrYI7wGiOYIhNcM1BgXt/BGjxOOHiu80QZ57bMXd2gotnQIb7Gkkm+H8CbPvJQREd5SaFWvLcJrYI/wGiCaIxBeM1BjXNzC+4c/ZfTcPwPttktOm27CDm+xpUN4iyWVfDuEN3nmpYyI8JZCq3ptEV4De4TXANEcgfCagRrj4hbey36T0dRPAx19RKuWW9Y48RqPQnjTW2CEN721iWaG8Ka7Ph2zQ3gNdUJ4DRDNEQivGagxLk7hnTs30LkXZNSjR16nDeU4QyllQ3hLoZVsW4Q3Wd6ljobwlkqsOu0RXgN3hNcA0RyB8JqBGuPiFN6X/h1qwh2h1l4rp5/sx3GGUsqG8JZCK9m2CG+yvEsdDeEtlVh12iO8Bu4IrwGiOQLhNQM1xsUpvPfcF+rJp0PttGNO390C4S2lbAhvKbSSbYvwJsu71NEQ3lKJVac9wmvgjvAaIJojEF4zUGNcnMJ79Zg6ffiRdNjPsvrGijxwopSyIbyl0Eq2LcKbLO9SR0N4SyVWnfYIr4E7wmuAaI5AeM1AjXFxCW9Tc6BzzssUHjQx/LRWZXjgRElVQ3hLwpVoY4Q3UdwlD4bwloysKh0QXgN2hNcA0RyB8JqBGuPiEt7X/xvo5t9mtNqqef30YL6wVmrJEN5SiSXXHuFNjnU5IyG85VBLvg/Ca2CO8BogmiMQXjNQY1xcwvvQI6EefTzUttvktN22nN8ttWQIb6nEkmuP8CbHupyREN5yqCXfB+E1MEd4DRDNEQivGagxLi7hve6GjN5+J9AhB2b1zdU5v1tqyRDeUokl1x7hTY51OSMhvOVQS74PwmtgjvAaIJojEF4zUGNcHMKbzUkjz6lTPi+dcVqr6uqME+4iUQhveguN8Ka3NtHMEN5016djdgivoU4IrwGiOQLhNQM1xsUhvNHObrTDu+IKeR1+KOd3yykXwlsOtWT6ILzJcC53FIS3XHLJ9kN4DbwRXgNEcwTCawZqjItDeKOzu9EZ3i22yGnnHTm/W065EN5yqCXTB+FNhnO5oyC85ZJLth/Ca+CN8BogmiMQXjNQY1wcwhvdnSG6S8P+++a07joIbznlQnjLoZZMH4Q3Gc7ljoLwlksu2X4Ir4E3wmuAaI5AeM1AjXFu4Y3O7Z59XkbNzYFOG5pVjx58Ya2cciG85VBLpg/CmwznckdBeMsll2w/hNfAG+E1QDRHILxmoMY4t/B++KF09dg6LbN0XscexfndckuF8JZLLv5+CG/8jCsZAeGthF5yfRFeA2uE1wDRHIHwmoEa49zC+8RToe69P9Qm/fPaY1eEt9xSIbzlkou/H8IbP+NKRkB4K6GXXF+E18Aa4TVANEcgvGagxji38N4+MdTLr4Taa2BW396I4wzllgrhLZdc/P0Q3vgZVzICwlsJveT6IrwG1givAaI5AuE1AzXGuYX33Asymjs30AnHZtVvSYS33FIhvOWSi78fwhs/40pGQHgroZdcX4TXwBrhNUA0RyC8ZqDGOKfwTpkaZt64fgAAIABJREFU6PIrM+rVK6+TT+Q4QyVlQngroRdvX4Q3Xr6VpiO8lRJMpj/Ca+CM8BogmiMQXjNQY5xTeJ/7Z6A//CmjDTbIa9+9EN5KyoTwVkIv3r4Ib7x8K01HeCslmEx/hNfAGeE1QDRHILxmoMY4p/BOujOjfz0faNedc9psU+6/W0mZEN5K6MXbF+GNl2+l6QhvpQST6Y/wGjgjvAaI5giE1wzUGOcU3kuvyOjTaYGOPLxVy3/NOMkuGIXwprfoCG96axPNDOFNd306ZofwGuqE8BogmiMQXjNQY5xLeBsbA51/cUYNDXkNP5XjDJWWCOGtlGB8/RHe+Ng6khFeB8X4MxBeA2OE1wDRHIHwmoEa41zC+9JLgSZMymitNXM6cH+OM1RaIoS3UoLx9Ud442PrSEZ4HRTjz0B4DYwRXgNEcwTCawZqjHMJ7133hnr6mVA/2CGnrbZEeCstEcJbKcH4+iO88bF1JCO8DorxZyC8BsYIrwGiOQLhNQM1xrmE96pr6vTRx9KhP81q5W9w/91KS4TwVkowvv4Ib3xsHckIr4Ni/BkIr4ExwmuAaI5AeM1AjXEO4W1qDnTOeRmFoTT8tFZlQuMEu2gUwpvewiO86a1NNDOEN9316ZgdwmuoE8JrgGiOQHjNQI1xDuH9z2uBbr0to1VXyWvwIXxhzVEehNdBMZ4MhDcerq5UhNdFMt4chNfAF+E1QDRHILxmoMY4h/A+8FCox/8e6ntb57T99zm/6ygPwuugGE8GwhsPV1cqwusiGW8Owmvgi/AaIJojEF4zUGOcQ3jHXp/Ru+8GOvgnWa3xTc7vOsqD8DooxpOB8MbD1ZWK8LpIxpuD8Br4IrwGiOYIhNcM1BhXqfC2tkqjfl2nfF46/ZSsujUgvI7yILwOivFkILzxcHWlIrwukvHmILwGvgivAaI5AuE1AzXGVSq8b74V6PqbMlrh63kNOYzzu67SILwukv4chNfP1JmI8DppxpeF8BrYIrwGiOYIhNcM1BhXqfD+9dFQD/8l1Bab57TzTpzfdZUG4XWR9OcgvH6mzkSE10kzviyE18AW4TVANEcgvGagxrhKhffGWzL63xuBfvyjnNZfF+F1lQbhdZH05yC8fqbORITXSTO+LITXwBbhNUA0RyC8ZqDGuEqENzq3+6tz6pTNSqcNzapHD87vukqD8LpI+nMQXj9TZyLC66QZXxbCa2CL8BogmiMQXjNQY1wlwvv+B4GuuTajpZfK67ijOb9rLIsQXidNbxbC6+XpTkN43UTjyUN4DVwRXgNEcwTCawZqjKtEeP/+RKj7HgzVf+O8Bu6G8BrLgvA6YZqzEF4zUHMcwmsGGlMcwmsAi/AaIJojEF4zUGNcJcI7fkKoV14NtdceWX37WxxnMJYF4XXCNGchvGag5jiE1ww0pjiE1wAW4TVANEcgvGagxrhKhPfcCzKaOzcoHGeIjjXw8hHgSIOPpTsJ4XUT9eYhvF6ecaUhvAayCK8BojkC4TUDNcaVK7yTpwS64qqMevXK6+QTOc5gLEkhCuF1E/XlIbw+lnEkIbxxUPVnIrwGpgivAaI5AuE1AzXGlSu8/3g20B/vzmj99fP68d4Ir7EkCK8bpjkP4TUDNcchvGagMcUhvAawCK8BojkC4TUDNcaVK7yP/CXQXx7NaJedc9p8U+6/aywJwuuGac5DeM1AzXEIrxloTHEIrwEswmuAaI5AeM1AjXHlCu95F2Y0e06gE47Nqt+SnN81lgThdcM05yG8ZqDmOITXDDSmOITXABbhNUA0RyC8ZqDGuHKE97XXA90yPqMVV8jr8EM5zmAsx/wozvDGQdWTifB6OMaVgvDGRdabi/AaeCK8BojmCITXDNQYV47wTvhdRi+9HOiHA3L6v804zmAsB8IbB0xzJsJrBmqOQ3jNQGOKQ3gNYBFeA0RzBMJrBmqMK1V4m5ulcy+oUy7X9jjh7t05zmAsB8IbB0xzJsJrBmqOQ3jNQGOKQ3gNYBFeA0RzBMJrBmqMK1V4n30u0J13ZbTWmjkduD+7u8ZSLBTFkYa4yFaei/BWzjDOBIQ3Trq+bIS3CJbP//t/uviaCbrxstO+sDXCWwTEhJsgvAkDL2G4UoV33I0ZvfV2oH33zmqD9dndLQF1SU0R3pJwJdoY4U0Ud8mDIbwlI6tKhy4lvC2tWY25+Y/6/b2PqTWb05qrraRfnTxYX19uqQL86P0rxk1Sc3Ordth6E51xwsHKZEIhvFVZmxUNivBWhC/WzqUI78xZgS68JKOGeun0k1sVZmKdWpcOR3jTW36EN721iWaG8Ka7Ph2z61LCO2PmbN125yP6yV47qFfPHrrqxjv13zff08VnHqW33v1Ig088T7dcMUzLLtNPp4y6Rhutt7oG7TsA4e0ca3mhWSK86S1aKcL7l8dCPfLnUP03zmvgbtydIc6qIrxx0q0sG+GtjF/cvRHeuAl78ruU8C6K7NX/vqPTzhmj348bpevG36NZjXN0/GH7FJpFn51x/jhNHHPmQsI7Z26TBh1/roYctLu226p/oS1HGjyL0ZmC8DpperNKEd5LLsto2oxAgw/JatVVOM7grcTCaQhvnHQry0Z4K+MXd2+EN27CnvwuLbzRbu+/X3tLvxo6WL+8YJz6b7iWBg7YqkC2qblFW+52lJ69f8x84b3+klN1zPDLtMUm6+vAvXecXwGE17MYnSkIr5OmN6tY4X3v/UBjrsuod++8fnF8VkHgnQdpCG9nWQMIb7orhfCmuz4ds+uywvvx5Gn66Qm/1rUXDtUKyy+jk0eOLuzYDvj+ZvMrt/62g/TSn6/XC6+8oUvHTtS6a6yiXD6vU48+YKHqzprb2jmq3YVm2aMho5bWnFpz7Aqmrey9etRp9txWLa4yv7szr8efkHbcTvrhD7DduOsY/YWiZ7c6Nc7jv2dxsy41vy4TqD4Tam4zx3pKZZdE+9496pIYhjEqJNAlhXf6jEb97KTzdeLh++q7m25QQDjiwuu14Tqra59dv1f4fePsudpmz2P13ANjCzu8R5x6sTJhqP0Hbq8jBw1cWHjntFRYBrq7CfTo1i682cVplXtk8hZHoFePes2e16L8V5Qmm5WGj5LmzZPOGBpoqaWo4+K4Vvp5EATq2T2jRv4CXylKe/+6TKj6ukBzmxBeO1xDYO8l6g0pRMRNoMsJb3RO97ChF+qnPx6gnbb9bDf3hgn3afKU6Rp65H4F5i+++qZGXDBOk64bWRDesy66QddedLL2/fkIXTryGG2w9mrza8ORhriXaen5HGkonVlSPYo50vDKf0KNvz3kUcJJFUUSZ3gThF3iUBxpKBFYws050pAw8DKH61LCO2fuPB1+8sWFuzQseHQhYvf+R1N0yLHn6KbLT2+/S8Norbn6Sjri4D0W+tLaY0+9oPOuHK87xp6l7t0aCtgR3jJXX4zdEN4Y4VYYXYzw3jYx1L9fCXmUcIWsS+mO8JZCK9m2CG+yvEsdDeEtlVh12ncp4Z3wp78Udmqj/7Av+Lr1N8O10Xrf1D0PP6WLRt+uuU1N2nrzjTRy6GA1NNR/7rZkUUZdXUbDjjsI4a3Oul3sqAjvYhFVrcHihHfevEC/vjBTOPLAo4STKxPCmxzrUkdCeEsllmx7hDdZ3uWO1qWEt1xIi+vHDu/iCCX/OcKbPPNiR1yc8D7zbKg/3R3yKOFigZraIbwmkDHEILwxQDVGIrxGmDFGIbwGuAivAaI5AuE1AzXGLU54x16f0bvv8ihhI/KiohDeojBVpRHCWxXsRQ+K8BaNqqoNEV4DfoTXANEcgfCagRrjvkp4p00PdMnlPErYiLvoKIS3aFSJN0R4E0de0oAIb0m4qtYY4TWgR3gNEM0RCK8ZqDHuq4T3z38NFf3YpH9ee+zKLZiM2BcbhfAuFlHVGiC8VUP/pQMHLc3SzKkKZs/S8ht/K30TZEafI4DwGhYFwmuAaI5AeM1AjXFfJbw8StgIusQohLdEYAk2R3gTgJ3LKWicoaBxuoIZ06TG6QpnTZdmTlMwa3rhh2a1vxf9vnle26Tq6rXkb/+cwAQZolICCG+lBLktmYGgPwLh9TN1JX6Z8L7zbqBrr+dRwi7OpeYgvKUSS649wlse62DubAWzphVEVjOnK2xsk9YOoQ1mTlfQOE2Fn2fPLGmQfF2D8n2WVL7v0lr6vLEl9aVxdQggvAbu7PAaIJojEF4zUGPclwnvn+4K9cxzobbZOqcdvp8zjkhUMQQQ3mIoVacNwtvOPZdVOHN6Yfe1IKmFnddP23ZhI5mNpHbWjM8+b20uvmBBoPwSvdsktle/ws/qFQltP6n3ksr1WlKKPuvd/nmPJeZnc4a3eMzVbInwGugjvAaI5giE1wzUGPdFwpvLSudelFHTvEAnHJtVvyV5lLAReVFRCG9RmKrSqJaFN5jT2H5cYNpnMhuJ7MxpbbuxhR/txwrmNJbEP9+th/KRtEbyWpDWvlKffoUfBXHt3a8gtblIZHv2iR43WFJ+R2OEtyxsiXdCeA3IEV4DRHMEwmsGaoz7IuF9+ZVQt0/kUcJGzCVHIbwlI0usQ2cS3qC1Zf6518IZ2FltIpsvSOuC52Hbz8tmS/hyaphRrlcfqXfbDmxBZntHO68dAtsuttF7fZdSdOwgiRfCmwTlysdAeCtnyKOFDQzdEQivm6gv74uE97e3h3r1PzxK2Ee59CSEt3RmSfWoqvDm84XzrYWjBNFua2P7zutCX+xaYBe2aW5JWHLR0YDoCEG049ousdFubL6wE7vAEYLovZ69pWDhJ6WWNFhMjRHemMCaYxFeA1B2eA0QzREIrxmoMW5R4eVRwka4FUQhvBXAi7mrW3g7bqn1mcTOKJyJ1cxP5x8h6LgjQeHLXLniz9Tn6+qUj44OtO++RuLatgu7VGFHNpLYXO/2z6P36upiphd/PMIbP2PHCAivgSLCa4BojkB4zUCNcYsK71PPhLr7Xh4lbERcVhTCWxa2RDotVnjLvaVWkbPPL9Fr/pnXgry278AWft1+BrbwXnTUYIleRabWTjOEt3PUEuE11AnhNUA0RyC8ZqDGuEWFd8y1Gb33AY8SNiIuKwrhLQtbbJ2CuXPabpk1a7oa5sxUt3kzNXvylC++pdacWVK++C96dtxSqyCr0W5rx1GC6NzrAscL5stsmIntOmshGOHtHFVEeA11QngNEM0RCK8ZqDFuQeF9651A427IqFu3vIadUsKXV4zzIaqNAMIb80rIZhXOmlb4QlfhfrDtt9WKvtRVeKjBgveJjT7LthQ/oQpuqVX8ILT8MgIIb+dYGwivoU4IrwGiOQLhNQM1xi0ovGOuy+i99wNt+X85DfhB8ecEjdMhqp0Awlv6UogeKzv/tlkdEhvdTmuh+8S2f6Fr7uySBljwllph337K9F1STT3abqvlvKVWSZOi8RcSQHg7x8JAeA11QngNEM0RCK8ZqDGuQ3iffzHQxEkZde+e10nH59Stofh/kjVOhyiEd/4aSOsttRZ7hpdVXFUCCG9V8Rc9OMJbNKovb4jwGiCaIxBeM1BjXCS8H06ZqwsuzaixMdAuO+e0+abs7hoRlxVVkzu8NXJLLYS3rCWdWCeENzHUFQ2E8FaEr60zwmuAaI5AeM1AjXGR8E74Y5MeeiTUssvkdNSQXLkPODLOiqjOIrxBdJ/XxukL3Be29m+phfCm+88nwpvu+nTMDuE11AnhNUA0RyC8ZqDGuJ4NPXTqmS1qaZV+NiirVVbmKIMRb9lRVRPeXPRlruihBjPaz8O2nXmNnswV3Se28IWu6NeFn6ercA/ZEl61cEsthLeEglehKcJbBehlDInwlgFt0S4IrwGiOQLhNQM1xt13f4P+/lROa62Z04H7c5TBiLaiKKfwBnNnL/QY2TZpna4g+kJXJLHRl7o6HjM7e1ZJ8+6Kt9RCeEtaIok3RngTR17WgAhvWdgW7oTwGiCaIxBeM1BT3IcfBbp6TEaZUDrmyKyWWordXRPaimO+UnjjvqVWzz7K9Wl/oEH7fWDVt1/bQw0KDzro1/bY2d59le+2RMXX2tkCEN50VwzhTXd9OmaH8BrqhPAaIJojEF4zUFMctyEzgTTELHpLrWDWDPVqnqnG6OEG7UcJ5t9yq4JbakVP4MpFj5rt049bapVZN4S3THAJdUN4EwJd4TAIb4UAo+4IrwGiOQLhNQM1xP371VC3TQi1RA/ppOOzqq9nd9eAdX5E0No8/2EGbedd287ARmdh5x8f6DgrGx0pyJbwoI8wo1yvPp/ttPaKdluXXEBg23ZnC7ux0dO66hqcl9blsxDedC8BhDfd9WGH11gfhNcI0xSF8JpAmmJyWenS32Q0fUagA/bJaL31m5TDd7+abi6nYPbMtrsRLHDmNXpSV9vTuma0n4ttPxMb3b2ghFeuxxLSgo+R7b2kei67jGbV92571GzvJT97zGzP3lIQlJBOUycBhNdJ05+F8PqZxpHIDq+BKsJrgGiOQHjNQCuM+9sToe5/sO02ZKOGddPk6XO7pPAGTXMUHR2IjgzMv/NA9ISuGW2Plw0apxUktnB3gtkzpXzxfyvI19UpHx0daN99zbc/kUt9llK+V7QbG+3Atn8evVdXt1BVnV9aq3C50H0RAghvupcEwpvu+rDDa6wPwmuEaYpCeE0gDTFz5wa6+LJQTc1B4TZkm3+7uz6ZVtvCW//wHQree/Oz+8V2PGo2OnZQwivfo6cicS18eat911V9l54vtm1f9OpXENp8z14lJH++KcJbEb5YOyO8seKtOBzhrRhhIgHs8BowI7wGiOYIhNcMtIK4P94d6h/Phlp7rZx+sl9OHY8WruUjDd0vOUnhay98jlqab6mF8FawyGPuivDGDLjCeIS3QoAJdUd4DaARXgNEcwTCawZaZtwnkwNdOTpTOP55/DFZLdk33yWEt+65xxTMmalc+1GC+edho3OzKX0hvCktjCSEN721iWaG8Ka7Ph2zQ3gNdUJ4DRDNEQivGWiZcdfdkNHb7wTa8v9yGvCDtodMdIUd3jJxVbUbwltV/F85OMKb3togvOmuzYKzQ3gNtUJ4DRDNEQivGWgZcf95LdCtt2XUvXteJx2fU7eGti9gIbxlwEygC8KbAOQyh0B4ywSXUDd2eBMCXeEwCG+FAKPuCK8BojkC4TUDLTEul5MuvzKjT6cF2m2XrDbd5LO7DSC8JcJMqDnCmxDoMoZBeMuAlmAXhDdB2BUMhfBWAK+jK8JrgGiOQHjNQEuMe/iRQH99PKPlls3r6CMWfsABwlsizISaI7wJgS5jGIS3DGgJdkF4E4RdwVAIbwXwEF4DvJgiEN6YwBYR+8JLge6YlCm0HHxIVquusvC9ZBHeIiBWoQnCWwXoRQ6J8BYJqkrNEN4qgS9xWIS3RGBf1JwdXgNEcwTCawZaZNwbb4a66dZQ0ZGGXXbOafNN276otuAL4S0SZsLNEN6EgZcwHMJbAqwqNEV4qwC9jCER3jKgLdoF4TVANEcgvGagRcS9/36gcTdl1NIibbVlTj/Y4fOyG8UgvEXArEIThLcK0IscEuEtElSVmiG8VQJf4rAIb4nAvqg5wmuAaI5AeM1AFxM3ZWqoa64L1DQv0EYb5rTPnl8suwhvsnUpZTSEtxRaybZFeJPlXepoCG+pxKrTHuE1cEd4DRDNEQivGehXxM2aFejqsaEaGwOt8c2cDtw/pzD88g7s8CZXm1JGQnhLoZVsW4Q3Wd6ljobwlkqsOu0RXgN3hNcA0RyB8JqBfkncvHlBYWd36tRQK66Q188GZVVX99VjI7zJ1KbUURDeUokl1x7hTY51OSMhvOVQS74PwmtgjvAaIJojEF4z0C+Ia22Voiepvf9BoKWXzunwn+ULD5lY3AvhXRyh6nyO8FaHezGjIrzFUKpeG4S3euxLGRnhLYXWl7RFeA0QzREIrxnoInHZnHTr+FD//V+oXr3yOuKwnHr3XrzsRjEIb7y1KTcd4S2XXPz9EN74GVcyAsJbCb3k+iK8BtYIrwGiOQLhNQNdIC6flybcEerlV0J1654v7Owus/SXf0lt0ZkgvPHVppJkhLcSevH2RXjj5VtpOsJbKcFk+iO8Bs4IrwGiOQLhNQNdIO6ue0M9/UyoTEY6dFBWK65Y3M5uRwTCG19tKklGeCuhF29fhDdevpWmI7yVEkymP8Jr4IzwGiCaIxBeM9D2uEcfD/XQI6GCQDrogFzhrgylvhDeUokl0x7hTYZzOaMgvOVQS64Pwpsc60pGQngrodfeF+E1QDRHILxmoJL+8lhGjz8eqLlFGrh7Tv2/XbrsRrNCeP21cSQivA6K8WQgvPFwdaUivC6S8eYgvAa+CK8BojkC4fUBnfxJoDt+H+rDj4NC6C4Dctp8s/JkF+H11cWdhPC6ifryEF4fyziSEN44qPozEV4DU4TXANEcgfBWDjS6E8Njj2f0l0cD5XLSssvmtd/eOS27XGlndhedCTu8ldcmjgSENw6qnkyE18MxrhSENy6y3lyE18AT4TVANEcgvJUBjXZ1J0wK9fEnQeG87lZbZrX9tnmFmcpy2eGtnF9cCQhvXGQrz0V4K2cYZwLCGyddXzbCa2CJ8BogmiMQ3vKALrqru1S/vPbdJ6cVvl7Zru6Cs2GHt7zaxN0L4Y2bcPn5CG/57JLoifAmQbnyMRDeyhkK4TVANEcgvKUDXXRX9/82y2nH7XOLfVRwqSMhvKUSS6Y9wpsM53JGQXjLoZZcH4Q3OdaVjITwVkKvvS/Ca4BojkB4iwe66K7ukkvmtc+eOa38Dd+uLju8xdejWi0R3mqRX/y4CO/iGVWzBcJbTfrFj43wFs/qS1sivAaI5giEd/FAm5ul518M9c/nA733XtsdGL7TP6cBP8ipoWHx/cttwQ5vueTi7Yfwxsu3knSEtxJ68fdFeONn7BgB4TVQRHgNEM0RCO+XA41uL/bkk4FefDlUa2tbu2hXd+Buea2+Wvm3Gyu2hAhvsaSSbYfwJsu7lNEQ3lJoJd8W4U2eeTkjIrzlUFukD8JrgGiOQHgXBtrSKr3wYqhnng30wQdtu7lhKK27dk7f2UQF0Y3uxpDEC+FNgnLpYyC8pTNLqgfCmxTp8sZBeMvjlnQvhLdC4nc9kNOcplbV1eWVyUiZOqkuI9W3/7rtvfwCvw8Kbevr1N4+r4b6CidB988RQHjbkEyZGuqpp6V//issPCEtevXrmy9I7iab5LREj3jO6X7VkkR40/kHFuFNZ12iWSG86a1NNDOEN9316ZgdwlthnQ49rt0iKsyJutdFEryAHEfiHAlz9H7h1+1S3SbLgerq2yQ7+qwg0fVSmAnn/76jbyYTCXdQkPGCXIefCXhdXfR+++/bx6gFAe+KwhsdT3j3vVBvvRPo7bfz+uDDQPPmfbZtu8F6kegmc2wB4TX8ByHhCIQ3YeAlDIfwlgCrCk0R3ipAL2NIhLcIaNNmzNJp54zRC6+8oaWW7KNfDR2s/huuWeh5571ZfTqjRa2tgVqzKpyJzLYGymbziv4ZOZtt+33HZ9HPbe9JLR2ftZ+jLGIqiTbpEPDCjnXHjnS7gEdyXZDsBqlb+4+GhrwauoXqVp8vfOmpviGvbt0CdWvIF6S9rU9bVuHX9W073XG8uoLwNs0LCnL71tt5vf1OqPfe//yZhKWWyus7/aX+G1dnN/eLassObxwrvvJMhLdyhnElILxxkfXkIrwejnGnILxFED7l7Gu04vLL6Oif7qWXXn1DJ551le666Vx1j0xPst2HN/on59aWSJbVLsht4rygHBc+a23/vF2c295rl+qWXOHX0XsLCveCAj7/8wX6V1PAo53o+vp2KW6X4flyHL1fLzVEu9z1HbIciXZYOCbSIc1R+0ieC23q81qqb72yuVbl1da/7bP4BLuIZVRWk6bmQI2NUuNsKZ9tE9r/vJ7Xa68Hmjwl/Fzm0kvn9I2VpFVWbvuxzNLxfwmt1AtDeEsllkx7hDcZzuWMgvCWQy25PghvcqwrGQnhXQy9XC6vLXc/Sn++41L16N4muMcMv1x7/3Abbbvlt63CW0kh4+gb3baqQ47bdqg/E/Do9y0tgZqapahdc1P0c6Dmllzhn9FbmlX4LPoRtWtpaRPw6OdI2KNfRzuU1Xh1CPZ8SV5IstuOhnTsSLftRi8g2NEudiTP7TvVHYI9X7YXEPMOwY6us6mlnU9zvsBL+YWvfW6T9MYb0qxZkdwGbT83BgVOX/WKnoC2ysp5rbqytPIqefVcIvkzuaXWEOEtlVgy7RHeZDiXMwrCWw615PogvMmxrmQkhHcx9D6ePE0HHnO2HrztwvktL75mgpbs20uD9/th4b1Zc3zneCspZmftGwlgQYab2wSvpSXfLshBQQ4XfS+S6UisI3nu2BUv/Nyan//7qE/0XnNLXq2Fn6PzrJ2TUN++efXtE6h3r7x69w7Ur29QENzVVm47UtLZXr161Gv2vBbl0+/mnQ1tRfMNgkA9u2fUODelZ6wqurrO3bkuE6q+LtDcpmznvpAanX3vJfjmeWcoLcK7mCq9/d7HOvr0S/Wnm86d3/KqG/6gaOf36MF7doYaM8dFCLTtOrftTLdJc5tgz38/2o1uaZPnz9q1/b7wXkG429pEfaId68/6tmVFn8+d2zZw4YxzN6l74Ueg7t1VOOMc/b5wxrmbtER3qUePQH37qCC3HT/3XILyQQACEIAABCBQKQGEdzEEP5kyXT8ecmbhSEPH67wrx2uZpfrqZ/uzw1vpAoyrf49uGbW05tSaZRsxLsbl5rLDWy65ePuxwxsv30rS2eGthF78fdnhjZ+xYwSEdzEU8/m8vrv70bpv/AXq06ttu23IKRfrR7tuq+237l/4PQ+ecCxFb0ZXuEuDl1hyaZzhTY51KSNxhrcUWsm25QxvsrxLHY0zvKUSq057hLcI7r+8YJyW7tdXxwztiywlAAAPWUlEQVRuu0vD0cMu0723nq+e0b9DI7xFEEy+CcKbPPNiR0R4iyWVbDuEN1nepYyG8JZCK/m2CG/yzMsZEeEtgtrMxjk6/Zyxeu6l19SnV0+dccLB+u6mG8zvyQ5vERATboLwJgy8hOEQ3hJgJdgU4U0QdolDIbwlAku4OcKbMPAyh0N4ywS3YDeE1wDRHIHwmoEa4xBeI0xjFMJrhGmOQnjNQM1xCK8ZaExxCK8BLMJrgGiOQHjNQI1xCK8RpjEK4TXCNEchvGag5jiE1ww0pjiE1wAW4TVANEcgvGagxjiE1wjTGIXwGmGaoxBeM1BzHMJrBhpTHMJrAIvwGiCaIxBeM1BjHMJrhGmMQniNMM1RCK8ZqDkO4TUDjSkO4TWARXgNEM0RCK8ZqDEO4TXCNEYhvEaY5iiE1wzUHIfwmoHGFIfwGsAivAaI5giE1wzUGIfwGmEaoxBeI0xzFMJrBmqOQ3jNQGOKQ3gNYBFeA0RzBMJrBmqMQ3iNMI1RCK8RpjkK4TUDNcchvGagMcUhvAawCK8BojkC4TUDNcYhvEaYxiiE1wjTHIXwmoGa4xBeM9CY4hBeA1iE1wDRHIHwmoEa4xBeI0xjFMJrhGmOQnjNQM1xCK8ZaExxCK8BLMJrgGiOQHjNQI1xCK8RpjEK4TXCNEchvGag5jiE1ww0pjiE1wAW4TVANEcgvGagxjiE1wjTGIXwGmGaoxBeM1BzHMJrBhpTHMIbE1hiIQABCEAAAhCAAATSQQDhTUcdmAUEIAABCEAAAhCAQEwEEN6YwBILAQhAAAIQgAAEIJAOAghvOurALCAAAQhAAAIQgAAEYiKA8JYBNpfL69e/+a3uefhJ1ddndPhBu2u/PbYrI4kulRCYNmOWTjtnjF545Q39fzv3HuVzncdx/FXuFIkuIsllN2m6TKfsca921+WEJVFRxpQoM8ZlGc24hTEzLmEmMu6MoYS0DRGpzHZQUVtEam2ho1DZVuQya8/nY2fWrOQ334/f+Jzj+f3HOb/5vT+fj8d7+L1+3+/n+73yivIaMSBakRF1zhhy62df6rmMRdqxc7fKlC6lqE4t9Ui7+1ympvYcAqH2Jm8Y82+qa9xo1a5RVcP6R+EbRoFPtv9Dg1NnaP93B3VTrepKHdxDV1W64hdnnL9ktRYue1NHjx1Xi3vu1p97dgrjyhj6lddzlD5rqY4dO6HfN75TQ/o+pmLFLj0DZnpWtpauWKeTJ6Ubq1fRqPjHValieQCLWODwkZ8VPyrDzpqeFFfEszNdYQUIvIUVk7Rk+Tobdicn99HhI0f1aGySxgzuqXq/rRFgNEqCCsQnZajqtZUV0629tmzfqX7PTlH2vGSVLlWywJDmQ6TG9dfqjlvq2A/5jj2Ga/rYAap9Y9WgU1N3DoFQe5M3zJxFK5Wz8WNVv+5qAm8Yf7tyc/+tVl3iNbjPY2pcP0Im0K7ftFWTR/c5Y1bTk/c+3KakQU+oYoXLw7gqhjYCX+7+RtH9UjU/PVFXVa5og9StN9dUVMcWBYA+2vqFRk6Yp8z0RJUtU8oG5IP/PGTDMUfRCezd971iEibq9nq1te/ADwTeoqMPPBOBNwBdz/jxerRDczW86xZbnbn4DX2z73sNePqhAKNREkTAnBFs0KaX3lo8UWVKnwq4sYPT9ECrJmrW4PZfHTI2cZLatmhkz6BwnH+BwvbGfNAPSsqwZ943bv6UwHv+W5I/4ifbdio5PUsLpgyxr5leNW3fWyvmp+ryy8rmv88E4+aPDNDiac/qigqXhXFFDJ0nMHPhCv3r0GH16d7BvrT9i10aMmaWXp42vADSug1/0+Ll7yhtZG/7+qq339Padz9UamIPMItQ4KfDP+uzv+/WsePHlbVkNYG3CO2DTkXgDSDXsvNAzRw/UNddW9lW52z8RFlLV2tqar8Ao1ESRODb/T+oS2ySVr84Lr/cbFswH87RD7U665DHj59Qyy7xykxLUJVrKgWZmppzCBSmNyZwdeubokExj2jX1/u0YdNWAm8Yf8OyV6+3Z3TNWdu846GnRigx7lFF3HRj/ms7du5RYsoM1Y+sa8+8X16urPr1eFCREb8J4+ou7qGHjp1lff/UopGFMNtIGrTupU2rphWAMf+Hdeubqt9F3qyIujX1wtxltp+1anDF6kL8Bq3/YKsWvLKGwHsh8As5J4G3kGDm7c0e6KMlM0bk75l6/6Pt9rLSvLSEAKNREkTgqz3f2stJr81Lzi+fMmeZPWMVE93urEOmzVwi8838mdjOQaalJgSBwvRm9ouvy+yD69WtnVa9/T6BNwRfl7e8nP22tn2+S0NPu/wd1SdFT3Vtq/p31M0f2nyJjxuSpuSE7vpDk7u07fOvFJM4UcszU+1ldI7zLzBw5FTd2yjS7pXOO+o1i9KWt2brkksuKTDhmzmbNWLCXJkz8a3/2ED9e3ZU8WLFzv+iGPGcAgTecxJ58wYCb4BWtOwcr4wx/VS96jW22lxOWvSXtZqa2j/AaJQEEdh34KA69RxutzTkHamTF6rylRX0+MO/fIb3pVfXanXOJr2Q3FclShQPMi01IQiE2pudu/ZqcMoMzU1LUInixQi8Idi6viV7zXrlbPjY3qiWdzzwxFAN7x9lzxbmHeZDfMrcV5WZ/r8v8dF9UxXXvYNuu7mW6zKo/wWBYeNmK+Kmmupwf1P700M/HVGTdr21+Y3pBd797vtb9PzsV5Qxpr9KlyyhlMkLZc76jhwYjesFECDwXgD0gFMSeAPA9UqYqAfvb5a/V3TWiytkLuNy1jAAZsCSkydPqmGbGK1cOFbl/7v3sGf8c7Yv9zWOPGPUV1e9q8XZ79gvKmXLlA44K2WhCITam4zM1zRjQbYuvfTUXegnTuQqNzdXNa6vomWzR4UyFe8ppMCnO77UsHFz8veFnsjNVaO2sVq1YKwqlC+XP5rZXmLuVTB7e/OOLjFJSozrorp1bijkrLw9FAFzk+D+Awfz7wUxT9MYNnaWls4cWaB8+Lg5uq1eLbVr2di+fuzYcTVpH6cN2VNCmYb3nGcBAu95Bg3jcATeALhmH5y58z/vKQ2de42yj4W581b2twXgDFxi9rxVqlhBsdGnntIQkzhJr2eNUbmypWUulZsPBbMnztzUMX/JGrvH2vyMI/wCZ+uNeSzcqEmZerprW3s2/vSDLQ3h74vZ8tO66zP2y3mju089pWHtXzdr1oR47dm7X1lL1yi+18N2IebpM22aN7RfIj/c8rkGjsrQ8nnJKlmyRPgXehHO8PU3B9S192i7Ne7UUxqmqk7NanrqsbZ2H/Xeb79Txzb3aM5LK2W+uIxO6G63Mbzxzgea9/IqzX8+8SJUu/B/ZQLvhe9BqCsg8IYq9X/vGz91kZatzLF7q6I6tfjVG6UCTkHZOQR+PHRYCaOna/OWHSp/WTn7WJ68J2c8/cwE3dsw0l4ebNo+Tt/98KNO3wbX8K4IbjIM42/Y2XpjbsRpZbcE9T/jsXAE3jA25LShzZ3lCcnTbYCqecN1Skl8UtWqXGVDrXndnNU1/6+Zs7yJKdPtn1dXrmj3/Z6+7aFoVntxzbLizY0aP/UlHTl6VI3r36qRA6LtFwxzFdHcSJiS8KTdvjA6PUsbNn1qn9F7jelNv666odqpLXYcRStA4C1ab5fZCLwuetQigAACCCCAAAIIeC9A4PW+RSwQAQQQQAABBBBAwEWAwOuiRy0CCCCAAAIIIICA9wIEXu9bxAIRQAABBBBAAAEEXAQIvC561CKAAAIIIIAAAgh4L0Dg9b5FLBABBBBAAAEEEEDARYDA66JHLQIIIIAAAggggID3AgRe71vEAhFAAAEEEEAAAQRcBAi8LnrUIoAAAggggAACCHgvQOD1vkUsEAEEEEAAAQQQQMBFgMDrokctAggggAACCCCAgPcCBF7vW8QCEUAAAQQQQAABBFwECLwuetQigAACCCCAAAIIeC9A4PW+RSwQAQQQQAABBBBAwEWAwOuiRy0CCCCAAAIIIICA9wIEXu9bxAIRQAABBBBAAAEEXAQIvC561CKAAAIIIIAAAgh4L0Dg9b5FLBABBBBAAAEEEEDARYDA66JHLQIIIIAAAggggID3AgRe71vEAhFAAAEEEEAAAQRcBAi8LnrUIoAAAggggAACCHgvQOD1vkUsEAEEEEAAAQQQQMBFgMDrokctAggggAACCCCAgPcCBF7vW8QCEUAAAQQQQAABBFwECLwuetQigAACCCCAAAIIeC9A4PW+RSwQAQQQQAABBBBAwEWAwOuiRy0CCCCAAAIIIICA9wIEXu9bxAIRQAABBBBAAAEEXAQIvC561CKAAAIIIIAAAgh4L0Dg9b5FLBABBBBAAAEEEEDARYDA66JHLQIIIIAAAggggID3AgRe71vEAhFAAAEEEEAAAQRcBAi8LnrUIoAAAggggAACCHgvQOD1vkUsEAEEEEAAAQQQQMBFgMDrokctAggggAACCCCAgPcCBF7vW8QCEUAAAQQQQAABBFwECLwuetQigAACCCCAAAIIeC9A4PW+RSwQAQQQQAABBBBAwEWAwOuiRy0CCCCAAAIIIICA9wIEXu9bxAIRQAABBBBAAAEEXAQIvC561CKAAAIIIIAAAgh4L0Dg9b5FLBABBBBAAAEEEEDARYDA66JHLQIIIIAAAggggID3AgRe71vEAhFAAAEEEEAAAQRcBAi8LnrUIoAAAggggAACCHgvQOD1vkUsEAEEEEAAAQQQQMBFgMDrokctAggggAACCCCAgPcCBF7vW8QCEUAAAQQQQAABBFwECLwuetQigAACCCCAAAIIeC9A4PW+RSwQAQQQQAABBBBAwEWAwOuiRy0CCCCAAAIIIICA9wIEXu9bxAIRQAABBBBAAAEEXAQIvC561CKAAAIIIIAAAgh4L0Dg9b5FLBABBBBAAAEEEEDARYDA66JHLQIIIIAAAggggID3AgRe71vEAhFAAAEEEEAAAQRcBAi8LnrUIoAAAggggAACCHgvQOD1vkUsEAEEEEAAAQQQQMBFgMDrokctAggggAACCCCAgPcCBF7vW8QCEUAAAQQQQAABBFwECLwuetQigAACCCCAAAIIeC9A4PW+RSwQAQQQQAABBBBAwEWAwOuiRy0CCCCAAAIIIICA9wIEXu9bxAIRQAABBBBAAAEEXAQIvC561CKAAAIIIIAAAgh4L0Dg9b5FLBABBBBAAAEEEEDARYDA66JHLQIIIIAAAggggID3Av8B4hKuQBsEKqYAAAAASUVORK5CYII="
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig = go.Figure()\n",
    "fig.add_trace(go.Scatter(x=grid, y=values))\n",
    "fig.add_trace(go.Scatter(x=grid, y=dummy_values))\n",
    "fig.show(renderer=\"png\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "21c56b35",
   "metadata": {},
   "source": [
    "We can see that the dummy function is lower than the highest achieved value of `jennrich_sampson` but higher than the start values. It is also rather flat. Fortunately, that is all we need. "
   ]
  }
 ],
 "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.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}