24141253_LAB3_Lagrange_Interpolation.ipynb

{"cells":[{"cell_type":"markdown","metadata":{"id":"-H-_6qE4TzPX"},"source":["Make sure you remove `raise NotImplementedError()` and fill in any place that says `# YOUR CODE HERE`, as well as your `NAME`, `ID`, and `SECTION` below:"]},{"cell_type":"code","execution_count":1,"metadata":{"id":"cZBZsJ7uTzPZ","executionInfo":{"status":"ok","timestamp":1711043413147,"user_tz":420,"elapsed":8,"user":{"displayName":"AHMED RIASAT","userId":"08082187387545029877"}}},"outputs":[],"source":["NAME = \"Ahmed Riasat\"\n","ID = \"24141253\"\n","SECTION = \"6\""]},{"cell_type":"markdown","metadata":{"id":"6jRU48H7TzPa"},"source":["---"]},{"cell_type":"markdown","metadata":{"deletable":false,"editable":false,"nbgrader":{"cell_type":"markdown","checksum":"472ee35503d72c0c10415479acd0c38f","grade":false,"grade_id":"title","locked":true,"schema_version":3,"solution":false,"task":false},"id":"JUqq2ULKTzPb"},"source":["# CSE330 Lab: Polynomial Interpolation using Lagrange Form\n","---"]},{"cell_type":"markdown","metadata":{"deletable":false,"editable":false,"nbgrader":{"cell_type":"markdown","checksum":"4a90b5ab8ffd9fe85f9eac720e5f3254","grade":false,"grade_id":"intro1","locked":true,"schema_version":3,"solution":false,"task":false},"id":"tqv6INGMTzPb"},"source":["### Importing necessary libraries"]},{"cell_type":"code","execution_count":2,"metadata":{"deletable":false,"editable":false,"nbgrader":{"cell_type":"code","checksum":"61b07e7eb90c9dd5d430e8e51f42846e","grade":false,"grade_id":"import","locked":true,"schema_version":3,"solution":false,"task":false},"id":"BQv_kayFTzPb","executionInfo":{"status":"ok","timestamp":1711043413148,"user_tz":420,"elapsed":6,"user":{"displayName":"AHMED RIASAT","userId":"08082187387545029877"}}},"outputs":[],"source":["import numpy as np\n","import matplotlib.pyplot as plt"]},{"cell_type":"markdown","metadata":{"deletable":false,"editable":false,"nbgrader":{"cell_type":"markdown","checksum":"81f5a760fa7cea3e76a5396bc16a1c44","grade":false,"grade_id":"intro2","locked":true,"schema_version":3,"solution":false,"task":false},"id":"D6j56_qsTzPc"},"source":["### The Lagrange_Polynomial class\n","General form of an $n$ degree Lagrange polynomial:\n","\n","\\begin{equation}\n","p_n(x) = \\sum_{k=0}^{n} f(x_k)l_k(x) = \\sum_{k=0}^{n} y_kl_k(x),\\tag{1}\n","\\end{equation}\n","\n","where\n","\\begin{equation}\n","l_k(x) = \\prod_{j=0, j\\neq k}^{n} \\frac{x-x_j}{x_k-x_j}. \\tag{2}\n","\\end{equation}\n","\n","Note that the Lagrange method is more efficient than the matrix method because *we do not need to calculate any inverse matrices*.\n","\n","1. **The constructor `__init__(self, data_x, data_y)` is written for you.**\n","\n","     * Here, we check whether the input vectors (numpy arrays) are equal or not.\n","     * We store `data_x` and `data_y`\n","     * We calculate and store the degree of the polynomial.\n","$$\\$$\n","\n","2. **The `_repr__(self)` function has been written for you.**\n","\n","    * This is similar to the `toString()` method in Java. This returns a formatted string of the object whenever the object is printed.\n","$$\\$$\n","\n","3. **You have to implement the `l(self, k, x)` function.**\n","    * This function would take `k` and `x` as inputs and calculate the Lagrange basis using the Equation $(2)$.\n","$$\\$$\n","\n","4. **You have to implement the `__call__(self, x_arr)` function.**\n","    * This function makes an object of a class callable.\n","    * The function calculates the lagrange polynomial from a set of given nodes. `self.data_x` and `self.data_y` contains the coordinates of the given nodes of the original function. Using Equation $(1)$, you have to use `self.data_x`, `self.data_y`, and the `l(k, x_k, x)` function to find the interpolated output of the polynomial for all elements of `x_arr`.\n","`x_arr` is a numpy array containing points through which we want to plot our polynomial."]},{"cell_type":"code","execution_count":3,"metadata":{"deletable":false,"nbgrader":{"cell_type":"code","checksum":"ee5fe25841a28719cac4c1675aa3b244","grade":false,"grade_id":"lag_ans","locked":false,"schema_version":3,"solution":true,"task":false},"id":"ATV3QRfpTzPd","executionInfo":{"status":"ok","timestamp":1711043415184,"user_tz":420,"elapsed":3,"user":{"displayName":"AHMED RIASAT","userId":"08082187387545029877"}}},"outputs":[],"source":["class Lagrange_Polynomial:\n","    def __init__(self, data_x, data_y):\n","        '''\n","        First we need to check whether the input vectors (numpy arrays) are equal\n","        or not.\n","        assert (condition), \"msg\"\n","        this command checks if the condition is true or false. If true, the code\n","        runs normally. But if false, then the code returns an error message \"msg\"\n","        and stops execution\n","        '''\n","        assert len(data_x) == len(data_y), \"length of data_x and data_y must be equal\"\n","\n","        '''\n","        Lagrange polynomials do not use coefficeints a_i, rather the nodes\n","        (x_i, y_i). Hence, we just need to store these inside the object\n","        '''\n","\n","        self.data_x = data_x\n","        self.data_y = data_y\n","\n","        self.degree = len(data_x) - 1\n","        # we assume that the inputs are numpy array, so we can perform\n","        # element wise operations\n","\n","    def __repr__(self):\n","        # method for string representation\n","        # you don't need to worry about the following code if you don't understand\n","        strL = f\"LagrangePolynomial of order {self.degree}\\n\"\n","        strL += \"p(x) = \"\n","        for i in range(len(self.data_y)):\n","            if self.data_y[i] == 0:\n","                continue\n","            elif self.data_y[i] >= 0:\n","                strL += f\"+ {self.data_y[i]}*l_{i}(x) \"\n","            else:\n","                strL += f\"- {-self.data_y[i]}*l_{i}(x) \"\n","\n","        return strL\n","\n","    def l(self, k, x):\n","\n","\n","        l_k = 1.0 # Initialization\n","\n","        # --------------------------------------------\n","        for j in range(self.degree+1):\n","          if j==k:\n","            continue\n","          else:\n","            numerator = x - self.data_x[j]\n","            denominator = self.data_x[k] - self.data_x[j]\n","            l_k *= numerator/denominator\n","\n","        return l_k\n","\n","\n","    def __call__(self, x_arr):\n","\n","        # initialize with zero\n","        p_x_arr  = np.zeros(len(x_arr))\n","\n","        # --------------------------------------------\n","        for k in range(len(x_arr)):\n","          for i in range(self.degree+1):\n","           p_x_arr[k] += self.data_y[i] * self.l(i,x_arr[k])\n","\n","\n","\n","        return p_x_arr"]},{"cell_type":"markdown","metadata":{"deletable":false,"editable":false,"nbgrader":{"cell_type":"markdown","checksum":"478a49eee2ad8430e336e36f9c3192a5","grade":false,"grade_id":"cell-5e302a53e7d98109","locked":true,"schema_version":3,"solution":false,"task":false},"id":"n6dM75WXTzPe"},"source":["### Calling the LagrangePolynomial object and plotting the polynomial."]},{"cell_type":"markdown","metadata":{"deletable":false,"editable":false,"nbgrader":{"cell_type":"markdown","checksum":"501f9a97334297861b6356643408c3a9","grade":false,"grade_id":"plot","locked":true,"schema_version":3,"solution":false,"task":false},"id":"l0NtndljTzPf"},"source":["\n","\n","*Note that in the plot the given nodes will be marked in red.*"]},{"cell_type":"code","source":["import numpy as np\n","\n"],"metadata":{"id":"DDdiwFRleBAI"},"execution_count":null,"outputs":[]},{"cell_type":"code","execution_count":6,"metadata":{"id":"rhysBDi7TzPf","colab":{"base_uri":"https://localhost:8080/","height":507},"executionInfo":{"status":"ok","timestamp":1711046831857,"user_tz":420,"elapsed":765,"user":{"displayName":"AHMED RIASAT","userId":"08082187387545029877"}},"outputId":"003ceac9-dc69-4034-d59d-fb4497976378"},"outputs":[{"output_type":"stream","name":"stdout","text":["LagrangePolynomial of order 6\n","p(x) = + 4.0*l_0(x) + 1.193*l_1(x) + 6.156*l_2(x) + 2.0*l_3(x) + 1.803*l_4(x) + 2.558*l_5(x) \n"]},{"output_type":"display_data","data":{"text/plain":["<Figure size 640x480 with 1 Axes>"],"image/png":"iVBORw0KGgoAAAANSUhEUgAAAj4AAAHHCAYAAAC/R1LgAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAABgCUlEQVR4nO3deVhU9f4H8PeZYd9m2BcBARcEU3FFVBLLUstyKbPypnbNSrM02+x3by51y+5NS69ttolli2lmZV1LS0tTwD0X3BAEkUXWYR2YmfP7Y5xRFBRwZs4s79fzzCOcOXPmw4DMm+8qiKIogoiIiMgByKQugIiIiMhSGHyIiIjIYTD4EBERkcNg8CEiIiKHweBDREREDoPBh4iIiBwGgw8RERE5DAYfIiIichgMPkREROQwGHyIiNooNTUVgiAgJydH6lJMQhAELFy4sM2P2759OwRBwPbt201eE5G5MPgQWSnDm+vevXulLsUuCIJgvMlkMoSFheH222/nmzaRg3GSugAiIku57bbbMHnyZIiiiOzsbLz77ru45ZZb8OOPP2LUqFFSlyeZuro6ODnx7YAcA3/Siei6ampq4OnpKXUZN6xr167429/+Zvx83Lhx6NmzJ5YtW+bQwcfNzU3qEogshl1dRDasoaEB8+fPR9++faFQKODp6Ynk5GRs27btqnNLS0vx0EMPwcfHB0qlElOmTMGhQ4cgCAJSU1ON502dOhVeXl7IysrCHXfcAW9vb0yaNAkAsGPHDkyYMAGRkZFwdXVFREQEnn76adTV1TV5LsM18vPzMXbsWHh5eSEwMBDPPvsstFptu+oCgOPHj+Pee++Fn58f3Nzc0K9fP3z//fftfv169OiBgIAAZGdnG4/99ttvSE5OhqenJ5RKJcaMGYPMzMxrXmfKlCkICAhAY2PjVffdfvvtiI2NNX4uCAJmzZqFjRs34qabboKrqyu6d++OzZs3X/XYAwcOYNSoUfDx8YGXlxduvfVWpKWlNTnH0CW6c+dOPPXUUwgMDIRSqcRjjz2GhoYGVFRUYPLkyfD19YWvry+ef/55iKLY5BpXjvE5e/YsZs6cidjYWLi7u8Pf3x8TJkywmzFN5NjY4kNkw1QqFT766CM88MADmD59OqqqqvDxxx9jxIgRyMjIQEJCAgBAp9PhrrvuQkZGBmbMmIFu3brhu+++w5QpU5q9rkajwYgRIzBkyBAsWbIEHh4eAIB169ahtrYWM2bMgL+/PzIyMrBixQqcO3cO69ata3INrVaLESNGIDExEUuWLMHWrVuxdOlSdOrUCTNmzGhzXUePHsXgwYPRoUMHzJs3D56envj6668xduxYfPPNNxg3blybX7/y8nKUl5ejc+fOAICtW7di1KhRiImJwcKFC1FXV4cVK1Zg8ODB2L9/P6Kiopq9zkMPPYRPP/0UP//8M0aPHm08XlhYiN9++w0LFixocv7OnTuxYcMGzJw5E97e3vjvf/+Le+65B7m5ufD39zd+vcnJyfDx8cHzzz8PZ2dnrFy5EikpKfj999+RmJjY5JpPPvkkQkJCsGjRIqSlpeGDDz6AUqnErl27EBkZiddeew0//fQT3njjDdx0002YPHlyi6/Lnj17sGvXLtx///0IDw9HTk4O3nvvPaSkpODYsWPGnwcimyQSkVVatWqVCEDcs2dPi+doNBpRrVY3OVZeXi4GBweLf//7343HvvnmGxGAuGzZMuMxrVYr3nLLLSIAcdWqVcbjU6ZMEQGI8+bNu+r5amtrrzq2ePFiURAE8ezZs1dd4+WXX25ybu/evcW+ffu2q65bb71V7NGjh1hfX288ptPpxEGDBoldunRp7uVpAoA4bdo08cKFC2JxcbGYnp4u3nrrrSIAcenSpaIoimJCQoIYFBQklpaWGh936NAhUSaTiZMnTzYeM3xvsrOzjTWHh4eLEydObPKcb775pigIgnjmzJkmdbi4uIinT59u8hwAxBUrVhiPjR07VnRxcRGzsrKMx86fPy96e3uLN99881W1jBgxQtTpdMbjSUlJoiAI4uOPP248ptFoxPDwcHHo0KFXvTYLFiwwft7c93n37t0iAPHTTz81Htu2bZsIQNy2bdtV5xNZK3Z1EdkwuVwOFxcXAPrWk7KyMmg0GvTr1w/79+83nrd582Y4Oztj+vTpxmMymQxPPPFEi9c2tMpczt3d3fhxTU0NSkpKMGjQIIiiiAMHDlx1/uOPP97k8+TkZJw5c6bNdZWVleG3337Dfffdh6qqKpSUlKCkpASlpaUYMWIETp06hfz8/Ba/FoOPP/4YgYGBCAoKQmJiIv7880/MnTsXc+bMQUFBAQ4ePIipU6fCz8/P+JiePXvitttuw08//dTidWUyGSZNmoTvv/8eVVVVxuOff/45Bg0ahOjo6CbnDx8+HJ06dWryHD4+PsbXRqvV4pdffsHYsWMRExNjPC80NBQPPvggdu7cCZVK1eSa06ZNgyAIxs8TExMhiiKmTZtmPCaXy9GvX78m34PmXP59bmxsRGlpKTp37gylUtnk54rIFjH4ENm41atXo2fPnnBzc4O/vz8CAwPx448/orKy0njO2bNnERoaelUXhaGL50pOTk4IDw+/6nhubq4xGBjG7QwdOhQAmjwfoB8wGxgY2OSYr68vysvL21zX6dOnIYoiXnrpJQQGBja5GbqRiouLm/1aLjdmzBhs2bIFW7duRXp6OkpKSrB06VLIZDKcPXsWAJqMxzGIi4tDSUkJampqWrz25MmTUVdXh2+//RYAcOLECezbtw8PPfTQVedGRkZedezy1+bChQuora1tsRadToe8vLxrXlOhUAAAIiIirjp++fegOXV1dZg/fz4iIiLg6uqKgIAABAYGoqKi4qrvM5Gt4RgfIhu2Zs0aTJ06FWPHjsVzzz2HoKAgyOVyLF68GFlZWe2+rqurK2Sypn8XabVa3HbbbSgrK8MLL7yAbt26wdPTE/n5+Zg6dSp0Ol2T8+Vyebuf/0qGaz/77LMYMWJEs+e0FOIuFx4ejuHDh5usrsvFx8ejb9++WLNmDSZPnow1a9bAxcUF991331XntvTaiFcMOm6Llq7Z3PHrPc+TTz6JVatWYc6cOUhKSoJCoYAgCLj//vuv+j4T2RoGHyIbtn79esTExGDDhg1NujmuHEzbsWNHbNu2DbW1tU1aV06fPt3q5zp8+DBOnjyJ1atXNxkYu2XLlnbX39q6DN09zs7OZgsuHTt2BKBvqbnS8ePHERAQcN0p/ZMnT8bcuXNRUFCAL774AnfeeSd8fX3bXEtgYCA8PDxarEUmk13VkmNK69evx5QpU7B06VLjsfr6elRUVJjtOYkshV1dRDbM8Nf85X/Bp6enY/fu3U3OGzFiBBobG/Hhhx8aj+l0Orzzzjs39FyiKGL58uXtqr0tdQUFBSElJQUrV65EQUHBVde5cOFCu2swCA0NRUJCAlavXt3kDf7IkSP45ZdfcMcdd1z3Gg888AAEQcDs2bNx5syZJmsGtYVcLsftt9+O7777rskU8qKiInzxxRcYMmQIfHx82nXt1j7/la1CK1asuGopAiJbxBYfIiv3ySefNLvGy+zZszF69Ghs2LAB48aNw5133ons7Gy8//77iI+PR3V1tfHcsWPHYsCAAXjmmWdw+vRpdOvWDd9//z3KysoAoElrUUu6deuGTp064dlnn0V+fj58fHzwzTffXHe8yLW0pa533nkHQ4YMQY8ePTB9+nTExMSgqKgIu3fvxrlz53Do0KF212HwxhtvYNSoUUhKSsK0adOM09kVCkWr9rIKDAzEyJEjsW7dOiiVStx5553truVf//oXtmzZgiFDhmDmzJlwcnLCypUroVar8Z///Kfd122N0aNH47PPPoNCoUB8fDx2796NrVu3GqfaE9kyBh8iK/fee+81e3zq1KmYOnUqCgsLsXLlSvz888+Ij4/HmjVrsG7duiZ7UMnlcvz444+YPXs2Vq9eDZlMhnHjxmHBggUYPHhwq1budXZ2xg8//ICnnnoKixcvhpubG8aNG4dZs2ahV69e7fra2lJXfHw89u7di0WLFiE1NRWlpaUICgpC7969MX/+/HY9/5WGDx+OzZs3Y8GCBZg/fz6cnZ0xdOhQ/Pvf/75qZlZLJk+ejE2bNuG+++6Dq6tru2vp3r07duzYgRdffBGLFy+GTqdDYmIi1qxZc9UaPqa2fPlyyOVyfP7556ivr8fgwYOxdevWFsdXEdkSQbyR0XREZNM2btyIcePGYefOnRg8eLDU5RhZa12t8d1332Hs2LH4448/kJycLHU5RHQFBh8iB1FXV9dkfRatVovbb78de/fuRWFhYZP7WFf7jR49GpmZmTh9+nSruhCJyLLY1UXkIJ588knU1dUhKSkJarUaGzZswK5du/Daa69JGi6sta62+uqrr/DXX3/hxx9/xPLlyxl6iKwUW3yIHMQXX3yBpUuX4vTp06ivr0fnzp0xY8YMzJo1i3WZgCAI8PLywsSJE/H+++/DyYl/VxJZIwYfIiIichhcx4eIiIgcBoMPEREROQx2Ql9Bp9Ph/Pnz8Pb25uBEIiIiGyGKIqqqqhAWFnbVXoOXY/C5wvnz5826Bw4RERGZT15eHsLDw1u8n8HnCt7e3gD0L5w598IhIiIi01GpVIiIiDC+j7eEwecKhu4tHx8fBh8iIiIbc71hKhzcTERERA6DwYeIiIgcBoMPEREROQwGHyIiInIYDD5ERETkMBh8iIiIyGEw+BAREZHDYPAhIiIih8HgQ0RERA6DwYeIiIgcBoMPEREROQwGHyIiInIYDD5E1KxGrQ6l1WqIoih1KUREJsPd2YkcVINGh8LKepwrr8W58rrL/tV/XKiqh04E4kJ9MKFvOMb27gA/TxepyyYiuiGCyD/nmlCpVFAoFKisrISPj4/U5RCZTG2DBit+O429OWU4V16HQlU92vK/31kuYHhcMO7rF4HkLgFwkrPBmIisR2vfv9niQ+QAMgtUmPXFfmRdqGly3NVJhg6+7gj39UC4r/vF26WPnWUy/PDXeXy9Nw9H8lX435FC/O9IIYK8XXFP33BM6BuOmEAvib4qIqK2Y4vPFdjiQ/ZEFEV8np6LlzcdQ4NGhyBvVzw7IhZdgrwQ7uuBAC8XCILQqmsdO6/Cun152HggH+W1jcbjfTv64r5+4bizZxi8XPm3FBFJo7Xv3ww+V2DwIXtRWduIeRv+wv+OFAIAhsUGYsmEXvD3cr2h6zZodPg1swjr9p3D9hPF0F38DeLuLMc9fTvgxVFx8GQAIiILY/BpJwYfsgf7zpbjqS8PIL+iDs5yAS+M7Ia/D46GTNa61p3WKlLVY8P+fKzbm4czJfputE6Bnnj/b33RJdjbpM9FRHQtDD7txOBDtkynE/H+H1lY+stJaHUiIv08sOKB3ugVoTTr84qiiJ2nS/DsukMoUqnh7izH6/f0wJiEDmZ9XiIiAwafdmLwIVtVXFWPuWsPYefpEgDAXb3C8Nq4m+Dt5myxGkqq1Zj91QH8eboUAPC3gZF4aXQ8XJ3kFquBiBxTa9+/OR+VyA78cfIC7li+AztPl8DNWYZ/39MD/70/waKhBwACvFzx6d8T8dQtnQEAa9JyMeH93cgrq7VoHURELWHwIbJhjVodXv/fcUz+JAMl1Q2IDfbGD7OGYGL/yFbP1jI1uUzA3Ntjserh/lB6OOOvc5UYvWInfs0skqQeIqLLMfgQ2bDn1/+F93/PAgBMSozEd7MGW82g4mGxQfjxqWT0ilCisq4R01bvxX82H4dGq5O6NCJyYAw+RDbqu4P5+PZAPuQyAW8/2BuvjusBN2frGkvTQemOdY8lYeqgKADAu9uz8LeP01FcVS9tYUTksBh8iGxQfkUd/rnxCABg1rDOGN0zTOKKWubiJMPCu7tjxQO94ekiR9qZMtz5351IP1MqdWlE5IAYfIhsjE4n4pmvD6KqXoOECCVmXRxIbO3u6hWG72YNQddgL1yoUuPBj9LxefpZqcsiIgfD4ENkYz7aeQZpZ8rg4SLHWxMT4GxDm4V2DvLCxicGY1zvDtDqRPzj2yP4MiNX6rKIyIHYzm9MIsKx8yq88fMJAMD80fGIDvCUuKK283Bxwpv39cIjQ6IBAP/37WGs33dO4qqIyFEw+BDZiPpGLeasPYBGrYjhccGY2D9C6pLaTRAE/OPOOExJ6ghRBJ5bfwjfHcyXuiwicgAMPkQ24j+bT+BkUTUCvFzx73t6SLZOj6kIgoCFd3fHAwMiIYrA3K8P4afDBVKXRUR2jsGHyAbsPFWCT/7MBgC8cW/PG95h3VoIgoBXx96Ee/uGQ6sT8dSXB/DL0UKpyyIiO8bgQ2TlKmob8My6gwD0e18N6xYkbUEmJpMJ+Pc9PTE2IQwanYgnvtiPbceLpS6LiOwUgw+RFRNF/cynIpUaMYGe+Mcd8VKXZBZymYAlE3rhzh6haNSKeGzNPuw4dUF/p1YLbN8OfPml/l+tVspSicjGOUldABG17NsD+fjxcAGcZAKWTUyAu4t1rcxsSk5yGZbdn4BGrQ6/HCvCI6v34vuwYsQu/idw7rJZX+HhwPLlwPjx0hVLRDaLLT5EViqvrBbzvzsKAJgzvAt6hiulLcgCnOUyvP1gH9zSLQgpR3eiyxNTIZ67Yqp7fj5w773Ahg3SFElENk0QRVGUughrolKpoFAoUFlZCR8fH6nLIQel1Yl44IM0ZOSUoV9HX6x9LAlymW3P4mqL+voG1IRFwLe8uPm/zgRB3/KTnQ3I7bcVjIhar7Xv32zxIbJCK//IQkZOGbxcnfDWxASHCj0A4Ja2C/4thR4AEEUgLw/YscOSZRGRHWDwIbIyR/Ir8eYvJwEAC+6KR4Sfh8QVSaCglev5tPY8IqKLGHyIrIh+deaD0OhEjLopBPf2DZe6JGmEhpr2PCKiixh8iKzIun3ncLq4GoHernhtnO2vztxuycn6MTwtfP2iIAAREfrziIjagMGHyEpodSI+/OMMAGDWsM7w9XSRuCIJyeX6KevAVeFHBwCiCN2bb3FgMxG1GYMPkZX435EC5JbVwtfDGff1s90NSE1m/Hhg/XqgQ4cmhwt9AvD42P/Dfzy7S1QYEdkyLmBIZAVEUcT7v2cBAKYMirLrhQrbZPx4YMwY/eytggIgNBTpXtH4ef0R4PcsdPT3wAMDIqWukohsCIMPkRXYlVWKI/kquDvLMSUpSupyrItcDqSkGD8dB+BshRrLtp7CPzceQbivO5K7BEpWHhHZFnZ1EVkBQ2vPxP4Rjj22p5Vm39oF43p3gFYnYuaa/ThRWCV1SURkIxh8iCR2JL8SO06VQC4TMG1ItNTl2ARBEPD6PT0wINoPVWoN/p66B8WqeqnLIiIbwOBDJLGVF2dy3dUz1DEXK2wnVyc5PnioL2ICPJFfUYdHPt2L2gaN1GURkZVj8CGSUG5pLX786zwA4LGhnSSuxvYoPVyw6uH+8PN0wV/nKjH7q4PQ6rj9IBG1jMGHSEIf7jgDnQgM7RqIuFBuitseHf098cFDfeHiJMOWY0VY/FOm1CURkRVj8CGSSGm1Gl/vzQMAPM7WnhvSL8oPSyb0AgB8tDMbn+3OkbYgIrJadhV8Fi5cCEEQmty6desmdVlEzVq9KwdqjQ69whUYGOMndTk27+5eYXhuRCwAYMH3R7HteLHEFRGRNbKr4AMA3bt3R0FBgfG2c+dOqUsiukqNWoPVu88C0Lf2OOyeXCY2M6UTJvQNh04EZn2xH8cLVVKXRERWxu6Cj5OTE0JCQoy3gIAAqUsiusraPXmorGtEdIAnbu8eInU5dkMQBLw6rgeSYvxR06DFtNS9uFCllrosIrIidhd8Tp06hbCwMMTExGDSpEnIzc295vlqtRoqlarJjcicGrU6fLwzGwAwPTkGchlbe0zJxUmG9/7WB9EXp7k/+tle1DdqpS6LiKyEXQWfxMREpKamYvPmzXjvvfeQnZ2N5ORkVFW1vKrr4sWLoVAojLeICG4OSeb1w6HzyK+oQ4CXK8b36XD9B1CbKT1c8PGUflC4O+NAbgWeX/8XRJHT3IkIEEQ7/m1QUVGBjh074s0338S0adOaPUetVkOtvtQUrlKpEBERgcrKSvj4cHoxmZYoihi5bAdOFFXh+ZGxmJnSWeqS7Nqu0yWY/EkGNDoRTw/vitnDu0hdEhGZiUqlgkKhuO77t121+FxJqVSia9euOH36dIvnuLq6wsfHp8mNyFy2n7iAE0VV8HJ1wqTEjlKXY/cGdQ7Av8beBAB4a+tJfH/ovMQVEZHU7Dr4VFdXIysrC6GhoVKXQgQAeO/iZqQPJkZC4e4scTWO4f4BkZierN8D7dl1h3Agt1ziiohISnYVfJ599ln8/vvvyMnJwa5duzBu3DjI5XI88MADUpdGhP255cjILoOzXMDfB3MzUkuaNyoOw+OC0KDRYfqn+5BfUSd1SUQkEbsKPufOncMDDzyA2NhY3HffffD390daWhoCAwOlLo0IKy+29oxN6IAQhZvE1TgWuUzA8vt7Iy7UByXVakxL3YNqNTc0JXJEdj24uT1aOziKqC2yLlRj+Ju/QxSBrXNvRucgb6lLckjnK+ow5p0/caFKjVu7BeGDyf24nACRneDgZiIr8uEfZyCKwPC4YIYeCYUp3fHh5H5wdZLh1+PF3NCUyAEx+BCZWbGqHhv25wMAZqTESFwNJUQosfS+Sxuafplx7UVOici+MPgQmdknf+agQatD/yhf9O3IzUitweieYZh7W1cAwEsbj2DX6RKJKyIiS2HwsQStFti+HfjyS/2/Wi6f7yg0Wh3W7c0DoN+egqzHk7d0xpiEMGh0Ih5fsw+ni6ulLomILIDBx9w2bACiooBhw4AHH9T/GxWlP052b1dWKUprGuDn6YJh3YKkLocuIwgC/n1PT/SJVEJVr8GUTzJQpKqXuiwiMjMGH3PasAG4917g3Lmmx/Pz9ccZfuyeYaXgO3qEwFnO/27Wxs1Zjo+m9EfMxQ1Np3ySgcq6RqnLIiIz4m9ic9FqgdmzgeZWCzAcmzOH3V52rL5Ri5+PFAIA7u7FzUitlZ+nC1b/fQACvV1xvLAKj37K3dyJ7BmDj7ns2HF1S8/lRBHIy9OfR3Zp+4kLqFJrEKpwQ7+OvlKXQ9cQ4eeB1Q8PgLerE9Kzy/D02oPQ6rjEGZE9YvAxl4IC055HNueHi91cd/UKg4yL5Fm9+DAfrJzcFy5yGf53pBCLfjgKru9KZH8YfMyltRujcgNVu1St1mBrZhEA4O5eYRJXQ601qFMA3pzYC4IAfLr7LN7dniV1SURkYgw+5pKcDISHA0ILf+kLAhARoT+P7M6WY4VQa3SIDvBE9zBufWJLRvcMw4LR8QCAN34+ga/35ElcERGZEoOPucjlwPLl+o+vDD+Gz5ct059HdueHQ/ouzLt6hUFoKfyS1Zo6OBozUjoBAF789jB+vdh6R0S2j8HHnMaPB9avBzpcMaMnPFx/fPx4aeoisyqvacAfJy8AYDeXLXt+RCzu6RMOrU7EE1/sx76z5VKXREQmwOBjbuPHAzk5+Pz1VDx113NY/eoqIDuboceO/e9IITQ6EfGhPugc5CV1OdROgiDg9Xt6ICU2EPWNOkxbvYerOxPZAQYfS5DLUTVwCL6PH4pDnRLYvWXnvj+k35D07gS29tg6Z7kM707qg14RSlTUNnJ1ZyI7wOBjIb4ezgCA8toGiSshcyqsrEd6dhkA/fgesn0eLk5YNZWrOxPZCwYfC1F6uAAAKvgL065t+us8RBHo19EXHZTuUpdDJnLl6s6PrN6DarVG6rKIqB0YfCzE1xB8ahl87NnlixaSfYnw80Dqw/3h7eqEPTnlmPRROir5/5nI5jD4WIiSXV12L6ekBofOVUImAHf04MKU9qh7mAJfTB8IpYczDuVV4P4P01BSrZa6LCJqAwYfCzEEn8q6Ru4BZKc2/aVv7RncOQCB3q4SV0Pm0iNcga8eHYgAL1dkFqgwceVuDngmsiEMPhaidNd3dYkioOI4H7v0Pbu5HEa3EB98/dhAhCrckHWhBvet3I1z5bVSl0VErcDgYyEuTjJ4uToBYHeXPTpeqMLJomq4yGUY0T1E6nLIAmICvfD1Y0mI8HPH2dJa3Pf+bmSX1EhdFhFdB4OPBV0a58MWH3vz/UF9a09KbCAU7s4SV0OWEuHngXWPDUKnQE+cr6zHfSt342RRldRlEdE1MPhY0KWZXWzxsSeiKOKHi+N7uGih4wlRuGHtY0noFuKNC1VqTFy5G0fyK6Uui6Sg1QLbtwNffqn/V6uVuiJqBoOPBbHFxz4dyKtAXlkdPFzkuLVbsNTlkAQCvFzx1aMD0StcgfLaRjzwYRr39nI0GzYAUVHAsGHAgw/q/42K0h8nq8LgY0Fs8bFPhm6u2+KD4e7C7UgcldLDBWseScSAKD9U1Wvw0Mfp2J1VKnVZZAkbNgD33gucO9f0eH6+/jjDj1VxkroAR2LYtoKLGNoPrU7Ej4cLAHAndgK83ZyR+vf+eOyzfdhxqgRTV2Vg5UN9kRIbJHVp1A7Vag0KK+tRpNLfClX1EEUgVOGGMKU7whTuCPFyhsvs2fopu1cSRUAQgDlzgDFjuE+jlWDwsSDFxRYfzuqyH+lnSnGhSg2FuzOSuwRKXQ5ZAQ8XJ3w4uR9mfbEfWzOLMf3TvXh1bA9M6BcOQRCkLo8uc6FKjf255SioqEOhSo3ii+GmUFWPYpW6VduSDMz9C19d2dJzOVEE8vKAHTuAlBTTFU/txuBjQWzxsT+GtXvu6BECFyf2HJOem7Mc7/2tL55eexCb/irA89/8hd1nSvGvsTfB05W/dqVS16BFenYpdp4qwc7TJTheeP0ZeN6uTghWuCHExw3BPm4AgILKOhRU1iO/og5B1a0by5Vz+DQ6Dh3K8GsF+D/QgnzZ4mNXGjQ6/O9IIQAuWkhXc5bL8N/7eyM+zAdLfzmJbw/k49C5Crw7qQ+6hfhIXZ5D0OpEHMmvxM7TJdh5qgT7zpajQatrck63EG9EB3gi2McNIQo3BPu46j++GHSuFVRFUYRqsxz44Y3r1jJv9wUU1m/H2N4dMDahA6ICPG/466P2YfCxIM7qsi9/nLyAyrpGBHm7IjHaX+pyyArJZAJmpnRG/yg/PPnFAZy5UIMxb/+JhXd3x/39I/jXvxnkV9Th9xMXsPP0BezKKr2qhb2D0h1DOgdgcJcADO7kD3+v9m8vIwgCFLffCoSH6wcyNzPORxQEVPgF4XB0T9SU1mLZ1lNYtvUU+kQqMa53B9zZMwx+ni7troHajsHHgjiry74Yurnu7BkKuYxvYNSy/lF++Gl2MuZ+fRDbT1zAixsOI+1MKV4d10O/ortWqx8DUlAAhIYCyckcCNsGlbWN+OlIAb7dn4+MnLIm93m7OmFgJ38kdwnAkM4BiA7wNG3glMuB5cv1s7cEoWn4EQQIAHw/eBcZd47AL8cKsWF/Pv48XYL9uRXYn1uBRT8cQ0psIMb1DsetcUFwc+b33dwEUWxuKLrjUqlUUCgUqKyshI+PaZujc0trcfMb2+DmLMPxV0aZ9NpkWbUNGvR9ZSvqGrX4duYg9I70lboksgE6nYgPdpzBGz+fgFYnIibAE6uV5xCx6MWmU6HDw/VvpuPHS1eslVNrtNh+4gK+3Z+P344XG7uwBAHoG+mL5C6BGNIlAL3CFXCSW2D83YYNwOzZTb+PERHAsmVXfR+LVfX4/tB5fHsgH0fPq4zHfT2c8cztsXhgQCT/mGqH1r5/M/hcwZzBR1XfiJ4LfwEAHH9lJJO9Dfvh0Hk8+eUBRPp54PfnUthlQW2yN6cMT355AD0zfsN7G1+DAKDJT5Dh52n9eoafy4iiiH1ny/HtgXxs+qsAlZdt+NwtxBvjenfA3QlhCFW4S1NgO1ruThZVYeOBfGw8kI/zlfUAgO5hPnh5THf07ehniartBoNPO5kz+IiiiM7/+B+0OhFpL96KEIWbSa9PljP9073YcqwITwzrhOdGdJO6HLJB5ao6iFFRUJYXN7+SrCDoW36ys+2326uVQeHMhWpsPJCPbw/mI6+szng82McVYxI6YFzvDogLte0B4xqtDmvSzmLplpOoqtdPox/fpwPmjeqGIG++V7RGa9+/OcbHggRBgNLdGaU1DSivbWDwsVGVdY34/cQFAMDdvTpIXA3ZKt/96UB5ccsn2Pv6L811DV3WxVdZ24gf/jqPb/afw4HcCuMpni5yjLwpFON6d0BSJ3+76RJyksswdXA0RvcKwxubT2Dt3jxs2J+PX44WYc7wLpgyKArOluiycwAMPham9LgUfMg2bT1WhAatDl2DvRAb4i11OWSrCgpMe54tMWzxcEWHg3hxi4ePZv8Hb3h2R4NGP25HLhOQ3CUA43p3wG3xwfBwsd+3rgAvV/z73p54IDESC747gkPnKvGvHzPx1Z48LLq7OwZ3DpC6RJvH+Ghhl2Z2cUq7rdp+Ut/ac3t8iMSVkE0LDW3VaQUeSvPWYWlarb6lp5lRFoIoQhRF3Pnxv6FpaERssDf+eWccdr94C1IfHoAxCR3sOvRcLiFCiW9nDsa/7+kBP08XnC6uxqSP0jHz833Ir6i7/gWoRQw+FqbkIoY2TasTseOUPvikxHKLCroBycn6rp0WBsbrAJz3DsDNuzV4dt0h5JbWWrY+c9mx4+rNPC8jAxBWVYJf+8uweU4yHkmOcdgxLjKZgIn9I7HtmRRMHRQFmQD8dLgQty7djhW/nkJ9o1bqEm0Sg4+FcdsK23boXAUqahvh7eaEhAil1OWQLTOs/wJcHX4EAYIgYMPkZ9EIGdbvO4dhS7fjhfV/Ia/MdgNQZW0j9uw+0qpzozVVnC15kcLDGQvv7o4fn0rGgGg/1DfqsHTLSYx5+0+cLq6Wujybw+BjYb4XV+gsr2GLjy0yDGpO7hJgmbVByL6NH6+fst7hikHy4eEQ1q/HrLdfwLczB+HmroHQ6kSs3ZuHYUu248UNh22muyOvrBaf7MzGAx+koc+/tmDp4evvjwWg1V2BjiQu1AdrHx2I5fcnIMDLFSeKqnD32zvx3cF8qUuzKY7RWWpFuG2FbTOM70npGiRxJWQ3xo8HxoxpcVp370hffPr3Adh3tgxvbTmFnadL8GVGLtbvy8PE/hEY1zscPcMVVjPjR6cTcTi/Elszi7DlWNFVG4Gq+g1E1ZYQeJUWQWhuNRXDNP7kZAtVbFsEQcCYBP2Mtqe+PIC0M2WY/dVBpJ0pw4K74rk+XCsw+FiYYXBzZR1bfGxNWU0D/jpXAQC4uSvH95AJyeXXnbLet6Mf1jySiIzsMry15SR2nynFmrRcrEnLhaeLHAOi/ZDUyR+DOgUgLtTHotO8a9Qa7Mkpw5ZjRdiaWYQildp4n1wmoH+UL26LD8FtccGI9PcAot5pcYsHAPrVju117SITCfJ2w+ePDMTyrSexYttpfJmRi0N5+k1wuQHqtTH4WJjSnS0+NufiImtndh9F4lkVVP0Gcg0mksyAaD98+ehApJ0pxae7c4wbcW47cQHbLnbFKtydMTDGD0kx/hjUOQBdgrxMMl6mvlGL08XVOFlUhZNF+n9PFFZd1e3m6SLH0NhA3BYfjGGxQcZJHUaGLr7m1vFpZosHap5cJmDu7bHoF+WHOWsP4liBCqNX7MS/7+mJO3uyq7AlXLn5CuZcuRkAdmeV4oEP0xAT6Infnkkx+fXJxJpZZE0VEAyfle/ylzNZBZ1ORGahCruzSrE7qxTp2WWoVmuanBPg5YKBMf64qYMCrk4yOMtlcJHL4OwkwFl+6XOXi/c5ywXIBAFny2px6mK4OVVcjbOlNdC18I4RpnDDsG5BuC0+GEmd/OHq1IoWG27OajKFlfV48sv92JNTDgCYktQR/3dnXOu+D3aCW1a0k7mDz/FCFUYu2wE/Txfsf+k2k1+fTKilRdYu7rjMfZTIGmm0OhzOr8TuM/ogtCenDPWNOpNd39fDGV2DvREb4o0uwd6IDfZG12Cvq1t1yOI0Wh2W/HIS7/+eBQDoGa7AOw/2QYSfh8SVWQaDTzuZO/gUqeqR+NqvkAnA6VfvgMxOllu3O1otEBXV8nojjrCPEtkFtUaLQ3mV2JVVgtzSWjTqRDRqdGjU6tCg1aHh4seNWtF4rFGrg0YrIlThhtgQb3QNvnQL8HLhNHMr99vxIsz9+hAqahvh4+aEJRN64fbu9r/gKvfqslKGWV06Eaiq10Bx8XOyMtdZZM3u91Eiu+HqpB/4PCCaO307ilu6BePHp5Ix64v9OJBbgUc/24dpQ6Ixb1Q3q5n9JyW+Ahbm6iSHh4u+hYCrN1sxR95HiYhsXgelO9Y+moRHhkQDAD7emY0HP0xDsape4sqkx+AjAV9uW2H9Wrt4GhdZIyIr5eIkwz9Hx2PlQ33h7eqEPTnluHPFTmRkl0ldmqQYfCRg6O6qqOOUdqt1nX2UIAhARAQXWSMiqzeiewi+f3IIYoO9caFKjQc+TMPHO7PhqEN8GXwkcGmHdrb4WK3L9lEScfU+SgC4yBoR2YzoAE98+8QgjEkIg1Yn4pVNx/DklwdQc8XSB46AwUcChgHN5TVs8bFq48dDXLcOxT7+TY+Hh3MqOxHZHA8XJyybmICFd8XDSSZg018FGPfunzhzwbE2OuWsLglc2qGdLT7WLjPpNox+7GMMKTiOj27vAJeIDlxkjYhsliAImDo4Gjd1UGDm5/txsqgad7/9J5ZM6IWRN9n/lHeALT6SuDS4mS0+1u73kxegk8nhdOswuDw0ST91naGHiGxcvyg/bHpqCAZE+6FarcHja/bh9f8dh0ZrusUurRWDjwSUnNVlM7afKAYApMRyU1Iisi/6jU4TjVPe3/89C5M/yUBptfo6j7RtDD4SuNTVxRYfa1ZV34h9Z/X73gzlbuxEZIec5fop728/2BseLnLsyirF6BU7jb/77BGDjwS4jo9t+PN0KTQ6EdEBnujo7yl1OUREZjO6Zxi+e2IwYgI9UVBZj3vf34UX1v+FEjts/WHwkYCSLT424feTFwCwtYeIHEOXYG9898RgjO/TAaIIrN2bh2FvbMdHO86g0Y7G/jD4SIDr+Fg/URTx+8XxPUM5voeIHIS3mzPevC8B38xIQo8OClSpNfjXj5kYtXwHdpy6IHV5JsHgIwFDi09NgxYNGvtJ0fbkdHE1zlfWw8VJhoHR/td/ABGRHenb0Q8bnxiM18f3gL+nC04XV+OhjzPw6Kd7kVta276LarXA9u3Al1/q/9VqTVlyq9ll8HnnnXcQFRUFNzc3JCYmIiMjQ+qSmvBxc4bs4uK/bPWxTttP6P+yGRjjD3cXTl8nIscjlwm4f0Akfns2BX8fHA25TMAvx4ow/K3fseTnE6htaMOqzxs2AFFRwLBhwIMP6v+NitIftzC7Cz5r167F3LlzsWDBAuzfvx+9evXCiBEjUFxcLHVpRjKZAIX7xdWbOc7HKnF8DxGRnsLdGfPvisfm2ckY0jkADRod3t52Grcu/R3fHzp//T2/NmwA7r0XOHeu6fH8fP1xC4cfQbSzXcoSExPRv39/vP322wAAnU6HiIgIPPnkk5g3b951H69SqaBQKFBZWQkfHx+z1XnLku04U1KDrx4diIEx7EqxJjVqDXq/vAUNWh1+fWYoOgV6SV0SEZFVEEURvxwrwiubjuFceR0AoHekEj07KODj7gwfN2co3J3h4+4EHzdn+LjI0C2pJ+Tn86/c9VBPEPTbAGVn3/DisK19/7arLSsaGhqwb98+vPjii8ZjMpkMw4cPx+7du5t9jFqthlp9abqeSqUye53A5TO72NVlbdLOlKJBq0O4rztiAjiNnYjIQBAEjOgegqFdA/HhH2fwzvbTOJBbgQO5Fc2ePzD3L3x1Pr/lC4oikJcH7NihXxnfAuwq+JSUlECr1SI4OLjJ8eDgYBw/frzZxyxevBiLFi2yRHlNcNsK62UY35MSGwhBaPZvFCIih+bmLMeTt3bBPX3D8b8jhaiobUBlXSNUdY1Q1WuMH3fJrmrdBQsKzFvwZewq+LTHiy++iLlz5xo/V6lUiIiIMPvzctsK6ySKIrafvDiNvWuQxNUQEVm3MKU7pl3c8qJZ20Xg68XXv1BoqOmKug67Cj4BAQGQy+UoKipqcryoqAghIc3vOuvq6gpXV1dLlNeEYduKSrb4WJXskhrkldXBWS5gUCeOvSIiuiHJyfoxPPn5+m6tKxnG+CQnW6wku5rV5eLigr59++LXX381HtPpdPj111+RlJQkYWVX8/Vki481Mszm6h/lB09Xu/q7gIjI8uRyYPly/cdXDh0wfL5s2Q0PbG4Luwo+ADB37lx8+OGHWL16NTIzMzFjxgzU1NTg4Ycflrq0Jjid3TpdPr6HiIhMYPx4YP16oEOHpsfDw/XHx4+3aDl29yftxIkTceHCBcyfPx+FhYVISEjA5s2brxrwLDVuW2F96hu1SDtTCoDje4iITGr8eGDMGP3srYIC/Zie5GSLtvQY2F3wAYBZs2Zh1qxZUpdxTYYxPmzxsR5pZ0qh1ugQqnBD12Cu3UNEZFJyucWmrF+L3XV12QolW3yszuWrNXMaOxGRfWLwkYivp2EBw8brL/dNFvE7x/cQEdk9Bh+JGMb4aHQiqtRt2OiNzCK3tBZnSmrgJBMwqHOA1OUQEZGZMPhIxM1ZDjdn/ctfUcNxPlL7/eKihX06+sLHzVniaoiIyFwYfCRknNlVx3E+UuNu7EREjoHBR0Jcy8c6NGp12J1lmMbO4ENEZM8YfCTEtXysw5H8StQ0aKFwd0Z8qI/U5RARkRkx+EjIMLOrvIbBR0ppZ8oAAInRfpDJOI2diMieMfhI6NIO7ezqkpJhtebEGG5KSkRk7xh8JGRYvZldXdLRaHXYm6Nv8RkY4ydxNUREZG4MPhLyZYuP5I6cV6GmQQsfNyd0C+H4HiIie8fgI6FLXV1s8ZFK+sVurgHR/pBzfA8Rkd1j8JGQoaurso4tPlIxjO9hNxcRkWNg8JEQW3ykpdHqsCenHAAwkAObiYgcAoOPhJSGwc3cskISxwpUqFZr4O3mhDiu30NE5BAYfCRkGNxcpdagUauTuBrHY+jmGhDlx/E9REQOgsFHQgp3ZwgX328rOLPL4tLPGKaxs5uLiMhRMPhISC4TjDuBcy0fy9LqRGRkM/gQETkaBh+JGWZ2cS0fy8osUKFKrYG3qxPiwzi+h4jIUTD4SIwzu6RhGN/TP5rje4iIHAmDj8S4bYU0jPtzRXP9HiIiR8LgIzHDzC4ObrYcju8hInJcDD4SU3CMj8VlFqigqtfAy9UJ3Tm+h4jIoTD4SOxSiw+7uiwl/WJrT78oXzjJ+V+AiMiR8Le+xC7N6mLwsZRL+3Oxm4uIyNEw+Ejs0qwudnVZgo7je4iIHBqDj8TY1WVZxwurUFnXCE8XOW7i+B4iIofD4CMxJQc3W5Shm6tflB/H9xAROSD+5peYr+elFh9RFCWuxv6lZ19cvyeG6/cQETkiBh+JGQY3N2pF1DZoJa7Gvul0onFGF8f3EBE5JgYfibk7y+FyscuFM7vM62RxFSpqG+HhIkePDgqpyyEiIgkw+EhMEATjOB+u3mxeaVn6bq6+HX3hzPE9REQOib/9rYAvNyq1iLQz7OYiInJ0DD5WgDO7zE+nE5GRYwg+HNhMROSoGHysANfyMb9TxdUoq2mAu7McPToopS6HiIgkwuBjBXw9L7b41LDFx1wM09j7RfnCxYk/9kREjorvAFZAyTE+ZmdYuDAxmt1cRESOjMHHCvgaZ3Ux+JiDKIpI58BmIiICg49VMLT4VNSxq8scThdXo7SmAW7OMvQMV0pdDhERSYjBxwoo3Tmry5wM3Vx9O3J8DxGRo+O7gBW4fL8uMr00wzYV0ezmIiJydAw+VsAwxqe8hsHH1PTjewwbkzL4EBE5OgYfK2AY46Oq10Cj1UlcjX3JulCNkuoGuDrJ0CuC+3MRETk6Bh8rYBjjAwCVHOBsUoZtKvpE+sLVSS5xNUREJDUGHyvgJJfB280JAAc4m5phYDOnsRMREdCO4DNlyhT88ccf5qjFoXHbCtMTRRHp2dyfi4iILmlz8KmsrMTw4cPRpUsXvPbaa8jPzzdHXQ7n0iKGbPExlTMlNbhQpYaLkwy9IpRSl0NERFagzcFn48aNyM/Px4wZM7B27VpERUVh1KhRWL9+PRob+abdXgpuW2Fy6cbxPUq4OXN8DxERtXOMT2BgIObOnYtDhw4hPT0dnTt3xkMPPYSwsDA8/fTTOHXqlKnrtHts8TE9ju8hIqIr3dDg5oKCAmzZsgVbtmyBXC7HHXfcgcOHDyM+Ph5vvfWWqWp0CL5s8TEpURQv25iUwYeIiPTaHHwaGxvxzTffYPTo0ejYsSPWrVuHOXPm4Pz581i9ejW2bt2Kr7/+Gi+//LI56rVbSg9uW2FKOaW1KL44vqd3pFLqcoiIyEo4tfUBoaGh0Ol0eOCBB5CRkYGEhISrzhk2bBiUSqUJynMcnNVlWobWnoQIju8hIqJL2hx83nrrLUyYMAFubm4tnqNUKpGdnX1DhTmaSy0+DD6msCdHP7A5MZrT2ImI6JI2B5+HHnrIHHU4vEstPuzqMoW9OeUAgP5RDD5ERHQJV262EhzcbDpFqnrkltVCJoDje4iIqAkGHyuh5HR2kzG09sSF+sDbzfk6ZxMRkSNh8LEShuCj1uhQ16CVuBrbZhjfw24uIiK6EoOPlfBydYKTTADA7q4btfesPvj0i/KVuBIiIrI2DD5WQhAEKDnO54ZVqzU4dl4FAOjXkS0+RETUFIOPFeG2FTfuQG45dCIQ4eeOEEXLSy4QEZFjYvCxIpzZdeP2GKaxs7WHiIiaYVfBJyoqCoIgNLm9/vrrUpfVaty24sbtzTGM72HwISKiq7V5AUNr9/LLL2P69OnGz729vSWspm2MixjWsMWnPRq1OhzIrQAA9OfAZiIiaobdBR9vb2+EhIRIXUa7KD0vjvGpY4tPexw7r0JdoxZKD2d0CvSSuhwiIrJCdtXVBQCvv/46/P390bt3b7zxxhvQaDTXPF+tVkOlUjW5SUXpzjE+N8Kwfk+/jr6QXVwagIiI6HJ21eLz1FNPoU+fPvDz88OuXbvw4osvoqCgAG+++WaLj1m8eDEWLVpkwSpbxlldN8awYjPH9xARUUusvsVn3rx5Vw1YvvJ2/PhxAMDcuXORkpKCnj174vHHH8fSpUuxYsUKqNXqFq//4osvorKy0njLy8uz1Jd2Fa7j036iKBoXLuT4HiIiaonVt/g888wzmDp16jXPiYmJafZ4YmIiNBoNcnJyEBsb2+w5rq6ucHV1vdEyTYItPu2XU1qLkuoGuDjJcFMHhdTlEBGRlbL64BMYGIjAwMB2PfbgwYOQyWQICgoycVXm4evJFp/2MozvSQhXwtVJLnE1RERkraw++LTW7t27kZ6ejmHDhsHb2xu7d+/G008/jb/97W/w9bWNrg/DOj6VdY3Q6kTIOUC31S6t32Mb32siIpKG3QQfV1dXfPXVV1i4cCHUajWio6Px9NNPY+7cuVKX1mqGWV2iCKjqGo0tQHR9hoHN3JGdiIiuxW6CT58+fZCWliZ1GTfExUkGL1cnVKs1KK9tYPBppZJqNc6U1EAQgD6RbPEhIqKWWf2sLkejcOcihm1laO2JDfaG4mJ3IRERUXMYfKyMr2H1Zg5wbjWO7yEiotZi8LEyxh3aa9ji01p7znJ8DxERtQ6Dj5XhIoZtU9ugwdH8SgBcsZmIiK6PwcfKcBHDtjmYVwGNTkSYwg0dlO5Sl0NERFaOwcfKsMWnbbg/FxERtQWDj5Vhi0/bGFZs5v5cRETUGgw+VsaXLT6tptHqsP8sW3yIiKj1GHysjJItPq12vLAKNQ1aeLs5oWuwt9TlEBGRDWDwsTKGMT5cx+f6DOv39O3oy33NiIioVRh8rIxhjE85W3yui+v3EBFRWzH4WBlDi09doxb1jVqJq7FeoiheWrG5Iwc2ExFR6zD4WBkfNydjtw3H+bTsXHkdilRqOMsF9IpQSl0OERHZCAYfKyMIApTuhu4ujvNpiWEae48OCrg5yyWuhoiIbAWDjxVSejD4XM+eHI7vISKitmPwsUK+xpld7OpqyaUd2Rl8iIio9Rh8rBC3rbi28poGnCquBqCfyk5ERNRaDD5WiIsYXtu+i9PYOwd5wc/TReJqiIjIljD4WKFL+3Wxxac5e85yfy4iImofBh8rdKmriy0+zTHuyN6R43uIiKhtGHysUICXPvhcqFJLXIn1qW/U4q9zFQA4o4uIiNqOwccKRfh6AADyymslrsT6/HWuEo1aEYHerojwc5e6HCIisjEMPlYowk8ffM6V1UGnEyWuxroYFi7sH+ULQeDGpERE1DYMPlYoVOEGuUxAg1aHoqp6qcuxKpf252I3FxERtR2DjxVyksvQQanvxsktZXeXgU4nYi93ZCciohvA4GOlIi92d+WWMfgYnCyuQlW9Bh4ucsSFektdDhER2SAGHytlGOeTV14ncSXWw7A/V59IXzjJ+aNLRERtx3cPK2Vo8clji4/Rpf25uHAhERG1D4OPlTJM1WZX1yV7uSM7ERHdIAYfK8UxPk3lV9Qhv6IOcpmAhAil1OUQEZGNYvCxUobgc6FKjboGrcTVSC8juxQAcFMHBTxdnSSuhoiIbBWDj5VSuDvD203/Bn+OKzgjI1s/vicxmt1cRETUfgw+VkoQBHZ3XSb9YvAZwPE9RER0Axh8rJhhzy5HDz4XqtQ4c6EGgsCBzUREdGMYfKxYpD+DD3Bpf67YYG8oPJwlroaIiGwZg48Vi+BaPgAuje8ZwPE9RER0gxh8rNilRQwde/VmBh8iIjIVBh8rdvngZlEUJa5GGpV1jcgsVAHgwGYiIrpxDD5WLEzpBkEA6hq1KKlukLocSew7WwZRBKIDPBHk4yZ1OUREZOMYfKyYq5McoRff7B11gDOnsRMRkSkx+Fg5wwBnR13EkON7iIjIlBh8rJxxnE+p4wWf2gYNDp+rBMDgQ0REpsHgY+UcefXmg7kV0OhEhCrcEO7rLnU5RERkBxh8rFyEAwef9Mu6uQRBkLgaIiKyBww+Vs6RFzHk+B4iIjI1Bh8rZ+jqKlDVo0Gjk7gay2nQ6LA/txwAd2QnIiLTYfCxcgFeLnB3lkMUgfwKx1nB+XB+BdQaHfw8XdAp0EvqcoiIyE4w+Fg5QRAccoDz5ev3cHwPERGZCoOPDYjw089ocqTgw/E9RERkDgw+NsDRBjhrdSL25ejH9zD4EBGRKTH42IBIBws+mQUqVKk18HJ1Qlyoj9TlEBGRHWHwsQGONsbH0M3VL8oXchnH9xARkekw+NiAy7etEEVR4mrMj+N7iIjIXBh8bEC4rz74VKk1qKxrlLga8xJFERk5+uDD9XuIiMjUGHxsgLuLHIHergDsv7sr60I1ymoa4OokQ48OSqnLISIiO8PgYyMuDXC270UMM7L1s7n6RPrCxYk/nkREZFp8Z7ERjjLAOSO7FADQn91cRERkBgw+NsIRdmkXRdG4YjPH9xARkTkw+NiICF/96s32vJbPufI6FFTWw0kmoHekUupyiIjIDjH42AhH6OoyTGPvEa6Ah4uTxNUQEZE9YvCxEZH++uBzvqIOGq1O4mrMg+v3EBGRudlM8Hn11VcxaNAgeHh4QKlUNntObm4u7rzzTnh4eCAoKAjPPfccNBqNZQs1k2BvN7jIZdDoRBRU1ktdjllw/R4iIjI3mwk+DQ0NmDBhAmbMmNHs/VqtFnfeeScaGhqwa9curF69GqmpqZg/f76FKzUPmUxAuJ/9jvMprqpHdkkNBAHo25HBh4iIzMNmgs+iRYvw9NNPo0ePHs3e/8svv+DYsWNYs2YNEhISMGrUKLzyyit455130NDQYOFqzSPC137H+ey5uH5PtxAfKNydJa6GiIjslc0En+vZvXs3evTogeDgYOOxESNGQKVS4ejRoy0+Tq1WQ6VSNblZK3se4GxYv4fdXEREZE52E3wKCwubhB4Axs8LCwtbfNzixYuhUCiMt4iICLPWeSOMqzeX29/qzekc2ExERBYgafCZN28eBEG45u348eNmreHFF19EZWWl8ZaXl2fW57sR9rqIYUVtA04UVQEA+kcx+BARkflIuljKM888g6lTp17znJiYmFZdKyQkBBkZGU2OFRUVGe9riaurK1xdXVv1HFK7tF+XfQWfvTnlEEUgJtDTuBkrERGROUgafAIDAxEYGGiSayUlJeHVV19FcXExgoKCAABbtmyBj48P4uPjTfIcUou4OKurrKYBVfWN8Hazj0HAnMZORESWYjNjfHJzc3Hw4EHk5uZCq9Xi4MGDOHjwIKqrqwEAt99+O+Lj4/HQQw/h0KFD+Pnnn/HPf/4TTzzxhM206FyPt5szfD30YceedmnnwoVERGQpNrMvwPz587F69Wrj57179wYAbNu2DSkpKZDL5di0aRNmzJiBpKQkeHp6YsqUKXj55ZelKtksIv08UF5bibzyWsSH+Uhdzg2rUWtwJL8SAMf3EBGR+dlM8ElNTUVqauo1z+nYsSN++uknyxQkkQg/Dxw6V2k343wO5FZAoxPRQemO8IvrFBEREZmLzXR1kZ69zewyrN/Dbi4iIrIEBh8bY2+LGHL9HiIisiQGHxtjT8FHrdHiQF4FAAYfIiKyDAYfG2MIPufK66DTiRJXc2P+OleJBo0OAV4uiAnwlLocIiJyAAw+NiZU4Qa5TECDRofiKrXU5dyQy6exC4IgcTVEROQIGHxsjJNchjClGwDb7+7aceoCACAx2l/iSoiIyFEw+NggexjnU1XfiL055QCAoV1Ns3o3ERHR9TD42CB7CD5/ni6FRiciyt8DURzfQ0REFsLgY4MMa/mcs+Hg8/vJYgBASmyQxJUQEZEjYfCxQbbe4iOKIraf0I/vGRrLbi4iIrIcBh8bFOFr28HnZFE1Cirr4eokQ1IMBzYTEZHlMPjYIEOLT3GVGnUNWomrabvtJ/TdXEmd/OHmLJe4GiIiciQMPjZI6eEMb1f9/rLnym2v1cfQzZXC2VxERGRhDD42SBAE4wDnPBsLPtVqDfae1S9cyIHNRERkaQw+Nso4wLnUtoLPn6dL0KjlNHYiIpIGg4+NivBzBwDkltVJXEnbGLu52NpDREQSYPCxUbY4pV0URfx+cWAzp7ETEZEUGHxslHGMjw0Fn1PF1TjPaexERCQhBh8bFXnZ4GZRFCWupnW2Hde39gyM4TR2IiKSBoOPjerg6w5BAGobtCitaZC6nFa5NL6H3VxERCQNBh8b5eokR4iPGwDbGOfDaexERGQNGHxsmC2N8zFMY+/o74FoTmMnIiKJMPjYsEgbCj6Gbq5hbO0hIiIJMfjYMFuZ0s5p7EREZC0YfGyYrQQfTmMnIiJrweBjwwyrN+dZ+erNht3YOY2diIikxuBjwwyDm89X1qFBo5O4mpZxGjsREVkLBh8bFujlCjdnGUQROF9hna0+1WoN9uRwGjsREVkHBh8bJgiC1Y/z2cVp7EREZEUYfGyctQef7ScvdnN1ZTcXERFJj8HHxoX7Wu9aPvpp7IbxPezmIiIi6TH42DhrbvE5XVyN/Io6uDjJMJDT2ImIyAow+Ni4y3dptzaG2VwDY/zh7sJp7EREJD0GHxsX6X+xxafU+oLPtovr93B8DxERWQsGHxsXcXGMj6peg8raRomruaTpNHYGHyIisg5OUhdAN8bdRY4AL1eUVKtxpqQavSN9pS4JAKexE9GN02q1aGy0nj/oSFrOzs6Qy2982ASDjx3oFa7Ar8eL8cfJEqsJPpdPYxcEQeJqiMiWiKKIwsJCVFRUSF0KWRmlUomQkJAbel9h8LEDI24Kwa/Hi7H5aCFmD+8idTmcxk5EN8QQeoKCguDh4cE/ngiiKKK2thbFxfqxo6Ghoe2+FoOPHRgeFwy5TEBmgQpnS2vQ0V/ariVOYyei9tJqtcbQ4+/P3x90ibu7fmPu4uJiBAUFtbvbi4Ob7YCfpwsSo/0AAD8fLZS4Gk5jJ6L2M4zp8fDwkLgSskaGn4sbGfvF4GMnRt4UAgDYfMQKgs9JTmMnohvD7i1qjil+Lhh87MSI7vrgsz+3AkWqesnqqFFrsCe7HACnsRORY0lJScGcOXOkLuOGCIKAjRs3SvLcCxcuREJCgtmfh8HHTgT7uKFPpBKAtN1du7JK0aDVIdKP09iJyLFs2LABr7zySqvOzcnJgSAIOHjwoHmLMjNLhRVTYvCxI9bQ3bXdsFpzLKexE5Fj8fPzg7e3t8Wfl2sdtQ2Djx0xdHelZ5ehrKbB4s8viqJxYDO7uYjI0Vze1RUVFYXXXnsNf//73+Ht7Y3IyEh88MEHxnOjo6MBAL1794YgCEhJSTHe99FHHyEuLg5ubm7o1q0b3n33XeN9hpaitWvXYujQoXBzc8Pnn3+O1NRUKJVKbNy4EV26dIGbmxtGjBiBvLy8JjW+99576NSpE1xcXBAbG4vPPvvsml/TCy+8gK5du8LDwwMxMTF46aWXjEErNTUVixYtwqFDhyAIAgRBQGpqKgCgoqICjzzyCAIDA+Hj44NbbrkFhw4danLt119/HcHBwfD29sa0adNQX2+ZYRoMPnako78n4kJ9oNWJ2JpZZPHnP5xfaZzGnhQTYPHnJyL7JIoiahs0ktxEUWx33UuXLkW/fv1w4MABzJw5EzNmzMCJEycAABkZGQCArVu3oqCgABs2bAAAfP7555g/fz5effVVZGZm4rXXXsNLL72E1atXN7n2vHnzMHv2bGRmZmLEiBEAgNraWrz66qv49NNP8eeff6KiogL333+/8THffvstZs+ejWeeeQZHjhzBY489hocffhjbtm1r8Wvw9vZGamoqjh07huXLl+PDDz/EW2+9BQCYOHEinnnmGXTv3h0FBQUoKCjAxIkTAQATJkxAcXEx/ve//2Hfvn3o06cPbr31VpSV6bcy+vrrr7Fw4UK89tpr2Lt3L0JDQ5sEPHPiOj52ZmT3EGQWqPDzkULc1y/Cos/90Y5sAMAdN4VwGjsRmUxdoxbx83+W5LmPvTwCHi7te6u84447MHPmTAD6lpO33noL27ZtQ2xsLAID9a3i/v7+CAkJMT5mwYIFWLp0KcaPHw9A3zJ07NgxrFy5ElOmTDGeN2fOHOM5Bo2NjXj77beRmJgIAFi9ejXi4uKQkZGBAQMGYMmSJZg6daqxprlz5yItLQ1LlizBsGHDmv0a/vnPfxo/joqKwrPPPouvvvoKzz//PNzd3eHl5QUnJ6cmX8POnTuRkZGB4uJiuLq6AgCWLFmCjRs3Yv369Xj00UexbNkyTJs2DdOmTQMA/Otf/8LWrVst0urDFh87Yxjns+NUCarVGos977nyWvx4uAAAMP3mGIs9LxGRterZs6fxY0EQEBISYlx5uDk1NTXIysrCtGnT4OXlZbz961//QlZWVpNz+/Xrd9XjnZyc0L9/f+Pn3bp1g1KpRGZmJgAgMzMTgwcPbvKYwYMHG+9vztq1azF48GCEhITAy8sL//znP5Gbm3vNr/vQoUOorq6Gv79/k68jOzvb+HVkZmYaA5pBUlLSNa9rKmzxsTNdg70QHeCJ7JIabDtejLt6hVnkeVf9mQOtTsSQzgHoHqawyHMSkWNwd5bj2MsjJHvu9nJ2dm7yuSAI0Ol0LZ5fXV0NAPjwww+vCgVXrlLs6Wn+WbO7d+/GpEmTsGjRIowYMQIKhQJfffUVli5des3HVVdXIzQ0FNu3b7/qPqVSaZ5i24DBx84IgoAR3UPw/u9Z2Hy00CLBp7KuEV9l6P8CYGsPEZmaIAjt7m6yVi4uLgD0W3QYBAcHIywsDGfOnMGkSZPafE2NRoO9e/diwIABAIATJ06goqICcXFxAIC4uDj8+eefTbrM/vzzT8THxzd7vV27dqFjx474xz/+YTx29uzZq76Oy78GAOjTpw8KCwvh5OSEqKioZq8dFxeH9PR0TJ482XgsLS2t9V/sDbCvnyQCoO/uev/3LGw7Xoz6Ri3cbuAvltb4MiMXNQ1axAZ74+YuHNRMRHQ9QUFBcHd3x+bNmxEeHg43NzcoFAosWrQITz31FBQKBUaOHAm1Wo29e/eivLwcc+fOveY1nZ2d8eSTT+K///0vnJycMGvWLAwcONAYhJ577jncd9996N27N4YPH44ffvgBGzZswNatW5u9XpcuXZCbm4uvvvoK/fv3x48//ohvv/22yTlRUVHIzs7GwYMHER4eDm9vbwwfPhxJSUkYO3Ys/vOf/6Br1644f/48fvzxR4wbNw79+vXD7NmzMXXqVPTr1w+DBw/G559/jqNHjyImxvx/PHOMjx3qFa5AqMINtQ1a7DxVYtbnatDosOpP/aDm6TfHcO0eIqJWcHJywn//+1+sXLkSYWFhGDNmDADgkUcewUcffYRVq1ahR48eGDp0KFJTU43T36/Fw8MDL7zwAh588EEMHjwYXl5eWLt2rfH+sWPHYvny5ViyZAm6d++OlStXYtWqVU2m0l/u7rvvxtNPP41Zs2YhISEBu3btwksvvdTknHvuuQcjR47EsGHDEBgYiC+//BKCIOCnn37CzTffjIcffhhdu3bF/fffj7NnzyI4OBiAfkbYSy+9hOeffx59+/bF2bNnMWPGjHa+mm0jiDcyV88OqVQqKBQKVFZWwsfHR+py2m3h90eRuisH9/YNx5IJvcz2PN/sO4dn1h1CsI8rdjx/C1ycmKWJqP3q6+uRnZ2N6OhouLm5SV2OzUhNTcWcOXNQUVEhdSlmda2fj9a+f/Ndyk4ZZndtzSxCo7blwXQ3QhRFfLjjDABg6qBohh4iIrJ6fKeyU/2j/ODv6YKK2kZkZJeZ5Tn+OFWC44VV8HSR48HESLM8BxERkSkx+NgpuUzAbfH6vlRz7d314R/61p6J/SOhcHe+ztlERGQuU6dOtftuLlNh8LFjIy52d/18tBA6nWmHch09X4mdp0sglwl4eHCUSa9NRERkLgw+dmxQJ394uzqhuEqNA3nlJr22cXuKHqGI8PMw6bWJiIjMhcHHjrk6yXFLXBAA03Z3na+oww+HzgMApidff4olERGRtWDwsXMju+u7uzYfLbyhXYYvl7orBxqdiIExfugZrjTJNYmIiCyBwcfODY0NhJuzDHlldThWoLrh66nqG/FFun57ike5PQUREdkYmwk+r776KgYNGgQPD48WNzkTBOGq21dffWXZQq2Mh4sThnYNBAD8bILurrUZeahWa9A5yAspXYNu+HpERESWZDPBp6GhARMmTLjuktarVq1CQUGB8TZ27FjLFGjFDIsZbj56Y8GnUavDJ4btKZKjIZNxewoiInNbuHAhEhISpC4DqampVrG7+o2ymU1KFy1aBED/wl+LUqlESEiIBSqyHbd0C4aTTMDJompkXahGp0Cvdl3nx78KUFBZjwAvV4xJ6GDiKomITEyrBXbsAAoKgNBQIDkZkJt302Z7NnHiRNxxxx1tekxKSgoSEhKwbNky8xTVDjbT4tNaTzzxBAICAjBgwAB88skn1x3Qq1aroVKpmtzsjcLdGYM663dN/7mdrT6iKOKDPwzbU3Q0+47vREQ3ZMMGICoKGDYMePBB/b9RUfrj1C7u7u4ICrL9IQ52FXxefvllfP3119iyZQvuuecezJw5EytWrLjmYxYvXgyFQmG8RUREWKhayzLM7mrvOJ9dWaU4VqCCu7MckxI7mrI0IiLT2rABuPde4Ny5psfz8/XHzRR+UlJS8NRTT+H555+Hn58fQkJCsHDhwibn5ObmYsyYMfDy8oKPjw/uu+8+FBUVNTnn9ddfR3BwMLy9vTFt2jTU19df9VwfffQR4uLi4Obmhm7duuHdd9+9bm2zZs3CrFmzoFAoEBAQgJdeeqlJ40B5eTkmT54MX19feHh4YNSoUTh16pTx/iu7ugxdcJ999hmioqKgUChw//33o6qqCoB+Nenff/8dy5cvN467zcnJQXl5OSZNmoTAwEC4u7ujS5cuWLVqVWtf5hsnSuiFF14QAVzzlpmZ2eQxq1atEhUKRauu/9JLL4nh4eHXPKe+vl6srKw03vLy8kQAYmVlZXu/LKtUrKoXo+ZtEju+sEk8V17b5sdP/jhd7PjCJnH+xsNmqI6ISK+urk48duyYWFdX174LaDSiGB4uikDzN0EQxYgI/XkmNnToUNHHx0dcuHChePLkSXH16tWiIAjiL7/8IoqiKGq1WjEhIUEcMmSIuHfvXjEtLU3s27evOHToUOM11q5dK7q6uoofffSRePz4cfEf//iH6O3tLfbq1ct4zpo1a8TQ0FDxm2++Ec+cOSN+8803op+fn5iamnrN2ry8vMTZs2eLx48fF9esWSN6eHiIH3zwgfGcu+++W4yLixP/+OMP8eDBg+KIESPEzp07iw0NDaIoXv3+u2DBAtHLy0scP368ePjwYfGPP/4QQ0JCxP/7v/8TRVEUKyoqxKSkJHH69OliQUGBWFBQIGo0GvGJJ54QExISxD179ojZ2dnili1bxO+//75Vr/G1fj4qKytb9f4tafApLi4WMzMzr3lTq9VNHtOW4LNp0yYRgFhfX9/qmlr7wtmiCe/tEju+sEn8ZOeZNj3ueIFK7PjCJjF63iYxp6TaTNUREZkg+Gzb1nLoufy2bZspyxZFUR8uhgwZ0uRY//79xRdeeEEURVH85ZdfRLlcLubm5hrvP3r0qAhAzMjIEEVRFJOSksSZM2c2uUZiYmKT4NOpUyfxiy++aHLOK6+8IiYlJV2ztri4OFGn0xmPvfDCC2JcXJwoiqJ48uRJEYD4559/Gu8vKSkR3d3dxa+//loUxeaDj4eHh6hSqYzHnnvuOTExMbHJ886ePbtJLXfddZf48MMPt1jrtZgi+Eg6uDkwMBCBgYFmu/7Bgwfh6+sLV1dXsz2HLRlxUwgycsqw+UghHh7c+hWXP9yhH9sz8qYQdPT3NFd5REQ3rqDAtOe1Uc+ePZt8HhoaiuLiYgBAZmYmIiIimgypiI+Ph1KpRGZmJvr374/MzEw8/vjjTa6RlJSEbdu2AQBqamqQlZWFadOmYfr06cZzNBoNFArFNWsbOHAgBOHSbNykpCQsXboUWq0WmZmZcHJyQmJiovF+f39/xMbGIjMzs8VrRkVFwdvbu9mvtyUzZszAPffcg/379+P222/H2LFjMWjQoGs+xpRsZlZXbm4uysrKkJubC61Wi4MHDwIAOnfuDC8vL/zwww8oKirCwIED4ebmhi1btuC1117Ds88+K23hVmRE92C8sukY9uSUoaRajQCv6wfCIlU9vjuYDwCYnswFC4nIyoWGmva8NnJ2dm7yuSAI0Ol0Jrt+dXU1AODDDz9sElIAQC7BjLX2fL2jRo3C2bNn8dNPP2HLli249dZb8cQTT2DJkiXmLNXIZoLP/PnzsXr1auPnvXv3BgBs27YNKSkpcHZ2xjvvvIOnn34aoiiic+fOePPNN5skYkcX7uuBHh0UOJxfiS/Sc5EY7YfKukbjTXXZx4ZbQWU9GrUi+kf5onekr9RfAhHRtSUnA+Hh+oHMzc3qFQT9/cnJFi8tLi4OeXl5yMvLM7b6HDt2DBUVFYiPjzeek56ejsmTJxsfl5aWZvw4ODgYYWFhOHPmDCZNmtSm509PT2/yeVpaGrp06QK5XI64uDhoNBqkp6cbW19KS0tx4sQJY23t4eLiAq1We9XxwMBATJkyBVOmTEFycjKee+45Bp8rpaamXnMNn5EjR2LkyJGWK8hGjbwpBIfzK/HmlpNtetwTwzqbqSIiIhOSy4Hly/WztwShafgxdPMsWybJej7Dhw9Hjx49MGnSJCxbtgwajQYzZ87E0KFD0a9fPwDA7NmzMXXqVPTr1w+DBw/G559/jqNHjyIm5lKL+6JFi/DUU09BoVBg5MiRUKvV2Lt3L8rLyzF37twWnz83Nxdz587FY489hv3792PFihVYunQpAKBLly4YM2YMpk+fjpUrV8Lb2xvz5s1Dhw4dMGbMmHZ/zVFRUUhPT0dOTg68vLzg5+eHhQsXom/fvujevTvUajU2bdqEuLi4dj9HW9lM8CHTuKdPONbtzYOqXgOFuzN83J2hMN6cLvvY2Xh/B6U7x/YQke0YPx5Yvx6YPbvplPbwcH3oGT9ekrIEQcB3332HJ598EjfffDNkMhlGjhzZZNmViRMnIisrC88//zzq6+txzz33YMaMGfj555+N5zzyyCPw8PDAG2+8geeeew6enp7o0aMH5syZc83nnzx5Murq6jBgwADI5XLMnj0bjz76qPH+VatWYfbs2Rg9ejQaGhpw880346effrqqO6stnn32WUyZMgXx8fGoq6tDdnY2XFxc8OKLLyInJwfu7u5ITk626PZSgiiaaMtuO6FSqaBQKFBZWQkfHx+pyyEicij19fXIzs5GdHQ03NzcbuxiXLnZyBpXUG6Pa/18tPb9my0+RERkn+RyICVF6irIytjVys1ERERE18IWHyIiIju3fft2qUuwGmzxISIiIofB4ENEREQOg8GHiIisDiccU3NM8XPB4ENERFbDsGZMbW2txJWQNTL8XNzI2kIc3ExERFZDLpdDqVQaN7r08PBosrEmOSZRFFFbW4vi4mIolcob2peMwYeIiKxKSEgIAFx3l29yPEql0vjz0V4MPkREZFUEQUBoaCiCgoLQ2NgodTlkJZydnU2yAz2DDxERWSW5XG6SNzqiy3FwMxERETkMBh8iIiJyGAw+RERE5DA4xucKhsWRVCqVxJUQERFRaxnet6+3yCGDzxWqqqoAABERERJXQkRERG1VVVUFhULR4v2CyHXBm9DpdDh//jy8vb0detEslUqFiIgI5OXlwcfHR+py7BJfY/Pi62t+fI3Ni69v24iiiKqqKoSFhUEma3kkD1t8riCTyRAeHi51GVbDx8eH/+HMjK+xefH1NT++xubF17f1rtXSY8DBzUREROQwGHyIiIjIYTD4ULNcXV2xYMECuLq6Sl2K3eJrbF58fc2Pr7F58fU1Dw5uJiIiIofBFh8iIiJyGAw+RERE5DAYfIiIiMhhMPgQERGRw2DwoVa5++67ERkZCTc3N4SGhuKhhx7C+fPnpS7LLuTk5GDatGmIjo6Gu7s7OnXqhAULFqChoUHq0uzGq6++ikGDBsHDwwNKpVLqcuzCO++8g6ioKLi5uSExMREZGRlSl2Q3/vjjD9x1110ICwuDIAjYuHGj1CXZFQYfapVhw4bh66+/xokTJ/DNN98gKysL9957r9Rl2YXjx49Dp9Nh5cqVOHr0KN566y28//77+L//+z+pS7MbDQ0NmDBhAmbMmCF1KXZh7dq1mDt3LhYsWID9+/ejV69eGDFiBIqLi6UuzS7U1NSgV69eeOedd6QuxS5xOju1y/fff4+xY8dCrVbD2dlZ6nLszhtvvIH33nsPZ86ckboUu5Kamoo5c+agoqJC6lJsWmJiIvr374+3334bgH6Pw4iICDz55JOYN2+exNXZF0EQ8O2332Ls2LFSl2I32OJDbVZWVobPP/8cgwYNYugxk8rKSvj5+UldBtFVGhoasG/fPgwfPtx4TCaTYfjw4di9e7eElRG1DoMPtdoLL7wAT09P+Pv7Izc3F999953UJdml06dPY8WKFXjsscekLoXoKiUlJdBqtQgODm5yPDg4GIWFhRJVRdR6DD4ObN68eRAE4Zq348ePG89/7rnncODAAfzyyy+Qy+WYPHky2FPasra+vgCQn5+PkSNHYsKECZg+fbpElduG9ry+REROUhdA0nnmmWcwderUa54TExNj/DggIAABAQHo2rUr4uLiEBERgbS0NCQlJZm5UtvU1tf3/PnzGDZsGAYNGoQPPvjAzNXZvra+vmQaAQEBkMvlKCoqanK8qKgIISEhElVF1HoMPg4sMDAQgYGB7XqsTqcDAKjValOWZFfa8vrm5+dj2LBh6Nu3L1atWgWZjI2x13MjP7/Ufi4uLujbty9+/fVX44BbnU6HX3/9FbNmzZK2OKJWYPCh60pPT8eePXswZMgQ+Pr6IisrCy+99BI6derE1h4TyM/PR0pKCjp27IglS5bgwoULxvv4F7Rp5ObmoqysDLm5udBqtTh48CAAoHPnzvDy8pK2OBs0d+5cTJkyBf369cOAAQOwbNky1NTU4OGHH5a6NLtQXV2N06dPGz/Pzs7GwYMH4efnh8jISAkrsxMi0XX89ddf4rBhw0Q/Pz/R1dVVjIqKEh9//HHx3LlzUpdmF1atWiUCaPZGpjFlypRmX99t27ZJXZrNWrFihRgZGSm6uLiIAwYMENPS0qQuyW5s27at2Z/XKVOmSF2aXeA6PkREROQwOJCAiIiIHAaDDxERETkMBh8iIiJyGAw+RERE5DAYfIiIiMhhMPgQERGRw2DwISIiIofB4ENEREQOg8GHiIiIHAaDDxERETkMBh8ismsXLlxASEgIXnvtNeOxXbt2wcXFBb/++quElRGRFLhXFxHZvZ9++gljx47Frl27EBsbi4SEBIwZMwZvvvmm1KURkYUx+BCRQ3jiiSewdetW9OvXD4cPH8aePXvg6uoqdVlEZGEMPkTkEOrq6nDTTTchLy8P+/btQ48ePaQuiYgkwDE+ROQQsrKycP78eeh0OuTk5EhdDhFJhC0+RGT3GhoaMGDAACQkJCA2NhbLli3D4cOHERQUJHVpRGRhDD5EZPeee+45rF+/HocOHYKXlxeGDh0KhUKBTZs2SV0aEVkYu7qIyK5t374dy5Ytw2effQYfHx/IZDJ89tln2LFjB9577z2pyyMiC2OLDxERETkMtvgQERGRw2DwISIiIofB4ENEREQOg8GHiIiIHAaDDxERETkMBh8iIiJyGAw+RERE5DAYfIiIiMhhMPgQERGRw2DwISIiIofB4ENEREQOg8GHiIiIHMb/A6ig0VlDZLxAAAAAAElFTkSuQmCC\n"},"metadata":{}}],"source":["data_x = np.array([-3.5, -2.231, -1.152, -0.5, 0.198, 0.785, 1.6])\n","data_y = np.array([4.0, 1.193, 6.156, 2.0, 1.803, 2.558, 0.0])\n","\n","p = Lagrange_Polynomial(data_x, data_y)\n","print(p)\n","\n","#generating 50 points from -2.55 to 1.189 in order to create a smooth line\n","x_arr = np.linspace(-3.5, 1.6, 40)\n","p_x_arr = p(x_arr)\n","\n","# plot to see if your implementation is correct\n","#google the functions to understand what each parameters mean, if not apparent\n","plt.plot(x_arr, p_x_arr)\n","plt.plot(data_x, data_y, 'ro')\n","plt.legend(['interpolated', 'node points'], loc = 'lower right')\n","plt.xlabel('x')\n","plt.ylabel('y')\n","plt.title('Lagrange Polynomial')\n","\n","plt.show()\n"]},{"cell_type":"markdown","source":["###Problem related Lagrange interpolation\n","\n","> Suppose, you have a function f(x) = 4x and three nodes (2, 8), (4, 16), (6, 24). Using Lagrange basis, print out the value of the interpolating polynomial at x = -1.5. Also, display the actual interpolation error at x = -1.5.\n","\n","Hint: Interpolation error = | f(-1.5) - p(-1.5) | where p is the interpolating polynomial.\n","\n","`You have to solve this problem using Lagrange_Polynomial class`"],"metadata":{"id":"-ESx-i4BGtBO"}},{"cell_type":"code","source":["#your code here\n","data_x = [2,4,6]\n","data_y = [8,16,24]\n","P = Lagrange_Polynomial(data_x,data_y)\n","print(P)\n","x_arr = np.array([-1.5])\n","p_x_arr = P(x_arr)\n","print(p_x_arr)\n","fx = 4*-1.5\n","#print(fx)\n","print()\n","error = fx - p_x_arr\n","print(f\"Interpolation ERROR at  {-1.5} : {abs(error)}\")"],"metadata":{"id":"g0jZzstDIndZ","colab":{"base_uri":"https://localhost:8080/"},"executionInfo":{"status":"ok","timestamp":1711043873604,"user_tz":420,"elapsed":7,"user":{"displayName":"AHMED RIASAT","userId":"08082187387545029877"}},"outputId":"af915cd4-6da6-4db8-8063-19f523c883f4"},"execution_count":4,"outputs":[{"output_type":"stream","name":"stdout","text":["LagrangePolynomial of order 2\n","p(x) = + 8*l_0(x) + 16*l_1(x) + 24*l_2(x) \n","[-6.]\n","\n","Interpolation ERROR at  -1.5 : [0.]\n"]}]},{"cell_type":"code","source":[],"metadata":{"id":"6jT2xQHiismZ"},"execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"codemirror_mode":{"name":"ipython","version":3},"file_extension":".py","mimetype":"text/x-python","name":"python","nbconvert_exporter":"python","pygments_lexer":"ipython3","version":"3.7.4"},"colab":{"provenance":[]}},"nbformat":4,"nbformat_minor":0}