neville.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAiMAAAGdCAYAAADAAnMpAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAABJxElEQVR4nO3deXgUVcL24d8pkkBYQgiECWGTHRd2FBdQwAVRxxEBhSjgwiYOgig47ygOM+q8AoIgoOz7IGAUg4iC8jKOKIIoIgrIIooggUASlpCQpc73Rw/5jKwNSVfSee7rykX6dHX3k5MmeVKnutpYay0iIiIiHnG8DiAiIiLFm8qIiIiIeEplRERERDylMiIiIiKeUhkRERERT6mMiIiIiKdURkRERMRTKiMiIiLiKZURERER8ZTKiIiIiHgqxOsA/khJSSE7Oztf7zM6OpqkpKR8vU85neY5cDTXgaF5DgzNc2AU1DyHhIRQoUKF82+X749cgLKzs8nKysq3+zPG5N6v3qKn4GieA0dzHRia58DQPAdGYZhnLdOIiIiIp1RGRERExFMqIyIiIuIplRERERHxVJE6gPVsrLUcP378og68SU9PJzMzswBSyW9pns+vZMmSlCxZ0usYIiIBFxRl5Pjx45QsWZKwsDC/bxsaGpqvr9CRM9M8n5u1lvT0dNLS0ihTpozXcUREAioolmmstRdVREQKC2MMpUuXzvfz6IiIFAVBUUZEgsWp1/uLiBQnKiMiIiLiKb/LSHJyMq+99hqPPPIIDzzwAE899RS7du06522+//57nnnmGeLi4hg4cCD//ve/LzZv0OjSpQvPP/+81zEuWdWqVfnwww8Lzf0EQrB870REAEI3bYL27X3/esSvA1iPHz/O8OHDufLKK/nrX/9KREQE+/fvP+cBdwcPHuTll1/m1ltvZeDAgXz33XdMnjyZyMhImjZteqn5i6xp06YRGhp6wdv/8ssvXHvttaxYsYKrrrqqAJMVrDFjxvDhhx/y0Ucf5RnfuHEj5cuX9yiVf4rr905EglP4W2/B6tWE165NZuPGnmTwq4wkJCRQsWJFBgwYkDtWuXLlc95m5cqVVK5cmZ49ewJQrVo1tm3bxvvvv1/oysimpE28uO5Fnmv1HE2imxToY13IGwcVlKysLL9+mQbC+Z5HhYmX3zsRkfxQYu9enORkMIZSS5cCUCohgRNdu4K1uFFR5FSrFrA8fi3TbNiwgdq1azN27Fh69+7NsGHD+Pjjj895mx07dtCoUaM8Y02aNGH79u1nvU1WVhYnTpzI/UhPT8+9zhhz2kd+id8Rz+f7P+ftHW/n232eze939bdq1YrXXnuNIUOGUL9+fa6++mrmz5+fe/21114LQIcOHahatSpdunTJvW7BggXcdNNN1K5dmxtvvJHZs2fnXvfLL79QtWpVEhIS6Ny5M7Vr1+add95h0aJFXH755Xz44YfccMMN1K5dm7i4OPbt25cn55w5c7j++uu57LLLaNOmDfHx8ef8ul566SVat25NnTp1uO666xg1alTuS3oXLVrE2LFj2bJlC1WrVqVq1aosWrQIyLtMcyrz8uXL6dKlC3Xq1OGWW25hw4YNeR7rX//6Fy1btqROnTo8+uijTJkyhcsvv/ys2X47F3fffTe1a9emffv2rF27Ns92a9eu5c4776RWrVo0a9aMf/7zn3le5ZKf37szOdNz3J+P/LgPfWieC8uH5rlgPv7QqhXRHTsSffvtOIcOAeAcPkz07bcT3bEjf2jVKt++fxfCrz0jBw8e5KOPPuLOO++kU6dO7Nq1i1mzZhESEkLbtm3PeJvU1NTTdr+XL18+9yRYZ3pJ7pIlS/L80qtVqxYjR44kOjr6jI+Rnp5+0X/p/3LsF5LTkzHGsPRHXztM+DGB7ld0x1pLVHgU1ctVv6j7PhdjDCVKlMjNbYxh6tSpPPPMMzz55JMsW7aM//mf/6FNmzbUrVuXFStW0KFDB+Lj42nQoAFhYWGEhoYSHx/PK6+8wv/+7//SqFEjNm/ezFNPPUW5cuXo1q0bISG+b/HLL7/MiBEjaNSoEaVKlWL16tWkp6czYcIEJk2aRFhYGM888wyPP/4477//PgDvv/8+f/vb33jhhRe46aabWLlyJUOGDKF69eq0bt0692v57dcRERHBhAkTiImJYevWrQwZMoSIiAgGDhxI586d2bFjB6tXr+att97K3f7UbU/dz6nMo0aNYsSIEdSuXZt//vOf/PnPf2bdunWEhISwbt06/vKXvzB8+HBuv/12PvnkE0aOHAlw1ufCqft98cUXeeGFF2jQoAGTJ0/moYceYsOGDURFRbF//3569uzJ/fffz6RJk9ixYwdPPfUU4eHhDBs2LF+/d2cSFhZGlSpVLuYplUdMTMwl34ecn+Y5MDTPBWD+fOxDD2GyszlVGcypE4eGhMDs2fnys+hC+VVGXNelTp06xMXFAb6SsGfPHj766KOzlpGL0alTJ+66667cy6faVVJS0hnPw5CZmXnRJ9RqMa/F/3+c/35LDqcf5pa3bskd39dn32m3u1TWWnJycnJzW2tp164dPXr0AKB///5MnjyZTz75hJo1a+YWunLlyhEVFQX49iCNGjWK559/ng4dOgAQGxvL1q1bmTNnDp07d86dr0cffTR3G/C9VXRWVhYvvPBC7nLZq6++yk033cT69etp1qwZkyZN4r777svN1Lt3b7788ksmTpxIq1atcu/rt1/HwIEDc8erVKlC//79SUhIYODAgYSEhBAeHo7jOLlfw6mv47f3cypzv379cp9XQ4YMoV27duzYsYO6desybdo02rVrR9++fQHo0aMH69ev5+OPPz7rc+HU/T700EPcfvvtgG9PzqpVq5g3bx4DBgxg+vTpVKlShRdeeAFjDJdddhlDhgzhn//8J4MGDcJxnHz73p1JZmYm+/fvP+N1F8IYQ0xMDImJiXrL9QKkeQ4MzXPBsVc2xPlTW2LePn11I2nZMrIbN4ZL+Fl0SkhIyFl3JOTZzp87rVChAtV+t4ZUrVo11q1bd9bbREZGcuTIkTxjR44cITw8/KwnKgsNDT3rX475/YR8/ZbXeWLVE2TbbCy++z71b4gJ4dW2r+br453LFVdckfu5MYbo6GgOHz581u1PnDjBTz/9xFNPPcXQoUNzx3NycihXrlyebZs0Of0YmJCQkDzH7dStW5fy5cuzY8cOmjVrxs6dO3nggQfy3Obqq69mxowZZ82UkJDAzJkz+fnnn0lLSyMnJ4eyZcuedftz+e2Sy6ljSg4dOkTdunXZtWsXHTt2zLN906ZNz7tsCNCixf8voCEhITRp0oQdO3YAsHPnTlq0aJFn9+LVV19NWloa+/fvp2rVqme8T3+/d+eSH89xa61+eAeA5jkwNM/5y/3yU+zciTgHfD+jrDEYa3P/hfz/XXs+fpWRBg0a8Ouvv+YZ+/XXX8/ZeurVq8fGjRvzjH377bfUr1/fn4cuMF3qd6F2udrcvuT2065bds8yGlVqdIZbFYxTywinGGNwXfes26elpQEwevRomjVrlue6EiVK5LkcHh6eTynPbsOGDQwcOJCnnnqKtm3bUq5cORISEpg6depF3d9v5+NUOTjXfHjJ3++diEig2cyT2EXTsf9ZAYDb4HJyfjhITtVqhD32GFlvvEGJX38lp2LFgGfz6wDWO++8kx07dvDOO++QmJjImjVrWLVqVZ7d/wsWLGDixIm5l2+77TYOHjzI/Pnz2bdvHytWrMg9QLCwObVMYyh8Z8E8tafot7/goqOjiYmJ4eeff6ZWrVp5PmrUqHHe+8zOzmbTb15XvnPnTo4cOUK9evUA356S3x80+uWXX+Ze/3sbNmygWrVqDBo0iCZNmlC7du3TDogNDQ3Nl1/SderU4Ztvvskz9vvLZ/P111/nfp6dnc23336b52v+6quv8vxV8OWXX1K2bNmLXj890/dORCSQbOJe3P8d6isixmDuvA87YhwH1n/J4eXLoV8/Di9fzoF163BjYwOez689I3Xr1uXpp59mwYIFvP3221SuXJlevXrRpk2b3G1SUlI49N8jc8G3e/0vf/kLc+bMYfny5VSsWJH+/fsXqpf1VixVkejwaGLLxNK9YXfe3PYmv6b9SsVSgW+HZ1OpUqXcA0+rVKlCyZIliYiI4KmnnmL48OFERETQtm1bMjMz+fbbb0lNTaVfv37nvM/Q0FCGDx/OP/7xD0JCQnj22Wdp3rx57l6Wxx57jP79+3PllVfSpk0bPvroIz744AMWLlx4xvs7VT4SEhJo0qQJq1at4oMPPsizTfXq1dmzZw/fffcdsbGxlClT5qLeqfaRRx7h3nvvZcqUKdx666189tlnrF69+oKO3p49eza1atWiXr16TJ06lSNHjtCtWzcAevXqxfTp03nuued4+OGH2bVrF2PGjKFv3744zsWdsPhs3zsRkUBwv1iNnf8GnMyAcuVxeg/BXPHfveklSsCpn5vGgEfvHO73u/a2aNEiz5r77z3++OOnjV155ZWMGjXK34cKmNiysazrvo4wJwxjDA82fJBMN5OSJQrP27mHhITwwgsv8Oqrr/LKK6/QqlUr4uPjiYuLIzw8nDfeeIMXX3yR0qVL07BhQ3r37n3e+wwPD2fAgAH8+c9/JjExkWuuuYYxY8bkXn/77bfz97//nSlTpvC3v/2N6tWrM3bsWK6//voz3t9tt91Gnz59ePbZZ8nMzOTmm29m8ODBjB07NnebO+64g+XLl3Pfffdx5MgRxo4dy/333+/3fFx99dW8/PLLjB07llGjRtG2bVv69OmT52XNZ/PXv/6VSZMm8f3333PZZZcxa9as3ANLq1Spwrx583jxxRe59dZbiYyMpHv37gwaNMjvjKec7XsnIlKQ7MmT2DcnYz9b5Rto0Ain91OYyKhz39ADxhaho4KSkpLO+CqEo0ePXvRfmsX1re0XLVrEiBEj2Lp1a0AeLxDzPHToUHbu3MmSJUvOeH1ROBPqpTyXwXesSpUqVdi/f78O+CtAmufA0DxfPLtvD+6UkbD/F9+yzF3dMHfdh3FKnLZtQc5zaGho/r+aRqQwmTx5Mm3atKF06dK55y755z//6XUsERHPWGuxn6/CLpgMmZlQvoJvb0hDb07zfqFURqTI2rhxI6+//jppaWnUqFGDf/zjH7nnwBERKW5sRjr2X5OxX6z2DVzRFOfRIZiISE9zXQiVkWLq/vvvv6hjNQqTKVOm+LV99erVT3t1j4hIMLB7d+NOGQWJ+8A4mD/FYTp2wVzkgfeBpjIiIiJSRFlrsZ+uwC6cDlmZEFkRp8/TmPpXeh3NLyojIiIiRZBNP4GdNwn75ae+gUYtcR4ejClX9E4doDIiIiJSxNg9u3zLMgf3Q4kSmE49MLfeU2SWZX4vaMqItdavtysWKWx0hlYROR9rLfbfy7GLZ0B2NkRF4/QdiqnT0OtolyQoykjJkiVJT0+ndOnSXkcRuSiu63Ls2DHKlCnjdRQRKaTsieO4cybC15/7Bppcg/PwIEyZcue+YREQNGUkLS2NI0eO+L13JCwsjMzMzAJKJqdons+vTJkyp73hnogIgN29A3fqKDh0AEqEYLr0wtx8d9CsCATNT76L+YtSZ/cLDM2ziMjFsdZiVy3Fxs+BnGyo9AecvsMwtc78hqVFVdCUERERkWBi047hzhoPm9b7Bppfj9Prz5jSZb0NVgBURkRERAoZu2sb7tTRkJwEISGY+x7FtL0jaJZlfk9lREREpJCwrov96F3sknmQkwOVq+D0G4apUcfraAVKZURERKQQsMeO4s4aB5s3AGCuboPp8TgmPPhfKaoyIiIi4jG7/Xvcaa9A6mEIDcN0641p0yFol2V+T2VERETEI9Z1sR/EY5cuANeFmKq+ZZlqtbyOFlAqIyIiIh6wR1NxZ4yFLd8AYK5th3mgP6ZUuLfBPKAyIiIiEmB227e408fAkRQIC8PE9cdcf3OxWZb5PZURERGRALFuDnbZYuyyRWBdqFIdp98zmKo1vI7mKZURERGRALCpyb5lmW3fAmBuuAXTvR+mZEmPk3lPZURERKSA2S0bcaePhWNHoGQpzAOP4VzXzutYhYbKiIiISAGxOTnYpW9iP3gLrIVql/neW6ZKNa+jFSoqIyIiIgXAphzGnTYadmwBwNx4O+b+RzFhWpb5PZURERGRfGY3f4U781U4fhRKhWN6PI5zzY1exyq0VEZERETyic3Oxr47H7viHd9Ajdq+k5hVjvU2WCGnMiIiIpIP7OEk37LMrm0AmHZ3Yro+jAkN8zhZ4acyIiIiconsN+twZ42HE8chvAxOr4GYFtd7HavIUBkRERG5SDY7C/v2XOzHCb6By+rh9B2KiY7xNlgRozIiIiJyEWxSou+ddndvB8Dccjemcy9MSKjHyYoelRERERE/2a8/x509AdLToHRZnIcHYZq28jpWkaUyIiIicoFsVhb2rZnY1e/7Bmo38C3LVKzsbbAiTmVERETkAtiDv+JOGQ17dgFgOtyLuedBTIh+lV4qzaCIiMh5uF+uwc6dABnpULYcziNPYhq19DpW0FAZEREROQubeRK7aAb2Px/6BupegdPnaUxUJW+DBRmVERERkTOwiXtxp4yCvT+BMZiOXTF3d8eUKOF1tKCjMiIiIvI77hersfPfgJMZUK48Tu8hmCuaeR0raKmMiIiI/Jc9eRL75hTsZx/7Bho0wun9FCYyyttgQU5lREREBLC/7vEty/y6x7csc9f9vg9HyzIFTWVERESKPfezVdgFb0BmJpSvgPPoEMzlTbyOVWyojIiISLFlM9KxCyZj1672DVzRFOfRJzERFbwNVsyojIiISLFk9/7kW5ZJ3AvGwfwpDtOxC8ZxvI5W7KiMiIhIsWKtxX66ErtwGmRlQmRF37lD6l/pdbRiS2VERESKDZt+AjtvEvbLT30DV7XwnU21XIS3wYo5lRERESkW7J5dvmWZg/vBcTD39sTceo+WZQoBlREREQlq1lrsvz/ALp4O2dkQVQmn7zBMnYZeR5P/UhkREZGgZU+k4c6dAF997htocg3Ow4MwZcp5G0zyUBkREZGgZHfvwJ06Cg4dgBIhmM69MLfcjTHG62jyOyojIiISVKy12FXvYeNnQ042VKyM028YplZ9r6PJWaiMiIhI0LBpx3BnvwbfrPMNNL8Op9dATOmy3gaTc1IZERGRoGB3bcOdOhqSkyAkBNP1EUy7O7UsUwSojIiISJFmXRf7UQJ2yVzIyYHoGJx+z2Bq1vE6mlwglRERESmy7LGjuLPGweYNAJir22B6PI4JL+1tMPGLyoiIiBRJdscW3GmvQMohCAnFdO+DadNByzJFkMqIiIgUKdZ1sR++jU34F7gu/KEqTv9hmGq1vI4mF0llREREigx7NBV3xquwZSMA5tq2mAcew5QK9ziZXAq/ysjixYuJj4/PMxYbG8u4cePOepv333+flStXcujQISIiImjVqhVxcXGEhYVdVGARESme7A+bcaeNgSPJEBaG6d4Pc8MtWpYJAn7vGalevTrDhw/Pveyc4w2G1qxZw4IFC3jssceoX78++/fv5/XXX8cYQ69evS4usYiIFCvWzcFdthj73kKwLlSp7nu1TNUaXkeTfOJ3GXEch8jIyAva9ocffqBBgwa0bt0agMqVK3PDDTewY8cOfx9WRESKoZzkQ7hjn8du+xYAc8PNvj0iJUt5nEzyk99lJDExkX79+hEaGkr9+vWJi4ujUqVKZ9y2QYMGfPrpp+zcuZO6dety4MABNm7cSJs2bS45uIiIBDd3yzckznwVm5oMJUthHngM57p2XseSAuBXGalXrx4DBgwgNjaWlJQU4uPjef755xkzZgzh4acfPNS6dWuOHj2au6yTk5PDrbfeyr333nvOx8nKyiIrKyv3sjEm9/7zc23w1H1pvbFgaZ4DR3MdGJrngmVzcnDfexP7/mKwFqpdRol+wzBVqnsdLSgVhuezsdbai71xWloaAwYMoFevXrRv3/6067///nvGjRtHt27dqFevHomJicyaNYubb76ZLl26nPV+f3+gbK1atRg5cuTFxhQRkSIi+9BBkkc/x8nvvgagzO2diOz7FI6WZYLaJb20t0yZMsTGxpKYmHjG6xctWsSNN97IzTffDECNGjXIyMhg6tSp3HvvvWc9+LVTp07cdddduZdPtbWkpCSys7MvJXIexhhiYmJITEzkEjqZnIfmOXA014GheS4Y7ndf4U4fC8ePQslwSvT6M1F/ul/zXMAK8vkcEhJCdHT0+be7lAfJyMggMTHxrMeAnDx58rTdPud69c0poaGhhIaGnvG6gnhCWmv1RA8AzXPgaK4DQ/OcP2x2NjbhX9gP3/YN1KiN03cYJqaq73rNc0B4Oc9+lZG5c+fSsmVLKlWqREpKCosXL8ZxnNxXy0ycOJGoqCji4uIAaNGiBe+//z61atXKXaZZtGgRLVq0uKBSIiIiwc0mJ/neaXfXNgBMuzt877YbqnNRFSd+lZHk5GTGjx/PsWPHiIiIoGHDhrz00ktEREQAcOjQoTx7Qjp37owxhoULF5KcnExERAQtWrSge/fu+ftViIhIkWM3rcedNR7SjkF4aZxeAzEtbvA6lnjgkg5gDbSkpKQ8r7K5VMYYqlSpwv79+7ULsABpngNHcx0YmudLY7OzsO/MxX6U4BuoWRen3zBMdEye7TTPgVGQ8xwaGlrwx4yIiIj4wyYl+t5pd/d2AMwtd2M698KEnPk4QSkeVEZERCQg7NdrcWe/BulpULoMzsODME2v9TqWFAIqIyIiUqBsVhY2fhb2/5b5Bmo3wOk7FFOxsrfBpNBQGRERkQJjD+73vVrm550AmA6dMPf0wITo14/8f3o2iIhIgXC/XIOdOwEy0qFsOZyHB2MaX+11LCmEVEZERCRf2axM7KLp2E8+9A3UvQKnz9OYqDO/qaqIyoiIiOQbm7gPd8oo2LsbjMF07IK5Ow5TooTX0aQQUxkREZF84X7xb+z81+FkBpQrj/PoEMyVzbyOJUWAyoiIiFwSe/IkduFU7JqPfAMNGuH0HoKJrOhtMCkyVEZEROSi2V/3+JZlft3jW5a5637fh6NlGblwKiMiInJR3M9WYRdMhsyTUL6Cb1nm8iZex5IiSGVERET8YjPSsQsmY9eu9g1c3sS3LBNRwdtgUmSpjIiIyAWze3/yLcsk7gXjYP4Uh+nYWcsycklURkRE5Lystdg1H2HfnApZmRAZ5Tt3SP2rvI4mQUBlREREzslmnMDOex27/j++gaua4zzyJKZceW+DSdBQGRERkbOye370Lcsc/BUcB9OpB+a2ThjH8TqaBBGVEREROY21FvvJB9hFMyA7C6Iq4fQZiql7udfRJAipjIiISB72RBp27kTsV5/5Bppcg/PQE5iyEd4Gk6ClMiIiIrnsTztwp46GpEQoUQLT+SHMLXdjjPE6mgQxlREREfEty6x6Dxs/G3KyoWJlnH7DMLXqex1NigGVERGRYs6mHced/Rp884VvoNm1vmWZ0mW9DSbFhsqIiEgxZn/8wbcsc/gghIRguj6CaXenlmUkoFRGRESKIeu62I8TsO/MhZwciI7xLcvUrOt1NCmGVEZERIoZe/wo7sxxsHkDAKZla0yPxzGly3gbTIotlRERkWLE7tyCO/UVSDkEIaGYbn0wN3bQsox4SmVERKQYsK6LXfEO9t354Lrwh6q+ZZnqtbyOJqIyIiIS7OzRVNyZr8L3GwEwrW7CPPgYplRpj5OJ+KiMiIgEMfvDd7jTXoEjyRAWhuneD3PDLVqWkUJFZUREJAhZNwe7/C3s0oVgXahS3bcsU7Wm19FETqMyIiISZOyRFNwZY2HrJgDM9Tdj4vphSpbyOJnImamMiIgEEbt1E+70MXA0FcJKYh54DOf69l7HEjknlRERkSBgc3KwyxZi318M1kLVmjj9nsFUqeZ1NJHzUhkRESnibOph30Gq278HwLS5zXf+kLCSHicTuTAqIyIiRZj97ivcGa/C8aNQMhzTYwBOq5u8jiXiF5UREZEiyObkYBPmYz942zdQvZZvWeYPsd4GE7kIKiMiIkWMTU7yLcvs3AqAaXsH5r5HMKFh3gYTuUgqIyIiRYjd9CXurHGQdgzCS+P0GohpcYPXsUQuicqIiEgRYLOzsEvmYVe+6xuoWdd3ErPoGE9zieQHlRERkULOHjqAO3U07N4OgLn5j5jOD2FCQz1OJpI/VEZERAoxu/EL3Nnj4UQalC6D89AgTLNrvY4lkq9URkRECiGblYV9ezZ21Xu+gVr1fcsyFSt7G0ykAKiMiIgUMvbgft+yzM87ATC3dcJ06oEJ0Y9sCU56ZouIFCJ2wxrcuRMh/QSUKYfzyGBM46u9jiVSoFRGREQKAZuViV08A/vvD3wDdS/H6fM0Jira22AiAaAyIiLiMZu4D3fKKNi7GwDTsQvmTw9gSpTwOJlIYKiMiIh4yF33CXbe63AyHcqVx3nkScxVzb2OJRJQKiMiIh6wJ09iF03DfrrSN9CgEU7vIZjIit4GE/GAyoiISIDZ/b/4lmX2/QzGYO68H/PH+zGOlmWkeFIZEREJIPfzVdh/TYbMkxARidP7KczlTbyOJeIplRERkQCwJzOw/5qMXft/voHLm+A8OgRTvoK3wUQKAZUREZECZvf9jDt5JCTuBeNg7u6OuaOLlmVE/ktlRESkgFhrsWs+wr45FbIyITIKp/fTmAZXeR1NpFBRGRERKQA24wR23hvY9Z/4Bq5q7nvZbrny3gYTKYRURkRE8pnd86PvvWUO7APHwdzTA9OhE8ZxvI4mUiipjIiI5BNrLfaTD7CLZkB2FlSohNP3aUzdK7yOJlKoqYyIiOQDeyINO28SdsMa30Djq3EeHoQpG+FtMJEiQGVEROQS2Z93+k5ilpQIJUpg7u2FufVPGGO8jiZSJPhVRhYvXkx8fHyesdjYWMaNG3fW26SlpfHmm2+yfv16jh8/TnR0NL169aJ5c733gogUbdZa3FXLsPEzITsbKlbG6TsUU7uB19FEihS/94xUr16d4cOH5152znFAVnZ2Ni+++CIREREMGTKEqKgoDh06ROnSpS8urYhIIeEeO4r7+v9iN671DTS9FuehJzBlynobTKQI8ruMOI5DZGTkBW37f//3fxw/fpwXXniBkBDfQ1WuXNnfhxQRKVTsjz+QOH0M9uB+CAnBdHkE0/5OLcuIXCS/y0hiYiL9+vUjNDSU+vXrExcXR6VKlc647VdffUW9evWYMWMGGzZsICIightuuIF77rnnnHtUsrKyyMrKyr1sjCE8PDz38/xy6r70A6RgaZ4DR3NdsKy12JXv4r4zB3JyIDqGEv2ewVxW1+toQUnP58AoDPNsrLX2QjfeuHEjGRkZxMbGkpKSQnx8PMnJyYwZMya3LPzW4MGDSUpKonXr1nTo0IHExESmT59Ox44d6dq161kf5/fHptSqVYuRI0f6+aWJiOSfnKOpJL/6dzLWfwpAeOtbiHriORwty4hcMr/KyO+lpaUxYMAAevXqRfv27U+7ftCgQWRmZjJp0qTcPSHLli1j6dKlTJ069az3e7Y9I0lJSWRnZ19s3NMYY4iJiSExMZFLmAY5D81z4GiuC4bdsYWcaaMh+RCEhOJ070Ps/Q9z4MABzXMB0vM5MApynkNCQoiOjj7/dpfyIGXKlCE2NpbExMQzXh8ZGUlISEieJZmqVauSmppKdnZ27nEkvxcaGkpoaOgZryuIJ6S1Vk/0ANA8B47mOn9Y18WueAf77nxwXagci9NvGE7NOhhjNM8BonkODC/n+ZLOTZyRkUFiYuJZD2ht0KABiYmJuK6bO7Z//34qVKhw1iIiIlIY2GNHcCf8A/vOXHBdzDU34Qwfi6lR2+toIkHHr0Ywd+5cWrZsSaVKlUhJSWHx4sU4jkPr1q0BmDhxIlFRUcTFxQFw2223sWLFCmbPns3tt99OYmIiS5YsoWPHjvn/lYiI5BO7/Tvcaa9AajKEhmG698W0vlUHUooUEL/KSHJyMuPHj+fYsWNERETQsGFDXnrpJSIifKc7PnToUJ7/rJUqVeLZZ59lzpw5DB06lKioKDp27Mg999yTr1+EiEh+sG4Odnk8dumbYF2oUh2n3zBM1ZpeRxMJapd0AGugJSUl5Tmw9VIZY6hSpQr79+/XemQB0jwHjub64tmjKbjTx8LWTQCY69pjHuiPKVnqtG01z4GheQ6Mgpzn0NDQgj+AVUQkGNitm3Cnj4GjqRBWEvNAf5zrb/Y6lkixoTIiIsWWdXOw7y3Cvr8IrIWqNX3LMlWqex1NpFhRGRGRYsmmHsadNga2fweAaXMb5v4+mJIlPU4mUvyojIhIsWO/+xp35qtw7AiUDMf0GIDT6iavY4kUWyojIlJs2JwcbMJ87Adv+waq1fIty8RU9TaYSDGnMiIixYJNTvKdO2TnVgBM2zsw9z2CCQ3zNpiIqIyISPCz336JO3McpB2D8NI4Pf+Madna61gi8l8qIyIStGx2NnbJPOzKJb6BmnVx+g7FVK7ibTARyUNlRESCkj18EHfKKNi9HQBz8x8xnR/CnOVNOEXEOyojIhJ07MYvcGePhxNpULoMzkODMM2u9TqWiJyFyoiIBA2blYV9ezZ21Xu+gVr1fcsylf7gbTAROSeVEREJCjYp0bcs8/NOAMxt92A69cCEaFlGpLBTGRGRIs9+9RnunAmQfgLKlMN5eDCmydVexxKRC6QyIiJFls3KxC6eif33ct9A3ctx+jyNiTr/u4SKSOGhMiIiRZI98CvulJHwy24ATMcumLvjMCH6sSZS1Oh/rYgUOe66T7DzXoeT6VA2AufRIZirmnsdS0QuksqIiBQZNvMkduE07KcrfQP1r8Lp8xQmsqK3wUTkkqiMiEiRYPfv9S3L7PsZjMHceR/mrm6YEiW8jiYil0hlREQKPffz/8P+6w3IPAkRkb5lmSuaeh1LRPKJyoiIFFr2ZAZ2wRTs56t8A5c38RWR8hW8DSYi+UplREQKJbtvj29ZZv8vYBzM3d0wd3TFOFqWEQk2KiMiUqhYa7FrPsIunAqZmVA+ynfukAZXeR1NRAqI43UAESneNiVtouuyrmxK2oTNOIGdMRY7d6KviFzZDOf5cSoiIkFOe0ZExFPxO+L5fP/nxG+cyVX/yYQD+8BxMPc8iOlwL8bR30wiwU5lREQCbu+xvSRnJGOMIWFXAgAJO9+l87F62JjyVLrvz1Rv1M7jlCISKCojIhJwrRa2yv3cWMBAckg2d1611Tf4xRfsa7TPm3AiEnDa/ykiATeh3QRCjO9VMdaQ598QE8KEdhM8SiYiXlAZEZGAstZyz89hJGy5/IzXL7tnGffWvTfAqUTES1qmEZGAsSeO486ZAF+vhdI5ABgMFpv7r4gUP9ozIiIBYXdvx/3HYF8RKRFC9B09iQ6PpnGlxrzc+mUaV2pMdHg0FUvpTe9EihvtGRGRAmWtxX6UgH1nDuTkQHQMTt+hVL2sHuty7ifMCcMYw4MNHyTTzaRkiZJeRxaRAFMZEZECY9OO4c4aD5vW+wZaXI/TcyCmdBmAPMXDGKMiIlJMqYyISIGwO7fiThsNyYcgJBRz/6OYmzpijPE6mogUMiojIpKvrOtiVyzBvjsPXBcqx+L0G4apUdvraCJSSKmMiEi+sceO4M4cB999BYC55iZMj8cwpUp7G0xECjWVERHJF3b7d7jTXoHUZAgNw3Tvi2l9q5ZlROS8VEZE5JJYNwe7PB679E2wLsRU8y3LVLvM62giUkSojIjIRbNHU3Cnj4WtmwAw17XHPNAfU7KUx8lEpChRGRGRi2K3bsKdMRaOpEBYScwD/XGuv9nrWCJSBKmMiIhfrJuDfW8R9v1FYC3E1vAty8TW8DqaiBRRKiMicsFs6mHfsswPmwEwbW7D3N8HU1InKxORi6cyIiIXxH6/0bcsc+wIlAzH9BiA0+omr2OJSBBQGRGRc7I5OdiEf2E/iPcNVKvlW5aJqeptMBEJGiojInJWNvmQ79whO7cAYNp2xNz3KCY0zNtgIhJUVEZE5Izs5g24M1+F48egVDim50Ccq1t7HUtEgpDKiIjkYbOzse/Ow65Y4huoWRen71BM5SreBhORoKUyIiK57OGDuFNHw48/AGDa34Xp8jAmNNTjZCISzFRGRAQA+80XuLNegxPHoXQZnF5PYJpf53UsESkGVEZEijmbnYWNn41d9Z5voFZ937JMpT94G0xEig2VEZFizCYl4k4ZBT/vBMDcdg+mUw9MiJZlRCRwVEZEiin71We4cyZA+gkoUw7n4UGYJtd4HUtEiiGVEZFixmZlYhfPxP57uW+gTkOcPkMxFaO9DSYixZbKiEgxYg/8ijt1FOz5EQBze2fMnx7AhOhHgYh4Rz+BRIoJd/1/sHMnwcl0KBuB8+iTmKtaeB1LRERlRCTY2cyT2EXTsf9Z4RuofyVO76cxFSp6G0xE5L9URkSCmN2/F3fKSNj3MxiDuaMr5o/dMSVKeB1NRCSXX2Vk8eLFxMfH5xmLjY1l3Lhx573tZ599xvjx42nZsiXDhg3zK6SI+M9duxr7rzfgZAZEROI8OgRzRVOvY4mInMbvPSPVq1dn+PDhuZcdxznvbQ4ePMi8efO4/PLL/X04EfGTm5FOzqzx2M8+9g00bIzT+ylM+QreBhMROQu/y4jjOERGRl7w9q7rMmHCBO677z62bt1KWlqavw8pIhfI7vuZA/8YhN3zIxgH88dumDu7Yhwty4hI4eV3GUlMTKRfv36EhoZSv3594uLiqFSp0lm3j4+PJyIigvbt27N169YLeoysrCyysrJyLxtjCA8Pz/08v5y6r/y8Tzmd5rngWWuxn32Mu2AKZJ6E8lE4fZ7GadjI62hBSc/pwNA8B0ZhmGdjrbUXuvHGjRvJyMggNjaWlJQU4uPjSU5OZsyYMbll4be2bdvGuHHjGDVqFBEREUyaNIm0tLTzHjPy+2NTatWqxciRI/34skSKDzf9BCkT/5cT//4AgFLNryXqqX9QIjLK42QiIhfGrz0jzZo1y/28Zs2a1KtXjwEDBrB27Vrat2+fZ9v09HQmTJhAv379iIiI8CtUp06duOuuu3Ivn2prSUlJZGdn+3Vf52KMISYmhsTERPzoZOInzXPBsb/sJmfySDiwDxwH554HqfTwnzlw8CB2/36v4wUtPacDQ/McGAU5zyEhIURHn//szpf00t4yZcoQGxtLYmLiadcdOHCApKSkPHs0Tn2R3bp1Y9y4ccTExJzxfkNDQwkNPfMbdRXEE9Jaqyd6AGie84+1FvufFdiF0yA7CyIr4vQdilP/SozjaK4DRPMcGJrnwPByni+pjGRkZJCYmEibNm1Ouy42NpZXXnklz9jChQvJyMjgoYceOudxJiJydjb9BHbeJOyXn/oGGrXEeXgwppx/eyBFRAoLv8rI3LlzadmyJZUqVSIlJYXFixfjOA6tW7cGYOLEiURFRREXF0dYWBg1atTIc/syZcoAnDYuIhfG/rzLdxKzpEQoUQLTqSfm1j9hLuAl9iIihZVfZSQ5OZnx48dz7NgxIiIiaNiwIS+99FLuMSGHDh3SUc8iBcBai139PvatmZCdDRUr4/R5GlOnodfRREQumV9lZPDgwee8fsSIEee8/vHHH/fn4UQEsCeO486ZAF+v9Q00bYXz0CBMmbLeBhMRySd6bxqRQszu3o47ZRQcPgglQjBdHsLc/EftgRSRoKIyIlIIWWuxHy/Fvj0HcrKh0h9w+g7D1KrndTQRkXynMiJSyNi0Y7izxsOm9b6BFtfj9ByIKV3G22AiIgVEZUSkELE7t+JOGw3JhyAkBHNfb0zbjlqWEZGgpjIiUghY18WuXIJdMg9cFypXwek3DFOjjtfRREQKnMqIiMfssSO4M8fBd18BYK65EdNjAKZUaW+DiYgEiMqIiIfs9u99yzKpyRAahunWB9PmNi3LiEixojIi4gHrutgP4rEJC8C6EFPNtyxT7TKvo4mIBJzKiEiA2aMpuDNehS3fAGCua4eJ648pFe5tMBERj6iMiASQ3boJd8ZYOJICYSUxcf1xbrjZ61giIp5SGREJAOvmYJctwi5bBNZCbA3fskys3jRSRERlRKSA2dRk3Olj4IfNAJjWt2K69cWULOlxMhGRwkFlRKQA2e83+pZljh2BkqUwDw7Aubat17FERAoVlRGRAmBzcrBLF2A/iPcty1S7zLcsE1PN62giIoWOyohIPrPJh3CnvwI7tgBgbrodc9+jmDAty4iInInKiEg+sps34M58FY4fg1LhmJ5/xrm6jdexREQKNZURkXxgs7Ox787DrljiG6hRB6ffUEzlWG+DiYgUASojIpfIHk7yndJ91zYATPu7MF0exoSGepxMRKRoUBkRuQT2m3W4s8bDieMQXgbnoYGY5td7HUtEpEhRGRG5CDY7C/v2HOzHS30Dl9XD6TsUEx3jbTARkSJIZUTETzYpEXfqaPhpBwDm1j9h7u2JCdGyjIjIxVAZEfGD/fpz3NkTID0NSpfFeWQwpsk1XscSESnSVEZELoDNysS+NQu7+n3fQJ2GOH2GYipGextMRCQIqIyInIc9+CvulFGw50cATId7Mfc8iAnRfx8Rkfygn6Yi5+Cu/w923iTISIeyETiPPIlp1MLrWCIiQUVlROQMbOZJ7KLp2P+s8A3Uu8K3LFOhorfBRESCkMqIyO/YxL2+ZZm9P4ExmDu6Yv7YHVOihNfRRESCksqIyG+4X6zGzn8DTmZAufI4vYdgrmjmdSwRkaCmMiIC2JMnsW9Oxn62yjfQoBFO76cwkVHeBhMRKQZURqTYs/v24E4ZCft/AeNg/tgNc2dXjKNlGRGRQFAZkWLLWov9fBV2wWTIzITyFXD6PI1p0MjraCIixYrKiBRLNiMd+6/J2C9W+wauaIbz6JOYiEhPc4mIFEcqI1Ls2L27fa+WSdznW5b5UxymYxeM43gdTUSkWFIZkWLDWov9dAX2zWmQnQWRFX3LMvWv9DqaiEixpjIixYJNP4GdNwn75ae+gUYtcR4ejCkX4W0wERFRGZHgZ/fs8i3LHNwPJUpgOvXA3HqPlmVERAoJlREJWtZa7L+XYxfPgOxsiIrG6TsUU6eh19FEROQ3VEYkKNkTx3HnTISvP/cNNLkG5+FBmDLlvA0mIiKnURmRoGN378CdOgoOHYASIZguD2Fu/iPGGK+jiYjIGaiMSNCw1mJXLcXGz4GcbKj0B5y+wzC16nkdTUREzkFlRIKCTTuGO2s8bFrvG2h+PU6vP2NKl/U2mIiInJfKiBR5dtc23KmjITkJQkIw9z2KaXuHlmVERIoIlREpsqzrYj96F7tkHuTkQOUqOP2GYWrU8TqaiIj4QWVEiiR77CjurHGweQMA5uo2mB6PY8JLextMRET8pjIiRY7d/j3utFcg9TCEhmG69ca06aBlGRGRIkplRIoM67rYD+KxSxeA60JMVd+yTLVaXkcTEZFLoDIiRYI9moo7Yyxs+QYAc207zAP9MaXCvQ0mIiKXTGVECj277Vvc6WPgSAqEhWHi+mOuv1nLMiIiQUJlRAot6+Zgly3GLlsE1oUq1XH6P4OJreF1NBERyUcqI1Io2dRk37LMtm8BMDfcguneD1OypMfJREQkv6mMSKFjt2zEnT4Wjh2BkqUwDz6Gc207r2OJiEgBURmRQsPm5GCXvon94C2wFqpd5ntvmSrVvI4mIiIFSGVECgWbfAh3+iuwYwsA5sbbMfc/ignTsoyISLBTGRHP2c1f4c4cC8ePQalwTI/Hca650etYIiISICoj4hmbnY19dz52xTu+gRp1cPoNxVSO9TaYiIgElMqIeMIeTsKdNhp2bQPAtLsT0/URTGiox8lERCTQ/CojixcvJj4+Ps9YbGws48aNO+P2H3/8Mf/5z3/45ZdfAKhduzbdu3enbt26F5dWgoL9Zh3urPFw4jiEl8HpNRDT4nqvY4mIiEf83jNSvXp1hg8fnnvZcZyzbrtlyxZuuOEGGjRoQGhoKAkJCbz44ouMHTuWqKioi0ssRZbNzsKNn4P9OME3cFk9nL5DMdEx3gYTERFP+V1GHMchMjLygrZ94okn8lzu378/69atY/Pmzdx0003+PrQUYdmJ+8gZ+RfYvR0Ac8ufMJ17YkK0LCMiUtz5XUYSExPp168foaGh1K9fn7i4OCpVqnRBtz158iTZ2dmULVv2nNtlZWWRlZWVe9kYQ3h4eO7n+eXUfek9TgqW/XotiXNeg7TjULosziODcZq28jpWUNJzOjA0z4GheQ6MwjDPxlprL3TjjRs3kpGRQWxsLCkpKcTHx5OcnMyYMWNyy8K5TJ8+nU2bNjFmzBjCwsLOut3vj02pVasWI0eOvNCYUkjYrExSp4/j+LLFAIQ1bEzFZ14ipHIVj5OJiEhh4lcZ+b20tDQGDBhAr169aN++/Tm3fffdd0lISGDEiBHUrFnznNuebc9IUlIS2dnZFxv3NMYYYmJiSExM5BKmQc7AHviVnCmjYM8uAMp17kn6bfdCiRIeJwtuek4HhuY5MDTPgVGQ8xwSEkJ0dPT5t7uUBylTpgyxsbEkJiaec7ulS5fy7rvvMnz48PMWEYDQ0FBCz/ISz4J4Qlpr9UTPR+6Xa7BzJ0BGOpQth/PIECI7/JH0/fs1zwGi53RgaJ4DQ/McGF7O8yWVkYyMDBITE2nTps1Zt0lISOCdd97h2WefpU6dOpfycFLI2cyT2EUzsP/50DdQ9wqcPk/jVDx/KxYRkeLLrzIyd+5cWrZsSaVKlUhJSWHx4sU4jkPr1q0BmDhxIlFRUcTFxQG+pZnFixfzxBNPULlyZVJTUwEoVaoUpUqVyt+vRDxlE/fiThkFe38CYzAdu2Lu7o7RsoyIiJyHX2UkOTmZ8ePHc+zYMSIiImjYsCEvvfQSERERABw6dCjP0bgfffQR2dnZjB07Ns/9dOnShfvuuy8f4kth4H6xGjv/DTiZAeXK4/QegrmimdexRESkiPCrjAwePPic148YMSLP5UmTJvmbR4oQe/Ik9s0p2M8+9g00aITT+ylMpE5oJyIiF07vTSMXxf66x7cs8+se37LMXd0wd92HcbQsIyIi/lEZEb+5n63CLngDMjOhfAXf3pCGjb2OJSIiRZTKiFwwm5GOXTAZu3a1b+CKpjiPPomJqOBtMBERKdJURuSC2L0/+ZZlEveCcTB/isN07II5xxslioiIXAiVETknay3205XYhdMgKxMiK+L0eRpT/0qvo4mISJBQGZGzsuknsPMmYb/81DdwVQucR57ElIvwNpiIiAQVlRE5I7tnl29Z5uB+cBzMvT0xt96jZRkREcl3KiOSh7UW++8PsIunQ3Y2REXj9B2KqdPQ62giIhKkVEYklz2Rhjt3Anz1uW+gyTU4Dw/ClCnnbTAREQlqKiMCgN29A3fqKDh0AEqEYLr0wtx8d57T+4uIiBQElZFizlqLXfUeNn425GRDxco4/Z7B1KrndTQRESkmVEaKMZt2DHf2a/DNOt9A8+tweg3ElC7rbTARESlWVEaKKbtrG+7U0ZCcBCEhmPsexbS9Q8syIiIScCojxYx1XexHCdglcyEnB6JjfMsyNet4HU1ERIoplZFixB47ijtrHGzeAIC5ug2mx+OY8NLeBhMRkWJNZaSYsDu24E57BVIOQUgopnsfTJsOWpYRERHPqYwEOeu62A/fxib8C1wXYqri9BuGqVbL62giIiKAykhQs0dTcWe8Cls2AmCubYt54DFMqXCPk4mIiPx/KiNByv6wGXfaGDiSDGFhmLj+mOtv1rKMiIgUOiojQca6Odj338K+txCsC1Wq+14tU7WG19FERETOSGUkiNgjKbjTx8C2bwEwN9yM6d4PU7KUx8lERETOTmUkSNgt3/iKyLEjULIU5oHHcK5r53UsERGR81IZKeJsTg72vTexy98Ca6FqTd+yTJVqXkcTERG5ICojRZhNOYw7/RXY/j0A5sYOmPt7Y8JKepxMRETkwqmMFFH2u698L9s9fhRKhmN6Po5zzY1exxIREfGbykgRY7OzsQn/wn74tm+gRm2cvsMwf4j1NpiIiMhFUhkpQuzhJNxpo2HXNgBMuzswXR/BhIZ5nExEROTiqYwUEXbTetxZ4yHtGISXwek1ENPieq9jiYiIXDKVkULOZmdh35mL/SjBN3BZPZy+QzHRMd4GExERyScqI4WYTUr0vdPu7u0AmFvuxnTuhQkJ9TiZiIhI/lEZKaTs12txZ78G6WlQugzOw4MwTa/1OpaIiEi+UxkpZGxWFjZ+Fvb/lvkGajfwLctUrOxtMBERkQKiMlKI2IP7caeOhp93AmA6dMLc0wMTom+TiIgEL/2WKyTcL9dg506AjHQoWw7n4cGYxld7HUtERKTAqYx4zGZlYhdNx37yoW+g7hU4fZ7GRFXyNpiIiEiAqIx4yCbuw50yCvbuBmMwHbtg7o7DlCjhdTQREZGAURnxiPvFv7HzX4eTGVCuPM6jQzBXNvM6loiISMCpjASYPXkSu3Aqds1HvoEGjXB6P4WJjPI2mIiIiEdURgLI/rrHtyzz6x7fssxd9/s+HC3LiIhI8aUyEiDuZ6uwCyZD5kkoX8G3LHN5E69jiYiIeE5lpIDZjHTsgsnYtat9A5c3wek9BBNRwdtgIiIihYTKSAGye3/yLcsk7gXjYP4U53vFjON4HU1ERKTQUBkpANZa7JqPsG9OhaxMiIzynTuk/lVeRxMRESl0VEbymc04gZ33Onb9f3wDVzXHeeRJTLny3gYTEREppFRG8pHd86NvWebgr+A4mE49MLd10rKMiIjIOaiM5ANrLfbfH2AXz4DsLIiqhNNnKKbu5V5HExERKfRURi6RPZGGO3cCfPW5b6DJNTgPPYEpG+FtMBERkSJCZeQS2J924E4dDUmJUCIE07kX5pa7McZ4HU1ERKTIUBm5CNZa7Kr3sPGzIScbKlbG6TcMU6u+19FERESKHJURP9m047izX4NvvvANNL8Op9dATOmy3gYTEREpolRG/GB//MG3LHP4IISEYLo+gml3p5ZlRERELoHKyAWwrov9OAH7zlzIyYHoGN+yTM26XkcTEREp8lRGzsMeP4o7cxxs3gCAadka0/PPmPDS3gYTEREJEioj52B3bsGd+gqkHIKQUEy3PpgbO2hZRkREJB+pjJyBdV3sinew784H14U/VPUty1Sv5XU0ERGRoKMy8jv2aCruzFfh+40AmFY3YR58DFNKyzIiIiIFwa8ysnjxYuLj4/OMxcbGMm7cuLPeZu3atSxatIikpCRiYmJ44IEHaN68+UWFzW+hmzbBAw8QOmwYmY0bY3/4DnfaK3AkGcLCMN37YW64RcsyIiIiBcjvPSPVq1dn+PDhuZedc7wJ3A8//MD48eOJi4ujefPmrFmzhtGjRzNy5Ehq1KhxcYnzUfhbb8Hq1ZSqVYuMn7dg31sI1oUq1XH6PYOp6n1GERGRYOd3GXEch8jIyAvadvny5TRt2pS7774bgG7durF582Y+/PBD+vbt6+9D54sSe/fiJCeDMZRauhSA8IULOP5DLcBir7sJd8AwTMlSnuQTEREpbvwuI4mJifTr14/Q0FDq169PXFwclSpVOuO227dv56677soz1qRJE7788stzPkZWVhZZWVm5l40xhIeH535+Kf7QqlXu5/a//zonMoj5bKvvwmfb2D/kb5f0GJLXqe+ZlrsKnuY6MDTPgaF5DozCMM9+lZF69eoxYMAAYmNjSUlJIT4+nueff54xY8bkloXfSk1NpXz58nnGypcvT2pq6jkfZ8mSJXmOTalVqxYjR44kOjran7hnNn8+PPQQZGdzatpzpz8kBGbPpkqVKpf+OHKamJgYryMUG5rrwNA8B4bmOTC8nGe/ykizZs1yP69Zs2ZuOVm7di3t27fPt1CdOnXKs0flVFtLSkoiOzv70u68fXtCli0j+vbbT7sqadkyshs3hv37L+0xJA9jDDExMSQmJmKtPf8N5KJprgND8xwYmufAKMh5DgkJuaAdCZf00t4yZcoQGxtLYmLiGa+PjIzkyJEjecaOHDly3mNOQkNDCQ0NPeN1+TlR1hiMtbn/5vf9S17WWs1vgGiuA0PzHBia58Dwcp7P/lKYC5CRkUFiYuJZy0X9+vXZvHlznrFvv/2WevXqXcrDXrKcihXJiY4mq3FjmDyZrMaNyYmOJqdiRU9ziYiIFEd+7RmZO3cuLVu2pFKlSqSkpLB48WIcx6F169YATJw4kaioKOLi4gC44447GDFiBO+99x7Nmzfns88+Y9euXZ69kuYUNzaWA+vWYUqWpEpsLIf/+EfsyZNQsqSnuURERIojv8pIcnIy48eP59ixY0RERNCwYUNeeuklIiIiADh06FCeo3EbNGjAE088wcKFC3nzzTepUqUKQ4cOLRTnGKFkSTiV1RgVEREREY/4VUYGDx58zutHjBhx2th1113Hdddd58/DiIiISDFySceMiIiIiFwqlRERERHxlMqIiIiIeEplRERERDylMiIiIiKeUhkRERERT6mMiIiIiKdURkRERMRTKiMiIiLiqUt6195ACwkpmLgFdb+Sl+Y5cDTXgaF5DgzNc2AUxDxf6H0aq/dlFhEREQ8V62Wa9PR0nnnmGdLT072OEtQ0z4GjuQ4MzXNgaJ4DozDMc7EuI9Zadu/ejXYOFSzNc+BorgND8xwYmufAKAzzXKzLiIiIiHhPZUREREQ8VazLSGhoKF26dCE0NNTrKEFN8xw4muvA0DwHhuY5MArDPOvVNCIiIuKpYr1nRERERLynMiIiIiKeUhkRERERT6mMiIiIiKeC+oT/W7ZsYenSpezevZuUlBSefvpprrnmmnPe5vvvv2fu3Ln88ssvVKxYkc6dO9O2bdvABC6i/J3ndevWsXLlSn766Seys7OpVq0aXbt2pWnTpoELXQRdzPP5lG3btjFixAiqV6/O6NGjCzhp0XYx85yVlUV8fDyffvopqampVKhQgc6dO9O+ffsApS56LmaeP/30U5YuXcr+/fspXbo0TZs2pUePHpQrVy5AqYueJUuWsH79evbt20dYWBj169fnwQcfJDY29py3W7t2LYsWLSIpKYmYmBgeeOABmjdvXmA5g3rPyMmTJ7nssst49NFHL2j7gwcP8vLLL3PllVcyatQo7rzzTiZPnsw333xTsEGLOH/neevWrTRu3Jj/+Z//yZ3vkSNHsnv37gJOWrT5O8+npKWlMWnSJBo1alRAyYLLxczzq6++ynfffUf//v0ZN24cgwYNOu8P++LO33netm0bEydOpF27dowdO5YhQ4awa9cupkyZUsBJi7YtW7bQoUMHXnrpJZ577jlycnJ48cUXycjIOOttfvjhB8aPH0/79u0ZOXIkV199NaNHj2bPnj0FljOo94w0a9aMZs2aXfD2K1eupHLlyvTs2ROAatWqsW3bNt5//3391X4O/s7zQw89lOdyXFwcGzZs4KuvvqJWrVr5nC54+DvPp0ybNo0bbrgBx3H48ssvCyBZcPF3nr/55hu2bNnCxIkTKVu2LACVK1cuqHhBw9953r59O5UrV+aOO+4AfHN8yy23kJCQUFARg8Kzzz6b5/Ljjz9O7969+fHHH7niiivOeJvly5fTtGlT7r77bgC6devG5s2b+fDDD+nbt2+B5AzqPSP+2rFjx2l/PTZp0oTt27d7lKh4cF2X9PT03B/kkn9Wr17NgQMH6Nq1q9dRgtaGDRuoU6cOCQkJ9OvXj0GDBjF37lwyMzO9jhZU6tevz6FDh/j666+x1pKamsoXX3xxUQW9ODtx4gTAOX/ebt++/Yy/C3fs2FFguYJ6z4i/UlNTKV++fJ6x8uXLk56eTmZmJmFhYR4lC27vvfceGRkZXHfddV5HCSr79+9nwYIF/P3vf6dEiRJexwlaBw4cYNu2bYSGhjJ06FCOHj3KjBkzOH78OAMGDPA6XtBo2LAhTzzxBOPGjSMrK4ucnBxatGjh97Jlcea6LrNnz6ZBgwbUqFHjrNud7XdhampqgWXTnhHx1Jo1a4iPj+fJJ5887ckvF891XV577TW6du2qYxcK2KmTWD/xxBPUrVuX5s2b07NnTz755BPtHclHe/fuZfbs2XTp0oWXX36Zv/71ryQlJTFt2jSvoxUZM2bM4JdffmHw4MFeRzmN9oz8RmRkJEeOHMkzduTIEcLDw7VXpAB89tlnTJ48mSFDhtC4cWOv4wSV9PR0du3axe7du5k5cybg+6VpraVbt24899xzXHXVVR6nDA6RkZFERUVRunTp3LGqVatireXw4cNUqVLFw3TBY8mSJTRo0CD3OIaaNWtSqlQpnn/+ebp160aFChU8Tli4zZgxg6+//pq///3vVKxY8Zzbnu13YWRkZIHlUxn5jXr16rFx48Y8Y99++y3169f3KFHwWrNmDW+88QaDBw8u0JeLFVfh4eG88sorecZWrlzJd999x5AhQ3SAZT5q2LAhX3zxBRkZGZQqVQrwLZEZY877Q18u3MmTJ09bbnQc3859vcXa2VlrmTlzJuvXr2fEiBEX9H+/fv36bN68mTvvvDN37Ntvv6VevXoFljOol2kyMjL46aef+OmnnwDfS3d/+uknDh06BMCCBQuYOHFi7va33XYbBw8eZP78+ezbt48VK1awdu3aPN8QOZ2/87xmzRomTZpEz549qVevHqmpqaSmpuYeWCVn5s88O45DjRo18nxEREQQGhpKjRo1cn9pyun8fT63bt2acuXK8frrr7N37162bNnC/PnzadeunfaonoO/89yyZUvWr1/PypUrc4/TmTVrFnXr1iUqKsqLL6FImDFjBp9++imDBg0iPDw89+ftb5cQJ06cyIIFC3Iv33HHHWzatIn33nuPffv2sXjxYnbt2sXtt99eYDmDes/Irl27+Pvf/557ee7cuQDcdNNNPP7446SkpOQ+8cH3UrG//OUvzJkzh+XLl1OxYkX69++vl/Weh7/z/PHHH5OTk8OMGTOYMWNG7vip7eXM/J1nuTj+znOpUqV47rnnmDlzJn/5y18oV64c1113Hd26dQt49qLE33lu27Yt6enpfPjhh8ydO5cyZcpw5ZVX8uCDDwY8e1GycuVKAEaMGJFnfMCAAbkn9Dx06BDGmNzrGjRowBNPPMHChQt58803qVKlCkOHDj3nQa+Xyljt3xIREREPBfUyjYiIiBR+KiMiIiLiKZURERER8ZTKiIiIiHhKZUREREQ8pTIiIiIinlIZEREREU+pjIiIiIinVEZERETEUyojIiIi4imVEREREfGUyoiIiIh46v8BFPMW1jU9kcUAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "5.5\n",
      "5.5\n"
     ]
    }
   ],
   "source": [
    "from matplotlib import pyplot as plt\n",
    "from matplotlib import style\n",
    "import numpy as np\n",
    "from scipy.interpolate import make_interp_spline\n",
    "\n",
    "\n",
    "def neville(xlist, ylist, a):\n",
    "    n = len(xlist)\n",
    "    p = n * [0]\n",
    "    for i in range(n):\n",
    "        for j in range(n - i):\n",
    "            if i == 0:\n",
    "                p[j] = ylist[j]\n",
    "            else:\n",
    "                p[j] = ((a - xlist[j + i]) * p[j] + \\\n",
    "                        (xlist[j] - a) * p[j + 1]) / \\\n",
    "                       (xlist[j] - xlist[j + i])\n",
    "\n",
    "\n",
    "    if n <= 3:\n",
    "        x = np.array(xlist)\n",
    "        y = np.array(ylist)\n",
    "        style.use('ggplot')\n",
    "        plt.plot(x, y)\n",
    "        plt.plot(xlist, ylist, 'r*')\n",
    "        plt.plot(a, p[0], 'g*',label='interpolating point')\n",
    "        plt.legend()\n",
    "        plt.show()\n",
    "    else:\n",
    "        x = np.array(xlist)\n",
    "        y = np.array(ylist)\n",
    "        X_Y_Spline = make_interp_spline(x, y)\n",
    "        X_ = np.linspace(x.min(), x.max(), 500)\n",
    "        Y_ = X_Y_Spline(X_)\n",
    "        style.use('ggplot')\n",
    "        plt.plot(X_, Y_,label='lagrange polynomial curve')\n",
    "        plt.plot(xlist, ylist, 'r*')\n",
    "        plt.plot(a, p[0], 'g*',label='interpolating point')\n",
    "        plt.legend()\n",
    "        plt.show()\n",
    "\n",
    "    return p[0]\n",
    "\n",
    "\n",
    "print(neville([1,2],[5,6], 1.5))\n",
    "print(np.interp(1.5,[1,2,3,4], [5,6,7,9]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h1><center>IMPLEMENTATION OF NEVILLE'S ALGORITHM IN PYTHON</center></h1>\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h3><center>Abstract</center></h3>\n",
    "<center>This report is concerned with the implementation of neville's algorithm in python. neville's algorithm is the recursive algorithm that provides the numerical solution for polynomial interpolations. So, in this report we make the function of neville's algorithm to find the estimate of polynomial on arbitrary value of x.</center>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h2>Introduction:</h2>\n",
    "<h3>Nevilles's Algorithm:</h3>\n",
    "<body>It is appropriate to discuss about neville's algorithm in detail before going into the explanation of code. Neville's algorithm is a recurrence relation used for polynomial interpolation. When there is a N degree polynomial p(x) so it have N+1 data points from which it goes through (x,y) which is the abscissa and ordinate of an underlying function f(x). So, we have to estimate p(x) on any value of x by following method.\n",
    "so if we take the polynomial</body> \n",
    "\n",
    "$$P_{i,j}$$ \n",
    "<body>which is the degree of N=j-i and passes through the points</body>\n",
    "\n",
    "\n",
    "$$(x_{k},y_{k})$$\n",
    "<body>where k=i,i+1,....,j\n",
    "So, the polynomial recurrence relation which it should testify is as follows:</body>\n",
    "\n",
    "$$P_{i,i}(x)=Y_{i}$$\n",
    "<body>when i is greater than 0 and lesser than N</body>\n",
    "\n",
    "$$P_{i,j}(x)=\\frac{(x_{j}-x)P_{i,j-1}(x)+(x-x_{i})P_{i+1,j}(x)}{x_{j}-x{i}}$$\n",
    "<body>when i is greater than 0 and j is greater than i and lesser than N</body>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h3>Description:</h3>\n",
    "<body>Neville's Algorithm which is derived by eric harold neville in 1934 is based on the newton method for the polynomial interpolation which also involves recursion. In 1966, Lyness and Moler exhibited how to process the Maclaurin series of the last adding polynomial, which brings about mathematical approximations for the subsidiaries of the capability at the beginning. They did this by involving undetermined coefficients for the polynomials in Neville's calculation. The selection of points for evaluation is unrestricted in any capacity, despite the way that \"this cycle needs more operations than are expected in limited difference approach.\" They likewise show how their methodology might be utilized to take care of direct issues.</body>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h3>Importance:</h3>\n",
    "<body>nevilles's algorithm has immense importance than other polynomial interpolation methods that is this algorithm is numerically stable that is it uses the given data directly and there is no requirement to represent the polynomial in the basis of 1,t,... . Although, it is fast and dynamic and it has the simple structure and easy to update.</body>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h2>Mathematical Explanation:</h2>\n",
    "<body>If we consider a polynomial P(X) having 'n' number of nodes so the formula is as follows</body>\n",
    "\n",
    "$$P_{a,...,b}(x)=(\\frac{x-x_{j}}{x_{i}-x_{j}})P_{a,...,j-1,j+1,...b}(x) + (\\frac{x_{i}-x}{x_{i}-x{j}})P_{a,..i-1,i+1,..,b}(x)$$\n",
    "<body>So, by using the formula recursively we can generate an array of polynomials </body>\n",
    "\n",
    "$$P_{0}(x)$$\n",
    "$$P_{1}(x)P_{0,1}(x)$$\n",
    "$$P_{2}(x)P_{1,2}(x)P_{0,1,2}(x)$$\n",
    "$$.$$\n",
    "$$.$$\n",
    "$$.$$\n",
    "$$P_{n}(x)P_{n-1,n}(x)....P_{0,...,n}$$\n",
    "\n",
    "<body>So, we can display the polynomial array as matrix form</body>\n",
    "\n",
    "$$S_{ij}(x)$$\n",
    "<body>where i and j are the position of element in the matrix</body>\n",
    "\n",
    "$$S_{00}(x)$$\n",
    "$$S_{10}(x)S_{11}(x)$$\n",
    "$$S_{20}(x)S_{21}(x)S_{22}(x)$$\n",
    "$$.$$\n",
    "$$.$$\n",
    "$$.$$\n",
    "<body>So, for all j=0 P(x) is constant and equals to f(x) and all values oj j>0 the general recursive formula is as given which is also used in the code. So, if there is (n+1) data points given:</body>\n",
    "\n",
    "$$(x_{0},y_{0}),(x_{1},y_{1}),(x_{2},y_{2}),....,(x_{n},y_{n})$$\n",
    "<body>for an underlying function f and given number of x for which we would like to estimate of f(x). We are interested in computing degree n polynomial given:</body>\n",
    "\n",
    "$$P_{0,1,...n}(x)$$\n",
    "<body>The general recursive formula for this is as follow</body>\n",
    "\n",
    "$$$$\n",
    "\n",
    "$$P_{i-j,i-j+1,..,i}(x)=\\frac{(x-x_{i-j})P_{i-j+1,..,i}(x)+(x-x_{i})P_{i-j,...,i-1}(x)}{x_{i}-x_{i-j}}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h2>About the code:</h2>\n",
    "<body>In this code we define a function called neville which has three arguments. These three arguments is used in the code and given when the function is called. In this functon we use the formula of neville's interpolation algorithm which is made recursive by the loops used in this code. Moreover,we import three libraries numpy,scipy and matplotlib. These libraries are used in making the plot of the given data and expressing the coordinate of the interpolating algorithm.</body>\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h2>Code:</h2>\n",
    "<h3>Importing libraries:</h3>\n",
    "<h4>Numpy library:</h4>\n",
    "<body>The numpy library is used here to get the result of interpolation directly from it so we can compare it with the estimate given by neville's algorithm and can judge the accuracy of the algorithm at the certain point.</body>\n",
    "<h4>The scipy library:</h4>\n",
    "<body>The function make_interp_spline is called from the scipy library to make the nodes to form the smooth curve ehen nodes are greater than three. We take the list x and y as an argument in the function that we import and assign it to X_Y_Spline method. Furthermore,we use numpy linspace function and set the curve smoothness to 500 and assign it to X_ and we put X_ as an argument in X_Y_Spline method and assign it to y.</body>\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h4>The matplotlib library:</h4>\n",
    "<body>The pyplot and style functions are imported from this library. pyplot is used to plot the nodes on the graph style function is used here to set the style of the grid to 'ggplot'. we use a third argument 'r*' and 'g*' which represents the input nodes to be red and interpolating point a to be green coordinate. To implement the plotting we use if-else loop to set the control flow. If the input nodes are lesser than four so the plot is without converting it to curve else when the nodes are greater than four so it first converts it to the curve using scipy implementation and than plot it on the grid using matplotlib.</body>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h3>Definition of function:</h3>\n",
    "<body>At first we define a function neville which have three arguments that is xlist,ylist and a. xlist is the input list of x and it is the list of size n, ylist is the input list of y which is corresponding to x in a sequece which is of the size n and 'a' is the arbitrary value of 'x' for which the value of polynomial to be estimated.</body>\n",
    "<h3>Setting up variables:</h3>\n",
    "<body>we set two variables inside the neville function which is n and p. n is the len(xlist) in which the length of xlist input is assigned to n and p uses n and generate a list od zero's whose length is of size n. The example is as follows in which we take xlist as a hardcoded list which have five elements and a print function to print p:</body>\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[0, 0, 0, 0, 0]\n"
     ]
    }
   ],
   "source": [
    "xlist=[1,2,3,4,5]\n",
    "n = len(xlist)\n",
    "p = n*[0]\n",
    "print(p)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h3>The for loops:</h3>\n",
    "<body>In this code we use two nested for loops who are responsible for the recursion of the neville's algorithm. The outer for loop has the variable i and runs from 0 to n-1 which runs n times and the inner for loops has a variable j and runs from 0 to (n-i-1) which is n-i times.i and j represent the position of element in the matrix. i changes when row chnges and j changes when column changes. So, when i=0 the inner loop runs n times and when i increases the occurence of inner loop decreases.</body>\n",
    "<h3>The if-else statement:</h3>\n",
    "<body>We use the if-else statement in the inner loop. The if statement executes when i=0 that is the degree 0 polynomial which is directly equivalent to their corresponding f(x) value. The else statement executes when i > 0. The else statement contains the recursive formula which uses the previous value of p[j] already re-stored in the list p. When the else statement executes the list p upgrades and move one consecutive forward in every iteration and in the last iteration p[n-1] changes which is the final estimated value.</body>\n",
    "<h4>example:</h4>\n",
    "<body>We take an example of four data points as follows:\n",
    "\n",
    "$$X_{1}(1,5) , X_{2}(2,6) , X_{3}(3,7) , X_{4}(4,8)$$\n",
    "<body>and we want to find the estimate of the underlying function f(x) at a=3.5. So, we make two list and one 'a' input as arguments to the function as follows:</body>\n",
    "\n",
    "\n",
    "$$xlist=[1,2,3,4]$$\n",
    "$$ylist=[5,6,7,8]$$\n",
    "$$a=3.5$$\n",
    "<body>the code given below runs the above example to elaborate it.</body>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "7.5\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "7.5"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import numpy as np\n",
    "def neville(xlist, ylist, a):\n",
    "  \n",
    "    n = len(xlist)\n",
    "    p = n*[0]\n",
    "    for i in range(n):\n",
    "        for j in range(n-i):\n",
    "            if i == 0:\n",
    "                p[j] = ylist[j]\n",
    "            else:\n",
    "                p[j] = ((a-xlist[j+i])*p[j]+ \\\n",
    "                        (xlist[j]-a)*p[j+1])/ \\\n",
    "                        (xlist[j]-xlist[j+i])\n",
    "    return p[0]\n",
    "print(neville([1,2,3,4],[5,6,7,8],3.5))\n",
    "np.interp(3.5,[1,2,3,4],[5,6,7,8])\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h2>Conclusion:</h3>\n",
    "<body>In this report the implementation of neville's algorithm is demonstrated using python and explained mathematically using the recurrence formula for neville's algorithm. The code is demostrated using python in which we made the function and implement the recurrence formula in the code. From this we can conclude that neville's algorithm is better the recurrence relation of newton form as it is numerically stable.</body>"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3.9.13 64-bit",
   "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.9.13"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "e89691e4663d647bc11c8a9859f5cdd56224e794f653bdd6c406cd9780fe241b"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}