Demos/Ch3_PolynomialInterpolation.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "Ch3_PolynomialInterpolation.ipynb",
      "provenance": [],
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/stephenbeckr/numerical-analysis-class/blob/master/Demos/Ch3_PolynomialInterpolation.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "1ddQFqPKPg27",
        "colab_type": "text"
      },
      "source": [
        "# Polynomial Interpolation (Lagrange form/Barycentric)\n",
        "\n",
        "There are excellent python demos already online, from the CS version of this class and from the Driscoll and Braun textbook; see below for links.\n",
        "\n",
        "References:\n",
        "- [Berrut and Trefethen paper on Barycentric interpolation, SIAM Review 2004](https://people.maths.ox.ac.uk/trefethen/barycentric.pdf)\n",
        "- Jed Brown's [CSCI 3656 Interpolation jupyter notebook](https://github.com/cu-numcomp/numcomp-class/blob/master/Interpolation.ipynb). [![Open In Colab](https://colab.research.google.com/assets/colab-badge.svg)](https://colab.research.google.com/github/cu-numcomp/numcomp-class/blob/master/Interpolation.ipynb)\n",
        "\n",
        "- Driscoll and Braun's textbook on numerical analysis has code in Matlab, Python and Julia. Here are their [core routines for interpolation](https://github.com/tobydriscoll/fnc-extras/blob/master/python/FNC09.py) as well as their [jupyter notebook with demos](https://github.com/tobydriscoll/fnc-extras/blob/master/python/Chapter09.ipynb) [![Open In Colab](https://colab.research.google.com/assets/colab-badge.svg)](https://colab.research.google.com/github/tobydriscoll/fnc-extras/blob/master/python/Chapter09.ipynb) \n",
        "  - If you run their jupyter notebook, you may need to download some files, using something like:\n",
        "\n",
        "```\n",
        "!curl -O https://raw.githubusercontent.com/tobydriscoll/fnc-extras/master/python/FNC[01-13].py\n",
        "!curl -O https://raw.githubusercontent.com/tobydriscoll/fnc-extras/master/python/FNC.py\n",
        "```\n",
        "\n",
        "  - See also the [main page for their book](https://tobydriscoll.net/project/fnc/)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "T_kHXe3xEKwQ",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "# Download the files from Driscoll and Braun:\n",
        "!wget https://github.com/tobydriscoll/fnc-extras/raw/master/python/FNC09.py\n",
        "from FNC09 import polyinterp\n",
        "import numpy as np"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "hqdfZ8DVwIoH",
        "colab_type": "text"
      },
      "source": [
        "## If you want to see how to use the code"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "W-0kFQeDEXd2",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "?polyinterp"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "X-ZGWp2hwNKB",
        "colab_type": "text"
      },
      "source": [
        "## Make our own functions and interpolate"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "tnNNZIliEazW",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 34
        },
        "outputId": "5f1c13e5-489d-46df-d53c-901851149964"
      },
      "source": [
        "f = lambda x : 1/x \n",
        "x  = [1,2,4]\n",
        "x  = np.array(x)\n",
        "p  = polyinterp( x, f(x) )  # returns a function that we can then query\n",
        "p(3)\n",
        "\n",
        "x   = np.array( [2,6] )\n",
        "p   = polyinterp( x, f(x) )\n",
        "p(4)"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "array(0.33333333)"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 6
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "yAXNCiEBEq4P",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 286
        },
        "outputId": "565a1c44-90b4-4f6c-9986-0159dbebb8dd"
      },
      "source": [
        "import matplotlib.pyplot as plt\n",
        "\n",
        "f   = lambda x : 1/x \n",
        "x   = np.array( [1,2,4,5])\n",
        "p   = polyinterp( x, f(x) )\n",
        "\n",
        "xGrid = np.linspace(.1,10,100)\n",
        "plt.plot( xGrid, f(xGrid) )\n",
        "plt.plot( xGrid, p(xGrid) );\n",
        "plt.plot( x, f(x),'o',markerSize = 12 );\n",
        "plt.ylim( -2,4 );"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "(-2.0, 4.0)"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 17
        },
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXIAAAD8CAYAAABq6S8VAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3deZwU9Z3/8denj+mei3sYkAFRxANBJU5UNCB4x1tjFH/R4MYEf5uYNYnJmsTdXzaJ2cSN8ZdzzbJqQowajdFoiPEKEHDjEVQUARHBA0RggHDMMFd3f/eP6uESmBmneqq76v18PPox1dXVVZ/2ePe3v/X9VplzDhERKV2xoAsQEZGeUZCLiJQ4BbmISIlTkIuIlDgFuYhIiVOQi4iUON+C3MziZvaSmc3ya58iItI5P1vk1wFLfdyfiIh0gS9BbmZ1wDnA7X7sT0REui7h035+CPwzUL2vDcxsOjAdoLKy8tjDDz/clwO/s2k7Le1ZDq3d56FFRELhhRde2OCcq9lzfY+D3MzOBdY7514ws8n72s45NwOYAVBfX+8WLFjQ00MDcMMDrzBn2Xqev/E0X/YnIlKszOztva33o2vlJOB8M3sL+A1wipn92of9dkllKkFTa6a3DiciUnR6HOTOua855+qccyOBqcBs59wVPa6si6pScZrasuRyuviXiERTyY8jr0x5vUPb27MBVyIiEgxfg9w5N9c5d66f++xMR5Cre0VEoqrkW+RV+SBvVJCLSESVfJCrRS4iUReCII8DapGLSHSVfJBX7WiR62SniERTyQe5ulZEJOpKPsh1slNEoq7kg1wtchGJupIP8oqkd7JTQS4iUVXyQR6LGZVlcRp1slNEIqrkgxx04SwRibZQBHlVOkFjm4JcRKIpHEGuFrmIRFgogryyTEEuItEVjiBPJXSyU0QiKxRBXpWKq0UuIpEViiDXqBURibJQBHlVKqEp+iISWaEI8spUgtZMjkw2F3QpIiK9LhRB3nHhrG0tapWLSPSEIshrqlMArN/WGnAlIiK9LxRBPqRvGoB1W1sCrkREpPeFIshrqxXkIhJdoQjywX3UtSIi0RWKIE8n4/QtT6pFLiKRFIogB6jtk1KQi0gkhSjI06zbqq4VEYme0AT54Oo069UiF5EI6nGQm1nazJ43s5fNbLGZfdOPwrqrtk+K9dtayeVcEIcXEQmMHy3yVuAU59zRwDHAWWZ2gg/77ZbaPmkyOcem7W29fWgRkUD1OMidpzH/NJl/9HqzuDY/BFEnPEUkanzpIzezuJktBNYDTzrnnvNjv90xuI83KWi9TniKSMT4EuTOuaxz7higDjjOzMbuuY2ZTTezBWa2oKGhwY/D7qa2j2Z3ikg0+TpqxTm3GZgDnLWX12Y45+qdc/U1NTV+HhaAmqqOrhW1yEUkWvwYtVJjZv3yy+XA6cBrPd1vd5UlYgysLGPdNrXIRSRaEj7sYygw08zieF8M9zvnZvmw324b3EdjyUUkenoc5M65V4DxPtTSY7V9UqxVkItIxIRmZid4l7NVH7mIRE24grxvmg2Nrbp3p4hESriCvE8K52BDo2Z3ikh0hCvIdacgEYmgcAW5JgWJSASFLMjzk4J0yzcRiZBQBfnAqhQxQ2PJRSRSQhXk8ZhRU61bvolItIQqyEG3fBOR6AldkA+uTqtFLiKRErog77jlm4hIVIQwyNNsamqjNZMNuhQRkV4RwiD3hiDqTkEiEhWhC/KDBlUB8EZDYydbioiEQ+iC/NBaL8iXr9sWcCUiIr0jdEHer6KMwdUplq1Vi1xEoiF0QQ5waG01y9erRS4i0RDKIB9dW8XydY3kci7oUkRECi6UQX5YbTXN7VlW/7056FJERAoulEE+urYagNd1wlNEIiCUQd4xcmWZglxEIiCUQV6dTnJA37SGIIpIJIQyyMHrXlm2TkMQRST8Qhvkhw2pZkVDI1mNXBGRkAttkI8eXEVbJsfbG5uCLkVEpKBCG+SHDdHIFRGJhtAG+SGDvZErr6ufXERCLrRBXlGWYPiAcg1BFJHQ63GQm9lwM5tjZkvMbLGZXedHYX44rLZaQxBFJPT8aJFngOudc2OAE4DPmdkYH/bbY6Nrq1nZ0ERbJhd0KSIiBdPjIHfOveecezG/vA1YCgzr6X79cFhtNZmc480NGrkiIuHlax+5mY0ExgPP7eW16Wa2wMwWNDQ0+HnYfTp6eD8A/vbWpl45nohIEHwLcjOrAn4HfME5t3XP151zM5xz9c65+pqaGr8Ou18jB1YwpE+aZ1Zu7JXjiYgEwZcgN7MkXojf7Zx70I99+sHMmDBqIM+u2IhzmuEpIuHkx6gVA+4Aljrnbu15Sf6acPBANja1aTy5iISWHy3yk4ArgVPMbGH+cbYP+927XPdGoEwYNRCAZ1ZsKEQ1IiKB82PUytPOOXPOHeWcOyb/eNSP4t7n2Z/DPR+HXLbLbxk+oIK6/uXqJxeR0CqtmZ1lFfDGUzD/B91624SDB/Lcm5t0D08RCaXSCvLxV8K4S2Hud+Gtp7v8tgmjBrJ5eztL175vMI2ISMkrrSA3g3NvhQEHwwNXQ2PXxqPv7CdX94qIhE9pBTlAqho+PhNaNsODn+lSf/nQvuWMHFjBs+onF5EQKr0gBxgyFj56M6ycA3O/16W3TBg1kOdWbiKT1XVXRCRcSjPIAT40DY65Aub9Byz7U6ebTxg1iG2tGV5do35yEQmX0g1yMzjnFhh6NDx4DWxcsd/NJx4yiETMeHTRe71UoIhI7yjdIAdIlsOld0EsBvddAa37nr3Zv7KMkw+t4ZGFa3RDZhEJldIOcoD+B8Ild0LDa/Dg9P3O/Lxg/DDWbm3huTd10lNEwqP0gxxg1Clw5ndh2R9hznf2udnpR9RSWRbn4ZfW9GJxIiKFFY4gBzj+Gu8E6Pxb4JXf7nWT8rI4Z44dwqOL3qOlvevT/EVEill4gtwMzr4FDjwJHv4crHp+r5tdeMwwtrVmmPPa+l4uUESkMMIT5ACJMu/kZ99hcM9lex3JcuKogQyqSvH7he8GUKCIiP/CFeQAlQPhEw94LfS7L4Gm3S9fm4jHOO/oocx5rYEt29sDKlJExD/hC3KAgaPg8t/A1jVw71Ro277byxeNH0ZbNsdDL60OqEAREf+EM8gBhh8HF/83rF4Av70Ksjtb3+OG9aX+wP789/w3adeUfREpceENcoAx53tXS1z+OPz+szvGmJsZn50yinc3N/PwQg1FFJHSFu4gB6j/FJzyr7Dofnjsq5C/CfOUwwZz+JBqbpv7hm44ISIlLfxBDjDxephwLTz/X/Dnb4FzmBmfm3IIKxqaeGLJ2qArFBH5wKIR5GZwxk1w7FXw9K07Ln179rihjBxYwX/OXYFzapWLSGmKRpBD/mqJ/9+79O1fvgfzvk88Zvzfk0fxyuotzFu+ofN9iIgUoegEOXhXSTz/x3DUVJh9E8z7Phd9aBjDB5Tz7VlLaMtoBIuIlJ5oBTlALA4X/iccdRnMvonUX/6db513JG+sb+T2p1cGXZ2ISLclgi4gELE4XPhzSKRg/i1MmdDCmWMu5Md/Xs55Rx3A8AEVQVcoItJl0Qxy8LpZzv0RJNLwzE+5dezfOcHO4d/+sJjrzk7xy8W/ZP7q+bRmW0nFU0yqm8RVR17F2EFjMbOgqxcR2cGCGK1RX1/vFixY0OvH3SvnvP7y+bewfNBkLshVku7/OtlcOzl29pnHiJFKpJhcN5nvTPwOyVgywKJFJIrM7AXnXP2e66PXR74nMzj1X3Fn3cwMW0R59SLac627hThAjhzNmWbmrJrDjfNv1HBFESkavgS5md1pZuvN7FU/9heERaNO5C99+tHWyT+RlmwLc1fP5dUNJftRRSRk/GqR/xI4y6d9BWLm4pm05rp216DWTCszl8wscEUiIl3jS5A75+YBm/zYV1DmrZ73vu6UfcmRY97qeQWuSESka3qtj9zMppvZAjNb0NDQ0FuH7bLWbGu3tm/JtBSoEhGR7um1IHfOzXDO1Tvn6mtqanrrsF2Wiqe6tX06kS5QJSIi3aNRK3mT6iYR6+I/jhgxJtVNKnBFIiJdoyDPm3bkNFKJrrXKy8yYdsSVBa5IRKRr/Bp+eC/wDHCYma02s6v92G9vGjdoHJPrJpOO77/LJJmDKdu2Mvaxb8DW93qpOhGRffNr1Mrlzrmhzrmkc67OOXeHH/vtTWbGdyZ+hynDp1CeKH9fN0uMGOl4mlS2ntia88m+9Qzutgmw+PcBVSwi4lHXyi6SsSQ3T7qZO864g9MOPI3yRDmGUZ4o5/SRp/OLs37BnGm3s/nwKzmj+Sbeyg2G306DB6fD9pIefSkiJUzXWvkAnHPcPv9NbnnsVb5aMYursr/DKvrD2bfAmAu8af8iIj7TtVZ8ZGZ8ZtLB3HPNR5gRv4xz225iLQO91vl9V8DmVUGXKCIRoiDvgWMPHMCfrpvI4UefyEkbb+S/yqaRXf4U/Ow4mH8rZNqCLlFEIkBB3kP9Ksr4waVHc+enJvCr2IVMarqZhWUfgj9/E247EZY95l0qV0SkQBTkPjn50Bqe/NIkLjn1RKZuvZbPZG5gY1Mr3HsZ/Op8eO+VoEsUkZDSyc4CWLulhVueWMYjL77FlcnZXF/2EBWZLXDkxTDl6zBo9F7f55xj0YZFujuRiOzVvk52KsgLaEVDIz+b8wazFy5nevyPXJ18jDLXih01FSZeD4MO2bFte66dG+ffyJxVc2jLtunuRCLyPgryAL21oYk7nn6TOS8uZlru90xLPEWSDNkjLiAx6Uu4IeO4Yd4NzFk1h5bsvq+qmI6nmTJ8CjdPulktc5EIUpAXgS3N7fx2wSpmPfMyZ2x9kE/Gn6TKmvnr0GO5rmIzLbn2TvdRnijnjjPuYFzNuF6oWESKiYK8iDjneO7NTTzy3BIGLr2blbVzeboyQa4LrewYMU4feTq3nHxLL1QqIsVkX0GeCKKYqDMzTjh4ICccPJGm1glMvO8Ecq5rN7bQ3YlEZE8afhiwylSCjOvexCHdnUhEdqUgLwLdvTtRzMX5weNLmbtsPVuaO+9XF5FwU9dKEZhUN4mn3n6qSzd/jjnHaU1buPKvZzF7/jHckDuGNQM+zCHD6ziqri9jh/VlzAF9qCjTv1qRqND/7UVg2pHTmP/ufJozzZ1uW5YoZ9r4TzFg1ctc8sZTTG2fS64xxpKlo5i36Ah+kjuCF91oagbVcMTQPowZ2ofDaqs5tLaauv7lxGLFNWxRk6BEek6jVoqAc+6DjSPPtsO7L8CK2bgVc2DNi1guQ44Yq8sO5sXcwfxP80hezo1ihTuAZDLJqJoqDhlcxaiaKg6uqWTkwEoOGlRJZar3v9M1CUqkezT8sMh1hNrc1XNpzbS+L9TK4mVMGT5l/6HW1gSr/wZv/xVWPe+FfOtWALKxMtamR7HcRvJy21AWNNXyem4Y6+gPGDXVKQ4cUMGIARWMGFhBXf8Khvcvp25ABbXVKRJxf0+nfOAvL59r0K8BKSUK8hLgnOPVDa96wfLufFoyLaQT6d2CpVtyOdi4HNYshLWveI91i2H7xh2bZOLlbE7XsSZ+AG9lB/F66wCWNPdjTW4g77kBbKGSeCzGkD5pDuiXZmjfcob2SzO0T5ohfdPU5v8OqkqR7EbYv9LwCp9+4tNd6k4qxCQo/RqQUqQgl50aG6BhKTQsg00rYeMK2LTCuyFGdvfx7JlYmsbkAP4e609Dri9rM5WsaaugIVfFFqrY4irZ4ipptAri5X1JV/Wnuk8/BvSppKY6xaCqFIOqyhhUlWJgVRkDKssYUFHGDfO/0vUTvD5PgiqGXwMiH4QmBMlOVTXe46BJu6/P5aBpPWx+B7auga1rSGx9l36N6+nXuI6DGtfD9jdwuY2Yy75/v1lgi/doI0GTS7OdFC2ujGZSbKOMDS5JK2XMOXgTuVjXGhE5cvzlndm4FXOxsgpIpCFZAclyKKuAsiqId73VvGjDIuaunrvfEAdoybYwd/VcXt3wqi6JIEVNQS47xWJQPcR77Ic5By1boGUzNG+G5r9D67b8Yyu0NVLW1kSytZGK5m20NTeRaWkk29ZMrr0Z2pvIWPd+CbZm27C7Ltjn6xlL0h6vIJOsIpusIlfWB9J9obwfsYr+JKoGkayuoaxPDTPfeojWLk6qas20MnPJTF0SQYqagly6zwzK+3mP/vvZDEjlH3tK/frDnbaId5WwMu4e83Pamptob2ki27Yd17Yd2hqx9mbimSbS7c1UtjbTh2aqrZG+rKOPNdGPRiptZ5fRvAPryMW61p+fI8fst+cy65U1VJTFKU8mKC+L55fjpJNxyvPL8SIb2inRoSCXQHRrEhQxThl5Cp84+fL9bteWydHYmmFbSzvbWjJsac3wbkuGxtYM27c3kW3aiGvaQMvWr3er1kyulS33X8vfXB3L3HCW5A5kK5Xv264sHiOdjO0I93QiTjoZI5XsCH3vtY716WScVMJ7fbe/iRip/DZl+eVUomN557qO54mYqQ8/4hTkEohuTYKKlzFtzLTOt0vEGJDwTqjuz09+/a1u/RpIYUytXEC89c871jVV1LGp7xjWVx/Ju5VjeCd9OI25MlraszS3ZWnJZL3l9hwt7Vk2b2+jpT1HSyZLa/5vS3uWlvbOv8g6Y+Z9iZQlYjv/7rKcjO/8EkjGYyTjtmNdWbxjXYxkwnZ/nt9uz+VE3NsuETcSsRhlCe9vIr9NIrb7tsn8ax3LxTYpLQwU5BKIcYPGMblucpdHjnR76OV+dPfXwMkHnU582vdh21pv+Obal6l87xUq17zE8Pee4FiAWAIOGA8HngiHf8T7m6rqdP/OOVozOdqyXuC3tue85xkv7Nsy3vPW9ixt2dwur+983pbN7diuPb/csa5juTWTo6k1Q3vW0ZbfruO4mWzOW5/11hd6IJsZJGMx4jHbLfwTMSOx4wvCiMfy63c8974w4jEjGd/9+Z7vief3F+9Yb/nX4kbMvHWxmO2+7S6Pjm12PI/t3E9sj+1jtpdl847tbU/++N57K5Jx3+dlaPihBMaXSVAfQHfGsKfjae488859j1ppbPAmXq16Ft5+xlvOtUMsCcOPg1FT4NCzoHasl2AlIJvLh30uR3vGC/n2fMh3LGdyO9dldrzuyOR2Ps/k3I4viUwu/7djmz1e89Y7svn9ZnNuxzYd6zve623jyLn3b7Nj/S7rcs5bXyx++Q8fZvJhgz/Qews6jtzMzgJ+BMSB251z39vf9gpy6eD7JKguHrNg48jbm2HVc7BiDqyY7U3CAug73Av0MefDgSdBLO7DJ5HuyOV2fllknSPb8cXgdn5R5PLL2V23zS/v+tqu+8g6R67jea7ji6Rje3a+nn/tzCOHMHxAxQf6DAULcjOLA68DpwOrgb8BlzvnluzrPQpyCVqv/RrYtg6WPwGvP+YFe/t2qKyBI86Doy6D4ceXTEtdglfIIJ8A/Jtz7sz8868BOOe+u6/3KMilGPT6r4G2Jlj+JCz5Pbz+uBfq/UfCUVNh/Ceg3wh/jyehU8ggvwQ4yzn36fzzK4HjnXPX7rHddGA6wIgRI459++23e3RckZLW2givzYKX74WVf/HWjT4djr0KRp8JcY1DkPcLfIq+c24GMAO8FnlvHVekKKWq4Oip3mPzKnjxV/DSXfCb/+O1zI+7Bj50pTc7VaQTfoyBeRcYvsvzuvw6EemKfsPhlBvhC6/CpXdB3xHwxI1w6xh4/Ebvujci++FHkP8NGG1mB5lZGTAVeMSH/YpESzzhjWr5hz/CNfPg8HPg2dvgR0fDH66Dv78VdIVSpHoc5M65DHAt8DiwFLjfObe4p/sVibShR8PFM+CfXoTxV8LCe+Enx8KsL8IW/eCV3WlCkEgp2LoG5t0CL84Ei8Px02Hil70Ll0lk7Otkp7/zREWkMPocAOfeCp9/AcZ+DP76U/jxeHhuhnfvVok0BblIKek/Ei66Da75C9QeCX/6Ctx2Erw5P+jKJEAKcpFSNPRomPYHmHovZJph5rnwu097F/aSyFGQi5QqMzj8bPjc8zDpn2HJw/DT47wx6QGc+5LgKMhFSl2y3BuH/o/PwJCx8Mjn4a4LNVwxQhTkImEx6BCYNgvO+QGsXuD1nb94l1rnEaAgFwmTWAw+/Gn47LPejS4euRbuuwKaNgRdmRSQglwkjPoNh08+Amfc5F1G97YT4c15QVclBaIgFwmrWAxO/Dx8Zo538a1fXQBzvwe5bNCVic8U5CJhN2SsF+bjLoW53/VOhKqrJVQU5CJRkKqCi34O5/8U3nkOZkyGNQuDrkp8oiAXiQoz7xrnVz/ujWS580x4+b6gqxIfKMhFouaA8d4U/7oPw0PT4alvQi7X+fukaCnIRaKochBc+ZB3a7mnb4UHroL25qCrkg9INwYUiap4Es79IQwcDU/8C2xZDZffB1U1QVcm3aQWuUiUmcGJ18Jlv4Z1S7x+c03tLzkKchGBI86FTz4M2zfCHWfA2kVBVyTdoCAXEc+I4+FTj0MsAb84G955NuiKpIsU5CKy0+DD4eonoGow3HWxpvWXCAW5iOyubx1c9Sj0GwF3fxyWPxV0RdIJBbmIvF91LVz1Rxh0KNw7FZY9FnRFsh8KchHZu8qB3u3khoyF+6+EN9QyL1YKchHZt/J+3sShmsPgN5+AlXODrkj2QkEuIvtX3t+7tvmAUXDPVHj7r0FXJHtQkItI5yoGwLRHvBtW3HMZvPdK0BXJLhTkItI1HddnSfWBX18MG1cEXZHkKchFpOv61sEnfw8uB7+6ELauCboioYdBbmYfN7PFZpYzs3q/ihKRIjZoNFzxO2jeBHdfCi1bg64o8nraIn8VuBjQ9C+RKDlgPFw6E9Yvgd9Og2x70BVFWo+C3Dm31Dm3zK9iRKSEHHIanPdDWDEbZn3Bu+uQBELXIxeRD+5Dn4TNq2Def0D/g2DSl4OuKJI6DXIzewoYspeXbnTOPdzVA5nZdGA6wIgRI7pcoIgUuSlf965hPvvb3sShI84LuqLI6TTInXOn+XEg59wMYAZAfX29foOJhIUZnP8T2LQCHrwGrh4JQ8YFXVWkaPihiPRcMg1T74F0X2/2Z+P6oCuKlJ4OP7zIzFYDE4A/mtnj/pQlIiWneghcfo93l6H7NZKlN/V01MpDzrk651zKOVfrnDvTr8JEpAQdMB7O/zG881d48v8FXU1kaNSKiPjrqEvh3Rfg2f+EYcfCuEuCrij01EcuIv474yYYMQEe+TysWxx0NaGnIBcR/8WT8PFfehfYuv+T0Lot6IpCTUEuIoVRPQQuuQM2rYRZX9LMzwJSkItI4Yz8CEz+Giy6H166K+hqQktBLiKFNfF6OOhkePQr6i8vEAW5iBRWLA4fu93rL3/gU9DeHHRFoaMgF5HCqxoMF/0cGl7T+PICUJCLSO845FQ44bPw/Ax4/YmgqwkVBbmI9J5TvwGDj4SHPwuNDUFXExoKchHpPck0fOy/vdvDPfJ5DUn0iYJcRHpX7ZFw2jfg9T/BwnuCriYUFOQi0vuO/0cYcSI89lXYsjroakqeglxEel8sBhf+DHIZePhadbH0kIJcRIIx4GA4/Vuwcg4suDPoakqaglxEglN/NRw82RtbvnlV0NWULAW5iAQnFoPzfgwuB7O+qC6WD0hBLiLB6n8gnPKv8MaTsOi3QVdTkhTkIhK846+BYfXwpxugaUPQ1ZQcBbmIBC8Wh/N/4t2A4rGvBV1NyVGQi0hxqB0DE7/kXbt8xeygqykpCnIRKR4f+ZI3LPGPX4b2lqCrKRkKchEpHsk0nPMD2LQC/udHQVdTMhTkIlJcRp0CYz8G838AG1cEXU1JUJCLSPE5898hkYJHv6yx5V2gIBeR4lM9BE75F++k52uzgq6m6CnIRaQ41V/t3YTi8a/rPp+dUJCLSHGKJ+CjN8Pmd+B/fhx0NUWtR0FuZt83s9fM7BUze8jM+vlVmIgIB02EIy+Cp2/1Al32qqct8ieBsc65o4DXAU3JEhF/nf5twOCJfwm6kqLVoyB3zj3hnMvknz4L1PW8JBGRXfQbDhOvhyUPw5vzg66mKJnzaWiPmf0BuM859+t9vD4dmJ5/ehiwrBu7HwRE8Uo6UfzcUfzMEM3PHcXPDD373Ac652r2XNlpkJvZU8CQvbx0o3Pu4fw2NwL1wMXOr2+G3WtY4Jyr93u/xS6KnzuKnxmi+bmj+JmhMJ870dkGzrnT9ve6mV0FnAucWogQFxGR/es0yPfHzM4C/hk42Tm33Z+SRESkO3o6auWnQDXwpJktNLOf+1DT3swo0H6LXRQ/dxQ/M0Tzc0fxM0MBPrdvJztFRCQYmtkpIlLiFOQiIiWu6IPczM4ys2Vm9oaZfTXoegrNzIab2RwzW2Jmi83suqBr6i1mFjezl8wsMpe7M7N+ZvZA/lIXS81sQtA19QYz+2L+v+9XzexeM0sHXZPfzOxOM1tvZq/usm6AmT1pZsvzf/v7cayiDnIziwM/Az4KjAEuN7MxwVZVcBngeufcGOAE4HMR+MwdrgOWBl1EL/sR8Jhz7nDgaCLw+c1sGPBPQL1zbiwQB6YGW1VB/BI4a491XwX+7JwbDfw5/7zHijrIgeOAN5xzK51zbcBvgAsCrqmgnHPvOedezC9vw/sfe1iwVRWemdUB5wC3B11LbzGzvsAk4A4A51ybc25zsFX1mgRQbmYJoAJYE3A9vnPOzQM27bH6AmBmfnkmcKEfxyr2IB8GrNrl+WoiEGodzGwkMB54LthKesUP8eYk5IIupBcdBDQAv8h3Kd1uZpVBF1Vozrl3gVuAd4D3gC3OuSeCrarX1Drn3ssvrwVq/dhpsQd5ZJlZFfA74AvOua1B11NIZnYusN4590LQtfSyBPAh4Dbn3HigCZ9+ahezfL/wBXhfZAcAlWZ2RbBV9b78THhfxn8Xe5C/Cwzf5Xldfl2omVkSL8Tvds49GHQ9veAk4Hwzewuv++wUM9vrxddCZjWw2jnX8YvrAbxgD7vTgDedcw3OuXbgQeDEgGvqLevMbChA/u96P3Za7EH+N2C0mR1kZmV4J0QeCbimgjIzw+szXeqcuzXoenqDc+5rzrk65/l058IAAADLSURBVNxIvH/Hs51zoW+hOefWAqvM7LD8qlOBJQGW1FveAU4ws4r8f++nEoGTvHmPANPyy9OAh/3YaY+utVJozrmMmV0LPI53ZvtO59zigMsqtJOAK4FFZrYwv+7rzrlHA6xJCufzwN35hspK4B8CrqfgnHPPmdkDwIt4o7ReIoTT9c3sXmAyMMjMVgPfAL4H3G9mVwNvA5f6cixN0RcRKW3F3rUiIiKdUJCLiJQ4BbmISIlTkIuIlDgFuYhIiVOQi4iUOAW5iEiJ+1/ru1J6yJzzSQAAAABJRU5ErkJggg==\n",
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ]
          },
          "metadata": {
            "tags": [],
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "BoZefw6iwRD1",
        "colab_type": "text"
      },
      "source": [
        "## Take a look at the code for Barycentric interpolation"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "AvxsNHkGu3Z0",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "# Or, see the code from Driscoll and Bruan:\n",
        "import numpy as np\n",
        "def polyinterp(t,y):\n",
        "    \"\"\"\n",
        "    polyinterp(t,y)\n",
        "    Return a callable polynomial interpolant through the points in vectors `t`,`y`. Uses\n",
        "    the barycentric interpolation formula.\n",
        "    \"\"\"\n",
        "    n = len(t)-1\n",
        "    C = (t[-1]-t[0]) / 4       # scaling factor to ensure stability\n",
        "    tc = t/C\n",
        "\n",
        "    # Adding one node at a time, compute inverses of the weights.\n",
        "    omega = np.ones(n+1)\n",
        "    for m in range(n):\n",
        "        d = tc[:m+1] - tc[m+1]          # vector of node differences\n",
        "        omega[:m+1] = omega[:m+1]*d     # update previous\n",
        "        omega[m+1] = np.prod( -d )         # compute the new one\n",
        "    w = 1 / omega                   # go from inverses to weights\n",
        "\n",
        "    def p(x):\n",
        "        # Compute interpolant.\n",
        "        z = np.where(x==t)[0]\n",
        "        if len(z)>0:       # avoid dividing by zero\n",
        "            # Apply L'Hopital's Rule exactly.\n",
        "             f = y[z[0]]\n",
        "        else:\n",
        "            terms = w / (x - t)\n",
        "            f = np.sum(y*terms) / np.sum(terms)\n",
        "        return f\n",
        "        \n",
        "    return np.vectorize(p)  # converts to the right data type"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "ehhOcFwqu-jJ",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        ""
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "_FX1xAnyI0tL",
        "colab_type": "text"
      },
      "source": [
        "# Using Scipy and Numpy libraries\n",
        "\n",
        "- Task: experiment with interpolating a Gaussian\n",
        "\n",
        "We'll explore these libraries:\n",
        "1. `numpy`  ([`numpy.interp`](https://numpy.org/doc/stable/reference/generated/numpy.interp.html))\n",
        "2. `scipy.interpolate` [tutorial](https://docs.scipy.org/doc/scipy/reference/tutorial/interpolate.html) and [full documentation](https://docs.scipy.org/doc/scipy/reference/interpolate.html#module-scipy.interpolate)\n",
        "3. [https://github.com/chebpy/chebpy](`chebpy`), a Python wrapper for [chebfun](https://www.chebfun.org/)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "jLLH21GMI5ua",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "import numpy as np\n",
        "import scipy.interpolate as scinterp\n",
        "import matplotlib.pyplot as plt\n",
        "\n",
        "# Make a function\n",
        "gaussian = lambda x: 1/np.sqrt(2*np.pi) * np.exp(-.5*x**2)\n",
        "\n",
        "n = 6  # use n points, so polynomial can be of degree n-1\n",
        "a, b = -5, 5  # the interval [a,b] that we're interested in\n",
        "nodes = np.linspace(a,b,n)\n",
        "xGrid = np.linspace(a,b,100)  # much finer grid, for plotting"
      ],
      "execution_count": 15,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "-2yI2dNALMR5",
        "colab_type": "text"
      },
      "source": [
        "## Run numpy's interpolation"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "dY_2WZZoJgIU",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 265
        },
        "outputId": "2c2a2f70-8b8d-4b09-9fee-e945ca5114bf"
      },
      "source": [
        "# Run numpy's interpolation.  What kind of interpolation is this?\n",
        "#    It's restricted to piecewise linear\n",
        "pGrid = np.interp(xGrid, nodes, gaussian(nodes) )\n",
        "\n",
        "plt.plot(xGrid,gaussian(xGrid),'-',label='Original function');\n",
        "plt.plot(xGrid,pGrid,'-',label='numpy.interp')\n",
        "plt.plot(nodes,gaussian(nodes),'o',label='Data',markersize = 10);"
      ],
      "execution_count": 16,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD4CAYAAADiry33AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3deXxU5b348c93luxhSwKyrwkQAUHDvglSxBbRWttqq+JPe6292uXS1qVa770uva1tsb0ttbXVVmsttRWvSKlWRfY1yCIBgRCWhDULJCHrLM/vjzOhISRkMpktk+/79corM+ec58x3IPnmmec85/mKMQallFKxyxbpAJRSSoWWJnqllIpxmuiVUirGaaJXSqkYp4leKaVinCPSATSVnp5uBg0aFOkwlFKqQ9m+fXuJMSajuX1Rl+gHDRpEbm5upMNQSqkORUSOtrRPh26UUirGaaJXSqkYp4leKaVinCZ6pZSKcX4lehGZJyL7RSRfRB65zHGfExEjIjmNtj3qa7dfRK4PRtBKKaX81+qsGxGxA0uATwFFwDYRWW6M2dvkuFTgm8CWRtuygduAK4E+wPsikmWM8QTvLSillLocf3r0E4B8Y0yBMaYeWArc1MxxTwE/AmobbbsJWGqMqTPGHAbyfedTSikVJv7Mo+8LFDZ6XgRMbHyAiFwN9DfG/F1Evtuk7eYmbfs2fQERuQ+4D2DAgAH+Ra5UmBwqPs/bu07g9VpLendJdHLbhAGkxEfdbShKNavdP6kiYgMWA3cHeg5jzAvACwA5OTm6QL6KCi6PlxfWFvDzDw5S7/YiYm03Bn6/4Qg/uGU0M7OavRFRqajiT6I/DvRv9Lyfb1uDVGAUsFqs34QrgOUissCPtkpFpYLi83z9zzvIO1HBZ0b35r8WXElGajwA24+e5eE3drPwpa3cek0/nr55FAlOe4QjVqpl/iT6bUCmiAzGStK3AV9q2GmMKQfSG56LyGrgO8aYXBGpAV4TkcVYF2Mzga3BC1+p4KuodXHPH7ZRXuPi13dczbxRvS/af83A7qz4+jR+seogSz48hMMm/PBzYyIUrVKtazXRG2PcIvIg8C5gB14yxuSJyJNArjFm+WXa5onI68BewA08oDNuVDTzeg3ffn0XRWdr+PN9kxg/qEezxyU47Xz3+hEALPnwEOMGdOOL4/X6kopOEm01Y3NycowuaqYi5Ver83n2nf18f342904b3OrxHq9h4Utb2XqkjDfun8Lofl3DEKVSlxKR7caYnOb26Z2xSvlsPFTCT97dz2fG9OaeqYP8amO3CT+/bSxpyXF87U/bKa92hTZIpQKgiV4pwO3x8sRbeQzokcSPPjcGaZhi44e0lHiWfPlqjp+rYcnq/BBGqVRgNNErBSzbcZz8M+d5eN6IgObHXz2gO58d15c/bDzCyfKaEESoVOA00atOr9bl4WfvHeCqfl2ZN+qKgM/zH3OywMDP3z8YxOiUaj9N9KrTe3XzUU6U1/LwvBFtGrJpqn+PJL48aQCv5xaSf+Z8ECNUqn000atOraLWxZIP85memc6UYemtN2jFA7OGkei089N/7g9CdEoFhyZ61am9uO4wZ6tdPOSbE99e6SnxfGX6EP6x5xS7i84F5ZxKtZcmetVp1bk9vLr5KHNG9grq/PevTB9MSryD3284ErRzKtUemuhVp/XOnlOUVtVz1+SBQT1vaoKTW67uy993n6T0fF1Qz61UIDTRq07rj5uOMigtiWlBGJtv6o5JA6n3ePnr9qKgn1upttJErzqlvScqyD16ljsmDcRmC3ymTUuyeqUycXAP/rTlKB5vdC0zojofTfSqU3p1y1HiHTZuvaZfyF7jzskDKSyrYe2B4pC9hlL+0ESvOp3KWhf/t+M4C67qQ7ekuJC9ztzsK8hIjeePm4+G7DWU8ocmetXpvLnjONX1Hu4M8kXYpuIcNm4f358P95+hsKw6pK+l1OVooledztKthYzu25Ux/bqF/LVunzgAAf6aW9jqsUqFiiZ61ankn6lk78kKbrn6khr1IdG7ayKTh6bx9u6TRFvtB9V5+JXoRWSeiOwXkXwReaSZ/feLyMcislNE1otItm/7IBGp8W3fKSK/DvYbUKotlu86iU3gM6N7t35wkNw4pg+HS6rYc7wibK+pVGOtJnoRsQNLgBuAbOD2hkTeyGvGmNHGmLHAs8DiRvsOGWPG+r7uD1bgSrWVMYa3d51g0pA0enZJCNvr3jCqN067sHzX8bC9plKN+dOjnwDkG2MKjDH1wFLgpsYHGGMad1WSAf2MqqJO3okKDpdUceNVfcL6ul2TnMzMymDF7pN4dU69igB/En1foPGVpCLftouIyAMicgirR/+NRrsGi8gOEVkjItObewERuU9EckUkt7hY5xyr0Fi+6wROu3BDO9acD9SNV/XhZHktuUfPhv21lQraxVhjzBJjzFDgYeBx3+aTwABjzDhgEfCaiHRppu0LxpgcY0xORkZGsEJS6gKv1xq2mZGZEdK58y2ZM7IXCU6bDt+oiPAn0R8H+jd63s+3rSVLgZsBjDF1xphS3+PtwCEgK7BQlQrc9mNnOVleG/ZhmwbJ8Q6uG9mLlR+fwu3xRiQG1Xn5k+i3AZkiMlhE4oDbgOWNDxCRzEZPPwMc9G3P8F3MRUSGAJlAQTACV6otlu88QYLTxqeye0UshgVX9aGsqp4Nh0ojFoPqnFqtgmyMcYvIg8C7gB14yRiTJyJPArnGmOXAgyIyB3ABZ4GFvuYzgCdFxAV4gfuNMWWheCNKtcTrNbyTd4rZI3qSHEDh72C5dngGqfEOVu4+ycwsHaJU4ePXT70xZiWwssm2Jxo9/mYL7d4A3mhPgEq1166icxRX1jE3O/wXYRuLd9iZMTyDDz45g9drQrJqplLN0TtjVcx7f99p7DZh1vCekQ6Fudm9KDlfxy4tM6jCSBO9innv7z3DhEE96JrkjHQoXJvVE7tNeH/f6UiHojoRTfQqph0rrWb/6UrmRPAibGNdk5xMGNSD9/eeiXQoqhPRRK9iWkPPec7IyA/bNJiT3Yv9pys5VqpLF6vw0ESvYtr7+06T1SuFgWnJkQ7lgoY/Ou/p8I0KE030KmaVV7vYcriMOSOjY9imwcC0ZLJ6pfD+Xk30Kjw00auYtfrAGTxeEzXj843NGdmLrUfKKK92RToU1Qloolcx6/19Z0hPiWNsGCpJtdWc7F54vIbVB/SirAo9TfQqJrk9XlbvP8PsET2j8saksf26kZ4Sx/v7NNGr0NNEr2LSzsJzVNa6uTYKbpJqjs0mzMjMYP3BYjy6Rr0KMU30KiatPVCM3SZMHZYe6VBaNHN4BmerXew5Xh7pUFSM00SvYtKaA8WM7d+NromRvxu2JdOGpSNixapUKGmiVzGnrKqe3cfLo36FyLSUeEb37aqJXoWcJnoVc9YdLMYYmBHliR5gRmYGO46d1WmWKqQ00auYs+ZAMd2TnIzu2zXSobRq5vAMvAY2HCqJdCgqhmmiVzHF6zWsPVDC9MwM7FE4rbKpcf27kZrgYK0O36gQ8ivRi8g8EdkvIvki8kgz++8XkY9FZKeIrBeR7Eb7HvW12y8i1wczeKWa2neqgpLzdR1i2AbAYbcxdWg6aw4UY4xOs1Sh0Wqi99V8XQLcAGQDtzdO5D6vGWNGG2PGAs8Ci31ts7FqzF4JzAN+1VBDVqlQWHvAGgKZkRm90yqbmjk8g5PlteSfOR/pUFSM8qdHPwHIN8YUGGPqgaXATY0PMMZUNHqaDDR0TW4Clhpj6owxh4F83/mUCok1B84wsncXenZJiHQofmv49KGzb1So+JPo+wKFjZ4X+bZdREQeEJFDWD36b7Sx7X0ikisiucXF+sOuAlNV52b70bPMyOo4vXmAvt0SGZqRrIlehUzQLsYaY5YYY4YCDwOPt7HtC8aYHGNMTkZGxxhbVdFn65EyXB7D9GEd72doemYG246UUef2RDoUFYP8SfTHgf6NnvfzbWvJUuDmANsqFbCN+SXEOWzkDOoe6VDabMrQNGpdXnYc06LhKvj8SfTbgEwRGSwicVgXV5c3PkBEMhs9/Qxw0Pd4OXCbiMSLyGAgE9ja/rCVutT6/FJyBnYnwdnxrvdPGpqGTWBDvs6nV8HXaqI3xriBB4F3gX3A68aYPBF5UkQW+A57UETyRGQnsAhY6GubB7wO7AXeAR4wxuhnUxV0Jefr2HeyIqoXMbucLglOrurfTRO9CgmHPwcZY1YCK5tse6LR429epu0zwDOBBqiUPzYdKgXosIkeYOrQdJ5fc4jKWhepCdG7GJvqePTOWBUTNuSXkJrg6BDLHrRk6rB0PF7DloKySIeiYowmehUT1ueXMHlIWodY9qAlVw/sRoLTxnodvlFBpoledXjHSqspOlvToYdtAOIddsYP6sFGXeBMBZkmetXhNfSAO3qiB6sYyYHT5zlTURvpUFQM0USvOrwNh0ro1SWeoRnJkQ6l3Rr+WOmyxSqYNNGrDs3rNWzML2HqsHREOu74fIPs3l3oluRk/cHSSIeiYogmetWhfXKqkrPVLqYM7fjDNgA2mzBlaBqbC0p12WIVNJroVYfWcOFyytC0CEcSPJOHpHH8XA3HyqojHYqKEZroVYe2uaCUQWlJ9OmWGOlQgmay79NJw01gSrWXJnrVYbk9XrYUlF1IjLFiaEYyGanxbNREr4JEE73qsPJOVFBZ52ZyDA3bAIhY4/QbD+k4vQoOTfSqw2ro8U4eEluJHqz3VHK+jkPFWl5QtZ8metVhbSooJbNnChmp8ZEOJegaZhHp8I0KBk30qkOqd3vZdrgspmbbNNa/RyJ9uyXqBVkVFJroVYe0u+gcNS5PzI3PNxARJg9NY1NBKV6vjtOr9tFErzqkjYdKEYGJg2Mz0YM1Tn+u2sW+UxWRDkV1cH4lehGZJyL7RSRfRB5pZv8iEdkrIrtF5AMRGdhon0dEdvq+ljdtq1QgNh0qZeQVXeieHBfpUEKm4dOKDt+o9mo10YuIHVgC3ABkA7eLSHaTw3YAOcaYMcDfgGcb7asxxoz1fS1AqXaqdXnYfuxszI7PN+jTLZFBaUma6FW7+VNKcAKQb4wpABCRpcBNWHVgATDGfNjo+M3AHcEMUqnCikJe3vsyKwpWUOWqJm5IHEX2uRRWfI3+XfpHOryQKKwoJKXvcrbWrWPMy/UkOZOYP2Q+C7MXxux7VqHhz9BNX6Cw0fMi37aW3Av8o9HzBBHJFZHNInJzcw1E5D7fMbnFxcV+hKQ6k3VF67jl7Vt448AbVLmqAIPY69hS/A9uefsW1hWti3SIQdfwngtdH4KtDoOhylXFGwfeiNn3rEInqBdjReQOIAf4caPNA40xOcCXgJ+JyNCm7YwxLxhjcowxORkZGcEMSXVwhRWFLFqziFp3LW7jvmif27ipddeyaM0iCisKWzhDx9P4PXvxXLQvVt+zCi1/Ev1xoPHnxH6+bRcRkTnAY8ACY0xdw3ZjzHHf9wJgNTCuHfGqTublvS/j9rgve4zb4+aVva+EKaLQ64zvWYWWtLaWhog4gAPAdVgJfhvwJWNMXqNjxmFdhJ1njDnYaHt3oNoYUyci6cAm4CZjzF5akJOTY3Jzc9vxllQsmfTaJN9wzeWleA2bzsbGfPNJ3YUqP4qcpzhT2PSlTWGISHUEIrLdN3pyiVYvxhpj3CLyIPAuYAdeMsbkiciTQK4xZjnWUE0K8FdflZ9jvhk2I4HfiIgX69PDDy+X5JVqqtrl35rsVTaBITNDHE14VJev9es4f/4AKgX+zbrBGLMSWNlk2xONHs9pod1GYHR7AlSdW5Izya+EluxMgZuXhCGi0Evy81NMsrPj18hV4aF3xqqoNn/IfBxy+f6IQxzMHzI/TBGFXmd8zyq0NNGrqLYweyEOeytJz+7gruy7whRR6HXG96xCSxO9imr9u/Rn8YyfkmDA0WTigEMcJDgSWDxzcUzdQNS/S38Wz1xMgiPhkp69wxgSxB5z71mFliZ6FfWmlxezrOgEQ2szsZkEBCHFmcKtWbey7MZlTO83PdIhBt30ftNZduMybs26lRRnCoJgPPFc5+rOshMlTO82ItIhqg6k1emV4abTK9VFPG741UQ8NidZRY/xwOzhLPpUVqSjiohP/3wdw50nea74qzDxazDvB5EOSUWRy02v1B69im67/wKl+ezJfACPscVk2UB/TR6axt9PpOIe/UXY9juoOBHpkFQHoYleRS93Paz5IfQey1u1VxPvsDFuQLdIRxUxU4amUe/2smvo/WC8sPYnkQ5JdRCa6FX02vEKnDsGsx9n0+EycgZ1J8Fpj3RUETNhcA/sNmH16US4+i746GU4eyTSYakOQBO9ik6uGqvH2n8SZb1nsO9kxYWC2Z1VaoKT0X27WgXDZ3wHxA5rnm29oer0NNGr6JT7ElSehNmPs/lwGQCTOvH4fIPJQ9PYVXiOqvieMP4rsOvPUHKw9YaqU9NEr6JP3XlYtxiGXAuDp7PxUAnJcXbG9Osa6cgibsrQNNxew7YjZTDtP8CRCKv/J9JhqSiniV5Fn62/geoSmPU4YNVMnTC4B067/rjmDOyB0y5WecGUDJh0P+x5A07tiXRoKorpb46KLjXnYMPPIWse9B/P6YpaDhVXdfrx+QaJcXbGDejOpgJfHdkpX4f4rvChzqlXLdNEr6LLpiVQWw6zHrOe+gpjT47xQuBtMWVoGnuOl1Ne7YLE7lay3/93OL490qGpKKWJXkWPqlLY/CvIvhl6jwGsRN810Ul27y4RDi56TB6ShtfAlsO+Xv2k+yGxB6x6JrKBqailiV5Fjw0/A1c1zPrehU0bC0qYNKQHNj8qLnUWYwd0I8Fp+9fwTXyqdWH20AdwdGNkg1NRya9ELyLzRGS/iOSLyCPN7F8kIntFZLeIfCAiAxvtWygiB31fC4MZvIohladg629h9BcgYzgAhWXVFJbV6Ph8E/EOO+MH9WBjfum/No7/CqT0glVPQ5StX6Uir9VELyJ2YAlwA5AN3C4i2U0O2wHkGGPGYNWOfdbXtgfwn8BEYALwn746skpdbN1i8Lrg2ocvbNp4qASwxqTVxSYPTWP/6UqKK+usDXFJMP07cHQDFKyOaGwq+vjTo58A5BtjCowx9cBS4KbGBxhjPjTGNBT33Az08z2+HnjPGFNmjDkLvAfMC07oKmacK4Ttv4dxd0CPIRc2r88vpWdqPMN6pkQwuOg0bZj1KafhjyEA1yyErv1h1VPaq1cX8SfR9wUKGz0v8m1ryb3AP9rSVkTuE5FcEcktLi72IyQVU9b6buOf8d0Lm7xew8b8EqYOS8dXcF41cmWfrnRJcFw8fOOIh5kPWbNvDrwTueBU1AnqxVgRuQPIAX7clnbGmBeMMTnGmJyMjIxghqSiXekh2PEnyLkHuva7sHn/6UpKq+qZOkzH55tjtwlThqazPr+Ei2pKXHW79alo1TPg9UYuQBVV/En0x4HGNcv6+bZdRETmAI8BC4wxdW1pqzqx1T8EexxMW3TR5g351pDE1GE6Pt+SqZnpHD9Xw9HS6n9ttDvh2kfh9Mew9/8iF5yKKv4k+m1ApogMFpE44DZgeeMDRGQc8BusJH+m0a53gbki0t13EXaub5tScGYffPxXmHgfpPa6aNf6/BKGZiTTu2tihIKLfg3j9OvzSy7eMepzkDHCulvW445AZCratJrojTFu4EGsBL0PeN0YkyciT4rIAt9hPwZSgL+KyE4RWe5rWwY8hfXHYhvwpG+bUlYiikuBqd+6aHO928uWgjIdtmnFoLQk+nRNuPDp5wKb3bqzuPSg9YdUdXqO1g8BY8xKYGWTbU80ejznMm1fAl4KNEAVo07shH3LYeYjkNTjol07C89R4/Joom+FiDB1WDrv7TuNx2uwN76pbOSNcMUYa2XLUZ8DR1zkAlURp3fGqsj48BlI6AaT//2SXevzS7CJrj/vj2mZ6ZyrdrH3RMXFO0Rg9vfh3FHY+WpkglNRQxO9Cr9jW+DgP2HqNyHh0jXmN+SXMLpfN7omOiMQXMfSsNjbJeP0AJmfgn4TYM2PwVUb5shUNNFEr8Lvw6chOQMmfvWSXZW1LnYWnmOazrbxS8/UBIb3Sr10nB6sXv1134fKE9YNaarT0kSvwqtgDRxeC9O/DXHJl+zeUlCGx2uYquvb+G3qsHS2HSmj1uW5dOfgGdbXup9CfVX4g1NRQRO9Ch9jrLH51D5wzf9r9pC1B4tJdNq5ZpAuieSv6Znp1Lm9bD3cwoS2WY9DVTFsfSG8gamooYlehc/B96BwC8z8LjgTmj1k7YFiJg3pQbzDHubgOq6JQ3oQZ7ex9kALy4cMmAiZc2H9z6yiLqrT0USvwsMYa7GtbgNh3J3NHnKstJojpdXMyNJlMNoiKc7B+MHdWXvwMutEzX4cas/Bpl+FLzAVNTTRq/DY9zac2m3dnm9vfjbNGl+i0kTfdjMyMzhw+jwny2uaP6D3VTBygVWqsVrvWexsNNGr0PN6rLtg07NgzBdaPGztgWL6dktkSPqlF2nV5TX8cVx3oJnZNw1mfQ/qz1vF11Wnoolehd6eZVC8z+rN25ofe3d5vGw6VMqMrAxdljgAI65IpWdq/IVPRc3qORJGfx62/AYqT4cvOBVxmuhVaHlcsPoH0Gu0VfS7BR8dPcv5Ojczs3RaZSBEhOmZGaw/WILHe5miI9c+Ap56WL84fMGpiNNEr0Jr15+hrABmPwa2ln/c1h4sttZY1/VtAjYjK53yGhe7i861fFDaUBj3Zch9CcqLwheciihN9Cp03HWw5lnoew1kXb6C5NoDJVw9oBtdEnTZg0BNz8xAxPq3vKwZD1nf17apPpDqwDTRq9D56BUoL7Sm9l1m3L30fB17TpQzI1Nn27RHj+Q4xvTtevlplgDd+sM1d8OOV61PWyrmaaJXoVFfbfUYB06FIbMue6hVDk+nVQbDjKwMdhw7S3m16/IHTv822JzWJy4V8/xK9CIyT0T2i0i+iDzSzP4ZIvKRiLhF5NYm+zy+YiQXCpKoTmDb7+D86VZ78wCrPjlDWnIco/peupKlaptrh2fgNbTeq0+9Aib8G+z+CxTvD09wKmJaTfQiYgeWADcA2cDtIpLd5LBjwN3Aa82cosYYM9b3taCZ/SrW1FXC+udg6HUwcMplD3V7vKzeX8y1w3teXDhDBWRs/+50T3Ky6pMzrR889VvgTLLucVAxzZ8e/QQg3xhTYIypB5YCNzU+wBhzxBizG9Cy8wo2/xpqyqyZNq3YUXiO8hoX143sGYbAYp/dJswa3pPV+89cfpolQHIaTPp3q4j4yd3hCVBFhD+Jvi9Q2Oh5kW+bvxJEJFdENotIsxOpReQ+3zG5xcWtfORU0a3mLGz8BQz/tDXbphUf7DuDwyZMy9RplcEya0RPzla72Fl4tvWDJz9gFX/RXn1MC8fF2IHGmBzgS8DPRGRo0wOMMS8YY3KMMTkZGXpBrkPb+Auoq7CKU/th1SenmTC4h06rDKIZWRnYbcIH+/wYvknsBlO+AQf+AYXbQh+cigh/Ev1xoH+j5/182/xijDnu+14ArAbGtSE+1ZGcL7aGbUbdAleMavXwwrJqDpw+z+wROmwTTF0TnYwf1N2/cXqAifdDUrpV+UvFJH8S/TYgU0QGi0gccBvg1+wZEekuIvG+x+nAVGBvoMGqKLfhZ+Cusda08cOH+61EpIk++GaP6Mknpyo5fq6F1Swbi0+B6YugYDUcXhfy2FT4tZrojTFu4EHgXWAf8LoxJk9EnhSRBQAiMl5EioDPA78RkTxf85FArojsAj4EfmiM0UQfiypOWFMqr7od0jP9avLBvjMMTk9mSEZKiIPrfGaP6AXgf68+5x5I7W1VADOtXMRVHY7Dn4OMMSuBlU22PdHo8TasIZ2m7TYCo9sZo+oI1v7EWo545sN+HV5d72ZTQSl3ThoY4sA6p6EZyQxMS2LVvtP+/Rs7E2HGd+HviyD/A8icE/ogVdjonbGq/c4esZY7uPpO6O5f4t6QX0q928t1OmwTEiLC7BE92XiolJr6ZoqGN2fcndBtgFUJTHv1MUUTvWq/NT8GsVk9Qj+9t/cUqfEOcgb1CGFgndt1I3pR5/ayrrW7ZBs44mDmI3ByJ3zy99AGp8JKE71qn5J82PUajP8KdOnjVxO3x8t7e09z3ciexDn0RzBUJg7pQddEJ+/knfK/0ZgvQtowa6zeq/c/xgr9LVPts/oH4EiEaf/hd5OtR8o4W+3i+iuvCGFgymm3MWdkL97fexqXx8+kbXdYs6bO7IW8ZaENUIWNJnoVuFN7YM8bMOl+SPH/Rrd395wi3mFj5nC9OS7U5o26gopaN5sLSv1vdOUt0PNK625Zjzt0wamw0USvArf6fyC+K0z5ut9NvF7Du3mnmZmVQVKcX5O+VDtMz0wnKc7OO3vaMHxjs1nrFJUdsiqEqQ5PE70KzPGP4JMV1lopid39brar6BynKmqZN0qHbcIhwWln1vCevJt3uvVFzhob/mnoM85ar95dH7oAVVhooleB+fAZSOwBk77Wpmbv5J3CYROu893Qo0Lv+lFXUHK+jh3H/FjkrIGIVUug/Bh89HLoglNhoYletd3RTZD/Pkz7FiR08buZMYZ395xi8tA0uibpImbhMmt4BnF2W9uGb8CqJzBgsnUznMuPpRRU1NJEr9rGGOuGmpReMP7f2tR0/+lKjpRW67BNmKUmOJk6LI138k5h2nIjlAjM/j6cPwXbXgxdgCrkNNGrtilYDUc3wPTvQFxSm5q+s+cUIvCpbB22Cbd5o66g6GwNeScq2tZwkK/m7/rFVuUw1SFpolf+MwZWPQ1d+sE1C9vY1LB81wnGD+pBz9SEEAWoWjI3+wocNuHtXSfa3nj241Bdai1BrTokTfTKfwfeheO5MPO74IhvU9O8ExUUFFdx89i2FCdTwdI9OY6ZWRks33UCb1tm3wD0y4GsG6yiMjVtuKCrooYmeuUfr9fqzXcfDGO/3Obmb+08jtMu3KDj8xGzYGwfTpbXsu1IWdsbz/oe1JXDxl8GPzAVcprolX/2vQWnP7Z+4e1tmzHj9Rre3nWSmVkZdE+OCyKVaoMAABd5SURBVFGAqjWfyu5FotPOW4EM3/QeA1d+FjY/D1UlwQ9OhZQmetU6r8e6HT5jBIz6XJubbz1SxqmKWhbosE1EJcU5mHtlL1Z+fJJ6dwALll37qFVBbP1zwQ9OhZRfiV5E5onIfhHJF5FHmtk/Q0Q+EhG3iNzaZN9CETno+2rbFTwVHXa/DiUHrN68zd7m5m/tPEFSnJ05I3Xt+Ui7aWwfzlW7WJ/v59LFjWUMt1a33PY7qDgZ/OBUyLSa6EXEDiwBbgCygdtFJLvJYceAu4HXmrTtAfwnMBGYAPyniPh/v7yKPI8L1vwQrhgDI25sc/N6t5eVH59kbnYvXdsmCkwblkG3JCdv7Qxg+AZg5kPgdcPaHwc3MBVS/vToJwD5xpgCY0w9sBS4qfEBxpgjxpjdQNPPg9cD7xljyowxZ4H3gHlBiFuFy45XrQpSsx+3Frtqo7UHiimvcXGTDttEhTiHjU+P7s0/805TXR/AypQ9hsC4O6yKYmePBj9AFRL+/Ob2BQobPS/ybfOHX21F5D4RyRWR3OLiAD5SqtBw1Vo9t34TIHNuQKd4c8dxuic5mZaZHuTgVKBuuqoPNS4P77alIEljMx6yKoqtfTa4gamQiYqLscaYF4wxOcaYnIwMXaM8amz/A1Qct3rzIm1uXlZVzz/3nuKz4/rhtEfFj5oCxg/qwcC0JP6yrbD1g5vTtS/k3AM7/2xVGFNRz5/fvuNA/0bP+/m2+aM9bVUk1VfBup/AoOkwZGZAp1j2UREuj+GL4/u3frAKG5tN+EJOfzYXlHG4pCqwk0xfZN00t/p/ghucCgl/Ev02IFNEBotIHHAbsNzP878LzBWR7r6LsHN921S02/oCVBVbi1oFwBjD0m2FjBvQjeFXpAY5ONVen7+mH3abBN6rT+kJE79qVRg7nRfc4FTQtZrojTFu4EGsBL0PeN0YkyciT4rIAgARGS8iRcDngd+ISJ6vbRnwFNYfi23Ak75tKprVlsOGn1vj8gMmBnSKj46dJf/MeW7T3nxU6tklgdkjevK37UX+15Ntaso3ID7VusdCRTW/Bk6NMSuNMVnGmKHGmGd8254wxiz3Pd5mjOlnjEk2xqQZY65s1PYlY8ww39fvQ/M2VFBtft5a02TW9wI+xdKthSTH2Zk/pk8QA1PBdNv4/pScr2PVJ2cCO0FSD6vC2Ccr4MSO4AangkqvkKmLVZfBpiUw8karlFwAKmtdrNh9khuv6kNyvM6dj1YzszLo1SU+8OEbgEn/bpWSXPVM8AJTQaeJXl1sw8+tdcdnPRbwKd7edZIal0cvwkY5h93G56/pz+r9ZzhZHmAFqYQuMPVbkP8eHNsc3ABV0GiiV/9SeRq2/AZGfx56jgzoFMYYXt18lOG9Uhnbv1uQA1TB9oWc/hjgz1uOBX6SCf8GyT3hg6esmgUq6miiV/+y/jnw1MO1lyxn5LfNBWXsPVnB3VMHIQHMvVfhNSAtietG9ORPW45R6/IEdpK4ZJj+bTi6Hg6vCW6AKig00StL+XHIfRHGfgnShgZ8mhfXH6ZHchyfHadLHnQU90wbTGlVPW/tbMctLtfcDV36aq8+SmmiV5a1P7Z+QWc+FPApjpRU8cEnp/nyxAEkONu+yqWKjMlD0hjZuwsvrj/ctuLhjTkTYMZ3rQpkB/8Z3ABVu2miV1B2GHb8EXL+H3QbEPBpfr/hMA6bcOekgUEMToWaiHDvtMEcOH2e9fntKCoy7g7oPghWPWVVJFNRQxO9gjU/ApvDGmcNUHmNi79uL+LGq/rQs4sW/+5obryqN+kp8by4/nDgJ7E7reIkpz6Gff7ePK/CQRN9Z1e8H3b/xZo5kRp4PdelW49RXe/h3mmDgxicCpd4h527Jg9k9f5i8s9UBn6i0Z+H9OHW3bLeAC/uqqDTRN/Zrf4fcCZZc6EDVOvy8NKGw0weksaVfboGMTgVTl+eOIB4h43nVxcEfhKb3bqjumQ/fPy34AWn2kUTfWd26mPIexMmfQ2SA18v/rUtxzhdUcc3rssMYnAq3NJS4rlj0kDe3FEU+KqWACMXwBWjYfUPrAplKuI00Xdmq56BhK4w+cGAT1Hr8vD8mkNMGtKDyUPTghicioT7Zw4lzmHjFx8cDPwkNpu16unZI7DzT0GLTQVOE31nVZQLB/5hrUCYGPgdrK9uPkpxZR3/MScriMGpSMlIjeeuyYP4v53HOVR8PvATZc6FfuNhzbNWpTIVUZroO6tVT0FSGky8P+BTVNe7+fWaQ0wdlsbEIdqbjxX3zRhCvMPO/7anVy9iVSarOA4fvRy84FRANNF3RofXQcFqmLYI4lMCPs2rm49Scr5ee/MxJj0lnrumDGT5rhMcPN2OGTiDZ1oVytb+xKpYpiJGE31nYwysehpSe8P4ewM+zbnqep5ffYjpmenkDOoRxABVNPjqjKEkOe386J39gZ+koVdfdQa2/jZ4wak28yvRi8g8EdkvIvkicsmKVyISLyJ/8e3fIiKDfNsHiUiNiOz0ff06uOGrNsv/AAo3w4zvgDMx4NM8994BymtcfO/Tga1yqaJbj+Q4Hpg9jPf3nWbtgeLATzRgEgz7FGz4GdRWBC9A1SatJnoRsQNLgBuAbOB2Eclucti9wFljzDDgOeBHjfYdMsaM9X0FPiCs2s8Ya2y+6wAYd1fAp/nkVAV/3HyUOyYNZGTvLkEMUEWTe6cNZmBaEv/9dl7g5QYBZj9mVSzb/HzwglNt4k+PfgKQb4wpMMbUA0uBm5occxPQcMXlb8B1omvURp9PVsDJnXDtw+CIC+gUxhj+a3keXRKdLPqUjs3HsniHnSfmZ3OouIqXNx4J/ER9xsGI+bDpl1YFMxV2/iT6vkDjWmNFvm3NHuMrJl4ONEzDGCwiO0RkjYhMb+4FROQ+EckVkdzi4nZ8TFQt83qsefNpw2DMbQGfZuXHp9hcUMZ35g6nW1JgfyxUxzF7RE+uHZ7Bz98/SHFlXeAnmvWYVbls4/8GLzjlt1BfjD0JDDDGjAMWAa+JyCWf9Y0xLxhjcowxORkZGSEOqZPKexOK91mLTtkDq+NaUevimb/vZWTvLtw+IfBVLlXHISJ8f342tW4PP1i5L/AT9cqG0bdaFczOB1iMXAXMn0R/HGhc/LOfb1uzx4iIA+gKlBpj6owxpQDGmO3AIUA/74ebx20tMtUzG668JeDTPPn2Xk5V1PKDz47CbtORuc5iaEYKX7t2GG/uOM47e04GfqKZj4C7zqpkpsLKn0S/DcgUkcEiEgfcBjRdg3Q5sND3+FZglTHGiEiG72IuIjIEyATasWKSCsjupVB2yPr4bAvsQ9y7eaf42/YiHpg1jHEDugc5QBXtvj57GKP6duF7b+4JfAgnfRiMvR22vWhVNFNh0+pvvW/M/UHgXWAf8LoxJk9EnhSRBb7DXgTSRCQfa4imYQrmDGC3iOzEukh7vzFGr8aEk7sOVv/Id0HsMwGdouR8Hd9b9jFX9unC12frwmWdkdNu47kvjOV8nZtHl+0OvBLVzIfBeK2KZips/BqsNcasBFY22fZEo8e1wOebafcG8EY7Y1Tt8dErUH4MbnzOuoGljYwxPLrsYyrr3Pz5i2OJc+g9dp1VZq9UHrp+OE//fR+v5xbyxfEBXKfpNsCqL7v99zD1m9BD6xeEg/7WxjJXjXX7+YDJMPS6gE7xq9WHeG/vaR66fjhZvVKDHKDqaO6ZOpgpQ9N44q08dhedC+wk079tVTRb82xwg1Mt0kQfy7a9COdPWbehB9Cbf3/vaX7yz/0suKqPVo5SANhswi9uH0d6Sjz3vbKdMxUBrEzZpTeM/4p17aj4QPCDVJfQRB+r6s7D+sUw5FoYNK3NzQ+eruRbf9nJlX268KPPjUHvf1MN0lLi+e1dOZTXuLj/1e3UuQMoGTjtP8CRaFU4UyGniT5WbXkeqkth1uNtblpcWce/vZJLgtPOC3fmkBhnD0GAqiPL7tOFn37hKj46do5H3vgYr7eNF2eT063KZnnLrEpnKqQ00ceimrOw4ReQdQP0H9+mpqXn6/jSbzdzuqKO39x5DX26Bb7wmYptnx7dm+/MzeLNHcd57P8CSPZTvm5VOPvwB6EJUF2giT4WbVoCdeVWkeY2KKuq58u/20Lh2Wpeuns81wzU+fLq8h6YNYwHZw3jz1sLeWL5nrZNu0zsZiX7/SuhaHvoglSa6GNOVYm1SmD2zdB7jN/NzlTWcsfvtnC4pIrf3TVe678qv4gI356bxVdnDuHVzcf4/lt78LSlZz/xfqvS2YdPhy5IpYk+5qx/DlzVberN7z1Rwc2/3MDhkipeuCuHaZnpIQxQxRoR4ZF5Iy4k+6+8vI3KWpd/jeNTrQuzh1bBkQ2hDbQT00QfSypOwrbfwegvQMZwv5q8m3eKW3+9Ea+Bv94/mZlZuqicajsR4dEbRvL0zaNYe7CEzz2/kWOl1f41Hv8VSLnCqnwW6B236rI00ceSdT8Fr9tab74VtS4PT63Yy/2vbiezVyrLH5zKqL5dwxCkimV3TBrIH++ZwOmKOm785Xre2nm89XF7Z6JV8ezYRjj0QXgC7WQ00ceKc8dg+x9g3B3QY8hlD91x7Cyf/t91vLj+MHdMHMhf7ptEzy4J4YlTxbwpw9J564GpDMlI5ptLd/Lvf/qI0vOtLIR29V1W5bNVz2ivPgQ00ceKNT+y7n6d8d0WDzlTUcujyz7mc89vpLbew6v3TuSpm0eR4NR58iq4BqUn87f7p/DwvBF8sO8Mcxav4aX1h6l3t1CS0BEPMx+CEx9Zs3BUUEnAq9CFSE5OjsnNzY10GFGvsKKQl/e+zIqCFVS7qknyepifNIiFN/ya/l36X3Tsuep6Xlp/mN+uO4zL4+WOSQNZNDeLLgnOCEWvOpP9pyr577fz2HiolP49EvnO3OF8ZnRvHPYm/UyPm8Lnx/NygrAiKc76uXYmMX/IfBZmL7zk51pdTES2G2Nymt2nib7jWVe0jkVrFuH2uHEb94XtDrHjsDtZPHMx0/tN58DpSn6/4Qhv7iii1uVl/pjefPf64QxMS45g9KozMsaw5kAxP/zHJ3xyqpI+XRO4c/Igbhvfn+7JVknKdUXrWPThN3F76nE3WnLDIQ4cdseFn2vVvJhP9Jf0bmO4F1BYUcgtb99CrbvlxaQcEk+f84+Td8xJvMPGzWP7cvfUQYzsfUkVR6XCyuM1fLDvNH/YeISNh0qJc9iYNTyDScPhVwe+Rq2n5Z/rBEcCy25cFpO/08HIXzGd6Fvu3cZmL+DpzU/zxoE3LnqvTRljI7V+GncPX8Tnc/rTI1mLeKvo88mpCpZuLWTlxycpT/oLzu5bEWlhDB/rd/rWrFt5bNJjYYwytIKZv9qd6EVkHvBzwA78zhjzwyb744FXgGuAUuCLxpgjvn2PAvcCHuAbxph3L/dabUn0/vRuw94L8LisG5ZctdZ3d621LryrBtw11vaGbe7a5o+9qI312OuqwV1XzYxubqr8qNea4kxh05c2heENK9U+Xq9h4muTqPW0Pu8+2cCm2m6IMxEcCdbUTGeitRKmM6HJ96aPL9fG92UL38SEYOevyyX6VitM+Wq+LgE+BRQB20RkuTFmb6PD7gXOGmOGichtwI+AL4pINlaN2SuBPsD7IpJljAlgXdNLvbz3Zdyelnu2AG6Pm1c+fpHHxj54+aTaXEK+8LimmWOrm0naNRDgW/OKE7c9HrctnjqJp444aoyTKq+TCreTSk8itXSlqnuRX+erclUFFIdS4WazCXWeGr+OrQZWHbfTxVFDiqOCZHGRIC4STB1OU4fDW4fdW4ct0BRjj2v0ByCh5T8OLf2h8OuPj3Xcy3l/8C9/7X2l3Z9i/CklOAHIN8YUAIjIUuAmoHGivwn4L9/jvwG/FGsB85uApcaYOuCwr6bsBCAoXc0VBSsuO4QB4DZuVhz4K4+917bK8x5s1BF3IenW4fsucdQRRy3x1JFKLfHU4qTWxFFDPDW+x1XmX4m6xjipJY5aE3fR8bVYbWqJw9topmucw0Z6chw9UuPomZpA764J9OmWSL/uiSTuuYUaf3o+Tr3gqjqOJGeSX52TOHsSu6b/hhPltZwsr6Gksp7SqjrKqur51xI7BiceEqgngToSpJ4EXL7n9b7n9STb6kkSF0niItFWT6LUk+h1keCqt74a2lNPPJXEU0q8qSfetz3O1BPv2x+IFQP74bZdfoa727hZUbAiLIm+L1DY6HkRMLGlY4wxbhEpB9J82zc3adu36QuIyH3AfQADBvhfh7La5d8t1lVi4/WMb1Av8dTb4nHZ4nFJHC6Jp96WgEvifN/jLzz2igMRoWGQRATEitX32PpuF8H6vxLstobngtMmJNhsZNgEp92G0yE4bTbinTbiHTbiHXYS4+wk+b5SE5x0SXCSmuAgKc7eYqGPj2tvbHWM3iEO5g+Z7/e/o1KRNn/IfL9+rj+buYBFky5d3sPrNVTWuamocVFR66KqzkN1vZvqeg+1Lg91bi91vu9ur6He7cXl8eIxBo/HUOk1lBuDx2vwGoMx4DUGr7Hu3zLGYGj83XpdAxivF6epv/CJIs73OM5be2F7nKfRY9/2KrnsKPYFwfh07ldx8FAzxrwAvADWGL2/7fztBSTHpfCFu58KPMAosjB7IW8degu3+zK/EHYHd2XfFcaolGqf9v5c22xC10QnXRM7zr0hP3xtvX/5Kwifzv25M/Y40PhKQD/ftmaPEREH0BXroqw/bQM2f8h8HHL5v1Wx1rvt36U/i2cuJsGRcMl7d4iDBEcCi2cujrkpaCq2dcaf63DmL38S/TYgU0QGi0gc1sXV5U2OWQ4s9D2+FVhlrOk8y4HbRCReRAYDmcDWdkftszB7IQ57K/9QMdi7nd5vOstuXMatWbeS4kxBEFKcKdyadSvLblwWU9NJVefR2X6uw5m//J1e+WngZ1jTK18yxjwjIk8CucaY5SKSAPwRGAeUAbc1unj7GHAP4Aa+ZYz5x+VeS+fRK6U6i6iaRx9Ogd4Z+8reV1hRsIIqVxXJzmTmD5nPXdl3xdRHPaVU7AlW/or5RK+UUp3d5RK9LlOslFIxThO9UkrFOE30SikV46JujF5EioGjkY4jAOlASaSDCDN9z52DvueOYaAxJqO5HVGX6DsqEclt6UJIrNL33Dnoe+74dOhGKaVinCZ6pZSKcZrog+eFSAcQAfqeOwd9zx2cjtErpVSM0x69UkrFOE30SikV4zTRh4CIfFtEjIikRzqWUBORH4vIJyKyW0TeFJFukY4pFERknojsF5F8EXkk0vGEmoj0F5EPRWSviOSJyDcjHVO4iIhdRHaIyIpIxxIsmuiDTET6A3OBY5GOJUzeA0YZY8YAB4BHIxxP0ImIHVgC3ABkA7f7Ct/HMjfwbWNMNjAJeKATvOcG3wT2RTqIYNJEH3zPAQ9hlZOMecaYfxpzYSHtzVhVxGLNBCDfGFNgjKkHlmIVvo9ZxpiTxpiPfI8rsRLfJfWeY42I9AM+A/wu0rEEkyb6IBKRm4DjxphdkY4lQu4BLltYpoPqCxQ2et5skftYJSKDsIoKbYlsJGHxM6yOmjfSgQRTVBQH70hE5H3gimZ2PQZ8D2vYJqZc7j0bY97yHfMY1sf9P4UzNhVaIpICvIFVHa4i0vGEkojMB84YY7aLyLWRjieYNNG3kTFmTnPbRWQ0MBjYJSJgDWF8JCITjDGnwhhi0LX0nhuIyN3AfOA6E5s3ZoS0yH20EhEnVpL/kzFmWaTjCYOpwAJf6dQEoIuIvGqMuSPCcbWb3jAVIiJyBMgxxnS0FfDaRETmAYuBmcaY4kjHEwoi4sC60HwdVoLfBnzJGJMX0cBCSKzeystAmTHmW5GOJ9x8PfrvGGPmRzqWYNAxetVevwRSgfdEZKeI/DrSAQWb72Lzg8C7WBclX4/lJO8zFbgTmO37f93p6+mqDkh79EopFeO0R6+UUjFOE71SSsU4TfRKKRXjNNErpVSM00SvlFIxThO9UkrFOE30SikV4/4/WFWzQOVawO0AAAAASUVORK5CYII=\n",
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ]
          },
          "metadata": {
            "tags": [],
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "dUlHMnJLLOrq",
        "colab_type": "text"
      },
      "source": [
        "## Run Scipy's interpolation\n",
        "We now have many more options. We're not just linear, and we can control the degree of the polynomial"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "1CuAqJmWKktw",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 265
        },
        "outputId": "8ed30a40-bc28-4c39-8334-07d755b6e25d"
      },
      "source": [
        "pLinear = scinterp.interp1d(nodes,gaussian(nodes)) # this returns a function\n",
        "pQuad = scinterp.interp1d(nodes,gaussian(nodes),kind='quadratic') # this returns a function\n",
        "plt.plot(xGrid,gaussian(xGrid),'-',label='Original function');\n",
        "plt.plot(xGrid,pLinear(xGrid),'-',label='scipy.interp1d, linear')\n",
        "plt.plot(xGrid,pQuad(xGrid),'-',label='scipy.interp1d, quadratic')\n",
        "plt.plot(nodes,gaussian(nodes),'o',label='Data',markersize = 10);"
      ],
      "execution_count": 19,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3dd3hUZdrH8e8zJZMeSCEhJBA6hC6hioAUQVGwIqgI2N11V9e2ro3Vta6r6Cv2AuqqKLoKKoo0AUWR0GsgBAgJkEIgPZn2vH8MaNRAZmAmk0zuz3XlIjNzzpn7aPKbk+c8RWmtEUII0fgZ/F2AEEII75BAF0KIACGBLoQQAUICXQghAoQEuhBCBAiTv944NjZWp6Sk+OvthRCiUVq3bl2h1jquttf8FugpKSmkp6f76+2FEKJRUkrtP9lr0uQihBABQgJdCCEChAS6EEIECAl0IYQIEG4FulJqrFIqQymVqZS67xTbXaaU0kqpNO+VKIQQwh11BrpSygi8BJwPpAKTlVKptWwXAdwOrPF2kUIIIermzhV6fyBTa52ltbYCc4EJtWz3L+BpoMqL9QkhhHCTO/3QWwEHajzOAQbU3EApdRaQrLX+Sil1jxfrE6JeOJ2ahVsPsetw6S/PpSZGMaZbPEopP1YmhPvOeGCRUsoAPAdMc2Pbm4CbAFq3bn2mby2EV2w7WMzD87exbv9RAJSCE8sEDG4fwyPju9ExPsKPFQrhHncCPRdIrvE46fhzJ0QA3YHvjl/JJAALlFLjtda/GQqqtX4deB0gLS1NVtYQfqW15qlvdvLGyiyahwbxzOU9ueysJAwGhcOp+eDnbP6zKIPzX1jFX0Z05PZRHf1dshCn5E6grwU6KqXa4gryScBVJ17UWhcDsSceK6W+A+7+fZgL0dC8tDyT11ZkcWVaMveP60pUiPmX14wGxZSBbRjXoyWPfLGNmUt2ERViYtrZbf1YsRCnVmega63tSqnbgEWAEXhba71NKfUokK61XuDrIoXwti83H+Q/3+7ikj6teOqyHidtJ48OC+K5ib2psDp49MvttI4JZUSX+HquVgj3KH+tKZqWlqZlci7hDxuyjzLp9Z/omRTFf28YgMVkrHOfCqudK179kX2F5Xxy62C6toysh0qF+COl1Dqtda1jfWSkqGhSSqps3PzeOuIjg3ltSppbYQ4QGmTiran9iAg2c9N76VTZHD6uVAjPSaCLJuXFpbspKKtm1lV9iA4L8mjfhKhgnp/UmwNFlby2IstHFQpx+iTQRZOxp6CM2T/sY2LfZHomNTutYwxsF8O4Hi15ZUUmuccqvVyhEGdGAl00GY99uZ0Qs5G7x3Q+o+P844IuaA1Pfb3TS5UJ4R0S6KJJWL4zn+UZBfx1ZEfiIixndKyk5qHcPKw9X2w6yM97i7xUoRBnTgJdBDy7w8m/vtxOu7gwpg5O8coxbx3WnsSoYB75Yhv+6ikmxO9JoIuA9/XWw2QVlnPvmC4EmbzzIx8SZOTO8zqz7WAJ3+0q8MoxhThTEugioGmteX1lFm1jwzgv1bsDgsb3SiQhMpg3VkqPF9EwSKCLgLZmbxFbcou54Zy2GAzenTUxyGRg+tkprN5zhK25xV49thCnQwJdBLQ3V2URHRbEZWcl+eT4k/q3JizIyJur5Cpd+J8EughYmfllLNmRz5SBbQg2uzci1FNRIWYm9W/NF5sPcVD6pQs/k0AXAeut77OwmAxMGdTGp+8z/ewUAGb/sNen7yNEXSTQRUAqKrfy6fpcLuubRGz4mfU7r0tS81DG9WjJhz8foLza7tP3EuJUJNBFQJq/MRer3cm1Pr46P+HaQW0oq7azcMuhenk/IWojgS4C0rz0HHq0iqJLQv1Mc9u3TXPaxYYxb11OvbyfELWRQBcBZ9vBYrYfKuGKNN/0bKmNUorL+ibx894i9h8pr7f3FaImCXQRcOal5xBkNDC+V2K9vu+lZ7XCoOATuUoXfiKBLgKK1e5k/sZcRqfG0yzUs/nOz1TLqBCGdIzj03U5OJ0yv4uofxLoIqAs25nH0Qobl9djc0tNV/RN4mBxFav3HPHL+4umTQJdBJR56TnER1oY2jHOL+8/OjWeyGAT89Yd8Mv7i6ZNAl0EjILSar7bVcAlfZIwenneFncFm42M753IN1sPU1Jl80sNoumSQBcB4+uth3A4NZee1cqvdVx6VhLVdidLtuf5tQ7R9Eigi4Dx1eZDdGwRTqf4CL/W0Se5GYlRwTLISNQ7CXQREPJLq/h5XxEX9Gjp71JQSnF+j5as3FUozS6iXkmgi4DwzdbDaA3jevo/0MFVh9UhzS6ifkmgi4DQUJpbTpBmF+EPEuii0WtIzS0nSLOL8AcJdNHoNbTmlhOk2UXUNwl00eg1tOaWE6TZRdQ3CXTRqDXE5pYTpNlF1DcJdNGoLdmej9Y0yEAHuKBHAlaHk+8yCvxdimgCJNBFo7Z4+2FaR4fSKT7c36XUqndyc2LDg1gs7eiiHkigi0arvNrOD3uOMKprPEr5Z+6WuhgNihFdWvBdRj5Wu9Pf5YgAJ4EuGq1Vuwux2p2MSm3h71JOaVTXeEqr7KzdV+TvUkSAk0AXjdaSHXlEBpvolxLt71JOaUjHWCwmgzS7CJ+TQBeNksOpWbYznxFdWmA2Nuwf49AgE+d0jGXJjjy0lpWMhO807N8EIU5iQ/ZRisqtjEqN93cpbhnVNZ6co5Vk5JX6uxQRwNwKdKXUWKVUhlIqUyl1Xy2v36KU2qKU2qiU+l4pler9UoX41eLteZiNiqGd/LMykadGdHW18y/eJs0uwnfqDHSllBF4CTgfSAUm1xLYH2ite2itewP/Bp7zeqVC1LB4Rx4D28UQGWz2dyluaRERTO/kZizZIYEufMedK/T+QKbWOktrbQXmAhNqbqC1LqnxMAyQhkLhM3sKysgqKGdU18bR3HLC6NR4NuUUk1dS5e9SRIByJ9BbATVXvM05/txvKKX+rJTag+sK/a+1HUgpdZNSKl0plV5QICPnxOlZvjMfgJFdG3Z3xd87Ue+y4/UL4W1euymqtX5Ja90e+Dvw4Em2eV1rnaa1TouLaxxtn6LhWZ6RT6f4cJKah/q7FI90jo8gMSqY7zIk0IVvuBPouUByjcdJx587mbnAxWdSlBAnU1Zt5+e9RZzbuXFdnYNrsq7hXVrw/fEBUUJ4mzuBvhboqJRqq5QKAiYBC2puoJTqWOPhOGC390oU4lerMwuxOTTDG2GgA5zbuQXlVgfp+2XUqPC+OgNda20HbgMWATuAj7XW25RSjyqlxh/f7Dal1Dal1EbgTmCqzyoWTdryjALCLSbSUpr7u5TTMrh9DEFGg8y+KHzC5M5GWuuFwMLfPfdwje9v93JdQvyB1prvMvIZ0iG2wY8OPZkwi4n+baNZvjOf+y/o6u9yRIBpnL8VoknKyCvlUHEV53Zp3DfUh3eOY3d+GQeKKvxdiggwbl2hC9EQLN/paqbwdvt5ha2CQ+WHOFR+iEp75S/Ph5nCSAhPoGVYS0JMIV57v3O7tOCxr3bw3a4Cpgxs47XjCiGBLhqN5Rn5pLaMJD4y+IyOs+fYHtYcWsOG/A2sz19PfkXd3QgTwhLo06IPfVv0ZUDLAaREpZz2+7eLDaN1dCjf7cyXQBdeJYEuGoWSKhvr9h/llmHtTmv/w+WHWbh3IV9lfcWuo7sAiA+Np298Xzo170TLsJa0DGtJeNCvKx+VWktdV+5lh8g4mkH64XS+3vs1AKkxqYxrO44L2l1AbEisR7UopTi3cxwfp+dQZXMQbDae1jkJ8XsS6KJR+H53IQ6n9rj/+Z5je3hry1ss3LsQh3bQM64n/+j/D4YlDyMxLNGjlY601uSU5rD8wHK+2vsVz6Q/w8z1M5nQfgLTu0+nTaT7V9vDu7TgnR/3s2ZvEcMayQRjouGTQBeNwoqMAiKDTfRObubW9tkl2cxcN5Ml2UsIMYVwVderuLLzlR6F7u8ppUiOTObabtdybbdryTqWxYc7P+R/u//HZ5mfMSZlDHecdQeJ4Yl1HmtQuxgsJgMrdxVIoAuvkUAXDZ7WmhW7ChjSMRZTHd0VK+2VvLH5DeZsm4PZYObmnjdzdderaR7s/X7r7Zq144GBD3Bzr5t5d/u7fLjjQ5ZnL+eGHjcwrfs0LEbLSfcNNhvp3zaaFbsKeMjrlYmmSrotigZvd34Zh0uq6ryS/fnQz0z4fAJvbHmDMSlj+PKSL7mtz20+CfOaYkNiubPvnSy4eAHnJJ3DrI2zuHT+pWzM33jK/YZ1iiMzv4zcY5Wn3E4Id0mgiwZvxfFRlSdbzMLmsDFz3Uxu+PYGLEYLc8bO4clzniQutH6bMlqGt+S54c/x2ujXcGgH076ZxiubXsHutNe6/YkPqJW7ZNSo8A4JdNHgrdxdQKf4cFpG/bEveG5ZLtd8fQ1vb32byzpdxkcXfkTf+L5+qPJXgxMHM++ieZzf9nxe3vgy1y26jrzyPy5s0aFFOIlRwRLowmsk0EWDVmG1syariKEd/3i1vT5vPZO/nMyB0gM8P/x5ZgyaQai5YUypGxEUwZPnPMlT5zxFRlEGV311FdsKt/1mG6UUwzrH8f3uQmwOmX1RnDkJdNGgrckqwupwMqzzbwP988zPuf7b64m0RPLBBR8wss1IP1V4auPajePd89/FaDAy7ZtpfLPvm9+8PrRjHKXVdjYeOOanCkUgkUAXDdqKXQUEmw30S4kGXD1eXt30Kg/98BB94/vy/gXvn9GozfrQObozH477kC7RXbhnxT28t/29X14b3CEWo0FJs4vwCgl00aCt3FXAwHYxBJuNaK2ZuW4mL218iYvaXcQro14hyhLl7xLdEhMSw1tj3mJ0m9H8e+2/eW3Ta2itiQox0ye5GSsk0IUXSKCLButAUQVZheUM7RiHUzt5fM3jzN42mys7X8ljQx7DbDD7u0SPBBmD+PfQfzO+/XhmbZzF8+ufR2vNsE5xbMkt5khZtb9LFI2cDCwSDYo1O5sjs2dTsuALHOUVfGIKIpwLeaHIxkdHFzK9+3T+dtbfPBqy35CYDCb+dfa/CDGF8PbWt9FoRkWcR+XGTzl0zsPkV1ViCA0lcvxFxEyfTlDr1v4uWTQiSmvtlzdOS0vT6enpfnlv0TCVrVxJzu13oG02sP/ad9tpVFgNmi13jGXKdc812jCvSWvN42seJ2Phh/x9vgFtdWLSNXq6mEwos5mkF54nfOhQ/xUqGhyl1DqtdVptr0mTi2gQrNnZrjCvrPxNmAMYHJpgG/Sf9R22Awf8VKF3KaW4u+U13PO5wlht/22YA9jt6MpKcm6/A2t2tn+KFI2OBLpoEI7Mnu26Mj8FbbNzZM479VSR7x2d8w5m56n/2tA2W0Cds/AtaUMXDULJgi/+cGX+B3Y7JZ99QsvzAmN2wpLPPnXvnBcsoOXDMoWXqJsEumgQnBXura/prKyGb/7u42rqh7OqJVD3/QBnebnvixEBQQJdNAiG0FC3gssQFg737q2HinzP8MVonOV1f5AZwsLqoRoRCCTQRYMQOf4ijs375NRNECYTkRMmQGh0/RXmQ5Hjx7t3zuPH119RolGTm6KiQYiZPh1HHUtrKrOZmGlT66egehAzfTrKfOrBUYF2zsK3JNBFg7ApKI9nJmhsJs0fOn6YTKiQEJJeeD6gBtoEtW5N0gvPo0JCwPTbP5YdBk2VGeyP/i2gzln4lgS68LvCykLuXXkvRe3MtJpg5duUfjhCQkEpDOHhNJs4kXbzPw/IATbhQ4fSbv7nNJs4EUN4OFopKkwWwjvaeOZ6C3dV/5fi6mJ/lykaCRkpKvzKqZ3csvgWNuSl80F2Nofb3sWUzT358R8jal3QItBl5pcx6rkVfNZzDaZ9r3BtUiuGJg3j+XOfD4gRsuLMyUhR0WB9sOMDfjz0I/eWOegY1Y7Xy4eedHWipqB9XBitmoXwpm0sPUISuKPazLIDy/g883N/lyYaAQl04TeZRzOZuW4mw8LacHnePqpGPsqafSUM79zC36X5jVKKoZ3iWJFVin3kDKbk7qJ/WGue+vkpDpQGxrQHwnck0IVfWB1W7lt1H+HmUP6ZtQXVfgSr6Y3V4ax1ubmmZFinOMqq7awPPxdDUn8ez96NURm4f9X9J11wWgiQQBd+8tLGl8g4msEjlnbEVpbAeY+zYlchIWYjaSnN/V2eXw3uEIPJoFixuwDGPEFCSR4PRHRjY8FGZm+d7e/yRAMmgS7q3bYj25izbQ6XJo1k+OYv4KypEJ/Kil0FDGrvWp2oKYsMNnNWm+auVYyS+0H3yxm36QvOSxzCK5teIetYlr9LFA2UBLqoVzanjRk/zCAmOIa7DueAORTOfYB9heXsO1LBsE5Nu7nlhGGd4tiaW0JBaTWM+icA/zhWQag5lBmrZ+D8/XS7QiCBLurZnK1zyDiawQMpE4jcvRiG3gXhcazc7VpTUwLd5cR/h1W7C6BZMgz6M7HbPufe9lewsWAjc3fO9XOFoiGSQBf1Zm/xXl7d9CqjW49i5LqPIao1DLgVcC0G3SYmlJRYmYgKILVlJLHhQaw8sXj0kL9BWAsu2vwlZycO5vn1z3Ow7KB/ixQNjgS6qBdaax758RGCTcHcH9IR8rbC6EfAHEy13cHqPUeafO+WmgwGxdCOcazcXYjTqcESASMeRB34mYejBwDw2E+P4a+BgaJhcivQlVJjlVIZSqlMpdR9tbx+p1Jqu1Jqs1JqqVKqjfdLFY3Zl1lfsi5vHX/rdSuxq2ZC8gDodgkA6/YdpcLqkOaW3xnaKY6icitbDx4f+t/nGojvTuKqmfy5x82syl3F8gPL/VukaFDqDHSllBF4CTgfSAUmK6VSf7fZBiBNa90T+AT4t7cLFY1XqbWUZ9OfpUdsDy49mAXl+TDmCTg+lH3FrgLMRsWg9jF+rrRhOadjLErBdxnHm10MRjjvMTiWzVXFxXRo1oGnf36aSnulfwsVDYY7V+j9gUytdZbW2grMBSbU3EBrvVxrfWKm/p+AJO+WKRqzlze+TFFVEQ+kXo/hp5egxxWQ9OtUFN9lFNAvJZowi0zPX1NMuIWeraJc3RdPaH8udBqLedVz3N/zzxwsP8hbW97yX5GiQXEn0FsBNccc5xx/7mSuB76u7QWl1E1KqXSlVHpBQUFtm4gAk1GUwYc7P+TyTpfTbf2HridHzvjl9dxjlWTklTKiS9Md7n8qwzu3YEP2UY6WW3998rzHwF5Jv61fckHbC5i9dTbZJdn+K1I0GF69KaqUugZIA56p7XWt9eta6zStdVpcnLSXBjqtNU/+/CQRQRH8NWEobJkHg/7s6oZ33HcZ+QBNev6WUzm3Swucml+6dQIQ2xHSrof173BXynjMRjNPr33af0WKBsOdQM8Fkms8Tjr+3G8opUYBDwDjtdbV3ilPNGZLs5eyLm8df+l9G82WPQHh8a7udzUs31lAcnQI7eOku2JteraKIiYs6Nd29BOG3weWCFqseJZbet7CypyVrD642j9FigbDnUBfC3RUSrVVSgUBk4AFNTdQSvUBXsMV5vneL1M0NjaHjZnrZtI+qj2X2gxwYA2MeNDV/e64aruDHzILGd6phcz1fRIGg2JYpzhW7CrA4azRRTE0Gob9HfYs5SpzAq3CW/Gf9P/gcDr8V6zwuzoDXWttB24DFgE7gI+11tuUUo8qpU6sXvsMEA7MU0ptVEotOMnhRBMxN2Mu2aXZ3N3nr5iWPgrx3aH31b/Z5ue9RVTaHJzbRZrfTmV4lxYUlVvZnHPsty/0uxGi2xG0ZAZ39Pkru4/uZv6e+f4pUjQIbnUr0FovBBb+7rmHa3w/yst1iUasuLqYVze9yuDEwQw5sAmOZcO1813d7mpYvrOAIJOBQe1i/VRp4zC0YywGBcszCujTusZMlKYgGP0ofHQNYwpz+W9cL17c8CJjU8YSag71X8HCb2SkqPC6Vze9SpmtjLtSp8PKZ6HTWGg3/A/bLc/IZ1C7GEKCmvbsinVpFhrEWa2bs3xnLa2ZXS6ENkNQ3z3JPT1vpbCykLe3vl3/RYoGQQJdeFVuWS5zM+ZySYdL6LR+LtgrXd3sfmdvYTl7C8s5t7M0t7jj3C4t2JJbTH5p1W9fUArGPA4VRfTa/g1jU8by7vZ3Kaws9E+hwq8k0IVXvbzxZYzKyK2tRsH6d1zd62I7/mE76a7omeHHP/hW/L63C0Bib+h9Fax5ldvaXozVYeX1za/Xc4WiIZBAF16z59gevsz6ksmdJxG/8llXj5bhf5j6B4BlO/NpGxsmsyu6KbVlJPGRFpZnnKQT2YiHwGCizU+vcUnHS5i3ax65ZX/oXSwCnAS68JpZG2YRYgrhurD2sGeZq1tdaPQftiurtrMmq4iRMjrUbUopRnRpwcpdhVjttSxuEdkSzr4Dts/n5uh+GDDwysZX6r9Q4VcS6MIrthZuZUn2EqZ2nULzZU9CdDtXt7parNpVgNXhZFRqfD1X2biN6hrv+jDce6T2DQbfBhGJJKz4N5O7TOKLrC9kubomRgJdeMX/rf8/mlmaMaVKQ2GGqzudKajWbRfvyCMqxExam6a9GLSnzu4QS7DZwJLtebVvEBQGo2bAwQ1cb4gl2BjMrI2z6rdI4VcS6OKMrc9bz4+HfuT6zlcRvvJZaDPE1Z2uFg6nZvnOfM7tHIfJKD9+ngg2GxnSIY4lO/JPvrBFj4mQ2IfmK57h2i6TWbx/MRlFGfVbqPAb+Y0SZ+zVTa8SHRzNxLz9UHnU1Y3uJEP512cf5WiFTZpbTtOori3IPVbJzsOltW9gMLjmmi/J5ZriUsLN4by2+bX6LVL4jQS6OCMb8zfy46Efmdb2IkLXvuXqPpfY+6TbL9meh8mgGCqrE52WEV1dN5JP2uwC0GYwdB1P1I8vc3W7CSzev5jdR3fXU4XCnyTQxRl5bfNrNLc058qs9WAwubrPncKSHXkMbBdDZLC5nioMLC0igumV3IwltY0arWn0I+C0M+VgFqGmUOmX3kRIoIvTtqVgC9/nfs+1icMJ3fmVq9tcZMuTbr+3sJw9BeWM6irdFc/E6K4t2HTgGPklVSffKLodDLiZqM0fc1XySBbtW8SeY3vqr0jhFxLo4rS9tvk1ooKimLzjO4hIhMF/OeX2S3e4mglGdpX28zNx4r/fsrqu0s+5G0KjuTZrA8GmYLlKbwIk0MVp2Vm0kxU5K5gS3ZuwQ5tc3eWCTj3D37fb8+gcH0FytMwEeCa6JETQqlkIi0/Vjg4Q0gyG/4Pm+1YzKa4/3+z7RpaqC3AS6OK0vL31bcJMoUzethQS+7i6y51CYVk16fuKGNNNrs7PlFKK87rFsyqzkLJq+6k37jsdYjtzbcYPmJSJOdvm1EuNwj8k0IXHDpQeYNG+RUwMTSGy5KCrm5zh1D9Ki7fn4dQwtvvJ29iF+87v3hKr3Vn7lLo1GU0w5nFij2QxIaID8zPny0yMAUwCXXjsnW3vYFRGrtn5PXQd7+omV4evtx6mdXQoXVtG1LmtqFvfNs2JDbfwzdbDdW/cYRS0H8G03WuwO+38d/t/fV+g8AsJdOGRI5VH+Dzzc8YbY2hht7q6x9WhuMLG6sxCzu+eIGuHeonR4Gp2WZ6RT5WtjnVElYLzHqd1RTGjLfF8lPERpdaTDEwSjZoEuvDI+zvex+qwMi1rAwy42dU9rg5LduRhd2rGdk+ohwqbjvO7J1BhdbBiVy1zpP9efCqcNZXr9m2mzFbGvF3zfF+gqHcS6MJtFbYK5mbMZaQOISUo0tUtzg3fbDtMQmQwvZKa+bjCpmVguxiiQswscqfZBeDcB0h1mhiownhv+3tYHVbfFijqnQS6cNtnmZ9Rai1l+qEsGP4PV7e4OpRX21m5q4Cx3RMwGKS5xZvMRgOjU+NZvCOv9jnSfy88DobexXUHsyisLGTh3oV17yMaFQl04RaH08F/t79Hbzv0jGjr6g7nhuUZ+VTbndLc4iNjuyVQWmVn9R43e64MuJWBwfF0dBp4d9u7J5+1UTRKEujCLd/lfEdOWS5TjhS4Fn02mtza75uth4kJC6Jfyh9XLhJnbkjHWMKCjCza5mazizkYNeoRphwpYPex3aw5vMa3BYp6JYEu3PLeltkk2p2MSBgIHUe7tU+F1c6ynfmc1y0eozS3+ESw2cjIrvF8s/UwNocbzS4A3S7hgubdiHZq3ts627cFinolgS7qtO3INtYVbuKqkhJMY5446Vznv7dkRz4VVgcX9Ur0cYVN20W9EjlaYeP73W42uyiFZcyTXFlcwsqDq9lbvNe3BYp6I4Eu6vTe+pcJdTq5tN0EV/c3Ny3YeJD4SAsD2sb4sDoxrFMcUSFmFmw66P5OSWlMbHUuZq15f6MsgBEoJNDFKeWV57Ho4EourbASMXKG2/sVV9hYsSufC3smSnOLjwWZDJzfPYFvtx2m0lrHIKMaYkc/xrjyKubv+5ri6mIfVijqiwS6OKWP1/wHh9ZclXqNq9ubm77eegibQzOhtzS31IfxvRIptzpYurOOGRhrapbMNe0mUIWTT9P/z3fFiXojgS5OqtpWwSfZixhmUyQP+btH+y7YdJCUmFB6tIryUXWipgHtYmgRYWHBRg+aXYDOI/5JP6uTubs/xe6w+ag6UV8k0MVJfbPynxQpzTXdpoI52O398kqq+DHrCON7t5K5W+qJ0aC4sGci32UUUFzpQTBbIri6w6UcUg6+W/207woU9UICXdRKV5Xw/t6v6KBN9B94p0f7frn5EFq7mgFE/ZnQOxGrw+n+VADHDT/nIRKdivczPgLbKZa1Ew2eBLqo1YZlD7LDbOCqrlej6pjr/PcWbMylW2IkHVqE+6g6UZueSVG0iQll/qZcj/YzmoKY3H4C6WbIWPm4j6oT9UECXfzRsQO8v/9rIjEy7qw/ebTr7rxSNuUUc3HvVj4qTpyMUoqLe7di9Z4j5B6r9GjfSwbcTQiKDzI+gjI3Zm8UDZIEuviDw4vvZ2mIhcs6XEyo2bP1P+ety8FkUFxylgS6P1zeNwmt4dN1OR7tF2WJ4sLWo/kqJIijy/7pm+KEz0mgi9/KSeeDgyvRysCVvW70aFebw8n/1ucwsqJjgb4AAB5/SURBVGsLYsMtPipQnEpydChnd4jh4/QDOJ2eTbx1dZ8/UW1QfJz1BeRt91GFwpck0MWvtKZs0X3Mi4zgvNYjaRXu2VX2sp35FJZZmZiW7KMChTsmpiWTc7SSn7KOeLRf+2btGZIwgA8iI6hedL+PqhO+JIEufrXtM/5XvJMyg2Jqj+s93n1e+gFaRFgY1sn9AUjC+8Z0SyAi2MTH6Qc83ndqzxsoMhr4Kv9n2L3YB9UJX3Ir0JVSY5VSGUqpTKXUfbW8PlQptV4pZVdKXe79MoXP2aqwL5nBf5tH07dFX7rHdvdo9/ySKpZnFHBZ3yRMRrlO8Kdgs5GLe7fi662HPeuTDgxIGECX5p15JzoG56L7wWH3UZXCF+r8zVNKGYGXgPOBVGCyUur3MzRlA9OAD7xdoKgna15hsa2QQwbN1G5TPd790/W5OJyaK/om+aA44amJaclU252eTdiFq6fMtd2mkmWE78uzYf0c3xQofMKdS6n+QKbWOktrbQXmAhNqbqC13qe13gy4OSGzaFDK8tErn2VOfBIpkSkMSx7m0e5aa+alH6BfSnPaxUnf84age6tIuiRE8PFaz5tdxrYdS4vQFrwTnwzLn4DKYz6oUPiCO4HeCqj5U5Fz/DmPKaVuUkqlK6XSCwqkr2uDsfwJ0o0OtusqpqROwaA8azJZvecIWYXlXNmvtY8KFJ5SSjG5f2u25Baz8YBngWw2mLmm6zX8rKrZbi+FVc/6qErhbfXa2Km1fl1rnaa1TouLkxtnDULedlj/Dm+17kJ0cDTj24/3+BCzf9hHTFgQF/Zs6YMCxem69KxWhFtMvLN6n8f7Xt7pcsLN4byV0h3WvApFsghGY+BOoOcCNfuhJR1/TjR2WsO3D7AtrBk/WAuZkjqFYJP7k3ABHCiqYOnOPCb3b02w2eijQsXpiAg2c3nfJL7cfJD8Us/maIkIimByl8ksth1hb5AFlrg/F77wH3cCfS3QUSnVVikVBEwCFvi2LFEvMpfAnmW81bYnEeYIrux8pceHeO+n/RiU4uqB0tzSEF07qA02h+bDNZ63pV/d9WosRgtvt0+D7fNh/2ofVCi8qc5A11rbgduARcAO4GOt9Tal1KNKqfEASql+Sqkc4ArgNaXUNl8WLbzAYYNFD5AV25YlZXuZ1GUSEUERHh2iwmpn7s/ZjO2eQMuoEB8VKs5Eu7hwhnWK4/01+7HaPeuzEBMSw2WdLuPLin0cikqERfeDU/o9NGRutaFrrRdqrTtprdtrrR8//tzDWusFx79fq7VO0lqHaa1jtNbdfFm08IJ1c6Awg7fb9sZitHBN6jUeH+LzDQcpqbIzbXCK18sT3jNtcAr5pdV8s82zaXUBpqZOBRTvdBoEBzfAlnneL1B4jYwAaYoqj8F3T3IwZRBfHdnE5Z0uJzo42qNDaK15Z/U+UltGktamuY8KFd4wrFMcKTGhzPnB8xubLcNbcmH7C/n06BaOJPaCpY+AtcIHVQpvkEBvilY9CxVFvJ3cBRSnNZBoxa4CMvJKmX52iqxK1MAZDIqpg1NYn32MtfuKPN7/uu7XUe2o5p0O/aAkF358yQdVCm+QQG9qivbCmlc52PMyPj24gks7XEpCWIJHh9BaM2tZJolRwUyQec8bhSv7JRMdFsSsZZke79s2qi3ntz2fuYdWcaTL+fD9TCg55IMqxZmSQG9qlswAg5k3YmJRKG7s6dkUuQBr9haRvv8otwxvT5BJfoQag9AgE9cPacuKXQVsySn2eP9bet1CtaOaOUmdwGmDZY/5oEpxpuS3sSnZvxq2zyd3wA18vv9bLu3o+dU5wEvLM4kNt8g0uY3MlEFtiAg2MWv5bo/3bRvVlnFtxzF339cUpk2Fje/DwY0+qFKcCQn0psLpdHU7i2zFG8FODMrADT1u8PgwGw8cY9XuQm4a2lYGEjUykcFmpg9OYdG2PHbllXq8/829bsbmtDG7WTMIjYZvH3QNThMNhgR6U7FlHhzcQM6QvzI/6ysu73T5aV2dz1qWSbNQM1cPaOODIoWvTT+7LaFBRl5e7nlbepvINoxrN46P98yncMhfYd8qyFjogyrF6ZJAbwqsFa7uZol9eKX6AAZl4PrTWMBiS04xS3bkMX1wW8IsJh8UKnyteVgQVw9ozYJNB8nML/N4/1t63oLNaeMNQznEdnZdpdutPqhUnA4J9Kbgx1lQksuus//MF1lfcHXXq2kR2sKjQ2iteWLhDqLDgpg+JMU3dYp6ccuw9oQGmXjq650e75scmcylHS/l492fcGDYnVCUBWvf9EGV4nRIoAe6kkOubmapE/i/wysIDwo/ravz5Rn5/Jh1hDtGdSQy2OyDQkV9iQm3cOvw9izZkccaD9cdBVePF5MyMevoBmg/AlY8DRWe928X3ieBHuiWPQZOO+v7XMmKnBVc1/06oixRHh3C7nDyxMKdtIsNY3J/mYQrEFw/pC0to4J5YuEOnE7Pbmy2CG3BNanXsHDvQnYOvBGqS1yhLvxOAj2QHdwIG99H97+Zmbvn0iKkBVd3vdrjw3ycnkNmfhn3ju2CWdYLDQjBZiN3ndeZTTnFfLHZs2XqAKZ3n05kUCTP718AZ011NbsUet4dUniX/HYGKq1dN6xCo1neYSAbCzZyS+9bCDF5NitiWbWd5xbvol9Kc8Z0i/dRscIfLunTitSWkTyzKIMqm8OjfSODIrmxx438kPsDa7qdD6YQ+PYhH1Uq3CWBHqgyFsK+VdiG/52Zm18jJTKFSzpc4vFhnv02gyPl1TwwLlXmbAkwRoPiwXFdyTlaeVrdGCd3nUxCWALPbnsLx5A7YNfXkLXCB5UKd0mgByK71XV1HtuZuaHB7CvZxz397sFk8Kyr4Ybso8xZvY8pA9vQO7mZj4oV/jS4QyyX9mnFy9/tYefhEo/2tRgt/O2sv7GjaAcLWrSGqNaw6AFwena1L7xHAj0QrX0TirI4OuIfvLLlNQYnDuacVud4dAir3cl9n24hITKYe8Z09lGhoiF48MJUIkPM3PfpFhwe3iA9v+359IrrxQubXqZsxH2QtwU2fuCjSkVdJNADTUURrHgK2o/kpZJtVNgquLffvR43l7y+cg8ZeaX8a0J3IqSbYkCLDgvi4QtT2XjgGO/9uM+jfZVS3Nf/Po5UHeFNewEk9Ydl/4Jqz6cWEGdOAj3QrHgaqkvZPegm5u2ax8TOE2nfrL1Hh8jML+X/lmYyrmdLRqXKjdCmYELvRIZ1iuPfizI4UOTZAhbdY7szvv143t3+LjnD7oSyPPjhBR9VKk5FAj2QFO6GtW+iz5rK03s+Icwcxp96/cmjQ1RaHfz5/Q1EBJuYcVGqjwoVDY1Siscv6Y5RKW77YL3H64/+tc9fMRlMPHPgG+hxBax+EY55vjC1ODMS6IHk24fAFMLC9gNZc3gNt/e5nWbBnt3M/OeCbezKL2Xmlb1pERHso0JFQ5TUPJR/X96TTTnFHk8LEB8Wz809b2bZgWWs6DbW9eTSR3xQpTgVCfRAkfUd7PqakrNv45ktr9I9pjuXd7rco0P8b30OH6Uf4M/DOzC0U5xv6hQN2vk9WjJtcApv/7CXRR4uKn1t6rW0j2rPk9vepHLAza4ZPnPW+ahSURsJ9EDgdLi6izVrzYtmK0erj/LQoIcwGtyfr3xXXikPfLaV/m2juWNURx8WKxq6f1zQhR6torh73ib2FZa7vZ/ZaObBgQ+SW5bL61ERENbCNQe/zJlebyTQA8HG9yFvK9sG3cRHuz9hUudJpMa43/59uLiKaW//THiwiRcn98Ekw/ubNIvJyEtXnYXRoJg2+2eOlFW7vW9aQhrj249nTsYHZA2+FQ78BNs/92G1oib5zW3sqkth6b+wJffnkfyVxITEcFuf29zevbjSxrTZP1NSZWfO9H7ER0q7uYDWMaG8NTWNQ8VVXPdOOhVWu9v73pV2F6GmUB4p2YgzvhssfhhsVT6sVpwggd7Yff88lOczu9NgdhTt5IEBDxARFOHWrtV2Bze/l86egjJem9KXbomezcIoAlvfNtG8OLkPW3KOcdsHG7A73Ov5Eh0czT397mF9/gY+7H4eHMuGNa/4uFoBEuiN27ED8OMsMlPH8eq+LxiTMoZRbUa5tWul1cHN763jp6wi/nNFL87uEOvjYkVjdF63BP51cXeW7cznLx9uoNru3rD+Ce0nMKTVEF7IXsiBjiNg5bNQVuDjaoUEemO29BHswEOWasLN4dw/4H63diuutDHlrTWs3FXAU5f2YELvVr6tUzRqVw9ow0MXpvL11sPc8E465dV1N78opZgxaAYGZWBGZDBOeyV890Q9VNu0SaA3VjnpsGUe73QbydZju7h/4P1EB0fXuVtBaTWTX/+JTTnHmHXVWUySBSuEG64f0pZnLu/JD5mFXPPWGo5V1L2OaEJYAvek3cPaoq3M6zYK1s2B/B2+L7YJk0BvjLSGRfezo1kCs8p2MKr1KMa0GVPnbun7irjwxVXsLSznzan9uKBHy3ooVgSKK9KSefnqvmzLLeHCF79nc86xOve5tOOlDE4czLOVWewNjXJ1rxU+I4HeGG37jMqcn/l7QgLRlmgeHvTwKSff0lrzxsosrnz9J0LMRj69dTDDZOCQOA1juyfw8S2D0Bouf+VH/vvTfvQp+pkrpXh08KNYTMH8Pbkttj1LYfeSeqy4aZFAb2xsVbBkBs+0as++6qM8fs7jNA9uftLNs49UMG32Wh5fuIPRXeNZ8JchpCZG1mPBItD0Tm7Gl38ZwuAOMTz4+VZufm8dB49VnnT7+LB4Hhn8CDuqC3kxsQ18+wA43O8GKdwngd7YrHmFpdYC5pmtTOs+jYEtB9a6mdXu5KXlmYyeuYJ1+4/y6IRuvHLNWUTKVLjCC5qHBfH21H7cf0EXVu4uYPRzK3hzVdZJuzaOaD2CiZ0mMtui+bF0L6yfU78FNxHqVH8u+VJaWppOT0/3y3s3FtbsbI7Mnk3Jgi9wVlRgCA3BmFzM3cPDsHTown/P/y9m428D2uZw8tmGXGYtyyS7qILzuycw46JuJETJgCHhGweKKnh4/laWZxTQLi6M20d25MKeiRgNv20GrLRXMunLSRQf28eHWWWgJ1GycNHxn+1QIsdfRMz06QS1lhv1p6KUWqe1Tqv1NQn0hqls5Upybr8DbbOB/dc/T+0GjcOoiHrmMdqOveyX50urbHy+8SCvr9zDgaJKerSK4s7zOnFu5xb+KF80MVprvt2ex8zFu9h5uJR2cWHcPLQdF/ZMJMzy69KHmUczeer5y/nLp1UEOQ1Qc4UkkwllNpP0wvOEDx3qh7NoHAIi0P94tRq4n+jW7GyyJlyMrjx5u6QKCaHNZ5+x2RnOJ+ty+HLzISptDnolRXH7qI6c27mFLOos6p3Tqfl2+2FeWJrJjkMlhFtMXNQrkcv7JtEnuRn2nAPsvuhCDNW2kx5DhYTQbv7nAfl77Y0Ma/SBfrKr1UD9RD/0yCMcm/fJb8/1dxwGI0s7DGJm6sWEBhkZ3yuRSf1b0yspSoJc+J3WmnX7j/Lhzwf4astBqmxOYsMtPLD7CzqnL0U5TjHi1GSi2cSJtHz4ofor2Me8mWFnHOhKqbHAC4AReFNr/dTvXrcA7wJ9gSPAlVrrfac6pruB7u7VaoP/RHc6wV712y/bie+rwV7p+tdWScaUR3BW1j1wo9oSQs57XzK8cwvCa/xZK0RDUlJlY/nOfL7dnsd1T1xHqL3u2RsNIRY6z30STMHHvyxgDnH9awr59bHRAoaG3bfD2xl2qkCvMwWUUkbgJWA0kAOsVUot0Fpvr7HZ9cBRrXUHpdQk4Gngyjorc8OR2bNdn2qnoG02jsx5x71PdIe91hD99XHNoK0RuL/ZppbHJw3oKtf3jroD+gRnZUug7qtsi7WKC3smun1cIfwhMtjMhN6tmNC7FTv+5d7vgbOyCj6c5N4bGC2u0DcfD/5fPgROPHeKxyf7oPjl8Sn2cXO9Aa9n2Cm4c1nXH8jUWmcBKKXmAhOAmoE+Afjn8e8/AWYppZT2QntOyYIvTtn0AIDdzrFP5lIV+QMmZzUmpxWjswqT0/rLY5OzGqOuxqjdm1zoZGwqCJvBgk1ZsKsgrIbj/2LGqixUqzCqaUaVDqJKm6nCRIXBTBlmyh1GSuxGKnUQ1dpMFUFUY6Ya17ZOo4WI8HDuDHqeIKsbVzFhYWd0LkLUN0NoKM7yuhfNqDJbuNT+OAZHNcHKhgUrwdgIxopF2QhRNiJNdsKMdsIMdsKwEWy3E+ywElxtxYKNIGwEUUSQtmHW1Zid1Zi1FZO2/vL9mXAoI3ZDMA4VhN1gwW5w/eswWH557DAEYftkD9S1RqvdTsmCBfUS6K2Amqu95gADTraN1tqulCoGYoDCmhsppW4CbgJo7WbziLPCzRXIrU5KivJ+DUkdRdXxsKzGTJX+NTyrCaKKIGyYjm9/4jUL1ZixEkS1cj1vU0FYlQUrZuzKhEEZMKIwKIVBgUkZMBoUJoPCZFSYDAbMJgOWGl8hQSZCzUZCgoyEW0zEB5sIDzbRPDSI5qFBRIcFER9pISrEjFKKQ9X762xDx2Qicvx49/7bCNFARI6/yK2f7YQrLufTh/7M0Qob+aVVFJVbOVpu42iFlbJqO2VVdoqqbOTaHFTanFRa7VTbna4vmwObQ2N3OrE5NA6n68uunDg0aDQOpXFqJyZtx6yrXR8A2koQVizH/w3GigUrQZz4QDnxvevxLx8yquY2NizKioVygjmKBRvYHLjzF7c7H3R1qdeGV63168Dr4GpDd2cfdz/RDeHh9Hp07ZkV2EBEXHs1BZ9+zKmGACmzmZhpU+utJiG8IWb6dIo/n48+RaArk4mYaVNRShEd5rrgacwyvkpzL8O88Be3O3cTcoHkGo+Tjj9X6zZKKRMQhevm6BmLHH8RmOr43Amgq1Wbw8bfs2byzMXgtJj/eO4mEyokhKQXnm/YN4GFqEVQ69YkvfA8KiTkDz/bTqOBKjMsurEnxqTAmdI5cvxFOOta1tFLGeZOoK8FOiql2iqlgoBJwILfbbMAOHG5eDmwzBvt5+D6RFfmUw9XdxhVQFyt2pw2/r7q76zMWclFV8+g4xdf0mziRAzh4aAUhvBwmk2cSLv5nwdUN03RtIQPHUq7+Z//4Wc7+spJ7Jx5E2+Gr+fh1Q/jcJ7Z/a6GYungMKyGU7ehe+svbne7LV4API+r2+LbWuvHlVKPAula6wVKqWDgPaAPUARMOnET9WS81Q/dZnDyzMUw7Iq/cUOPG9w6XkNUZa/i7hV3syJnBfek3cO13a71d0lC+MUrm17h5Y0vMzZlLE8MeeIP01s0FlprXtzwIm9seYMby/py3hubXE1NPuyH7lYbutZ6IbDwd889XOP7KuAKt6o5DSc+0Y/MeYeSBQtwlpdjCAsjcvx4oq69msTs13hh/QtU2Cr4S5+/NLqBNaXWUm5behsb8jfw0MCHmNh5or9LEsJvbu11KyHGEJ5d9yyl1lKeG/4coeZQf5flEad28u+1/+b9He9zWcfL+PPAh3BclFtrhsVMm+q15tNGMVK0Lg6ng3/99C8+3f0p49uPZ8agGQQZG8eNlLzyPG5bdhuZxzJ5csiTjG071t8lCdEgfLb7M/754z/pEduDF859gZiQGH+X5JYKWwUP/vAgi/cv5pqu13Bvv3u9epF5qiv0hj3Eyk1Gg5EZg2bwp15/YsGeBdzw7Q0cqfTKPVmf2pi/kUlfTSK7JJsXR7woYS5EDZd0vIRnhz3LzqKdTP5qMtuPbK97Jz/LK89j2jfTWLJ/CXen3e31MK9LQAQ6uFZGubX3rTwz7Bm2H9nOVV9dxdbCrf4u66T+t/t/TF80nRBTCO9f8D5DWg3xd0lCNDij2ozinfPfQaOZ+vVUvsr6yt8lndT6vPVM/moy+0v28+KIF5nabWq9N/8GTKCfMDZlLHPGzsGJkykLp/DWlrdw6jpGadWj4upi7lt1HzNWz6B/Qn8+HPchHZp38HdZQjRY3WK6MXfcXFJjUrlv1X38c/U/Kbed+SAcb7E77by08SWmL5qOxWjhvQveY1jyML/UEhBt6LUpri7mkR8fYfH+xQxIGMAjZz9Cq3D/9m1dfXA1D/3wEEWVRdzU6yZu6nETRjfngxCiqbM5bLy08SVmb5tNy7CWPD7kcfrG9/VrTftL9vPA9w+wqWAT49uP5/4B9xNm9u2UHI1++tzTpbXm88zPefLnJ3FqJzf2uJFp3adhMVp8+r6/l1eexwvrX+CLrC9oG9WWJ4c8SbfYbvVagxCBYkP+Bu5fdT+5Zblc2vFSbutzG7EhsfVaQ6W9kjc2v8GcbXOwGC08OPBBxrUbVy/v3WQD/YTD5Yd5Zu0zfLv/W1pHtOYvff7C6DajfX51XGGr4N3t7/L21rexO+1cm3ott/S6hWCTLAcnxJkot5Uza8Ms5u6ci8Vk4eaeNzO5y2Sf/27ZnXa+3vs1szbM4mD5QS5qdxF3pt1Zrx8oTT7QT1h9cDVP//w0WcVZpESmcEOPG7ig7QVeH7hQWFnIBzs+4KOMjyixljC6zWj+1vdvJEck172zEMJtWcVZPJv+LCtzVhIdHM2kLpOY1HkSzYObe/V9qh3VLNizgLe3vE1OWQ6dm3fmvv73kZZQa676lAR6DQ6ng6XZS3l98+tkHM2gmaUZY1PGMr79eLrHdj/tu9JV9iq+z/2eb/Z9w/Ls5dicNka0HsG0btPo3aK3l89CCFHT2sNrmb11NqtyVxFsDGZkm5GMaTOGs1udfdpjUrTWbCrYxII9C/hm3zeUWkvpEduDG3vcyLDkYRiUf/qUSKDXQmvNDwd/YH7mfJZlL8PqtBIXEkdafBppCWl0je5KUkQSzSzN/hDyDqeDoqoisoqz2FSwic0Fm1l7eC0V9gqig6M5r815XN31alKiUvxzckI0UZlHM3l/5/ss3r+Y4upiIswR9E3oS6+4XvSK60XbqLZEB0f/IYy11hRVFZFblsv2I9tZe3gt6XnpFFUVEWIKYWTrkVzc4WL6J/T3+0h0CfQ6lFhLWLp/KT8e+pH0w+kUVBb88lqYOYxmlmauedCVkQp7BUcqj+CosVBG26i29I3vy3ltzqNfQj9MBlkOTgh/sjltrDm0hsX7F7Mubx37S/b/8ppJmYgNjSXYGIxDO3BqJ0VVRVTaf10iLj40nn4J/RiUOIiRrUf6vOeKJyTQPaC1Jrs0mz3H9pBblktOaQ6l1lIc2oFDO7AYLcSHxtMitAXJEcl0j+1OlCXK32ULIU7haNVRthRuIac0h4LKAvIr8ql2VGNURozKSJQliqSIJFqFt6J9s/YkhSf5/Ur8ZM54cq6mRClFm8g2tIls4+9ShBBe0jy4OUOTAn/K6YAbKSqEEE2VBLoQQgQICXQhhAgQEuhCCBEgJNCFECJASKALIUSAkEAXQogAIYEuhBABwm8jRZVSBcD+OjdseGKBQn8XUc+a2jk3tfMFOefGpI3WOq62F/wW6I2VUir9ZMNuA1VTO+emdr4g5xwopMlFCCEChAS6EEIECAl0z73u7wL8oKmdc1M7X5BzDgjShi6EEAFCrtCFECJASKALIUSAkEA/A0qpu5RSWikV6+9afEkp9YxSaqdSarNS6jOlVDN/1+QrSqmxSqkMpVSmUuo+f9fja0qpZKXUcqXUdqXUNqXU7f6uqb4opYxKqQ1KqS/9XYu3SKCfJqVUMnAekO3vWurBYqC71ronsAv4h5/r8QmllBF4CTgfSAUmK6VS/VuVz9mBu7TWqcBA4M9N4JxPuB3Y4e8ivEkC/fTNBO4FAv6ustb6W621/fjDn4Akf9bjQ/2BTK11ltbaCswFJvi5Jp/SWh/SWq8//n0proBr5d+qfE8plQSMA970dy3eJIF+GpRSE4BcrfUmf9fiB9cBX/u7CB9pBRyo8TiHJhBuJyilUoA+wBr/VlIvnsd1Qeb0dyHeJItEn4RSagmQUMtLDwD342puCRinOl+t9fzj2zyA60/09+uzNuF7Sqlw4FPgDq11ib/r8SWl1IVAvtZ6nVJquL/r8SYJ9JPQWo+q7XmlVA+gLbBJKQWu5of1Sqn+WuvD9ViiV53sfE9QSk0DLgRG6sAdvJALJNd4nHT8uYCmlDLjCvP3tdb/83c99eBsYLxS6gIgGIhUSv1Xa32Nn+s6YzKw6AwppfYBaVrrxjhrm1uUUmOB54BhWusCf9fjK0opE66bviNxBfla4Cqt9Ta/FuZDynVV8g5QpLW+w9/11LfjV+h3a60v9Hct3iBt6MIds4AIYLFSaqNS6lV/F+QLx2/83gYswnVz8ONADvPjzgamACOO/7/dePzKVTRCcoUuhBABQq7QhRAiQEigCyFEgJBAF0KIACGBLoQQAUICXQghAoQEuhBCBAgJdCGECBD/D08Ng1XbsGRnAAAAAElFTkSuQmCC\n",
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ]
          },
          "metadata": {
            "tags": [],
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "xaUHnQ2_MLtz",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 265
        },
        "outputId": "311fff5a-277b-4225-9488-b27b3ef0291c"
      },
      "source": [
        "pGrid = scinterp.barycentric_interpolate(nodes,gaussian(nodes),xGrid)\n",
        "plt.plot(xGrid,gaussian(xGrid),'-',label='Original function');\n",
        "plt.plot(xGrid,pGrid,'-',label='Barycentric interpolation')\n",
        "plt.plot(nodes,gaussian(nodes),'o',label='Data',markersize = 10);"
      ],
      "execution_count": 21,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3dd3xUVfrH8c+ZmVRCaAkJpNAjhCYQmhgB6Uix4NrrWn7uWlF37W11rcuu6+q66Lq2dVEUFRGlCYIKQigSElqoAUIINUAIycyc3x83kQApk2l3ZvK8Xy8lmblz5xmYfHPm3FOU1hohhBDBz2J2AUIIIbxDAl0IIUKEBLoQQoQICXQhhAgREuhCCBEibGY9cVxcnG7btq1ZTy+EEEFp5cqV+7XW8dXdZ1qgt23blqysLLOeXgghgpJSakdN90mXixBChAgJdCGECBES6EIIESIk0IUQIkS4FOhKqdFKqY1KqTyl1EO1HHeZUkorpTK8V6IQQghX1BnoSikr8DowBkgHrlJKpVdzXGPgHuBnbxcphBCibq600PsBeVrrrVrrMmAaMLGa4/4EvAiUerE+IYQQLnJlHHoSkF/l+11A/6oHKKV6Ayla66+VUg96sT4h/MLh1Mxau4ct+479eluXVrGM6pqIxaJMrEwI13k8sUgpZQGmADe6cOxtwG0Aqampnj61EB7TWvP9piJe+GYDG/YeBUApqNwmoEdyEx4e04WBHVqYWKUQrnEl0HcDKVW+T664rVJjoBuwSCkFkAjMVEpN0FqfNhVUaz0VmAqQkZEhO2sIU5WWO7jrf6uZl1tIavNo/nF1Ly7q3gqlFE6n5os1u3llzkauemsZE3q25i+/6UmYVQaGicDlSqCvADoppdphBPmVwNWVd2qtjwBxld8rpRYBD5wZ5kIEErvDyT3TjDD/4+jO/Pb8doTbToW1xaK4tHcyY7u34p+LtvDqgs3YLIpXLu8pXTAiYNUZ6Fpru1LqTmAOYAXe0VrnKKWeAbK01jN9XaQQ3qS15uEZ2czJKeTJ8encNKhdjcdGhlm5b0QaNoviL/M2ERsVxpPj06n4NCpEQHGpD11rPRuYfcZtT9Rw7BDPyxLCd174ZgPTV+7i7mGdag3zqu68sCOHSsp558dttGgUzl3DOvm4SiHqz7TVFoUww495+/nX4q1c0z+V+4a7HspKKR67qAsHj59kyvxNnN8pjl6pzXxYqRD1J1d4RINRWu7gkc+zaRfXiMfH1b/bxGJRPHtJdxJjI3l4RjblDqePKhXCPRLoosF4dcFmdhwo4blLuhEZZnXrHDERNp6Z2I0Ne48ydfFWL1cohGck0EWDsL6gmKmLt3J5n2TO6xBX9wNqMSI9gbHdE3l1wWa27T/upQqF8JwEugh5TqfmoRnZNI0K45GxXbxyzqfGdyXCZuHRz7PRWqZUiMAggS5C3tzcvfySf5iHx3ahWaNwr5yzZWwkD4w8h5+2HODHvANeOacQnpJAFyFNa82rC/JoH9eIS3olefXcV/ZLITE2klcXbJJWuggIEugipM1fv4/1BcX8fmhHrF6e4Rlhs3LHkA6s2H6IZVsPevXcQrhDAl2ELK01f1+wmdTm0Uw8t7VPnuOKvim0bBzB3xds9sn5hagPCXQRshZtLCJ79xHuHNoRm48W1YoMs3L74A4s3XqA5duklS7MJYEuQpLRd76ZpKZRXNLbu33nZ7q6XypxMeG89p200oW5JNBFSFqx/RBr8g9zx5AOPl/yNircyq2Z7VmyeT+5e4p9+lxC1EYCXYSkD5ftoHGkjct6J/vl+a7sm0pkmIUPf97hl+cTojoS6CLkFB09yTfrCpjUJ5mocPem+NdXk+gwxvdozRerd1NcWu6X5xTiTBLoIuR8kpVPuUNz7YA2fn3e6wa2oaTMweerdtd9sBA+IIEuQorDqfno552c16EFHeJj/PrcPZKb0iO5CR8u2yETjYQpJNBFSFm4YR+7D5/gOj+3zitdO6ANm/cd42cZwihMIIEuQsqHP+8gITaC4ekJpjz/+B6taRIVxofL5OKo8D8JdBEy8g+W8P2mIq7sm+rzoYo1iQq3MqlPMt+u20vR0ZOm1CAaLgl0ETJmrNqN1vCbvimm1nFl3xTsTs3MX/aYWodoeCTQRUjQWvP56l0MaN+cpKZRptbSKaEx3ZJi+Xz1LlPrEA2PBLoICavzD7P9QAmX9vLPRKK6XNIrmXW7i9lceNTsUkQDIoEuQsLnq3YTYbMwpnui2aUAMKFna6wWxYzVMiZd+I8Eugh6ZXYnX63dw4j0BBpHhpldDgDxjSPI7BTHl6t343TKmHThHxLoIugt2riPwyXlXOrjVRXr65JeSew5UsqybbJFnfAPCXQR9D5fvZu4mHAyO8WbXcppRqYnEhNhk6UAhN9IoIugdqSknAXr9zG+Z2vTxp7XJCrcyphuiczOLuBEmcPsckQDEFg/AULU07c5BZQ5nF7fANpbLumVxPEyB99t2Gd2KaIBkEAXQe3r7L2kNI+ie1ITs0upVv/2LYiLCWd2doHZpYgGQAJdBK3DJWX8lLefsd1boZQyu5xqWS2KUV0T+W7DPul2ET4ngS6C1tycQuxOzbjurc0upVYX9WjFiXIHCzdKt4vwLQl0EbS+zi4gpXkU3ZJizS6lVv3bGd0uX0u3i/AxCXQRlA6XlPFjgHe3VPq122W9dLsI33Ip0JVSo5VSG5VSeUqph6q5//+UUtlKqTVKqR+UUuneL1WIUyq7Wy7q3srsUlxyUXfpdhG+V2egK6WswOvAGCAduKqawP5Ia91da30u8BIwxeuVClFFZXdLoI5uOVO/ds2l20X4nCst9H5AntZ6q9a6DJgGTKx6gNa6uMq3jQBZvEL4TDB1t1SyWS3S7SJ8zpVATwLyq3y/q+K20yilfq+U2oLRQr+7uhMppW5TSmUppbKKiorcqVcI5uUGV3dLpcpul+83SbeL8A2vXRTVWr+ute4A/BF4rIZjpmqtM7TWGfHxgbXuhggec3MLad0kMmi6Wyr1a9ecptFhzM0pNLsUEaJcCfTdQNU9vZIrbqvJNOBiT4oSoiYlZXYWbypiZNfEoOluqWSzWhjWOYH56wspdzjNLkeEIFcCfQXQSSnVTikVDlwJzKx6gFKqU5VvLwI2e69EIU5ZvKmIk3YnI7smmF2KW0Z1TaC41M7PWw+aXYoIQXUGutbaDtwJzAHWA59orXOUUs8opSZUHHanUipHKbUGmAzc4LOKRYM2N6eQptFh9Gvb3OxS3JLZKZ7IMAtzc/eaXYoIQTZXDtJazwZmn3HbE1W+vsfLdQlxlnKHk/nrCxmRnogtwJbKdVVUuJXBafHMzSnkqfFdsViCq9tIBDaXAl2IQLB820GKS+3e727RGo7vhyP5UF5y6vbwGGiSDNEtwIv99SPTE5mTU0j27iP0TGnqtfMKIYEugsacnL1Ehlm4wNOdiQ7nw/YfIH8Z5C+HA1vAcbLm422R0KITpPSD1AHQ9nyIdX9BsGFdWmK1KObk7JVAF14lgS6CgtaauTmFXNApnqhwa/1PcHQv5HwO2Z/C7izjtogmkNIXOg43WuKxSRBZZaGv0mI4sstouRfmwNpPIOvfxn2pA6HbZZB+McTU7xdM0+hw+rdrztzcQv4wunP9X4sQNZBAF0Fh7a4j7C0u5cGu59TvgXvXwU9/h3WfgdMOid1h+FPQaSTEdwFLPfrinQ4j2DfPgezPYPYDMOcR6HEFnHc3xKe5fKpRXRN5cmYOW4qO0SE+pn6vSYgaSKCLoDAvtxCrRXFh55auPWDvOpj/FOTNg7BG0O826HMjxNfzF0JVFiu06mH8d8GDRrhnvQOrP4TVH0DncTDsSZeCfUR6Ak/OzGFuTiF3DJFAF94RnEMFRIMzf30hGW2a0axReO0HHiuCr+6Ff2UaXSsXPg6Tc2D0856FeXUSusJFf4H7cmDwH2HbYvjnQPjmISipfZx566ZRdG0dy4L1MmtUeI8Eugh4+QdL2LD3KCPSaxndorXRUn6tt9Fa7nc73L0aLngAopr5tsBGcTD0EbhrFfS+Hpb/y6gj+9NaHzYiPYGVOw9x4FgtF2SFqAcJdBHw5le0YmsM9KN74X9Xwpe/h8QecMdSGPOC74P8TDHxMO6vcPsSaN4BPvstTL8Rjh+o9vDhXRLQGr7bIIt1Ce+QQBcBb/76Qjq1jKFNi0Zn37lpLrwxALYugtEvwA1f1evipE8kdoOb58CwJ2D9LKO+bYvPOqxr61haNYn89ReWEJ6SQBcB7ciJcn7eepDhZ7bOnU5Y9CJ89BtjyOHtS2DAHfUbteJLVhtk3g+3LYKopvD+xfDTa0bXUAWlFMO7JLB4035Ky2WNdOG5AHn3C1G97zcVYXdqhnepEuilR2DaVbDoz8aQwZvnmt8qr0liN7j1O+h8Ecx9DD69CcqO/3r38PQETpQ7WLql+m4ZIepDAl0EtPm5hcTFhHNu5YzK4gL4z1jImw9jX4FL3oTwaHOLrEtEY/jN+zD8acj9Et6bYCw1AAxo35xG4VbmSbeL8AIJdBGwyh1OFm7cx4WdjanyFG2Ef4+AQ9vhmunQ71avrrHiU0rB+ffCFR9C4Tr490g4uI0Im5XB58SzYH0hTqfs3Cg8IxOLREDJL87nvdz3mLV1FsfLS9Cp4RyJHkn+pt2kzLgDrOFw02xo1dPsUt3T+SK4fib87woj1K+bQe/2Dhbs+5gBHz1KqeME0WHRjGs/jhvSbyAlNqXucwpRQWltTqsgIyNDZ2VlmfLcIjAt2bWEyd9Pxu6wY9f2X2+3YcXmtDOlxErm1V9Cs7bmFektRZvgg4tZosqYHNeEEw47Sp3axcimbNisNqYMnkJmcqaJhYpAo5RaqbXOqO4+6XIRASG/OJ/J30+m1F56WpgD2HFQalFMbhJGvtWNhbkCUXwa+Ze/zeRmkZQ6y04LcwC7tlNqL2Xy95PJL86v4SRCnE4CXQSE93Lfw+6w13qM3eng/dz3/VSR7723+zvsltp7Pe0Oe0i9ZuFbEugiIMzaOuuslvmZ7NrOrK2z/FSR7xmvufbx56H2moVvSaCLgFBSdaegWhwvP173QUGiIb5m4VsS6CIgRIe5Npa8UVg10/+DVEN8zcK3JNBFQBjXfhw2VXt/sk3ZGNd+nJ8q8r2G+JqFb0mgi4BwQ/r12LSz1mNsVhvXp1/vp4p874b0G7BZ6wh0S2i9ZuFbEugiIKSs/JApe/cSqWxoffrb0qZsRNoimTJ4SkhNtEmJTWHK4ClE2iLPaqnbUEQ6nUxJGBJSr1n4lgS6MF/2p/D9i2R2nsQferxD+aF+RFqjUShiwmKYlDaJGeNnhOQEm8zkTGaMn8GktEnEhMWgUGhHBINbT2BGVDcyf3wLNsw2u0wRJGSmqDDX3mx4ewQk9YbrvuCxWRv5bOVuVj8xgsiwEJlEVA/FpeX0fmYet2S256FhqfDuWDiwBW5dCHEdzS5PBACZKSoC04lD8PG1xs5Cl7+LtoYxP3cfF6TFNcgwB4iNDKN/++bGphfh0fCbD8AaBh9fAyePmV2eCHAS6MIcTifMuA2O7DaWlo1pSc6eYvYWlzIiPdHs6kw1oksCefuOsW3/cWiaApPegf2bYOadp22QIcSZJNCFORa/BJvnwpgXIaUvAHNzC7EoGHpOvMnFmWtYxWYe83Mr1khvPwSGPQk5n8OyN0yrSwQ+CXThf9sWw6IXoOdVkHHzrzfPzy2kT5tmtIiJMLE486U0j6ZzYuPTN70YdA90HgfznoBdK80rTgQ0CXThX8eK4LNboUVHY8ehig0qdh8+QW5B8elbzTVgI9ITyNp+kEPHy4wblIKJ/4DGrY1t7EqPmFugCEgS6MJ/nE744g7jYujl/4GImF/vWlDRGh1x5mbQDdSI9AScGhZu3HfqxqhmMOnfcGQXzLxb+tPFWSTQhf8sex3y5sGo5yCx+2l3zcstpH18I9rHx9Tw4IalW+smJMRGMC/3jL1GU/rBsMch9wtY+R9zihMBSwJd+EfBWpj/tNEP3PeW0+4qLi1n2dYDjJDull9ZLIphXRL4flMRpeVnLLF73j3Q4UL49hHYv9mcAkVAcinQlVKjlVIblVJ5SqmHqrl/slIqVym1Vim1QCnVxvuliqBVXmoMUYxuARNeO2tj58Wbiih3aIZLd8tpRnRJoKTMwbKtB06/w2KBiW9AWKTx9+ooN6dAEXDqDHSllBV4HRgDpANXKaXSzzhsNZChte4BfAq85O1CRRBb8AwUrYeLX4fo5mfdPS+3kOaNwumd2syE4gLXwA4tiA63MvfMbheA2FYw7m+wZxUsftn/xYmA5EoLvR+Qp7XeqrUuA6YBE6seoLVeqLWuXK1/GZDs3TJF0Nq6yOg773srdBx+1t1ldiffbdjHsM4tsVrU2Y9vwCLDrAxOi2d+biFOZzUXQLtebAz9XPwK5K/wf4Ei4LgS6ElA1V1qd1XcVpPfAt9Ud4dS6jalVJZSKquoqMj1KkVwKj0CX/wO4tJgxDPVHrJs6wGOltoZ1bVhzw6tyaiuiew7epI1uw5Xf8CYFyE2CT6/Hcpc2wFJhC6vXhRVSl0LZADVfgbUWk/VWmdorTPi4xv2bMAGYe5jcLQALn7TWJekGnNy9hIdbuX8TnF+Li44DO3cEptFMSdnb/UHRDYxurIOboGFz/m3OBFwXAn03UDVBZmTK247jVJqOPAoMEFrfdI75YmglbcAVr0P590NyX2qPcTp1MzLLWRwWnyDXYyrLk2iwhjYoQVzcwqpcWXUdhdAxm9h6euw82f/FigCiiuBvgLopJRqp5QKB64EZlY9QCnVC/gXRpjvq+YcoiEpLTYmvsSlwZCHazxsza7D7Dt6Urpb6jCyayLb9h8nb18tqy2OeBqapMCXv4PyE/4rTgSUOgNda20H7gTmAOuBT7TWOUqpZ5RSEyoOexmIAaYrpdYopWbWcDrREMx7HI7uOTW0rgZzcvZisyiGdm7px+KCz8iK4Zw1drsARDSGia/BgTz47lk/VSYCTe0bGlbQWs8GZp9x2xNVvj57+IJomLYtgZXvwnl3/bqKYnW01szNKWRghxY0iQrzX31BKCE2kl6pTZmTU8idF3aq+cD2Q6DPTcaKjN0uhaTqu7pE6JKZosJ7ykth1r3QrC0MeaTWQyvX+x4p3S0uGZmeSPbuI+w5XEd3yoinoVFLmHmPTDhqgCTQhfcsecX4yD/urzWOaqlU2X0wUmaHumRUV+PvaW5t3S5gjHoZ+zIUZhsXSUWDIoEuvKMwF374qzHRpcOFdR7+bc5eeqU2JSG25j52cUr7+Bg6tYzhm3V1BDpA+gRjzZxFz8PBrb4vTgQMCXThOacTvrrHaB2OrHss9M4DJazbXczYbq38UFzoGNO9FSu2H6ToqAujgse+DJYwmHWfLLPbgEigC8+tehd2LYdRf4ZGLeo8fPa6AgDGdJf+8/q4qHsrnLqO0S6VYlvD8CeNpRfWfebz2kRgkEAXnjlWBPOfgraZ0OMKlx4yO7uAnslNSG5Wez+7OF1aQgzt4xsxO7vAtQdk3Ayte8GcR2SHowZCAl14Zt7jxhoiF005a1nc6uQfLGHtriOM7S7dLfWllOKi7q1YtvUA+4+50O1isRoXqI8Xydj0BkICXbhv+w/wy/9g0N0Qn+bSQ76p6G6RQHfP2Ipul7k51SypW53WvYwNRVa8DXtW+7Y4YToJdOEeexl8fT80TYXMB1x+2NfZe+me1ISU5tLd4o7OiY1pF1ePbheAoY9CdJxxgdTpqPt4EbQk0IV7lv8LijbAmJfqHHNeadehEn7JPyytcw8opRjbPZGlWw9w8HiZaw+Kamrs47pntbFgmghZEuii/o7uhUUvQKdRcM4Ylx/2bcUY6rEyusUjY7u3wuHUdU8yqqr75ZB6Hix4GkoO+q44YSoJdFF/854ARxmMfr5eD5u1toCurWNp06KRjwprGNJbxdKmRTRf16fbRSljbHrpEblAGsIk0EX97PgJ1n5srHPeooPLD9t5oIQ1+YcZ37O1D4trGJRSjO/Rmh/z9rs2yahSYjdjK8Csd2DPGt8VKEwjgS5c57DD7D8Y625n3l+vh878xdgTRQLdOyae2xqnhq/X7qnfA4c+AtEtYPaDxgxfEVIk0IXrVr1nLPo08lmXL4SCsVTuF2v20K9dc5KaRvmwwIajU0Jj0lvF8sWaegZ6VFMY/pQxs3fdp74oTZhIAl245sQho++1bSakT6zXQ3MLisnbd4yJ50rr3JsmntuaNfmH2XHgeP0eeO410Opc41rIyVp2QRJBRwJduGbRi1B62LgQ6sKM0KpmrtmDzaJkMS4vq+y+mlnfVrrFYgw3PVoAP/7NB5UJs0igi7rt2wDLp0KfGyGxe70e6nRqZv6yh8Fp8TRrFO6b+hqo1k2j6NeuOV+s2V3zBtI1Se1vDGX88e9waIdvChR+J4Euaqc1zHkYImJg6GP1fviK7QcpOFLKBOlu8YmJ57ZmS9FxcguK6//g4U8b673Me9z7hQlTSKCL2m2aA1u+gyEPu7Q07pm+WLOHqDArI2RnIp8Y260VNoviy/p2uwA0SYLzJ0Pul8ZesCLoSaCLmjnKYe5j0KKTscBTPZWWO5idXcDIrglEh7u0H7mop2aNwhmcFs+Xa3bjcLqxkcV5d0JssrHErqzzEvQk0EXNVvwbDmw2hilaw+r98Hm5hRw5Uc6kPsk+KE5UmtQnmcLikyzeXFT/B4dFGRtL711rrJwpgpoEuqheyUFjT8r2QyBtlFun+CQrn6SmUZzXIc6rpYnTDeuSQPNG4UzPynfvBN0ug+S+sOBPMowxyEmgi+p9/xKcLDa2lavnMEWAPYdP8EPefi7rnYTVUv/HC9eF2yxMPLc183ILXV+BsSqlYNTzcGwv/Piq9wsUfiOBLs62Pw9WvAW9r4eErm6d4rOVu9AaJvVJ8XJxojqX90mh3KH5cs1u906Q0he6TYKf/g5Hdnm3OOE3EujibPOeAFuksTGCG5xOzfSVuxjYvgWpLWQjC39Ibx1L96QmTM/yIIyHP2UMU13wJ2+VJfxMAl2cbvsPsPFrOP8+iGnp1imWbz/IzoMlXJ4hF0P96fKMZHILilm3280NoZumwMDfwdppsl1dkJJAF6c4nTDnUWMY28Dfu32a6Vm7iImwMUam+vvVhJ6tCbdZ+HSlB6308ycb29XNfdxorYugIoEuTsmeDgVrYNgTxnA2Nxw5Uc7s7ALG92xFVLjVywWK2jSNDmdU10RmrNrFiTI3x5RHxsKQh2D7Etj4jXcLFD4ngS4M5SdgwTPGKnzdL3f7NJ+u3MWJcgfX9G/jxeKEq67pn0pxqZ2vfnFj5milPjdCXJqxJICj3Gu1Cd+TQBeGZW9A8S5jEpHFvbeF06n5cNkOeqc2pVtSEy8XKFzRv11z0hJieH/Z9vov2FXJGgYjnoEDebDyXa/WJ3xLAl3A8f2w5K9wzlhol+n2aX7I28+2/ce5fmBb79Um6kUpxXUD27JudzGr8w+7f6K00cba94tegFI3Fv4SpnAp0JVSo5VSG5VSeUqph6q5/wKl1CqllF0pNcn7ZQqf+v5FKC8xVt/zwPtLtxMXE86Y7oneqUu45dJeSTSOsPH+T9vdP4lSMPJPULJf1kwPInUGulLKCrwOjAHSgauUUulnHLYTuBH4yNsFCh/bn2dsGtznBohPc/s0+QdLWLBhH1f2TSXCJhdDzdQowsZlfZKZnb2X/cfqsYn0mVr3Mq6nLH0djrg5YUn4lSst9H5AntZ6q9a6DJgGnLYHmdZ6u9Z6LSC7zgabBU8Zk4iGPOzRaf77804UcHX/VK+UJTxz7YA2lDmcfLzCzfVdKl34OGgnLHzOO4UJn3Il0JOAqu+KXRW3iWC382dY/xUMusftSURgLJP7SVY+I9ITaC2bQAeEji1jGNSxBf9dtoNyhwftrGZtoP/tsOYj2JvtvQKFT/j1oqhS6jalVJZSKquoyI2lPoX3aG2sdR6T6NEkIoDpWfkcPF7GzYPaeak44Q03D2rHniOlzFrrwRBGgMz7IaopzHvSO4UJn3El0HcDVVdYSq64rd601lO11hla64z4+Hh3TiG8Zf1XsGs5DH0Ewhu5fRq7w8nUJVvpldqUfu2ae7FA4amh57QkLSGGNxdtdX8II0BUM7jgQdiywNi9SgQsVwJ9BdBJKdVOKRUOXAnM9G1Zwqcc5TD/KYjvDOde49Gpvs4uIP/gCe4Y3AHlxjK7wncsFsX/De7AxsKjLNy4z7OT9b0FmqYaC7c55VJZoKoz0LXWduBOYA6wHvhEa52jlHpGKTUBQCnVVym1C7gc+JdSKseXRQsPrXwXDm4xhila3d8aTmvNm99vpWPLGIZ3kT1DA9H4nq1JahrFm4u2enYiWwQMe9LoR8+e7p3ihNe51IeutZ6ttU7TWnfQWj9XcdsTWuuZFV+v0Fona60baa1baK3dW0Rb+F5psTFZpG2m2zsRVVq8eT/rC4q57YL2WGQTi4AUZrXw2/PbsXz7QVbuOOjZybpeaiwN8d2foLzUOwUKr5KZog3NT383JouMeNqtnYiq+ueiPBJjI7n4XBn0FMiu7JdC0+gw/ulpK91iMSYbHcmHn9/0TnHCqyTQG5LiAvjpH0ZLK6mPR6davu0gy7Ye5JbMdoTb5G0UyKLDbdx4Xlvmry90f630Su0ugE6jYMkUY99ZEVDkJ7EhWfRncNqN5XE9oLXm5TkbaNk4QlZVDBI3DWpHk6gw/jJ3o+cnG/4UlB2FJX/x/FzCqyTQG4p962H1h9DvVmju2XjxRRuLWLH9EHcP6yRrngeJJlFh3DGkAws3FrFiu4ct64R0Y3TU8qlwaLtX6hPeIYHeUMx/CsIbG+OJPeB0al6as5HU5tH8JkM2gA4mNwxsS8vGEbz07QbPxqWDMX9BWeG7Z71TnPAKCfSGYNsS2PQtZE6GaM8m/3ydXcD6gmImj0iTvvMgExVu5a5hnVix/RCLNnk4Uzu2tTHDOHu67D8aQOQnMtQ5ncbOM7HJxpocHih3OJkybxPnJDRmfM/WXipQ+NMVGSmkNo/m5W834nR62EofdI/sP/4j6CQAABaGSURBVBpgJNBDXc4MowV14WNu7xNa6b2ftrNt/3EeHHUOVhl3HpTCbRbuH5lGbkEx01d6uBJjZCwM/qOx/+jmud4pUHhEAj2U2U/CgqchsTv0uMKjUxUWl/K3+ZsZck48w7q4vzKjMN+Enq3p27YZL3yzgcMlZZ6dLOMmaN7BWBLAYfdOgcJtEuihbPlbcHgnjPiT2/uEVvrz7PWU2Z08Nb6rrNkS5JRSPDOxG8Wldl6e4+EwRmuYMYyxaAOs+a83yhMekEAPVScOweKXocMw6DDUo1Mt3XKAL9fs4f8Gt6dtnPsrM4rA0aVVLDcMbMtHy3eydpcHe48CdBkPKf1h4Z+h7Lh3ChRukUAPVYtfgdIjxu7tHih3OHniy3UkN4vid0M7eqk4EQjuHdGJuJgIHv9inWcXSJUyPgUe2ws/vea9AkW9SaCHooPbjEkf514Did08OtVr3+Wxed8xnhrflcgwmUQUSmIjw3jsoi78susI//5hm2cnS+0PXSbAj3+Ho3u9U6CoNwn0ULTgGbDY4MJHPTrNqp2HeH1hHpf2TmJ4uiyPG4om9GzNyPQEXp6zkQ17iz072fCnwFFmdL0IU0igh5r8FcZQxfPuMiZ/uOn4STuTP15DYmwkT02Q1ZBDlVKK5y/tTmxUGPdOW8NJu8P9k7XoYCwtsfoDKMz1XpHCZRLooaRyn9BGLeG8uz061XOz17PjYAl/+U1PYiPDvFSgCEQtYiJ4aVJ3Nuw9ypS5mzw72QUPQkRjYzKb8DsJ9FCy/ivIX2Z0tUTEuH2ab7IL+OjnndyW2Z4B7Vt4sUARqC7snMA1/VOZumQrizzZri66uRHqefMhb4H3ChQukUAPFfaTxuSO+M5w7rVunyZ3TzGTP/mFc1OaMnlkmhcLFIHu0Yu6cE5CY+7632q2Fh1z/0T9boOmbYwlAZwedOGIepNADxXL34JD22Dkc27vE3rg2ElufT+L2CgbU6/rQ4RNRrU0JNHhNt66PoMwq4Vb3s+iuLTcvRPZIowdsfblGP3pwm8k0EPB8QPw/UvQcTh0Gu7WKcrsTu747yqKjp1k6nUZtIyN9HKRIhikNI/mn9f0ZueBEu7532oc7o5PT78YUgcay+uWejh6RrhMAj0ULHoeyo4ZrXM3OJyayZ+sYfm2g7x0WQ96pjT1coEimPRv34KnJ3Zl4cYiHv08271JR0rBqOfgeBH88FfvFymqJYEe7PZtgKx3jEWSWnau98OdTs2Dn/7CrLUFPDymMxf3kg2fBVzTvw13X9iRaSvyeXJmjnsbYiT1gR5XwtLX4dAO7xcpziKBHuzmPgbhMTDk4Xo/1OnUPPJ5NjNW7eb+EWncPriDDwoUweq+EWncfkF7Pli2g2e/Xu9eqA97ApTF2DFL+JwEejDbNBfy5sEFD0CjuHo99KTdwQPTf2HainzuurAjdw3r5KMiRbBSSvHQmM7cNKgt//5hG49+sQ67w1m/kzRJMjbCyJkBO37yTaHiVxLowcpeBnMeMdai7v9/9XrooeNlXPf2cmasNlrmk0fI8ERRPaUUT4xL544hHfjo553c9O6K+o9+GXQPxCbBN3+UYYw+JoEerFa8BQc2w+jnwRbu8sPy9h3lkjd+ZM2uw/z9ql7cNayTrG8uaqWU4o+jO/PiZd1ZuuUAk/75EzsPlLh+gvBoY9XPvWtlzXQfk0APRseKYNGL0HEEpI1y6SFaaz5Yup1xr/1Acamd/93anwmyL6iohyv6pvLezf3Ye6SUMa8uZnpWvuv96t0ug5QBxsJxpUd8W2gDJoEejBY+C+XHYZRrq9rtKy7lpndX8PiXOfRr14Jv78mkT5vmPi5ShKJBHeOYfU8m3ZKa8OCna7njw1UcOHay7gcqBWNegOP7jTkTwick0IPN7lWw8j3odzvE1973XVru4PWFeQx9ZRFLtxzg6Qldee+mvjJpSHgkuVk0H906gIfGdGbBhkKGvLKIt5dspcxexwXT1r2g17Xw85tQ5OHWd6Jayq2hSF6QkZGhs7KyTHnuYJFfnM97ue8xa+ssSspLiA6LZlypgxsOHyHldysgskm1jyt3OPnqlz38Ze4mdh8+wYj0BB4Z24V2sn2c8LLNhUf509frWbypiHZxjZg8Io0x3RKxWWtoKx7fD6/1Jr9VN95LG8isbV+fem+3H8cN6TeQEpvi3xcRZJRSK7XWGdXeJ4EemJbsWsLk7ydjd9ix61O7qdu0xmYJZ8qFr5KZnHnaY46WljNteT7v/LiNgiOldGkVy+MXdeG8jvUb0ihEfS3cuI/nvl5P3r5jpDSP4reD2nF5RgqNIs5eV2jJgkeYvPNL7BYbdk616m3Khs1qY8rgKWe9t8UpIRHo1bZWQ/Q3en5xPpd+dSml9tIaj4m0RTJj/AwSGyXxw+b9zFi9m3m5eyktdzKgfXNuu6A9Q9JaYrHICBbhHw6nZl5uIVMXb2HVzsNEh1sZ3TWRS3oncV6HOKwWZby3Z15KqaPu93Yo/lx7I8OCPtBrbK2G6G/0Z5c9y2ebPjvttZ7JgpUENYTCbWM5cqKcptFhjOvRit9kpNAjWdZiEeZaueMQ07Py+Tq7gKOldpo3CiezUxxHG33MqkNzcNTy3rYpG5PSJvHoAM+2UAwk3swwjwNdKTUaeBWwAm9rrV844/4I4H2gD3AAuEJrvb22c7oa6PVprYbKb/QBHw3gePnxug90RjIq5m1GpCcw9JyWhNvkGrcILKXlDr7bsI95uYUs3lTEyaSHUda6R8XEhMWw9OqlfqjQ97ydYbUFep0LZyulrMDrwAhgF7BCKTVTa11108DfAoe01h2VUlcCLwJX1FmZC97LfQ+7o+bf5gB2h533c98Pmt/o5Q4nx0/aOVxSzuET5RwqKWP/0ZPsO3qSoqMnXQtzQFlO8srlPX1crRDuiwyzMrZ7K8Z2b4XTqen5wf0uPe5Y+XHeWJRHy8aRxMWE0yw6nKbRYTSJCiMmwlbzRdcA5M8Mc2UnhH5AntZ6K4BSahowEaga6BOBpyq+/hT4h1JKaS/058zaOqvWrgcAu7bz2aYvCTt82a+31fTEVSvSlUfpU8drrdEV3xt/VnyvjaOdWuNwGt87nBpHxZ92p8bh0NidTsocmjK7g3KH5qTdQWm5kxNlDk6UOzh20l7r8K7GETZUu0i0qvm3eaVGYTJqRQQPi0XRKCzapQaLdkbw0rc1D22MsFmIibARGWYlKtxKVJiVcJuFcKuFcJuFMKvCZrFgsyqslor/lPGnxaKwKLAohcKYCasUKCr/NIbNV86gVr/+zzimUk0TrM+8eca+mS5l2Kyts/wS6ElAfpXvdwH9azpGa21XSh0BWgD7qx6klLoNuA0gNTXVpQJLyl2bYlzmLOWtJVtPu02d9Vf76x1nfVn5D3rq64p/aE7946qKN4FFKawWKv489WaxVbyJwqyKcJuFyDALsZE2osKtxhsvzEpMhI1GETZiImw0jQ6raHWE07JxBHExEUSFW3l22Yo6+9Btysa49uNc+rsRIlCMaz/Opff2pC4Xc981o9h/tIyiY6XGp9mKT7THT9o5ftLO0ZN2SssdlJY7OFHmoMzhpMzupKTETnlF48ruMBpddoeuaIxpnNpomDkrGm/Oihac0YjTvzbm4FSDzvi6ihpajLqaOyI6nagx/Kty9ZN5bdzbq8xNWuupwFQw+tBdeUy0i7/RG4c3Yt1zYz0rMEDckH4DX275Eru9lje91cb16df7sSohPFef93Z0uI3UFjZSW0T7sULvG/DRMy5lmDc+cbvSEbUbqNpTn1xxW7XHKKVsQBOMi6MeG9d+HDZV+++dUGutpsSmMKXDVUQ6ndjO+JRhUzYibZFMGTwlZC4Ci4YjJTaFKYOnEGmLPOvn2qY1kZbwkHtv+zPDXAn0FUAnpVQ7pVQ4cCUw84xjZgI3VHw9CfjOG/3nYPxGt9Wx6XHItVZPHCbzh38xo7Qxk9ImERMWg0IRExbDpLRJzBg/I6SGaYqGJTM5kxnjZ5zx3m7EpDILMw6Xk9myt9klepU/M8zVYYtjgb9hDFt8R2v9nFLqGSBLaz1TKRUJfAD0Ag4CV1ZeRK2JV8aha7DZIpkyJLTGofPlncYyo7csgKTQenMLUaPtP8C7F8GA38No1xaeCxZLdi5i8nd3YleKqp1N3h6H7tLYH631bK11mta6g9b6uYrbntBaz6z4ulRrfbnWuqPWul9dYV5f1f5Gt0QwqfgoM3o/FFphvmUhrP4AzrtLwlw0LG3Ph4ybYdkbkL/C7Gq8KrOkhBm7CpjUsr9PP3EHxUzRapWfgL90hg5D4fJ3vVaXqU4eg38OBEsY3PEjhEWZXZEQ/lVaDG8MhIgYuO17CAuRlUE/nASF6+DebLCGeXQqj1voASksCs69BtZ/BcUFZlfjHXMfhcP5MPF1CXPRMEXGwvhXoWiDse5/KDi4FfLmQ+/rPQ7zugRvoAP0u9UYMLr0H2ZX4rlNc2Dlu0ZXS5uBZlcjhHk6DYc+N8FP/zD61YPdD38Da7jRneRjwR3ozdtB90mQ9Q4c98ooSXMcP2BcCG3ZFS58zOxqhDDfyGehWVv4/A6jGyZYHdkNaz4yNvZonOjzpwvuQAc4fzKUl8DP/zS7EvdoDbPuhROH4NJ/gS3C7IqEMF9EDFw6FYp3wbcPm12N+356DbQTBt3jl6cL/kBv2Rm6jIefpwbn5rOrP4D1M+HCRyGxu9nVCBE4UvrB+ffBmg9h3Qyzq6m/Y0VGN2rPK6FZG788ZfAHOkDm/XDyCKx42+xK6qcwF2b/AdoNhvPuNrsaIQLPkIchuS/MvNu4uBhMlr0O9lLjl5KfhEagt+4FHYfD0jegzPMFbvyi7DhMvxEiGsNlb4PFanZFQgQeaxhMesf4+Zh+I9jrXks9IJw4BMvfhq4XQ1wnvz1taAQ6wAV/gJL9Rp9VMPj6Adi/CS57C2Jaml2NEIGraSpc/AYU/AJzHze7GtcsfgXKjkHmA3592tAJ9NT+0PVSY4jQ4fy6jzfTynfhl49g8B+g/RCTixEiCHS+CAb8Dpb/C7I/Nbua2u3fDD+/aYw7T+zm16cOnUAHGPEMoGHeE2ZXUrMdS43WeYdhMPiPZlcjRPAY/jSkDjSG+O5ZY3Y1Nfv2YQiLhgv9/2kitAK9aQoMuhdyZsCOn8yu5myH8+Hja42PkJP+Lf3mQtSHLRx+8wFEt4BpV8OxfWZXdLZNcyFvntFYi4n3+9OHVqCDMd4zNgm++SM4HWZXc0pZifEmdJTBVdMgqpnZFQkRfGLi4aqPoOQgfHxdYF0ktZfBnIehRUfod5spJYReoIdHG10ve9caq7YFAke5cYV+b7YxoiU+zeyKhAherXrCxa9D/jL4/P/AWfMevX615BU4kAejnjc+TZgg9AIdoNtl0HkczH/auDJuJqcTZt4Fm+fARX+BtFHm1iNEKOh2mdFwy5kB3/7x9N3fzbBjKSx+GXpeBWkjTSsjNANdKZjwGjSKg89uMXds+vwn4Jf/wdBHoe9vzatDiFAz6B5jMbvlU40wNcuJwzDjNuPa2FgT6yBUAx0gujlc8qYxhGjOI/5/fq1h4fPGuPh+t8EFD/q/BiFC3fBnoOfVsPA5c+agaA1fT4bi3XDZv42JgiYK3UAHY4z3oLsrxn1P89/zag1zH4PvXzDWbB/9ovGpQQjhXRYLTPg7pF9s/MwtfN6/3S8r3oZ1n8HQhyG52j0n/Kr2nUtDwdDHYM9q+OJ3ENkUzhnt2+dzOmDWfbDqPeh3O4x+wXjTCSF8o3J5gJkxRiOq7Jix/K6vG1HZn8LsByFttLHqawAI/aSxhcOVH0GrHjD9Bt+OTy89AtOuMcI8834Y86KEuRD+YLEa18363W5sePPpzb69drZ5Pnx+uzHR6fJ3A2ZOScNIm4jGcM2n0CQFProC8pd7/zn2rYepQ41JBWNehmFPSDeLEP5ksRiNqOFPQe4X8PZwOLDF+8+z9Xv45Dpo2QWunhZQ20U2jEAHY8TL9V8YF0v/M9bY5cgbfW1OJ6z6AN4aZnzUu2EW9DdnUoEQDZ5SxnK1134GRwuMRtba6d75Wdfa2Bbvg0uMxuG1MyCyiefn9aKGE+gATZLh1oXQ7gKjn3vmnVB+wv3zFebCu2ON87TqaexSLvuBCmG+DhcaP49xHWHGLfD+RGPEm7tOHoPPfmts5H7OGLhlfkCukqq0SQPyMzIydFZWlinPjdMBi543xq42bg0X3A+9rnd9dldhjjH2ddUHxm/oEc8Yo1mkv1yIwOJ0GJ/GFzxjbDbR50ZjGLGra5SXn4Cs/8APU6DkgLHg1vn3mdqdqpRaqbWudkhNwwz0Stt/hO/+BDuXGh+hzr0G2mUaO6RU3dvTYYeiDbBrhXFle8cPYI2Ac682/oEbtTDvNQgh6na00Aj1tR+DsxzaD4Vul0JSBsSfc/pFzfITxnW27Utg9YdG1027wcbPekpf815DBQn02mgNWxYYC9LvXAZosEUaO3QrCyirMWmgvMQ4vmkq9L0Fel1n9McLIYLHsX3GKLSs/xg/1wDhMdC4FWiH0aI/WmAsoqcs0GaQsXJiu0xz665CAt1VJw4bwxq3/wDHi079AzdOhKQ+xn/N28voFSGCndMJB7fA7pWwK8vY7ayyAdc4AdpmGkMSI2PNrvQsEuhCCBEiagt0uYonhBAhQgJdCCFChAS6EEKECAl0IYQIERLoQggRIiTQhRAiREigCyFEiJBAF0KIEGHaxCKlVBGww5Qn90wcsN/sIvysob3mhvZ6QV5zMGmjtY6v7g7TAj1YKaWyapqlFaoa2mtuaK8X5DWHCulyEUKIECGBLoQQIUICvf6mml2ACRraa25orxfkNYcE6UMXQogQIS10IYQIERLoQggRIiTQPaCUul8ppZVScWbX4ktKqZeVUhuUUmuVUp8rpZqaXZOvKKVGK6U2KqXylFIPmV2PrymlUpRSC5VSuUqpHKXUPWbX5C9KKatSarVSapbZtXiLBLqblFIpwEhgp9m1+ME8oJvWugewCXjY5Hp8QillBV4HxgDpwFVKqXRzq/I5O3C/1jodGAD8vgG85kr3AOvNLsKbJNDd91fgD0DIX1XWWs/VWtsrvl0GJJtZjw/1A/K01lu11mXANGCiyTX5lNa6QGu9quLroxgBl2RuVb6nlEoGLgLeNrsWb5JAd4NSaiKwW2v9i9m1mOBm4Buzi/CRJCC/yve7aADhVkkp1RboBfxsbiV+8TeMBpnT7EK8yWZ2AYFKKTUfSKzmrkeBRzC6W0JGba9Xa/1lxTGPYnxE/68/axO+p5SKAT4D7tVaF5tdjy8ppcYB+7TWK5VSQ8yux5sk0GugtR5e3e1Kqe5AO+AXpRQY3Q+rlFL9tNZ7/ViiV9X0eisppW4ExgHDdOhOXtgNpFT5PrnitpCmlArDCPP/aq1nmF2PHwwCJiilxgKRQKxS6kOt9bUm1+UxmVjkIaXUdiBDax2Mq7a5RCk1GpgCDNZaF5ldj68opWwYF32HYQT5CuBqrXWOqYX5kDJaJe8BB7XW95pdj79VtNAf0FqPM7sWb5A+dOGKfwCNgXlKqTVKqTfNLsgXKi783gnMwbg4+Ekoh3mFQcB1wIUV/7ZrKlquIghJC10IIUKEtNCFECJESKALIUSIkEAXQogQIYEuhBAhQgJdCCFChAS6EEKECAl0IYQIEf8Pf4dxS5jWoL4AAAAASUVORK5CYII=\n",
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ]
          },
          "metadata": {
            "tags": [],
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "BierK9JtNzez",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        ""
      ],
      "execution_count": null,
      "outputs": []
    }
  ]
}