{ "cells": [ { "cell_type": "markdown", "id": "468aa069", "metadata": {}, "source": [ "# Tutorial 2 : using the `qmat.nodes` module\n", "\n", "📜 _The_ `NodeGenerator` _class from_ `qmat.nodes` _allows to generate sets of quadrature nodes associated to various types of orthogonal polynomials._\n", "_It is based on the book of W. Gautschi : [Orthogonal Polynomials: Computation and Approximation](https://doi.org/10.1093/oso/9780198506720.001.0001)._\n", "\n", "Gauss quadrature approximate integrals by a given quadrature rule on $M$ points :\n", "\n", "$$\n", "\\int_{-1}^{1} f(t)\\omega(t)dt \\simeq \\sum_{m=1}^{M} \\omega^M_m f(\\tau^M_m)\n", "$$\n", "\n", "where $f(t)$ is a function of interest, $\\omega(t)$ a weight function and $\\tau^M_m, \\omega^M_m$ are the **quadrature nodes and weights** \n", "associated to the Gauss quadrature on $M$ points. \n", "In particular, the quadrature nodes $\\tau^M_m$ are the roots a given polynomial belonging to \n", "an orthogonal basis with respect to the scalar product :\n", "\n", "$$\n", "\\langle p,q \\rangle = \\int_{-1}^{1} p(t)q(t) \\omega(t)dt\n", "$$\n", "\n", "This polynomial basis is solely determined by the weights function $\\omega(t)$, and classical polynomial basis exist already in the literature :\n", "\n", "- $\\omega(t) = 1$ : Legendre polynomials\n", "- $\\omega(t) = \\frac{1}{\\sqrt{1-t^2}}$ : Chebyshev polynomials of the 1st kind\n", "- $\\omega(t) = \\sqrt{1-t^2}$ : Chebyshev polynomials of the 2nd kind\n", "- $\\omega(t) = \\frac{\\sqrt{1+t^2}}{\\sqrt{1-t^2}}$ : Chebyshev polynomials of the 3rd kind\n", "- $\\omega(t) = \\frac{\\sqrt{1-t^2}}{\\sqrt{1+t^2}}$ : Chebyshev polynomials of the 4th kind" ] }, { "cell_type": "markdown", "id": "99c66e96", "metadata": {}, "source": [ "## Node generation\n", "\n", "The type of polynomial defines a **node type**, and the roots of the $M^{th}$ degree polynomial of this type are then the quadrature nodes $\\tau^M_m$\n", "and can be generated, _e.g_ for $M=4$ :" ] }, { "cell_type": "code", "execution_count": 1, "id": "b71b9640", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[-0.86113631 -0.33998104 0.33998104 0.86113631]\n" ] } ], "source": [ "from qmat.nodes import NodesGenerator\n", "\n", "gen = NodesGenerator(nodeType=\"LEGENDRE\", quadType=\"GAUSS\")\n", "nodes = gen.getNodes(nNodes=4)\n", "print(nodes)" ] }, { "cell_type": "markdown", "id": "fd22885d", "metadata": {}, "source": [ "💡 Note that those nodes symmetrically distributed, which is not necessarily the case for other types of nodes _e.g_ for the Chebyshev polynomials of the fourth kind :" ] }, { "cell_type": "code", "execution_count": 2, "id": "a4cd7d0a", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[-0.93969262 -0.5 0.17364818 0.76604444]\n" ] } ], "source": [ "gen = NodesGenerator(nodeType=\"CHEBY-4\", quadType=\"GAUSS\")\n", "nodes = gen.getNodes(nNodes=4)\n", "print(nodes)" ] }, { "cell_type": "markdown", "id": "644b2a74", "metadata": {}, "source": [ "Different types of nodes are available, checkout `qmat.nodes.NODE_TYPES` for the current list, in particular :\n", "\n", "- `LEGENDRE` : for Legendre polynomials\n", "- `CHEBY-1` : for the Chebyshev polynomials of the first kind\n", "- `CHEBY-2` : ...\n", "\n", "💡 You may noticed that those nodes are always **strictly included in** $[-1,1]$, hence usually named _Gauss points_.\n", "But four specific **quadrature types** can be considered for each node type (_i.e_ for each polynomial basis) :\n", "\n", "- `GAUSS` : nodes do not include $-1$ or $1$,\n", "- `LOBATTO` : nodes include $-1$ and $1$,\n", "- `RADAU-LEFT` : nodes include $-1$ (usually called Radau-I),\n", "- `RADAU-RIGHT` : nodes include $1$ (usually called Radau-II).\n", "\n", "The quadrature type is selected when instantiating the node generator, as for the node type :" ] }, { "cell_type": "code", "execution_count": 3, "id": "0b92caf1", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[-1. -0.4472136 0.4472136 1. ]\n" ] } ], "source": [ "gen = NodesGenerator(nodeType=\"LEGENDRE\", quadType=\"LOBATTO\") # usually called Gauss-Lobatto in the literature\n", "nodes = gen.getNodes(nNodes=4)\n", "print(nodes)" ] }, { "cell_type": "code", "execution_count": 4, "id": "0a713599", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[-1. -0.42720716 0.36290645 0.92144357]\n" ] } ], "source": [ "gen = NodesGenerator(nodeType=\"CHEBY-3\", quadType=\"RADAU-LEFT\")\n", "nodes = gen.getNodes(nNodes=4)\n", "print(nodes)" ] }, { "cell_type": "markdown", "id": "c720e329", "metadata": {}, "source": [ "> 📣 We use the naming convention `RADAU-RIGHT` / `RADAU-LEFT` as it is more explicit than the usual one in the literature, \n", "> and also because Radau-I and Radau-II nodes are usually associated to the Legendre polynomials." ] }, { "cell_type": "markdown", "id": "bb59a6cb", "metadata": {}, "source": [ "## Orthogonal polynomials\n", "\n", "The node computation process relies on the **three term recurrence coefficients** associated to each polynomial basis :\n", "\n", "$$\n", "\\begin{gather}\n", "\\forall j \\in \\mathbb{N}, \\quad\n", "\\pi_{j+1}(t) = (t-\\alpha_j)\\pi_{j}(t) - \\beta_j \\pi_{j-1}(t), \\\\\n", "\\pi_{-1}(t) = 0, \\quad \\pi_0(t) = 1.\n", "\\end{gather}\n", "$$\n", "\n", "Those coefficients are know analytically for each polynomial basis (_i.e_ node type), \n", "and are used to generate the tri-diagonal Jacobi matrix for the weight function $\\omega$\n", "\n", "$$\n", "J_\\infty^{\\omega} = \\begin{bmatrix}\n", "\\alpha_0 & \\sqrt{\\beta_1} & & & \\\\\n", "\\sqrt{\\beta_1} & \\alpha_1 & \\sqrt{\\beta_2} & & \\\\\n", " & \\sqrt{\\beta_2} & \\alpha_2 & \\sqrt{\\beta_2} & \\\\\n", " & & \\ddots & \\ddots & \\ddots\n", "\\end{bmatrix}.\n", "$$\n", "\n", "Computing the eigenvalues of the leading principal sub-matrix of size $M$ allows to retrieve the `GAUSS` nodes,\n", "and some small modifications of this sub-matrix allow to retrieve the `LOBATTO`, `RADAU-LEFT` and `RADAU-RIGHT` nodes.\n", "\n", "But one can also use the orthogonal coefficients to evaluate the orthogonal polynomial of any degree, _e.g_ for $M=5$ :" ] }, { "cell_type": "code", "execution_count": 5, "id": "0fe4c282", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjgAAAGwCAYAAACkfh/eAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjcsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvTLEjVAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAjFxJREFUeJztnQWYHFXWhr9xd/d4ZuKuJCGeQAJBg4VFwsKySJB/WRb3ZZfFHQLBCRCCBiIQd3eXyWR8JuMu/T/n3q4en8xMurukz/s8narpVHffrqquOvfId5xMJpMJDMMwDMMwBsJZ7QEwDMMwDMNYGzZwGIZhGIYxHGzgMAzDMAxjONjAYRiGYRjGcLCBwzAMwzCM4WADh2EYhmEYw8EGDsMwDMMwhsMVDkhtbS3S0tLg5+cHJycntYfDMAzDMEwbIOm+oqIiREdHw9m5dR+NQxo4ZNzExcWpPQyGYRiGYTpASkoKYmNjW93GIQ0c8twoO8jf31/t4TAMwzAM0wYKCwuFg0K5j7eGQxo4SliKjBs2cBiGYRhGX7QlvYSTjBmGYRiGMRxs4DAMwzAMYzjYwGEYhmEYxnCwgcMwDMMwjOFgA4dhGIZhGMPBBg7DMAzDMIaDDRyGYRiGYQwHGzgMwzAMwxgONnAYhmEYhjEcbOAwDMMwDGM42MBhGIZhGMZwsIHDMAzDMIzhcMhmm7aiptaEsyWVKK6oRudQH7WHwzAOR3VNLc6WVgImwM3FGf5ebnBxPndTPoZhrEd5VQ1S88sQ6uOBAG83qAUbOFbkYHohZryxDuF+HtjyyCRrvjXDMM1MKLadOos1R7Ox50wBDqYXIbekAiZT3TZuLk6ICfRCUpQ/hnYKxtgeYegW7sv7kmFsyKGMIsx6az2iAzyx4eGJUAs2cKxIqK+HWOaWVKK21gRnnjkyjNU5lVOCzzYlY/HOVOExbYyTE0A+m1oTUFVjwqncUvH4bV+G+H8ydq4YFIPZQ+Pg56ne7JJhjEpucYVYhpjviWrBBo4VCfZxt8wsC8qqEGT+m2GY8+dYVjH+t+ywxVAhArzcML5nGIZ2Dkaf6ADEBHkhyNtdhKUoXJVRWI7TuaXYmZKPTSdysfF4rvC0PvtrIV774yhuHJmA28d1hT8bOgxjNWiSX/+eqBZs4FgRd1dnccEl44Zc5WzgMMz5U1hehf/8fghfbUkRkweCjJobR3bCBd1DRa5Nc7i6OCM2yFs8RnULxd/Hd0NeSSV+3ZuOj9efxPHsEry18jgWbk3BA1N6YvaQOPa6MowVyC2WBk6ILxs4hoIOKBk4OcWV6Bau9mgYRt/8eSgT//p+n/DEEJOSIvCPaT3RI8KvQ+9Hk44bRiTgumHxWH4wEy/+fggnskvw8Pd78eOuVLx0VX9hEDEM03HOlphDVCp7cLhM3MpQ1jiRY45BMgzTfqpqavHUz/txy4JtwrjpFOKNr24bgQ//MqTDxk19KD9uau9ILJ03Fo/N6AUvNxdsOnEW015di593p/EhYxgreHCCzfdDteAQlZUJ9XNvcIAZhmkfWYXluPOLHdiWnCf+vvWCznhwSk94ubtYfVdSeIvef2JiOB74dje2J+fh7q92Yl9aAf4xNZFLzBmmA2SbJ/ihKoeo2INjZULMFquSRc4wTNs5liXLS8m48fNwxftzBksPiw2Mm/p0CvXBN7ePxB3juoq/31t9AnM/2YrSymqbfi7DGJFc8wQ/1E9dDw4bOFZGSarKaaZ8lWGYliHvyZXvbkRaQTm6hPrgp7svwJTekXbbZVR59c/piXjj2oHwdHPGysPZmDN/CwpKq/iwMUw7oCKb+ikbasEGjpVR6v5zitiDwzBthUq4r/9wE/JLqzAgLhDf/W2UamrgM/tH44u5I+Dv6SqMrtnvb+ScOoZpI6QBV+fB4RCVoQgze3AUHQCGYVqH1IhvWbAV5VW1GNcjDF/eNlx1/YzBCUH45o6RQpWcVFnZk8MwbZd1qDbLOaj9O2YPjo08OJyDwzDnZndKPm76mHJdajCmeyjemzMY3u7aqH1IjPTHwttHCoVyEgf8y8dbRJ85hmFaRqkg9vN0hYerbXPnzgUbOFZGqfvnKiqGaZ3k3BLcvGCrMBpGdAnG+3OGwNNN3QtiYyhM9sXc4Qj0dsOulHzc/tk2UcLOMEzzkAYcEaZymwaCDRwbeXCKKqpFR1WGYZqSX1opjBvqJdUnxh/z/zLU5pVSHaVnpB8+vWUYfNxdsP5YLh5dvA+m+h09GYbRnIoxwQaOlaHERHezdDzn4TBMUyqra3H7Z9uFgjB1G/7oL0Ph46GNsFRL9IsNxBvXDQT1z124LQXvrzmh9pAYRtMhqlD24BgPJycni+XKeTgM05Snf9mPzSfPwtfDFR/dPBTh/p662E0TEiOEJg/x798P4Y+DmWoPiWE03EncXe2hsAfHplo4LPbHMA34fscZfL7pNJycIPRmKJFXT9w0qhNuGBEPilDdt3AXUs6Wqj0khtEU2UqJOHtwjK1mrCRbMQwDHMooxL8W7xW74p4J3TE+MVyXHtrHZ/QWWj2F5dWipQTn2jFMcx4cTjI2JIrlypVUDCOhSqm/fU7GQC3G9gjDPRO763bXuLs6463rB4nKqr2pBXjmlwNqD4lhNEOuWQMuVGUNHIKTjG2A0mCMc3AYRvL0z/txMkcmFb86e4Dum1jGBHqJ70Ghti82n8bS/RlqD4lhtJVk7MceHEPCOTgMUwfd/L/ZdkYYA6/MHqC6uqm1uLBnOP46totYf/j7vcjm9iwMA0uZuAZ+5+zBsWVHcW7XwDg4WUXl4uZPkDEwvEsIjMT9k3sgMdJP6Pn8c9Ee1sdhHJryqhqL2jd7cAyKcmA5yZhxZEgM76Hv9oibf1KUvzAGjAZJ0b96zQChffXHoSx8vTVF7SExjOrhKXcXZ/hpQNuKPTg2bdfAHcUZx2XxzlSsPJwtknJfu2aA6n1pbAWVuv/f1J5i/blfDyK9oEztITGM6irGVHGoNmzg2LKKqqRStI5nGEeDvDZKddG9E7ujR4QfjMytF3TGoPhA4Z5//Mf9HKpiHJIcDakYE2zg2AAlibKm1oSCsipbfATDaJpnfzmAvNIqkZ+iJOIaGWdnJ/z7in5wc3HC8gOZ+H0fV1UxjkeuhvpQEWzg2AByyQd4uYn13BIOUzGOxZoj2fh+Z6qomnrh8r5wM/dmMzrkpfrbuK5i/Ymf9vPkhnE4stmD42il4qxmzDhWFcUjP8iqqb+M7ISB8UFwJO4c3w1dwnyQVVSB//x+SO3hMIxdYQ+OgxBqadfAHhzGcaAu2ylnyxAV4IkHzYm3joSnmwuev6yvWP9yy2nsSy1Qe0gMYzeUiEUY5+AYm1A/pZKKPTiMY5CaX4a3Vx0T6/+6KEl0C3dERnQJwSX9o0VDzid/4oRjxnHI0VAnccIxguNqiv2xB4dxEF5YclD0mhrWORgz+kXBkXn4okR4ublgW3IeftqdpvZwGMbOKsYOVEX19ttvo3PnzvD09MTgwYOxdu3aFrdNT0/Hddddh549e8LZ2Rnz5s1rdrtFixahV69e8PDwEMvFixdDkzk4rGbMOACbTuTilz3poBZTT87srQkNDDWJCvDC38fLhOMXlhxCaaVUd2UYI5PjaEnGCxcuFEbKI488gp07d2LMmDGYPn06Tp8+3ez2FRUVCAsLE9v379+/2W02btyI2bNnY86cOdi9e7dYXn311di8eTO0gtIqPof70zAGh+QQnvpZat5cNzwevaL91R6SJpg7pgvigr2QUViOt1ceV3s4DGPz6wDpX9VvOG14A+fll1/Grbfeirlz5yIpKQmvvvoq4uLi8M477zS7fadOnfDaa6/hxhtvREBAQLPb0HtMnjwZDz/8MBITE8Vy4sSJ4nmtEKZ0FGcPDmNwFm0/g4PphUIa4YHJjpdY3FrC8aMX9xLrH6w9wQrHjKHJL62EomurlYa6NjVwKisrsX37dkyZMqXB8/T3hg0bOvy+5MFp/J5Tp05t8T3JK1RYWNjgYS8PDufgMEYvC395+RGxfveEbgjSyIVNK0zpFYFhnYJRUV2LV5cfVXs4DGMzFEmUIG83uGpE+8qmo8jJyUFNTQ0iIiIaPE9/Z2R0XOmTXtue93zhhReEN0h5kAfJfv2ouIqKMS4LNpwSIZiYQC/cMCJB7eFoDspFemh6olj/dnsKjmYWqT0khrEJuRrLvyHsYmY1TjikLsPnm4TYnvekEFZBQYHlkZKSYjcPTlFFtZjlMozRKCitwtsrZVn4fZN7iJAM05TBCUGY1jtSuO9f/P0w7yLG0CrGIRrJv7G5gRMaGgoXF5cmnpWsrKwmHpj2EBkZ2a73pEorf3//Bg9b4+/pKlrGEyz2xxiRt1cfQ2F5NXpG+OGygTFqD0fT/N+0nnBxdsKKg5nYeuqs2sNhGBuqGDuIB8fd3V2UhS9fvrzB8/T3qFGjOvy+I0eObPKey5YtO6/3tDbkTVIyybldA2M00vLL8PH6U2L9oeny5s20TNcwX1w9RIbG//3bIe42zhiOXI2pGBM2lxq9//77RRn3kCFDhGHy/vvvixLxO+64wxI+Sk1Nxaeffmp5za5du8SyuLgY2dnZ4m8ylkjvhrj33nsxduxYvPjii7j00kvx448/YsWKFVi3bh20RJifB9IKyrlUnDEcb648hsrqWpFAO75nuNrD0QXzJnXH4p1nsD05D6uPZONC3m+MgcgpUkT+3B3HwCG9mtzcXDz99NNCxK9Pnz5YsmQJEhJkQiI911gTZ+DAgZZ1qsL68ssvxfanTskZI3lqvv76azz66KN47LHH0LVrV6G3M3z4cGjNwKkfm2QYI3AmrxTfbpN5bA9M6eHwon5tJcLfE3NGJOCDtSfxyoqjGNcjjPcdYzgPTogjeXCIO++8UzyaY8GCBU2eo4Thc3HllVeKh5axGDgs9scYiLdXHUdVjQmjuoZgeJcQtYejK/46tis+25SM3Sn5WHUkm71fjGHINt/nlPueFtBGsbpBUWKRbOAwRmqoqXhv7p3YXe3h6A66+JMXh3h1xVHOxWEMQzYbOI4Fe3AYo0Fl4ey9OX8vjqebs/TiHM620pFhGPUwmUyWVAz24DgInIPDGM178w17b6zsxTnCXhxG9xSUVYmJj5b6UBEcorIhiqIjh6gYI/DOKum9GdmFc2+s5sU5U8BeHEb3ZJvDU4HebvBw1Y7gJxs4dgpRtSVxmmG0SlZhOb7Zekas3zuJc2+scW24Ybj04ryzijuNMwbJv/HVToIxwQaOHTw4ZVU1KKnkdg2Mfpm//iQqa2oxtFMQRnDllFWYO6YL3FycsOXUWWxPZnVjRr9kazD/hmADx4b4eLjCx1266zhMxeiVwvIqfLlJalXdMa6r2sMxDJEBnpYWF++sOqH2cBimw2QVsoHjkHAlFaN3vth0WjSN7RHhy7otNsjFoR7B1KPqCHcaZ/TuwfFlD45DwQYOo2fKq2rw0fqTlpuxM/ecsirdwn0xtVekWH9vNXtxGH2SrUENHIJDVDZGOeDcUZzRI4t3poqLV1SAJy7pH632cAzJHRfKsN+Pu1JFKT7D6I1sNnAcE1YzZvRKTa0JH6yRXoVbL+gMd1eeD9mCAXGBovS+utaE+Wult4xh9EQ2GziOCYeoGL2y/EAGTuSUIMDLDdcOi1d7OIbmb2YvztdbTwvRNIbRE9lcReWYsJoxo1eo6zVBqrtUEcjYjjHdQ9Ezwg+llTX4Zqvs9cUweqCqphZnSyrFOicZOxjswWH0yJ4z+dienCd0Wm4cJQXpGNvh5OSEWy7oJNYXbDiF6ppa3t2MLsgtlsaNi7MTgry106aB4KC6jeF2DYwe+Xj9KbGc2S8a4X6eag/HIbh0QAyCfdxFovGyA5lqD4dh2pV/Qz2otFZlyQaOHauoamu5XQOjj7YMv+xJE+s3j+6s9nAcBk83F1w/XOY6fbSOk40ZfZBdXK7JEnGCDRwbE+IjDzpVSORz8iCjAz7flCyaag5JCELf2AC1h+NQUL4ThQW3Jedhd0q+2sNhmDZ7cLTo6WUDx8ZQaW2Qt5tY53YNjB6E/b7YLNsysPfG/oT7e4qwIPGxWWCRYbRMtkYbbRJs4NgBTjRm9MLPu9OQW1KJ6ABPTO0dofZwHBLFsPxlTzoyC6X7n2G0SrZGNXAINnDsWirOFytGu5hMJnxkTi6eM7ITXF348qAGFBYc1ilYhLU/25isyhgYpq1ksYHj2LCaMaMHNp88i4PphfB0c8a1w+LUHo5Dc9PoThbhv8pqLhlntEs2GziODYeoGD2geAsuGxiLQI3pWTgak3tFINzPAznFlVi6P0Pt4TCM7lSMCfZB27VUXAoiMYzWyCoqt9xIqZKHURc3F2dcY26P8dkmDlMx2iWbk4wdG/bgMFrn221nRM7HoPhA9Ir2V3s4DCDChKQOu+XkWRzJLOJ9wmiOkopq0V6EYA+OgxLmK/UBuEyc0WrX8C/NpeHXD2fvjVaICvDCpKRwizYRw2iNbLP3xtvdRZP96jhEZQdC/dwbxCoZRkusPpIl2gMEervh4n5Rag+HqcecETLZ+PsdqWK2zDBaIlvD+TcEGzh2rKKijqvUeZVhtMTnm6T35spBsaJdAKMdRnUNQedQHxRXVOOHXalqD4dhdJN/Q7CBYweowyrF0ut3XmUYLXAmrxQrD2eJ9evMfZAY7UDNC5X+VFTlRlpFDKMVsjVcIk6wgWOPnezsJDqtKtUqDKMVvt6SArpnju4Wgi5hvmoPh2mGKwfHwsPVGYcyirDjdB7vI0YzZLOBw9RvRMaJxoxWIAG5r7emiPUbOLlYs5Am0cz+sj+V0ieMYbRANoeoGIJEu4jMQk40ZrTB8gOZyCmuEOfmpF7cd0rLXGvWxFmyNx2F5VVqD4dhGkQkOETl4FCXYIKb5zFagdoAELOHxglhOUa7kD5R93BflFfV4qddaWoPh2Ea9KEK9+ccHIcmwnwCKCcEw6idXLzuWI5Yv3oI953SOk5OTsIQJRaaw4oMozaZ5oiEkoKhNXjaZicizB6crEJOMma0oVysJBfHBXurPRymDVw2MAZuLk7Ym1qA/WkFvM8YVamuqUVuSUWD+5vWYAPHzh6cTK6iYlSmttaE77afEevsvdEPIb4emNIrUqx/w14cRmVyiivFJIkkUEJ8tNmclw0cO6G48DjJmFGb9cdzhHKxv6crpvaWN0xGH1xtDlMt3pmK8irZA4hh1CDTHI2gIgWSQtEibODYCSUJK7e4Qrj2GEYtvtkmvTezBsawcrHOuKBbKGICvVBYXm3p/s4wqho4/toMTxFs4NiJEB8P4cqrNQG5JaxmzKhDfmml5cbI4Sn9QdcQEv4jONmYUZNMc8FMhEZVjAk2cOx4YVL6dXCpOKMWP+5KEwJ/vaL80ScmgA+EDrlqSCycnIANx3ORnFui9nAYByXL4sFhA4epn2jMYn+MSiiz/quHSC8Aoz9ig7wxpnuYWP9mG5eMM+qQZb6PRWi0RJxgD44dCbMkGnOpOGN/9qUW4EB6IdxdnHHpgBg+BDpmtlm76PsdqaihuDfD2JlMc0WwVkvECTZwHFXsr6YKKMkBirOpbljt0TB2QJntT+kdgSCNlnUybWNiUriogksvKMemE7mtb0y1vOUF8rdemA5UFMnnGMYaIn/NhaiqyoHsw/JcUxFXVT/dwVBd7K/0LLD7K2D/D0D6LqDGnOzs4gHEDgF6TAMG3gB4B6szPsZmUEnxDztTxTonF+sfTzcX0YCTmm8u2n4Go7uFyv8gwyVjD3BqPXBmK5C5Hyg4A1Q1ytVx9QQC44HIvkD0QKDbJCAskSSTVfk+jP7IKmzFg5N9EHj/QsA3EnjwMNSCDRxVcnDsbODU1gBbPgBWPS9ncg1wAmoqgOT18rHyeWDE34AxDwAevvYdJ2Mz/jyUJUqLowI8626GjK65fFCsMHB+25eBZ8ccgPfeL4CDPwH5rXQcd3IBTDVAdTmQc0Q+9i0Clj0KBMQD/WcDA+cAQQn2/CqMzqisJhXjypYNnAI5mUKAuqFwNnBUabhZYV+vzaK5wPE/5N9hScDQW4GuE4CgToCpFshLBk6sBHZ8AmTsBda9DOz5BrhyPhA/wn5jZWwG5WoQlHtDFX2M/hkUF4ArA4/ikpLv4P3+3rr/cPMGOo8FYodK7wz9zv2jpdeGPDSVJUBJNpBzTHp7kjcAp9YCBaeBNf8F1rwEJM0Axv0TiOyj5ldkNEp2sbyHUeuQIG+3phsUmhvC0nmnImzg2BFSfLRrDk5RJvDJTCDnMODqBUx9Fhh8M+DsUm8jFyC0m3wMnQscXgL8/k85C/z4ImDKM8DIv9tnvIxNIHHJVYezxPrlgzi52BCc2QanFU/ipfK14idcCyc495wODLgO6DoRcG+lv5i7j3yQ4dN9EjDmfqCyFDjyG7DjU+DEKuDgz/LR+3JgyrOqz8QZjZaI+3mKRrBNKJRiovBXt1qTk4ztiOLKowZlVbZWMy7JBT6ZIY0b/xhg7gppwDQwbhpBJ2rixcDfNgB9r5au7KX/ApY+wonIOuaXPemorjWhT4w/ekT4qT0c5nwoywN+uhv4cKLwuphc3PFx9VSMrXgVqdM/ApJmtm7ctAS9ps8VwI0/AnduloYNha/3fw+8ORTY8AZQU83Hjjl3grGGPDhs4NiRYG93uDo7iTzAHLOLz2YVUt/+RcbXybi56Zf2uZo9/IDL3wcmPyP/3vgm8Ot9XHmhU743JxdfNpC1b3TNgR+lsUFeFmLA9XC6eweWxt+HM6YwSxL5eROeCFz1MXDHWiBuuExQphwdmjC1lt/DOAxZSol4Sxo4GsnBYQPHnjvb2ckSprJpHs6yx2RM3d0XuGERENyl/e9B3pzR9wCz3gWcnIHtC4Dlj7GRozOOZxdjd0q+yLu5pL+6symmg1SVAT/PA765UebOhPYEbv4dmPU2EBiHKwZJw3XRjjMwWbP8myqs6HMueQNw9wNObwTeuUCGrhiHJtNSQdWSB8ds4NAEW0XYwFEp0dhmpeLHVwKb35Hrl38AhCed3/sNuBaY+bpcJzc1eXMY3bDYnFw8tnsowjTcM4ZpAUoE/mACsP1jGTK64D7pWUkYadlket8oeLm54ER2CXal5Ft3Vzo7A4NulJ8ZMwSoKAAW3gCs/i9PdhyYTEuIqhkPDumqFaXLdQ5RORYWD44tEo2pBPzHu+Q65dskXmSd9x00RyYaEssfB46usM77MjalttaExUp4yjzLZ3TEidUy1ybrAOATDsz5Hpj0JODa0FD19XDFtD6RDarlrE5wZ+CW34Hhd8i/Vz4rqzOrNSBayqjowfFs+p+lOWaNNSfALwpqwh4cI4n9rfq3zF4P6gxMesq67z3yLjmTo7Ly726RM0tG02w9dRap+WXiBjilV4Taw2Haw/ZPgM8vB8rzZbn3HeuktEMLKNVxP+1OQ0V1jW32tYsbMP1FYMargLMrsO874KtrZAUW41BkmyfoyoS92fCUb4Q8Z1SEDRyjiP1lHwG2vC/XZ7xsfZE+ysm56CWZdEhuakpiJjluRrMo3pvpfSKF8i2jAyiHZvV/gJ/vAWqrgT5XAn/5BfBr3UAd1TUUkf6eKCirwqrD2bYd45Cbgeu/k3o7x/80G2KNBUQZh/XgFGgjwZhgA8coYn/LqJS7GugxvdWZ3nlBrvGrPwO8Q4HMfcCKJ2zzOYxVWjP8uifdonjL6MS4oRDwyufk32P/AVzxIeB27maGIol8gEwi/2mXuUTXlnQdD8z5AfAIkMnHn10OVBTb/nMZ1amorkFeaVXLScYaKREn2MBRK0RlzRyc5I3A0WWAsxsw1XxxtBU0k5xlTmLe/C5wZKltP4/pECsOZqKoohoxgV4Y3pl7i2keSsz89QFggzmhf+oLwIRH2tUbSqmSE8e+XN6AbEr8cOCmnwGvICB1G/D1dezVdQCyzJNzd1dnBHi5tSLyxx4cx1UztmaIas1/5HLg9UBIV9icHlOA4X+T6z/+XbaDYDRZPXXpgGghT8Bo3HND6uHb5svETKpaHHlnu9+md7Q/uob5oKK6Fsv2Z8IuRPWXUhQkSXFytczPIx0uxvgaOP4eLagYpzmWgfP222+jc+fO8PT0xODBg7F27dpWt1+9erXYjrbv0qUL3n333Qb/v2DBArFjGz/Ky8t1pGZcKRqWnTdntsk4OCX9XXA/7Mbkp2T3YdLlIKVjRjPkl1Zi9RGZh3HZQPUvMsw5+PNZYMt7cv2yd4HBf+nQLqNrIPUaI37cbYcwlULMYODarwEXD+Dwr8CSB7mE3MBkmj04LYr8OVKIauHChZg3bx4eeeQR7Ny5E2PGjMH06dNx+nTzipgnT57ERRddJLaj7f/1r3/hnnvuwaJFixps5+/vj/T09AYPMoi0DjUmowZl9RuWnRfUHI/of419OwBTPg4JgNGMc/eXwDEuHdcK1F2aWjMkRvqhO7dm0DbrXgHWviTXL/6f/B2fB0qYav2xHEuli13oPAa4aoG8HpAo6CZzGJsxbh8q/xZ0tQrMIaqAWOMbOC+//DJuvfVWzJ07F0lJSXj11VcRFxeHd95p/gdA3pr4+HixHW1Pr7vlllvw0kvmi0C92UpkZGSDhx6gcVODMquEqXKPA0d+NwuA2dF7oxA3rE4Xg5RWK4rsPwamCUqSqZJ0ymiUnZ8DK56U65OfltpV50mnUB/0jwtETa0JS/aaxdbsBeluKXpZ1MPuMF2bGKORaSkR17bIn80NnMrKSmzfvh1Tpkxp8Dz9vWHDhmZfs3HjxibbT506Fdu2bUNVVV1st7i4GAkJCYiNjcWMGTOEt6clKioqUFhY2OBhiFLxrR/KZffJ9sm9aY4JjwKB8UBBiixvZVSFzqlNJ3PF+sx+6l9gmFZE/H6+V67T5GS0ed0KXGr24vy4y0aif60x8u/AIAqxmYBFtwLZh+0/BsamZOpE5M/mBk5OTg5qamoQEdFQw4H+zsjIaPY19Hxz21dXV4v3IxITE0Uezk8//YSvvvpKhKZGjx6No0ePNvueL7zwAgICAiwP8iCpSWSAPDEyCs7DwKGSzJ1fyPVht0M1SG+H9HEIckvnNH8MGPtApeGUszooPhBxwR3oKs3YHrrpL5xTp3Mz8XGrvv2MflGgvPIdp/NxOtfOInyUdEqhtk5jgMpi+T25fNxx+lDlp8glGTcqi/zZLcm4caY1NYRrNvu6le3rPz9ixAjccMMN6N+/v8jV+eabb9CjRw+88QblhDTl4YcfRkFBgeWRkmI+CCoR6e8llunn48HZ+40U3KNGmrbSvWkrPaYC3acAtVXAbw9xgqGKkJItwY01NUpxFvDFlfK3GzcCuPStdpWCt1Vri4T/iJ/32DHZWIFubFd+BPhGAjmHgV/u42uCgcgwT8yViXoDyJOvkfwbmxs4oaGhcHFxaeKtycrKauKlUaBcmua2d3V1RUhISLOvcXZ2xtChQ1v04Hh4eIik5PoPNYmyhgeHpNwJittTQzy1mfZvwMUdOP4HcPg3tUfjkNBsnZot0uz9on7qu4eZRlD59Dd/AfJPy3Yq13zZJhG/jqDkX/2wM9W6Hcbbim84cNXHgJOLnIyJZqGM3jGZTEg337eiAuREvVkDJ1DdKImCTe+M7u7uotx7+fLlDZ6nv0eNGtXsa0aOHNlk+2XLlmHIkCFwc3Nrcafv2rULUVH6uKhHmA0c5URpN1kHgfRdsjS832xoAsoBovg7QZoe3MbB7iiz9ZFdQ5pPAGTUZdljwOkNgIc/cN03gE/zEzZrMLV3JNxdnHE0qxiHMlRK/k8YJZuDEuTZzdyvzjgYq0HioaWVstcZtQZpMUQV4AAGDnH//ffjww8/xEcffYSDBw/ivvvuEyXid9xxhyV8dOONN1q2p+eTk5PF62h7et38+fPx4IMPWrZ56qmnsHTpUpw4cUIYNlSlRUvlPbWO4sHpcJLx7q/ksvtUwEe6ojXBmAdl7DU/2SxaxtiTnzk8pV32fgdsfqdO6yash00/jhRmxyeGifUf7dG6oSVG3Q30mCYTTxfdxt3HdU6GeVJO55eXu4tjh6iI2bNni5Lvp59+GgMGDMCaNWuwZMkSUQFFkH5NfU0cEgSk/1+1apXY/plnnsHrr7+OK664wrJNfn4+/vrXv4oycqq4Sk1NFe87bNgw6AHF8iUPTrvdx7U1wJ5v5PqAa6EpKOF4/L/q9HnK8tUekcNwJLNIzNRJY2lab314Mh2GjH3Aj3fJ9TEPAIkX2+VjFdE/Mnxra1UIUxGUX0R6WdS/Lms/8Ocz6oyDsaqBo0zSWw5RxUMLuNrjQ+68807xaA6qhmrMuHHjsGPHjhbf75VXXhEPvaKU15GSMTUtC/Zxb/uLT6yUOgPU/4USe7VG/+uADW/K5ML1r9a5qBm7aN+M6xGOAG/1qxcYM+WFwMIbgOoyWQww3n6q3xMSw+Hr4YrU/DLsTMnH4IQgdQ4L5eNc+ibw1TXy2kCeZxIGZIyVYOyIISqmKdSkLNRXGjXpBWXt20X7vpfLPldINWGt4eJaZ9RQ2XiBClocDgZ5AZXqqZn92XujGcg7SxVEeSflBf+K+YBzM259G+Hp5oLJvWQxh9JZXjV6Tq/Tx1l8hzT8GN2R3poHh4Rey/MdK0TFWFELh6owDv0q13vN0u6upYtZ/EiguhxY9bzaozE8e84U4PTZUnjVu6ExGoBy5fZ9JyuJyLjxtn9X94v7SoOXVI1VC1MpTH0eCOoku03/8ZS6Y2E6REZhWQOpk2ZbNHgGAJ7qViorsIGjEsoJktGeRONTa6WFTPFsqlDQKhR3J+l5YteXQPYRtUdkaBTvzaReEfB2t0vUmTkXOceAX82FEeMfBuKHq7LPxvQIhZ+Hq7jO7DidB9Vz9KhTuqLCnrxR3fEwHfbgRAa0IvIXoI38G4INHD1p4Rz4SS6TZtjV1d3hPlU9LwJMtcDqF9UejWGhWfkv5vJwFvfTCNVUMXQLUFUiFX3V6BNnxsPVBZN7S6/eL2qHqYgu44CBc+T6z/ewnIRuc3C8mv5nwWlNaeAQbOCoHKJqsxYOVU8d+kWuJ10CXXDhP+Vy3yIg65DaozEk25LzkFlYAT9PV4ztoSHJAEeGwi/pu2UhwOXvqz4ZodYNSpiKmnCqzpRnAN8IIOdIXSd1RhdkFLaSg2Px4Ggj/4ZgA0flUvE2e3BObwRKsgHPQKDzWOiCqP5A4gyZWMheHJugdIym3BuarTMqc2o9sPEtuU5tGDTQUfmCbmHCAM4qqsC2U2fVHo40/C76r1xf94oULmU0T1llDfJLq1pp03BGUxVUBBs4KqFYwG2uolLaH1ACrwaamLXbi7N/MV/IbBCe+m1feoNkUkZFqKnkD3+TBj2FYeykd9OWqk1SNiZ+NRvEqtPrUqDnxbLh6JL/415VOvLe+Li7iLwurbdpINjAUQnFAqbwQps4tkIutah90xqRfYGkmfKiv+rfao/GUFDSqAhPebjigu4cnlKdFU9IFW+awVLFkIaoC1NlaCNMRUx7AXD1lMUT+83yF4xmSTdPxqnVULPNsjWmgUOwgaOygVNcUY2icun2axFqzpd9SJabdh0P3THO7MU58AP3o7Eiymycw1Ma4MQqWRlEkHKvRspkFUZ3CxXy+jnFFdhyUgNhKiIooS4Be+mj0gPG6FPFuLpSCtASbOAwVM7r7+natjyco8vrKpMofq03IvtIlzSxhpMKrRae2psh1i/i8JS6kGid0ophyK2anIS4uThjmiVMpWJvqsaMvldq4xSlyfYujPZLxP2bqaCi40deehcPwEf2QNMC7MFREaXdfHpbDZzuk6Fbxv5fnRcn97jao9E9O1PyREycpPhJ64RRkWWPyPwDulEr+k8a5GJzmOr3fRmorqmFJnDzBKaZQ9eUnJ1zVO0RMS2Q2aYKqhjAWTtmhXZG4shqxq2J/VVXACdX6zP/pnEuDo2fdHHWv6b2aHTPr3uk92ZSUjhXT6nJidXAjk/l+qVvSzE7jTKyawiCvClMVYnNWglTEdRtnK4NtVXAUvv16mI6KvLXXAWV9vJvCDZwtC72l7weqCoF/KKAiD7QNUq8ndSNCzXkJtdx9RSHp1Skqgz4ZV5daKrTaGgZEabqE6kd0T8FSlid+gLg7AocXSqNRka7In/+rXhwNFRBRbCBo4Gu4q2GqCh5keg6UV4I9EzCSCB+lJypKVohTLuhztB0zlB4amwP7cS7HQ7KJzt7AvCNBCY9AT1wcV+py/P7vnTthKmI0G7AkFvk+rJHyYpXe0RMezw4VD1IBHaClmADRxMenFa0cE6ulUu9iPudizFmL862j4FSDbnJdcRv5uqpiUnhomM0owKZB4D1r8p1Eq2jBoM6YESXYIT4uCOvtAobjudCc9WWHv5Axh5gz0K1R8PUo7K6VlTgtZiDk3dKLikPTUOwgaPldg3lBUD6LrneeQwMQbdJMh+H+vRseV/t0egOk4nCU1w9pSrkXfj5XilSR2J1QudJH7i6OGOqOUylnEeawScEGPOAXP/zGaCyVO0RMWayiuQ9yt3FGcE+7mhCXnJd6b+GYANHA1VUSnZ6E5I3yKTckG6akHy3ChRmu+A+ub75Xda+aCe7UvKRml8m1ETHcXhKHbZ/BJzZArj7Su+NzkLH080GzvIDGhL9Uxh+h0xULUwFNr2t9miYJk02mxH5o0IYOl4Ee3CYxh4ccheXV9W0HJ6ijsRGotcsILgLUJZXV4HCtKv31MSkCA5PqUFhOrDiKbk+8QlZFqszRnQJMYv+VWKrFnpTNS4bp/2q9KkqyVF7RAzqa+C01IPKBLh5a0oDh2APjoqQ0J+XOYei2Uqqk2uMFZ5SoO7Ko+6R65veAWqq1R6RbsJTJLVPcPWUSiz9F1BRCMQMAYbeCj1C1VSkfq1o4miOPlcAUQOAymJp5DCa8uA0Ie+kXAYmaM6byQaOipCrr67pZiMDhxJwM/ca04ND9L8G8A4FCk4DB39UezS6YPeZAkt46sKe2popOQQ04aCeSU7OwIxXpKGuU5QwFRk4JDugKUgobsJjcn3LBywpoQHSW2vToNEEY4INHM0kGjeqpDq1Ti7DkgDfcBgONy9g6Fy5vuFN7ibcjvDUBA5P2Z+aKtn1WtG8ieoHPUPNWUlmgERGd53Jh+boNhGIHwnUVHALBw2QUVjWigdHmwnGBBs4KhMd2EK7BkowJjpdAMNCBg71LknbAZzeqPZodBCekgbOxX3l7JuxI5vfkw1vyes4Qf9qux6uLpiQGK7dMBWFOhQvDuXpnTWHQRhVSM0vb3C/agB7cJiWUE4YCj00IGWzXMaPMO7O8w2ToSrFi8O0yIH0QpzJK4OnmzPG9TCgR0/LFGUAq8z9kiY9qc+Gt62EqUgVmwxozUHK0F0nyHL81S+qPRqHJs18f4ppzsBRRP44RMU0JibQs8EJJCD9BxK7IuKGG3unjTR3YT68hJtwtsKy/ZliSaXhXu76zf3QJcseAyqLZGLxgOthFMb1DBMGc8rZMuxPK4QmmfCoXJLwX/ZhtUfjkFRU1yC7qOLcHhxKMtYYHKJSGeWEaWDgUMiGZi1+0UBALAxNWA+g+1RZZsjtG1pk2QFp4EzpxeEpu3JqPbD3G4qZSM0bDXVKPl+83V0tWkpL92swTEXEDAYSZ0g9sJXPqT0ah66g8nRzFs1aG0BSHyRIS3AODtNiiCqvrM5NrISn4oZpruzOJoy6q64JJ7dvaELK2VIcTC+Ei7OTJW+CsQMkX6AkFg++CYgZZLjdPr1PlDZVjesznnKenIADPwKZ+9UejcORap58072qicifkmDsEw64+0BrGGc6olOizWrGJZU1KCw368GkbHGM8JQClcFH9gOqy4Ct89UejWa9N8M6BSOoOZl0xjZsmw9k7Zc5NxMfN+RenpAUDjcXJxzLKsaxrCJokoheQK9L6xqcMnYlzZxgHNNqgrH2wlMEGzgqQ/kUSm8PEaYiL44lwdhBDByaFYy6W65Tf6qqVrqrOyDLzOGDKb2lOBtjB8iTuPJ5uU7VPN7Bhtzt/p5uuKBbqFj/zSwiqUnGmj1p+xcD2UfUHo1DkaZ4cMyTcb0kGBNs4GiAmHphKuQclXFNVy/p1XAUel8G+McAJVnA3m/VHo1mOFtSJ6evqM8ydoCqdsrzgfDeMjxlYHQRporsIxubUq7eWvbiqGLgBOqrRJxgA0cDRCuVVCT2p3hvKN7v0iihy8jQdx3217omnFosW1WBPw5mgoRme0f7IzbIW+3hOAY0ydj6oVyf+pyuFYvbAhnOlN9FUgSnczXcwXuc2YtDE6Dc42qPxgFzcDyb/qeSg6PBCiqCDRytaeFQl2IidigcjkE3Ss9V5j4geb3ao9EEXD2lxk5/VFYx9pgGdB0Po0N5XSO6BFs0cTRL9ECg+xRZUbXuZbVH4zCktaaBo/Sh4hwcpiWUE0ckc6XulE/GDnG8HUZ5DorwHzXhdHDKKmuw9mi2WOf8Gztx/E/gyO+Asysw5Vk4CtN6R2o/TEWM/Ydc7v66znvA2AyTyWRJMm4SoqL2JcoxCO6qyaPAHhwNoJw4OXn5QNYB85PGK0ltE8PvqBP+c/AL2Jqj2SivqkVcsBcSI/3UHo5jlIUvNbdhGHobENodjsLU3pEi139XSn7TvnhaIm4o0GW89LBxp3Gbk19ahbKqmub7UOWfBkw10uvuJ/O4tAYbOBoycHzyDsoThjQF/KPhkIQnygsYuaGposqBUdSLSdyvif4EY312fionGJ6BwDizp8BBCPf3xOD4IO32pqqPcmx2fg4UpKo9GofIvwn19YCnW6NcNCUPKriLZgUwtTkqB0NJ3oopPWR+YqBjCPy1xIi/yeWOz4CKYjgi1TW1+OOQYuBw9ZTNITXWP81KuRc+bNiy8NaYZu5NpVlVY4WEUUDCaKC2Ctj0ttqjcZD8G8+m/3nWbOCEdIFWYQNHA4T6eMDdxRl9nU/WGTiOTLfJMqZbUQDs/gqOyNZTecI9TBpJgxOM0dxR05CAXGkOENIdGHorHBEKUxFbTp5FXkklNM3oeXK5fYGU1WDsXyKee1zT+TcEGzgawNnZCVGBnujrdEI+4egGDrk7h98u1ze/B9TWwtFYdkDOoicmhsPVhX+mNuXsCSlNoJSFO5I8Qz3igr2RFOUvZAn+OJQFTdN9stQoqixm9XMbkmbuQ9WsgWPx4LCBw5yDzn5ANydzPNnRDRxiwHWAhz+Qe1RWtjhY5YIl/8Y8q2ZsyIongZpKmftFZcgOjBIOVdSzNQuF8EffK9fJOK3ScGK0QfpQNTsxINiDw5yLIR4pcHEyocgjAvDjnAt4+AEDb5A7Z7NjlYyT4BpdWLzcXDCmu5TRZ2xEylbZxJGaOZL3xpFz3+qFqaiCj2QKNE2fy4GAOKAkWzbqZeyngVNdKauolCRjjcK+b43Qx0m6+1I8eqg9FO0glI2dgGMrHKr/jOK9GdsjtGnlAmM9SC17+eN1HsOI3g6/d5Oi/BAb5CXkCcjI0TQUShx5l1zf8AZQq3GDzEgGTn6yrHR18wH8tOtlZgNHI3SqlDfwg87ajWfaneDOQM/pct2BSsZZvdhOkKDf6Q2Aqycw/l/2+lRNQ3IEJEtQ39DWNIPmAF7BUlFXeOIYa1FZXYusoorm2zTULxHXsNeTDRyNEF4sS8R3VGqzaZnqwn/kgi7Lh9FJOVuKg+mFojfQhMRwtYdjbFE/yr1RzrGAWLVHpBkU1WySKSC5Ak3j7lPXw279q9zDzopkFpYLJ6eHq7Oo5tRbiTjBBo4WqCyFV5Hsyrq+xEEF/lqi81ggvBdQVSKFvRzEezOsU7DoEcTYiF1fANmHAK8g4IL7eDfXY0hCkLihkUzBFnMne01DBg6p6abvBk6sUns0hkswjgn0aio0qoMScYINHC2QdRBOplpkm/xxqsIXheVVao9IO9APSykZ3/Ke4ePsSvUK956yIZWlwKoX5PrY/wO8Am35abqDZAlInkA3YSqfENmoV/HiMLbXwDmr/RJxgg0cLZC5VyyOOXVqcGIxZvpeLWfalLV/+DfD7pazJZXYap4xT2b1YttB6rdF6UBgPDB0rg0/SL8o8gTLD2QK2QLNM/LvgJOL9OBkyOspc36k5sn7UFTjHlRErvZLxAk2cLRAxj6xSPPsJpZnzrKB0wB3b2DwTXJdEWQzIH8czBQia72j/REb5K32cIxJSQ6wzjzLn/AY4Oqh9og0CckTkEwBhSn2pxVC8wQlAL0ukeubHEtWwlacMRs4Ta5FVeVAQYpcZw8Oc04y94tFgX9P84lVyjutMdTdmWZop9ZaDEKjwdVTdmDNf4HKIiCyH9DnSnt8oi4heYJxPcL0IfqnoJSM7/kGKNLJmDXMmXx5H4oLbhSiooo1mAB3P8BHniNahT04akPuX7OBUxXaq4HlzNQjIKZuhmZA4T8SVVtr1h3h/BsbcfZknaz/5Kc12wFZKyjnoWJ4a57YIUDccNmEc+uHao9G96ScbcGDk31YLkO7a7pEnOBfuNpQXgk1lXRxh1d0kngqhT04zTPiTrnc860MNRgIElUjcTWaLSVG+qk9HGPy5zPy5td1AtB1vNqj0TwkU0ByBYcyipCcWwJdXSPIkKVkcqZD1NSaLLmgTTw4OUflMkxGHLQMGzhqk2kOt4T1REyIv1hlD04LxA4FogcBNRXAto9hJJRqlclJkU1LMpnzJ3U7sG+RVMae9BTv0TYQ6O2OEV2CG5yfmidxhkweLzsL7Pla7dHolozCclTXmuDm4oRwv0ZJxjmKB0f7qvts4KiNkk8S0Ud081XE3phmoBv/iL/JdXJBUz8UA0BiaiSqRnB4ylYtGZ6Q6/1mA1H9bPEphsSiamzubq95XFyB4X+rSzau1bhQoUY5Y74HUYk4efEakGNum8MGDtPWEnEycJR+H4Xl1SgoYy2cZuk1C/CNBIozDCPNvvVUnhBVC/J2EyJrjJWhXmaUnO7iDkx4hHdvO1DkCrYl5yGnWMr2ax5q0uvhL2/EdOyZdqNEEeIa59+QwcghKqbdHpzIPvDxcEWIWb2WK6lawNUdGHprnZ6JHjQ6zoEyO56YFCFE1hgrQsKQiveGFG8pfMG0GZrB94sNED+zFXpJNvb0rxP+2/im2qPRJSnmPFBqvNqAwlSgqhRwdgWCtN9WiK+malJRbC65kx6c+icU5+G0wuCbARcPIG0HcGYr9AyJqCn5DVNY3M/67P4ayNoPeAYAYx6wwQcYH+W81E01FUHq507OwMnVLPx3Xho4Xs3n35DAH3Vz1zhs4KiJcrL4hAM+oWI1lvNwzo1vGND3KkOIeh1ILxRiap5uzhjTXduaErqjqgxY+ZxcJ+PGWybMMh1TNV53LAfFFdX62H3kqet1qSGuEWqQYs7BUfJCLVjCU9pPMCbYwFETRU8gPNHyFHtw2sgIc5dxysMpSIVeUbw3Y7uHwcvdRe3hGIvN70mXun8sMMzcz4xpN93DfdEpxBuV1bVYfVhqNekCFv6zvgfHooHDBo6Ft99+G507d4anpycGDx6MtWvXtrpzV69eLbaj7bt06YJ3320qz79o0SL06tULHh4eYrl48WLoDupmTITVGThKUhfn4JyDyL5AwgWAqQbY+gF0r15sniUzVqL0LLDuZblOicVuzfTTYdoEyRZM7a2zaiqChf86RFVNLdILWkgytlRQaV8Dxy4enIULF2LevHl45JFHsHPnTowZMwbTp0/H6dOnm93+5MmTuOiii8R2tP2//vUv3HPPPcKgUdi4cSNmz56NOXPmYPfu3WJ59dVXY/PmzdAVWYqBU3eyKBazoiLJtIJSMr59gS5FvcgNfDC9EFSFqXRvZqzEuleA8gIgvLcsDWfOC0W+4M9DWcKToxsU4T/SzaIeSsw5ySgoFz3x3F2dEerr0YKB0x16wOYGzssvv4xbb70Vc+fORVJSEl599VXExcXhnXeaj4uStyY+Pl5sR9vT62655Ra89NJLlm3o/yZPnoyHH34YiYmJYjlx4kTxvO49OMF1HhxddPFVk57TgcAEoCwP2PsN9Oq9GdY5GEHm6jnGChSckeEpYtKTgDOH/s6XgXFB4mZXVF6NTSdyoRtI+I9ClKU5ZqFHps0VVIFecK6vgUNe0RJziJJDVEBlZSW2b9+OKVOmNNiB9PeGDRua3bnknWm8/dSpU7Ft2zZUVVW1uk1L71lRUYHCwsIGD9WpLJFtGogw2aKBULRwSiprkFfKWjitQjcuKv0lNr2ru5JxpYmhIqbGWImVL0i1awphdp/Mu9UK0I1O0cRZrqdqKhL+GzZXrm/W3zVCDc4oPahaSjAmg9HDF3B0D05OTg5qamoQESF/GAr0d0ZG87Fcer657aurq8X7tbZNS+/5wgsvICAgwPIgD5LqiJPFBHiHAj4hDbr4hvtJtyDn4bRR1MvNB8g+KEtCdcLZkkpsPXVWrCs3DsYKZB4Adn8p1yc/pflmgHosFycDp5ZiGHph0F8AV08gYw9weqPao9E8Z1rSwFEiDjoJT9mtiqpxbx0KvbTWb6e57Rs/3573pBBWQUGB5ZGSkgIthqcUOA+nHXgFAgOu01056B8HM0Wcu1eUf9NSTOY8duzTgKkWSLpEJpkyVmNk1xD4uLuIPkV7Uwv0s2dJHqDf1XVeHKZVUlpSMc46KJfhvaAXbGrghIaGwsXFpYlnJSsrq4kHRiEyMrLZ7V1dXRESEtLqNi29J1Va+fv7N3hox8Bpmo1ePw+HaQPDzSXjR5YCucd1Vj3F3hurkbwBOPIb4OQCTDSrFzNWg7zLF/YM1181Vf1rxMFfgHwNTHD16MHJ2i+XEWzgCNzd3UW59/LlyxvsJ/p71KhRze7ckSNHNtl+2bJlGDJkCNzc3FrdpqX31LYGTl3+jYJiOSvJXsw5CO0GdKecLBOw5X3N766yyhqsPSqT9Tj/xpoNNR+X64P/Is8JxuooBrluuosrRPQGOo0xy0p8qPZo9KmBk3lALtmDU8f999+PDz/8EB999BEOHjyI++67T5SI33HHHZbw0Y03mvuGAOL55ORk8Tranl43f/58PPjgg5Zt7r33XmHQvPjiizh06JBYrlixQpSj6wbF3deMB4fF/s5jhrbzC6BcA0nkrbDmaDbKq2rFcU6K8lN7OMbg0C+ybYebNzDuIbVHY1jIg+Pq7ISjWcU4mVMCXcpK7PhEl7IS9qCiukaEIIkGofPiLFmJBqdm0yocNgeH9GqofPvpp5/GgAEDsGbNGixZsgQJCQni/9PT0xto4pAgIP3/qlWrxPbPPPMMXn/9dVxxxRWWbchT8/XXX+Pjjz9Gv379sGDBAqG3M3z4cOhGQj7vlFxv5mRRTixFLptpA10nSPGpyiJg5+ea3mV1vaciW81FY9pITTWw4im5PvLvgB9XpdmKAC83jOgiUwWW6y1M1WNaPVmJb9UejSZJzy8XzlAvNxdL42dBpjk8FdwZcNdPzqCrPT7kzjvvFI/mIOOkMePGjcOOHTtafc8rr7xSPHSJUkHlFQz4hLXqwTlXQjZjhvYRNdj79X5gy3tyXYP6J9U1tfjjEOffWJVdnwO5R+XvadQ91n1vptkwFfWlIkP9r2O76k9WYtkjMtmYOo7ztbUBSlpETJBXw/uODhOMCe5FpWb+DYWnmvmBRQV4CXXbiupaZBdX2H98eqX/NYBnoPSOUcKxBtmWnIf80ioEebthSEKQ2sPRPxRqIN0bYtw/AE8NFBAYnElJMg9n++k8ZBfp7PokZCW8gawDwKnWWwY5Ism50sBJaFzZaUkw7g09wQaOGtBssxU9AZLIJiOH4DBVO3D3kQmmxOZ3NB2empgUAVcX/vmdN3ScizNk9+ght5z/+zHnJDrQC/1iA0Qog+QOdCcr0f/aOnFQpgGnzWkR8SHeuk8wJvgKqwa5x+QypGXBpHizBa1Y1EwbGXqbLBM+uQbI2Kep3UbhRqW8VhFNY863oaa5PcuExwDXRn1zGLuI/ukOCl8Th5fU5UIyguTckqYenNraOlkTNnCYths4LZeyJpgt6FNs4LSPwDggaaZc3/S2pk7Gg+lFIq/K080ZY7o3zb1i2sna/wEVhUBEX6CPTvPxdMpkc3uRtcdyUFJRDV1BqQFUlCBkJT5QezTaDFGF+NQ9mXcSqCoFXDyA4C7QE+zBsTfk11XE6Fo1cOQJdtpsUTPtYNTdcrnnG6AwTTO7TvHejO0eBi937SVA6wrq46ZoHk2mhpp8KbMnPSJ8xSSMOouvOWJuwKgnhisl458BFcVqj0YzHubTzYWoKF9JMQypt5eO4KuCvSnKACqLZRglqNM5PTjJXCrefkiiP34UUFulKWl2S3l4by5jPm9WPg/UVAKdxwJdJ57/+zHtgipslDCVosqtK7pNkt6IigJgz9dqj0YT5BRXorSyRtS9NBD5U/JvdJZgTLCBo1Z4KigBcK2nM9AIzsGxkhdn28eaEP6jZPED6YWiOm5iopS7ZzoI5VbtNt+UJnFDTbVQDHVKNK6qqYWuII/fMHMuzub3ZJ6Jg3P6rIwWRAd4wcO1noeZmpQSbOAwba6gaiU8Vd+DQ12ni8qreMd2RNSLkrgpR2PHp6rvPyUZc1jnYATVF9Bi2s8fJOpnAnpfDsQM4j2oEoPig4QYXGF5NbacPKu/40BNet39gJwjwImVcHSSzfk3yuS6iYET1R96gz049qYN+TeEn6ebRUmSK6k6OEMbdVddl/EadY3EuuopDk+dFydWAUeXAc6uwIRHrXNwmA7h4uxk0cTRZTUVaSYNvL7Oi+PgJFsSjOsZOKT6TPluRGRf6A02cDRYQaWgJHopiV9MO+l3jVSKLjwD7F+s2u7LK6m0zHAnc3l4x6EwwjKzUTPkViBERyq6BkU5n5ftzxBJqrqDlI2pv9LRpXWTTwfldHMJxulm7w21uPDSnzApGzgaNnAULYJTXEnVMdw86+LsG16XFWwq8MehLNSagF5R/g0b2DHtY89CIGMv4OHPDTU1wgXdQ0XforSCcuxPUz/Xrd2Qkdx9ilxXqvLg6Bo4PoYITxFs4NgTCpMowlJtMXAspeLswekwQ2+V0ux0Y6TwhgrQ7Fbp4cOcR0uGP5+R62MeAHxkw0dGXTzdXDCuR1iD81x3jLhDLnd+oYmCBE2FqNJ3y2VUP+gRNnDsSV4yUFstb7h+Uefc3FIqzgZOx/EOBgbOkesb3oC9KauswZqjUieE82/OAxJtLEwFAuKA4eYbEqMJFMNdl+XiRJfxQGhPoLII2PUFHJHiimrkllS2HKKKGgA9wgaOKuGprm0SJqszcFjs77wYeSfg5Awc/8Pu7RvWHs1GeVWt0JVIivKz62cbhuIsYN0rcn3iEzL0yGiGCYnhIuH4UEaRPr3NJPyieHFIN6u2Bo5GsvkeQ02A/T3d5JOVJbLCjIhkDw5jxfyb+iGq9MJyVFQ73o/OapCgYq9LVfHiKLNa8t6QOBrTAVb9W4pjRg8E+lzBu1BjBHq7Y1in4AbVgrosSPAMlCkER36Ho3FaKRGv36IhkzqImwDfSMBPn+F19uDYkzZq4ChQmbiPu4vIjU05W2bbsTmK8N/eb2Wo0A5U19Raui1z9VQHyT4MbF8g16c8yy0ZNIruw1Tu3sCQm+tkJRyMZHMFVYMmmzrPvyHYwNGgBo4CzfgVi1pRmWQ6SMxgoMuFgKlGVlTZgW3JecgrrUKgtxuGdtJfiaUmWP6EPGY9LwY6XaD2aJgWUAz4bafOIre4Qp/7aehtsoXOqbV1uScOQnKrCcb6rKAi2MCxJ0oFVVDnNr/EUiqeo8PYttYY82Bdgz3qCWan3lMTEyPg6sI/tXZzcg1w5Dd505lM6sWMVokN8kbvaH8hh0CyCLokIAboPUuua6iHnT04bZ5AN+ginrpD1wnGBF917UV1BVBwRq4Ht8PACWWxP6tBHoC44UBNBbDxTdgSEj2zqBdzefh5ivrdAoR2t+4BYqyOUiWoGPa6ZMSddaHsIh1/j/P14FQUA9kH65oX6xQ2cOxFfopM2HLzkeq6bUQRXeJKKitASb6KF2frR0Cp7frnHEwvwpm8Mni6OWNs97Yfb8bM3m+ki5x6BV34T94tOkAx5NcdyxbyCLqEbuaxw2Sn+m0fwRGorK5FWn5ZwxyctJ2AqRbwjwX89Ntehg0ce5F3sq6ipx3VNJZScW7XYB26T5Ylj1UlNnVDK96bMd3D4OVerzMvc25o9rjiSbk+5n7AJ5T3mg5IjPQTcggki6BoP+mSEX+Ty23zgapyGJ2UvFIRWvR2d0GYn4d8MnWbXMYOhp5hA8denD3Z7vBUfQMn5WwpaugsZKzgxXlArpOBYyPlUsVNP4V7T7Uf0rwpSpeTASVkwGgeKoowRJgq6RLpuSjJBvYtgtE5mS3zbzqH+tRJWZwxGzgx+g1PEWzg2D3BuFO7XhYV4AV3V2dU1ZiQmsel4lYhaSYQ0h0oL5CzNCtDxuiB9EI4OwETzd2WmXb8ThStoinPsaifTsNUfxzKFDIJusTFFRh2W13JuB6biLaDkzl1Bo6F1O26z78h2MBRI0TVDkghtLM5s/1ETrEtRuZ4OLvI0Aex8S2gyrqG43KzFsjQTsEI9nG36nsbnmWPySTwzuOAxIvVHg3TToYkBAk13PzSKmw9laff/TfoRtlSJ3MvcGodjMwJs4HTRTFwCtOkB5WqF3VcIk6wgaPxEFV9y1qxtBkr0PcqIDBeuqEVITkrUVc9pd/kPNXKwg/+JNtqTPt3u3LVGG1AcgiK11Ix9HXbw67/tQ4h/HfSPHHuHObTMDwV3gtwr+fV0SFs4NgDcnF2QANHQTnxTphjpYwVcHEDLri/LueDulVbgbySSmw5KauzOP+mHdRUA78/LNeH3ApE9LLK8WDUE/0jQ5/kEnSL0tT18BLg7AkYP0Tla6gEY4INHHtQnAlUl8mZKXVDbieK65A9OFZmwPXSi0PHx0q5OCRyRrngSVH+iKsve860zo5PgMx9sh/Q+H/x3tIxJItA8ggkk0ByCbolrAfQbbKU99j8PoxISUU1Mgul8rSSCmGUBGOCDRx7hqcCYgHX9udkdDF7cNjAsTJ0LMY9JNfXvSrLk8+TZfvN4Smunmo7ZXnAn8/K9fGPyPAAo1tIFoHkEXTdfLNxyfjOz21WcakmJ83eG8oVDPB2k4K0ioETPwJ6hw0cuyYYtz88Vd91mJpfpl8BLS13EQ7uApTmAFvOb5ZGx0bR/2D14nZ2Cy87C4QlSdViRvcoBr6uy8WJrhOA0J5AZRGw8zMYvoIqbadM8icx2jb2TNQybODY04PTzgoqBbKuqWEjcSqX83CsXhI6zqyUS004z2OWtvZothA5iwn0Qq8of+uN0chQU0PFsJz2vDwejO6hRGOSSSC5BJJN0C2U6K54cTa9K3PFjGzgJK+Xy/iRhkjyZwPHHigJxh2ooFLgSiob0vdKILSHDJWch7rxMnPVCHlvLIJZTOv9pn59QErC95olZ8uMIaBJ2ZBOMtS44qDOvTj9rwG8Q4GC08CBH2BsA2ejXCaMghFgA0cHIar6J+CJbNbCsYkujtLvaMObQFl+u9+CRM3+MF/IFTVX5hzs+gI4swVw9wWmvcC7y2AYJkzl5lVXUbX+VUMJ/52or4FTWwOkbK7z4BgANnA0roGj0DXMt8EJyViZXpdJ3YeKgjol3XawLTkPeaVVIpQ4tFMQH55zQY1Olz8u18m49I/mfWYwFEN/y6mzQj5B1wy9VQr/ZewFTqyEETCZTDhpnjALKZLM/UBFoWxwG9kXRoANHFtTUSQTWM8jB4fgEJWNcXaWFTzEpreBwvR2vVyZpU5MjBBiZ8w5+OMpmVhMRqUyO2YMRXyIt2jAST30/jyUBV1DlX2kbkysfw1GIK+0CoXlMqeoE5WInzaHp+KGSa+2AeArsa3JS5ZLryDAM6DDb8Ol4naAWgPEDQeqSoFVL7RrJrTUXB6uiJwxrZCyFdj+iVy/+H9SdJExJIqat+7LxQlq/ErtC06sAtJ2wSgKxjGBXvB0cwGSN8j/SDBGeIpgA8fWFKTIZWDCeb2NsLAB0ePlrN7dvVqFEoMnPy3XqSQ0+3CbXkaVIlTCT+Jm43pI/Q+mBSjO/yspSJuA/tcZJpmRaT0PZ82RHJRX6VziIigB6HN5XcWlzjlRr4u4SPg/tVb+R8IFMAps4Nia/NNyGdh+BeP6kIVNlnZ9y5uxASRulThDVvaseLJNL1lqDk+RgiuJnDGtsPk9IGOP9GYqxiRjWHpH+4vrVllVDdYdNYfq9cyoe+Ry/+K66lgjVFBl7gVKc2XCv847iNeHDRy7GTjn58GpH6binlQ2ZuIT0hVNPWgUt21b1Iu5uea5w7V/PiPXJz0J+LK3y+iQXEL93lS6J6qflDOgCdDGt6FnTpg9OJ3IwKGwG5Ew2lAhYzZw7GXgdKAHVYul4lxJZfseNEpCIVX6tFIWejq3FIcyiuDi7IRJSeE2HpiOoX34yzyZ30QX0UE3qT0ixs5hqhUHs0TCse4Zfa9c7vgUKMmFXjlmrqDqFu5bZ+B0HQ8jwQaO3Tw48ef9VqyFY0cufBhw8wHObJXu6BZQZqXDOwcj0Lv9fcYcht1fA8f/BFw8gJmvy6o1xiEY2jkYAV5uIndwe3IedE/ncUBUf9lAect70CNVNbU4ZZ4odwt2rfNUd7kQRoKvMnZLMj5/D46ihXPc7FpkbIhfBDD6njovTmXzcvNK9RQ312yF4mxg6cN1mjeh+u9xw7QdNxdnTEwMbxDO1X0xwgX3yXVSPi8vgN5Izi1Fda0J3u4uiC7cDVSXA76RQFgijAQbOLakskQmblkpRNU9Qho4ZHmTBc7YwRUdEC+NVFIwbUR2UYUQ+CM4/6YVfn9ItsEg8bBRd9vueDGaRcnDWX4wU8gq6J6kS2QTTjJuzrNJrxocyyq2TJqdlPAUeW8M1mKGDRxbkm/23ngEAF6B5/12kf6e8PVwFZZ3MjfdtI9E+9Rn5fq6V5tUTVBrBrpW94sNQLS5wo1pxOHfgH2LACdn4JI3DJXAyLSdsT3C4O7qLDwHRzINUAVKQnhjH5TrG9+Sgq464liWHG93yr85tsKQ4SmCDRyd5N8oFQld6YQEcNQIFwm9zNQ6jwVqKoBljzb4Lw5PnQNKwPzJHOYbeRcQPdBmh4nRNj4erhjTLdQ4YSqi9+VAcFfpndw6H3r04PQLKJHtJ+AEdJuk9rCsDhs4toS6z1op/0ahu2LgmE9QxsaQy3bai7Js/ODPwHHZh6a4ohrrj8nw41QuD28KubZ+vQ8oyZJxfaUNBuOwTOmtlIvrvPmmgosrMOYBuU7961rI09NyBdXQyi3yidihhpRtYANHRx4cgg0cFYjoBQy7Ta7/9hBQXYFVh7NQWVMruvCKMkumIXu/Aw78CDi7Ape9B7h58h5ycCYmRYj5wt7UAqTll8EQ9LtaapxRv8HtH0MP1NaacDzLrIGTu0Y+2XMajAgbOPbIwbFCgrGCcjM9mqmvmK/uoeofnzAg57DIx1HUiyf3pou2sRLzzpvCNGCJeWY77iEgeoDaI2I0QKivB4YkBIn1FQeN4sVxA8bcX9eEs0r7hltaQZlQlvZ3qYD3mfXyyR7TYUTYwNGdB8fPIvZnCNEsvUDNUqe/KFZNa1/CqUM7xTqHp5oJTf34d1ldEj0IuMB88WeYetVUy8wTBENAPdX8Y4HiTGD7Amido+b0hlkBx+BEuYVUKRqeBCPCBo5ONHAUYoKo86szKqtrkXJWPzFfwyQVdp8Cp5pKPGp6DxG+bhgQe/7VcYbrNUWCfq6eMjRFeQoMY2ZyL9ldfNOJXBSUVhljv7i6A2PNHss1L2m+ouq42cCZ5rqjLjxlUC80Gzi2oqpcWvRW6kOlQC0BFME/TjS2M3QRuPh/qHD2wnDnQ/hn5FY4OxvzwtAh0nYByx+T65OfkS0vGKaRGnuPCF8hdbHycJZx9s3AOUBwF5mLs+ldaL2CyhXVGFRqDk8lzYRRYQPHVhSckUuS+6fwhhWx5OGYtQwY+1HrH4e3cI1Yn5n5DlCYzrufoFnrdzcDNZWyG7uSlM0wjZhi9uIYovlm/VwcpVJww+tA6Vlo2cC5wHkfPKsLAZ9w2RvOoLCBYyvyk+vyb6zs/lMqqY6xFo7d2ZmShzdLJ2IvusK1qkjmmxhBmfW8G2neD5w9IRPqSdDPoC5vxnrl4qsOZ6O8qsZYIeyIvkBFIbDuFWgRk8kkPP8XO2+ST/S6VIoWGhQ2cHSUf6PQzZxozCEq+0PJkbVwxs+dH5d5Jsf/ALbpS+TL6uz6Etj7jdQKuuJDwDtY7RExGqZvTIBQZS+trMHG4/rtxt0EaiA78XG5Tu0bqJpQY+QUV6KsrBRTXbbJJ3pfBiPDBo6OKqga96Q6nl0sNA0Y+81+FPXiAYOGA5OelP+x7DEg97hjHob0PcCv5gTL8f8C4keoPSJG45CsgqWaykhhKqL7ZCB+pGxeuVpWXWotPDXGeQ/8nUplc00aq4FhA0dHGjgKCcHecHNxEjMg0jRg7AN5zE7lloqeOuN6hAHDbgc6jwOqSoHv/wrUVDvWoaA8g4XXA9VlQLfJXBLOtDtMtfxAprHkLig0q0x8dnwKZOyDljiSWYQrXNbKP/pcLr1OBsbY305NClPlMiDW6m/t6uIsqhEIDlPZj6X75GyTeupQbx1xcZj1tmymmroNWPUCHAYy5iipmDyVQZ2BKz4w/MWSsR7DO4fAz9NVhEx2peQZa9eSF7PXLMBUC/z+T03l6KWcOY1JztvlHwNvgNHhK5Ktq6j8Y2zy9orgHysa24+lZnd6A3E/MmBnvCzX174EHDV35jU6fzwFnFglqwSv+dLqlYKMsSEv6ITEcOOJ/ilMfhpw8QBOrQUO/QKtEH36Z7g71SA/sDcQ0RtGx6YGTl5eHubMmYOAgADxoPX8/Pxz5jk8+eSTiI6OhpeXFy688ELs37+/wTb0HMVx6z+uuUaW7mqC2lqgyFw+HGAbA6dHhDRwDmVwqbg9IFHFfamFINmbiUnywmyh75XAkFvl+ve31Rm3Rk4qplJYYtZbslcXw7STujycTHHdNxRBCcCou+X6skelLprKmGprcUHR72K9ou91cARsauBcd9112LVrF37//XfxoHUyclrjP//5D15++WW8+eab2Lp1KyIjIzF58mQUFTW8kd92221IT0+3PN577z1oBhJ7Ij0QakHvF2WTj0iMkgbOYTZw7MLv5vAUudZDfD2abjD1eSCqP1B2FviW9GAMotLamBOrgZ/MF27qpGzwKgzGdlAem7uLM07mlIiCCcNxwX0ykTfvFLDpLbVHg6yD69HD6TQqTG4IHs4Gznlx8OBBYdR8+OGHGDlypHh88MEH+OWXX3D48OFmX0NW/KuvvopHHnkEl19+Ofr06YNPPvkEpaWl+PLLLxts6+3tLYwf5UEeopaoqKhAYWFhg4dNUWbwvhFSAMoGJEX6i+XRzGJU19Ta5DOYOpbskx65i/rWC0/Vh7plX/WJzMc5swVY8n+air1bhaxDwMI5QG010OcKYPyjao+I0TF+nm4Y1S1ErCvNaw2Fhy8w+Sm5vvq/0tBREdPmd+RQPMbCzdcxpBxs5sHZuHGjMDqGDx9ueW7EiBHiuQ0bNjT7mpMnTyIjIwNTpkyxPOfh4YFx48Y1ec0XX3yB0NBQ9O7dGw8++GATD099XnjhBUuYjB5xcdavbGo+wdg24SkiNsgL3u4uqKypFTMgxnakF5Rh5+l8USDRanPN4M7A5e9Lz932j2VfJqNQlAF8cRVQUSBLSy99m5OKGSuqGhvQwCH6zQY6jZGVhiSGqdakpyAVYSkyPLUrWkPpHHo1cMhQCQ9vlKsAiOfo/1p6DRERIWOzCvR3/ddcf/31+Oqrr7Bq1So89thjWLRokfD4tMTDDz+MgoICyyMlxVzCbSsUgScbJRgT1AOpZyTn4dizempIQhDC/T1b35ga11GCoXjhw8ZIOi7JBT69FCg4DQR3lUnF5LFimPNkUq9wMXHYnZIvJhKGg77cjFdlwjGJgu79Tp1xbP0QLqYabKpNgm+nQXAU2m3gUAJw4wTfxo9t26RKIq03F4Zq7vn6NP7/xq+h/JtJkyaJEBYlF3/33XdYsWIFduwwd0dtBHmB/P39Gzz0XEGlkGgxcGwccnNwlpgNnGl92phPRcmFVIJJZaJUSk1ieHqlLB/4bBaQfQjwiwbmfM9KxYzVCPfzxOD4oAZ5boYjtBsw7v/kOpWN27tPVXkBsO0jsfpx9TTLfcMRaLeBc9ddd4n8mtYeZHhQXkxmZlO3Y3Z2dhMPjQK9hmjs4cnKymrxNcSgQYPg5uaGo0ePQhPYIURFJJrzcDjR2HZkF1Vg6yl5QZrWp5XwVJOu468ACRfIvjSfXQbkHIMuG2h+eTWQsQfwCQP+8hMQ1EntUTEGY3pfOXH4ba9BDRxi1L1AeC9ZgPLLPPuGqja/B5Tn45gpBstrB6On+b7hCLTbwKG8l8TExFYfnp6eIqmYwkFbtmyxvHbz5s3iuVGjRjX73p07dxZGzvLlyy3PVVZWYvXq1S2+hqAy8qqqKkRF2aZiqd0UpNrVg3MwnUvFbQW1ZqBrUf+4QMQEerX9ha7uwLVfApH95EWNQjyKurUeoFnmp7OAlM2AZyAw5wcgtLvao2IMyHTzxGFr8llkFapfTm0T6HpAoqDOrsCBH6XUgr28NxvfFKuvVV0GH093RAc4TnjZZjk4SUlJmDZtmggnbdq0STxofcaMGejZs6dlOzKIFi9eLNYpDDVv3jw8//zz4rl9+/bhpptuEhVTVHJOHD9+HE8//bQIg506dQpLlizBVVddhYEDB2L0aI20fVdycGygYtycByc1vwyF5QYtS1YZxW2uXITbhWcAMGcxENIdKDwDfHqJPoycokxgwQypzkwCfvQdIvuoPSrGoEQHemFgfKCYSPxu7vVmSKIHyn5txG//AM6esP1nbnhDGDlFfl3wa+0I9IzwO2eKiJGwqQ4OVTr17dtXVEXRo1+/fvjss88abEMl4+TVUfjHP/4hjJw777wTQ4YMQWpqKpYtWwY/P+mtcHd3xx9//IGpU6cKQ+mee+4R7005OC4uGmj7XlsDFClJxtE2/agAbzdEma3xI6yHY3XySiqx8URuxw0cwicUuPEH2XSVLmgfT9d2Y04aG40xa7+UObhpCRDjOEmJjDpcZM5v+3WPWSDVqIyeBySMBiqLge9uta0AYN4pYL0U5FwReRtq4WwpTHEUXG355sHBwfj8889b3aaxgiVZl5TITI/moBJvCllpluIsqRPi5CxFnmwMnbDpBeVC0XhIJ8fQNrAXSiPAXlH+SAiRvb86BHnybv5NhqlyjwEfTZPJupF9oSlOrQMW3gCU5QEB8dIwC+mq9qgYB2B630g8t+Qgtpw6K/LewvyaEdM0As4uwGXvAu+OAdJ2AEseAC55U+btWZtljwI1FaJM/ZfKIZRR6FAJxgT3orJVeIoUjF1saj82CFNxJZUK4n4dMXIi+gIlWcD8qcChX6EJaJKxfYHMuSHjJmYwMHcFGzeM3YgN8kb/2ADjh6kI8uZe9bGcBO/8HNjygfU/4+AvwMGf5WdMfxEHzR5+R0owJtjAsTaUa2GHBGOFJHPLhkOcaGxVCsqqsP5YTvvKw8+Fbzhw089A53FAVQnw9fXAmpdk7zK1KC+U/bN+vheorQJ6Xw7c9Cvg13LVIsPYgoss1VQGD1MRXScAk56qKx0/tMR6712SI3/PxKh7kOfbHWkF5Q3uF44CGzg2q6Cybf6NghJTpVJxwzWsU5E/DmaiqsaE7uG+6Bbua703pqTdGxYBQ28j1wnw5zPA55cDhSpc1FO2AO+PA/Z+Czi5ABMfB678CHBrR7UYw1iJ6eaJxKYTucgtrjD+fiW9rAHXA6Ya4NubgJNrrZMD+v1fZeVmeG+R1Lw/TeqkJYR4i/YYjgQbODbTwLFtBZVCl1BfuLk4oaiiGmfyDKgEqhK/KdVT5lmlVaH+ZBe/BMx8HXD1Ak6sBN4ZBez51j76GOS1+fVBYP4UmfjsT+GzJbJ5pgNVWDDaIj7EG31i/FFrMmhvqsbQb42uAT0vlrkyX10DHF95fu+5/HGpmOzmLdvGuHpgf5os4ukd7VjhKYINHFsZOHYKUbm71mXGKycyc34UV1Rj9ZHs86ueaguD/wLcvlomG1MX8u/nAp/MBDIP2ObzqMP51g+BN4cAWynub5IzyDvWAvEjbPOZDNORMJU5/83wUJ4meU27jJeVVdTvraPtHKihp1nzRmjuREpphwPp0oPTO7rlhtRGhQ0cW4WobKxiXJ8+5hN3Xyq3bLAGKw9lobK6Fp1DfWxfdRDWE5j7h+zM7eoJnForvTnf/AXI2Gedz6gsAbbOB94cCvz6AFCcCQR3AW78SV4Ivbn6jtFWufiG47k4W1IJh4D6ul23EOh9mcyDW3QrsOT/2l5CXlMNLH0EWPms/Jtye+i9zCghql7swWH05sEhesdIA2dvKntwrIGixUGtGewiiuXqIXvV/H0LkHSJ9Kwc+AF4d7QMI+34VDa8bG8snmL6ZNC8nAT8ej+Qd1K2XLjoJeDOzUCXcbb6RgzTITqF+ghZBpJnWH7A4NVUja8BV8yXOjnElvfl7//wb62HrbMOSa+v4rmZ+ARwwby6VnKVNTiRXeywISrb1zE7EmRJF6Xb3cDpYz5x96UWtKmZKdN6eGrl4SyxPqOfnVt/BCUAsz8DMvcDa/4LHPhJtkqgB+4BogcAcSOA8CQguDPgFQy4ewPVldK9nZ8M5J4AzmyVrynPr3tv8tgMu102AfWwYtI0w1gZkmWgsMqvezMwe2i84+xf0siZ/BTQ6QLghzulZhbl5YT2BPpdDcQNB/wigaoyeY2gMvAjZADVAm4+wKy3GnhuiIMZhSKniXSFqLGpo8EGjjUh1z+dbNRvhEqC7URSlD9cnJ2QW1KJzMIKRDpQrxFrs+JAJiqqa9HFPJNUhYjewFULgKIMYPfXssopcx+QtlM+2gr1kEqcIS96VJbqzBFpRvtQYv9Ly45gw7Ec5JdWItDbHQ5F98nA3duBtS8BWz4Ecg7LasuWoN/41OflBKkR+5XwlFrXMpVhA8eaUKIolQG7+0lr3E54urmgW5gvDmcWCS8OGzgd55c9UqhxRv9o9T1hNFsjdzM9qIz85GogY6+cvZGgJInyUX4NubfdfWTlXmACENVPenpoSRVbDKMjuob5itw3UmcnNfGrhsTB4fD0ByY/DYx5ENi3CDj+J5CxR4aqScaBjBlq+TDgOpnH1wIHHLiCimADx5pQNcxDp4Bq+2s49I7xlwZOWgEm9WKRto5QUFplqZ6aae/w1LnwjwL6XyMfDOMAmjhk4CzZm+6YBk59Q2fIzfLRAfanOW4FFcE+a1tAM2o7w5VU58/SAxlC3I867naPcCzFT4bREhf3k/IMa4/miKa3TPupqqkVRqIje3DYwDEIfcyVVKyF03F+MVdP2T25mGGYBnQL9xO5hdW1JuP3prIRx7OLhdyFr4cr4oO94YiwgWMQFI0D6iye4wgy51aGNDeU3lOUf8MwjLrM7C8nGj/tMjcwZtrFvtS6BGNnZ8esrGUDxyCQlU6VP/Xjrkzb+X1fhtDeIFcuCfwxDKMuM/vJicamk7nIKmyj6B1jYXeKlInoF+uY+TcEGzgGQhH8o0oqpoPVU+aLKsMw6hIX7I2B8YFC504JHzNtZ88ZaeD0jwt02N3GBo6BUAT/9p5hA6c9ZBdViA7GBOffMIx2uMQcLv7ZPAFh2kZFdY2lB1X/WDZwGAOgWOq7zZY70zaosR+pfdL+o1kjwzDa4OJ+UaD0kZ2n85FytlTt4eiGQ+lFoiI0yNsNccFecFTYg2Mg+sYEiIsBJRpncsy6zfyyO12b2jcM4+BQe4ERXULEOntx2o4yye0XG6i+YKmKsIFjIHw8XNHDrN9CMx7m3GQUlGNr8lmxflFfNnAYRqthKq6maju7U2SagiPn3xBs4BiMAeYTepc5g55pnZ93p4kkxqGdghAd6LiuXIbRKtP6RMLNxUmI1h3NlMJ1TNs8OP0duIKKYAPHoAaOUiLItM7inalieekA+3V/Zxim7VCzzbHdwywTEqZ1isqrhMifEqJyZNjAMRgD4gMtJYKk68K0zJHMIlFpQLPDizk8xTCa5ZIB5jCV8Ljyda019qYWCK90TKAXwvzs3zZIS7CBYzC6h/vBx90FJZU1OJYlrXimeX4we28u7BmOIB933k0Mo1EmJUXA080Zp3JLLQq9TPPsMcuE9I9z7PAUwQaOwXBxdkJfc9x1V0qe2sPRLLW1JvxoloCfxeEphtF8AcXEpAix/uMuOTFhzqVgHOjwu4gNHAOiZM7vMmfSM03ZlpyH1Pwy+IkLZzjvIobROJcq1VS70zj83gIUvqNrGzHQwSuoCDZwDIhyYnMlVcv8YJ4FUoWGp5uLnY4MwzAdRYSSvd2QVVRhaYzLNORMXplQZqe8wv5s4LCBY0QGxAWJ5eGMQpRWVqs9HM1RWV2LX829bS4byNVTDKMH3F2dMdPsxfl+xxm1h6NJtpu9N72jA3jixh4cYxIZ4InoAE/RfmAXC/41YdXhLBSUVSHC3wPDzSqpDMNon8sHxYrl0v2ZKK7gyVtjtplFSwcnyEmuo8MhKoMypFOwWG49xYnGLYWnSPuGkrIZhtEHJFzXJdQHZVU1+H1fhtrD0Rzbk2WC8RA2cARs4BiUoZ0VA0da9IyksLwKKw5mifVLzdoaDMPoA+qrdPkgGVbmMFVTgT9KSyDYgyNhA8egUOsBYsfpPFTX1Ko9HM3w2950kYPTI8IXvaL81R4OwzDtRFEd33giF2n5Zbz/zFD/QUpLoO7h4f6evF/YwDEuPcL94O/pitLKGqHWy0i+2y6TE2cNjHHoLrsMo1figr0xvHOwUOtVws1MXYLxkATpvWfYg2NYnJ2dOA+nESdzSkROEqXdXGFOVmQYRn8ov9/vd6Ry64ZGBg6Hp+rgEJWBGWpONN7GeTiCRWbvzdgeYYhgFy7D6JbpfSPh4eos2tFw6wYI4cOdp9nAaQwbOA6Qh0OJxo7eoI4uAIvM2hlXDY5TezgMw5wHfp5umNI7Uqwrv2tHZn9ageg/SGkJPSL81B6OZmADx8BQTyoSx8oprhRN6hwZUj5NLyhHgJcbt2ZgGANwuVmkk1o3UOGAI7PxeK5YDuscwtIX9WADx8B4uLpggLnh2taTjl0u/q05PEWl4dyagWH0z5juoQj388DZkkqsOJgJR2bTCWngjOjCCcb1YQPH4AztLMNUm07KH4AjUlBahaX7pSgYh6cYxhi4ujjjqiEy2Xjh1hQ4KiQDogi6juzKyuz1YQPH4IzqGiqWG47lOmwezk97pAs7MdIPfWJY+4ZhjMLVQ2Q+3Zqj2Uh1UE2cfWmFom0Fhd+TIvn6Vh82cAwOlQxSHk5GYTlO5JTAEflum5zdXTk4lrVvGMZAJIT4YGSXEKGJ8635d+64+TfBQh6EqYMNHIND+SaD42WYasOxHDgaRzKLsPtMAVydnYS4H8MwxmL2UOnF+XbbGVEt6aj5N2ToMQ1hA8cBGN1NnvjrjzleHs5XW06L5YTEcIT6eqg9HIZhrMy0PpGiPJpCVFQt6UhUifwbWUAygg2cJrCB4wCM6hZq6d3iSDOc8qoai7jfdcPj1R4OwzA28lJfZvbOLnSwMNWeMwWiHU+gt5vIMWQawgaOA9AvJgB+Hq4oKKvCgTTH6Uv16550FJZXIzbIC2O7h6k9HIZhbMTV5jDVsv0ZomzcUVh3NMcSnuL8m6awgeMg5ZTDzfoI6487jgv3S3N46tph8fzjZxgD0zs6AH1jAlBVY8LinY7TgJOqx5T2M0xT2MBxsHJxR4lRH84oEs3nKLlY0cpgGMb4XpwvNyc7hCQG6Xsp/afYwGkeNnAchNHmPJwtJ8+K3BSjQxc5YnKvCIT7eao9HIZhbMysAdHwcXfB8ewSS+m0kSFvPKVUdg3zQUygl9rD0SRs4DgIPSJ8ERXgiYrqWpFsbGRKK6vxvdlNzcnFDOM4DTgvGySTjT/dKCc4RmbNEQ5PnQvXc27BGAInJyeMTwzHl5tPY9WhLIzvGQ6j8svudBSVVyM+2BujzaE5pu3U1NSgqqqKd5kGcHNzg4uLi9rD0A1zRnTC55tOY/nBTKQXlCEqwJieDQrBsYFzbtjAcSDIqCEDZ+XhbDxpMhlS1Zd++F+Yw1OcXNz+fZeRkYH8/HybHBumYwQGBiIyMtKQv1dr0zPSTyj6Uij+q82ncf+UnjAix7OLkVZQLlTqR3Rmgb+WYAPHgRjVNQTuLs44fbZUtG3oGuYLo7EzJV8oF9MPn5OL24di3ISHh8Pb25tvqBowOEtLS5GVlSX+joqKUntIuuDGkQnCwPlySwrumtBdXAuMxuojslhkWKdgeLmzh68l2MBxIHw8XEW5+NqjOVh5KMuQBs6C9afE8pL+0axc3M6wlGLchITwjFAreHnJEAsZOXRsOFx1bqb2jkSYnweyiyqwdH8GZvaPhtFYcSBTLC/syeXhrWE805ZplQvNuTcrD8tZoZHIKCjHkr3pYv2mUZ3UHo6uUHJuyHPDaAvlmHBeVNtwc3EW4WniMwMmG1N5+BZzewaqEmVahg0cB2O82eInF25xRTWMBOXeVNeahNu2T0yA2sPRJZznoT34mLSf64bFw8XZSRgCRlNvp8kptdyhyljqps60DBs4DkaXMF90DvURip+rDOTFIW0fSqAmbh7N3huGcWQiAzxFE05i/rqTMBLLzeEp9t6obODk5eVhzpw5CAgIEA9aP1eFxvfff4+pU6ciNDRUzFx27drVZJuKigrcfffdYhsfHx9ccsklOHNGNlVk2hajJn7bl2GY3fXz7jTkllQiOsCTf/gMw2DuBZ3FXvhpdyoyC8sNsUcqqmuw2qx/M7mXvI4zKhk41113nTBQfv/9d/GgdTJyWqOkpASjR4/Gv//97xa3mTdvHhYvXoyvv/4a69atQ3FxMWbMmCESJZlzM908s6FEYyOoGlO1ycfm5OI5IzuJ3lsMwzg2A+ODMLRTkPBWf7JBXh/0zqYTMrUg3M9DNFFmWsdmd4KDBw8Ko+bDDz/EyJEjxeODDz7AL7/8gsOHD7f4OjKAHn/8cUyaNKnZ/y8oKMD8+fPxv//9T2wzcOBAfP7559i7dy9WrFjR7GvI41NYWNjg4cj0iw0Qno7SyhpRUaV31h/LxYH0Qni6OeMacz8ahlGT3NxcUfV06lTrN9Yrr7wSL7/8st3G5WjMHdNFLD/flIwSA+QcUrd0YmJSBDcQVtPA2bhxowhLDR8+3PLciBEjxHMbNmzo8Ptu375dVBNMmTLF8lx0dDT69OnT4vu+8MILljAZPeLiHPsmSKG/qWYvzm/7ZNWRnnl39XGxvGZoPIJ83NUeDsOIa87MmTPRqVOnBp7nWbNmNdg7NJl77rnnHH7SZSsmJUWgU4g3Csur8e22FF2fmdU1tfjdnFag5BcxKhk4JBpGM5jG0HP0f+fzvu7u7ggKCmrwfERERIvv+/DDDwvPj/JISdH3iW4NpveJsugpVFbXQq/sPVOAdcdyRMXEreaYO8OoSVlZmfAyz507t8HzW7duxbBhwxo8169fP2EEffHFF3YepWNQ/7rw0fpTovpIr1APQcozDPZxF6KtjA0MnCeffFJ4AFp7bNu2rcXyRsqXsEXZY2vv6+HhAX9//wYPR2dwQhBCfd3FzGbDcf2Gqd5dI703M/tFIS6YNVwclcTExBavR6+//rpdx/Lbb7/B1dVVhOUJ8jjTpIw8zI888ogYU33PNhVJfPXVV3YdoyNx5eA4BHq7CQV3Ev7Tc489JYeStH4YGygZ33XXXbjmmmta3YZmJHv27EFmpixnq092drbwtnQU6slSWVkpKrTqe3FI6XPUqFEdfl9HnNmQF+ezTcn4aVeaRQBQT5zKKcFvZmG/28d1VXs4hoMmDWUqJaF7ubm0ayJERQe9evXCH3/8IZZkVHTp0kUYDhQqsidr1qzBkCFDLH+T+jAVQ5BRQ4UWdP3z9PS0/D95dSikRbmCNBljrAu1MrhxRAJe//MY3lp5TBgIetMWIi+7kk4wo5/xlJk1Y+BQaTY9zgXNXigctGXLFotbdvPmzeK58zFEBg8eLDrsLl++HFdffbV4Lj09Hfv27cN//vOfDr+vIzJrYIwwcH7fn4FnK6vh7a6vzh3vrz0B8jiTXHlSFHvlrA0ZN70eXwo1OPD01HadjxSeJq8JVWCSkUCGRHV1NcaMGSP+pqKEBx98EIcOHULPnrIB41//+ldRBFFUVCTkJtrDZZddhlWrVmHixIn47rvvGvwfJRZTXqCCs7Mz0tLSRAuM/v37N3mvmJgYYdzQd0hISGjXOJi2cfPozkIPZ39aIf48lCWSdPXEumPZwttO1VPUTJRpGzbzcyUlJWHatGm47bbbsGnTJvGgdSrnVi4wimuZZl8KZ8+eFRenAwcOiL+p4or+VvJrKEn41ltvxQMPPCBmazt37sQNN9yAvn37tlh5xTTPoPhAJIR4i2oqvbluqS3Dd9ul9tEd7L1xeKiKskePHhYPCF0zwsLCLN5imgDRNUKp4Dx58qQIpZOXp73GDXHPPffg008/bTEHp76HhqDrVHPGTf1+U9RYk7ENVHxww0hpPJInh7yTeoK87MRFfaOE951pGzadslPiHF0IlIonijW/+eabDbahCw55dRR++ukn3HzzzZa/lXDYE088IfJ/iFdeeUXM1siDQxcTmkUtWLCAG9G1E3LTzhoQg9f+OIrFO9Nw2cBY6IW3Vx0TbltqyzCcZzQ2CxORJ0Wtz24PFBInA0aBDBxK4FUgA4euF4qB88wzz4iKJqrK7Ajjx48XHpzmIA83hdDrQ+NpycChSR1BBhljO24b00Xo4exOyceaozkY10Mf+7uwvEp42YlLB3B4SjMGTnBwsNCoaY3GlvRNN90kHq1Bs6M33nhDPJjz47KB0sBZdzQbWYXlCPdvOPPUImn5Zfh6i6yEmze5u+7i6XqB9qtewpZk4NQvwSaDYtCgQZZrDIWNyHtM14yjR48iNTVVhLtJXsLaKNpcjT1MFNZqDjK+YmNj2xT6ZzpOqK8HrhuWgI/Wn8QbfxzF2O5SLV8PycXlVbXoHu6LAXGBag9HV3AqtoPTKdQHA+MDRS7LD7tSoQcoUbCyplZ4bkZ15ZuCo1NbW4v9+/c38NicOHHCks9C6/Hx8SJsfvz4cTz99NNCf4YMi8YGjmL0NH5QDk1boVYzNJ76XhwaIxlh9D71PdbE2rVrG+h6Mbbj9nFd4O7qjG3Jedh4PFcXu/obs37P1UPidGGQaQk2cBhcNVgKH361JQW1GteJOJNXavnB3ze5h9rDYTQAGS2Uv1LfwKFwEIW0qaJJMWQoP4dUzEllmJKRmzNwKGRFzzd+1E8aPhcUKqMqqm+++cby3LPPPouFCxeKhGIysBTKy8tFDiLlJzK2J8LfE9ea1c5fXHpY87k4RzKLsCslH67OTqIohGkfbOAwIq7r6+GKkzkl2KDxWc2bfx4TvWVI6GpEFxa7YoDu3buLGxV5aRR+/vln4UEZO3ZsA0PmoYcewksvvSTWjxw50qDgwZo89thjeO2114TnhqBCCAqL0TipokuBBAGpfJxU3hn78PcJ3USOF+XiaL244putcjI3PjEcYX4sIdBe2MBh4OPhKnJxiC82J2t6NqN4b+5n7w3TRihcpBg4s2fPFjo5FCqinBcS4OsIFIa66qqrsGTJEpE/QyrF9bnoootw++23C6OmNUjygnMJ7Uu4n6dF3fg/Sw+LFghapLSy2nK9mz3EsdsLdRR9ZBAyNuf6EfFCE2fZgUxkFpYLV67W+Pdvh0Su0NTeERjSibUgmLbx5ZdfNnmOQk5k+HSUpUvPrQ907733nnMb0uJh7M9fx3URk7kT2SVCbuKaYXXeP63ww840oX0TH+wtPDhM+2EPDiNIjPTHkIQg0avli82nNbdXNhzLEQJdFIt+aFqi2sNhGEbH+Hu64e/ju4n1V1ccRVmlOordLUGhzAUbTor1G0cmsPZNB2EDh7Fw02jZ+fizjac09YOnxOfnlhwU69cPj0eXMF+1h8QwjM65YUQCYgK9kFFYjndWy552WmqseSSzGN7uLriKw1Mdhg0cxsK03pHCHZpXWmWJ/WoBciGTxLqfhyvumdhd7eEwDGMAPN1c8MjFSWL93dXHcTpXO0rSH6w5IZZXDIpFgJeb2sPRLWzgMBZcXZxx2xiZfPfB2hOaSL7LK6nEC79J7w0ZNyG+XEnAMIx1oMabo7uFCFX0Z36V7YHUZl9qAVYezgZ1ZFCSoZmOwQYO04ArB8ch2McdZ/LK8Mse2b1WTV78/ZDwKCVG+llCaAzDMNaAhPOenNlb5PYtP5CJlYezVN+xb/x5VCwv6R8thFiZjsMGDtMAL3cX3GI2JF5dcQRVKnpxtifn4WuzDsSzs/rAzYVPV4ZhrEv3CD/cNEpe8x7/cR9KKqpV28WHM4qwdH8mSLD4rgkyCZrpOHzHYJpw0+jOCPFxx6ncUny7TXbstjflVTV4aNEesX71kFguC2cYxmbMm9xDJBynnC3Df5fKhqxqoHz2RX2i0C3cT7VxGAU2cJgmkKqxUkL52h9HhLFhb15efgTHsopFg7yHp8tEQIZhGFtd8164XHajX7DhFDafsL+iO/XGWnEwU5SE3z+F29BYAzZwmBaF/2hGk1lYgfdWy4x+e7H11FmR5Ez8+/K+CPLpmNoswzBMWxnbIwzXmPtU/WPRHruGqkgK43mzFMZ1w+LRlaUwrAIbOEyzeLi64KHpUlDvrVXHkJxbYpc9VVBWhfu/2QXqgXfV4FhM6hXBR4hhGLvwr4uTEB3gieTcUjz6wz67NeMkWY69qQXCk3TvJJbCsBZs4DAtMrNflKWE8smf9tv8x06zmAe+2S3i4LFBXnhsZi8+OoxdoU7j4eHhOHXqVKvbXXnllXj55ZftNi7GfgrHr107UISJFu9MtYseWFZRucV7M29SdxGWZ6wDGzhMqyWUT19K1UtOQpfh+x2tNw48X95bc0LEoN1dnfHuDYPFxYZh7MkLL7yAmTNnolOnOkmCefPmYdasWQ22e/zxx/Hcc8+hsLCQD5DBGNop2NLM9/Ef9wtdGltBk0aaPFLPqb4xAZZqLsY6sIHDtArFgudNUn7s+2wWqlq2PwP/XXpIrJMuRZ+YAD4yjF0pKyvD/PnzMXfu3AbPU6fwYcOGNXiuX79+wgj64osv+CgZkL+N64pxPcJQUV2LWz/ZivSCMpt8DnmIluzNEB6jf1/RV4itMtaD9yZzTu4Y1xXDOgWjpLIG93y10+pVVTtP5+Ger3eKTuHXDovHtcNkoh/DtJXExEThcWzu8frrr7fpPX777Te4urpi5MiR4u+qqiq4u7tjw4YNeOSRR8R7DR8+3LL9JZdcgq+++ooPkgFxdnbC69cORLdwX1FoceuCbSi2ctIxad488ZPsaP/AlB7oHc2TOmvDBg5zTmh28co1A0RPlN1nCvCP7/ZYLR9nf1oBblmwFeVVtRjfMwzPXNpb3EgYDUDHuLJEnUc7z6/FixeL5R9//IH09HScPn1aGCvffvstbr/99ja9x5o1azBkyBDL3y4uLli3bp1Y37Vrl3jfpUuXWv6fvDpbtmxBRUVFu8bK6AO63n1801CE+rrjQHohbvpoi9WMnOyiCsz9VF73qHrrjrFdrfK+TENcG/3NMM1CJePvXD8IN360BT/tTkNUgCf+OV3OmjvK7pR88X5UOdU/NgBvXjeIXbRaoqoUeD5anc/+Vxrg3naZ+oyMDGHQjB49Gh4eHsIgqa6uxpgxY8Tf//vf//Dggw/i0KFD6Nmzp3jNX//6V3z44YcoKiqCj4+PSCyOjq77vs7OzkhLS0NISAj69+/f5DNjYmKEcUOfnZCQYKUvzmiJuGBvfHzTMFz34SZsS84TRs78m4aeVwPMovIqEfaiYoqEEG+8cnV/4TFirA97cJg2M6pbKJ67rI8lIfjZXw+KyqeOsGRvOq55f5MwbgbFB+KzucPh48H2NtMx9u7dix49eghjhiADJywsDBERUmZg37596Nu3Lw4flkqxJ0+exLZt29ClSxdh3Cg5OJ6eng3ed+fOnc0aN4SXl5dYlpZqpws1Y336xgbgi7nD4e/pKoycy99ej1M5JR1uHnzdB5ux50yB6Pm34OZh3EDYhvAdhWkXs4fGo7LGhMd+2If5606KpOP/XS3DV22hrLIGLy07LF5LjOkeinduGCz0HxiN4eYtPSlqfXY72LNnjzBgFMjAoURgBTJwrr76aouB88wzz4jKqO3bt1u2CQ0NRV5eXoP3pfdpycA5e/asWJIhxRibfrGB+OqvI0QuzvHsElzy5jo8M6uPaIjZVi82VWPd+cUOnD5bKoybT28Zhs7cTNOmsAeHaTdzRiTgf1f1F+XcKw5mYcorq/HjrtRWvTk1tSYR2pr22hqLcTP3gs4ixs3GjUahCzeFidR4tDP0SQZOfYOmvoFD+WIUfpoxY4YIUR09ehSpqaki9NSnj/RIEgMHDsSBAweaeIbqv299yGiKjY0VhhFjfCgJ+Ke7RmNgfKAo67736124ecHWc5aRU0jqP78fwuVvbxDGDYX7v7l9BFeK2gGeNjMd4orBsege4Yu7v9opVD/px079oy4dECMqriIDPFFrMuFMXim2nsrDz7vTcCZPllpS/s7zl/fF+J7hvPeZ86a2thb79+8X2jQKJ06cwGWXXWZZj4+PR1JSEo4fP46nn35abPv5559j7NixltdMnToVDz/8sPDiBAUFWd6bjCfKxaFQVkBAXaXL2rVrMWXKFD6CDkS4vye+uX0k3l11HK//eRSrDmeLx5CEIEzrEyk8PSG+7kIc9WROCdYezRHXPiU5eXKvCLx0ZX8EeLPGlz1wMtlLi1pDkDgXXagKCgrg7++v9nB0DZWMv7/mhOgdVVTeeoUBxbDnjumCm0d3gh+L+GmK8vJykZfSuXPnJnkoWoc8MpR/k5ycLAwZgsT6qALqxx9/FAbL999/j08++QSDBg1CZGQklixZIhKQ33777QahLSoRv+mmmyyVV2QEPfTQQ8LAuf/++0WysrK/KL+HqqpGjBhh0++n52NjZE5kF+P1P44Kz/S5UhG7hvngoWmJwsDhKlH73b/ZwGEDxyqUVlbj1z3pWHM0B/tTC0QZpKuLE8L8PNAnOgATksIxMTECXu4u1vlAxqoY+SZKisOkZ/N///d/WLhwoTBoevXqJQyUlJQU8X8KZPhQtRWFn6iKqiXeeustYTwtW7bM5uM38rExApmF5cJLs+nEWRzKKBSFE24uzsJTTZ6dyb0iRcsbNmzsb+BwiIqxCt7urrhqSJx4MIyWoPDVnDlzxPrs2bPFkjwylDtT37ghLrroIkuOTlxcy+eym5sb3njjDRuPnNEDEf6ewjNND0ZbsAeHPTgMw14CDcMeHIbpmAeHq6gYhmEYhjEcbOAwDMMwDGM42MBhGIZhGMZwsIHDMAzDMIzhYAOHYRgLDiiLpXn4mDBMx2ADh2EYUfZMcONI7aEcE+UYMQzTNlgHh2EYuLi4IDAwEFlZWWJveHt7szCZBjw3ZNzQMaFjQ8eIYZi2wwYOwzACamFAKEYOow3IuFGODcMwbYcNHIZhBCQlHxUVhfDwcFRVVfFe0QAUlmLPDcN0DDZwGIZpAN1Q+abKMIze4SRjhmEYhmEMBxs4DMMwDMMYDjZwGIZhGIYxHK6OLJxFXUkZhmEYhtEHyn27LQKYDmngFBUViWVcXJzaQ2EYhmEYpgP38YCAgFa3cTI5oA54bW0t0tLS4OfnZ3UxM7IuyXBKSUmBv78/jIbRv58jfEf+fvqHj6G+Mfrxs+V3JJOFjJvo6Gg4O7eeZeOQHhzaKbGxsTb9DDqgRj1xHeH7OcJ35O+nf/gY6hujHz9bfcdzeW4UOMmYYRiGYRjDwQYOwzAMwzCGgw0cK+Ph4YEnnnhCLI2I0b+fI3xH/n76h4+hvjH68dPKd3TIJGOGYRiGYYwNe3AYhmEYhjEcbOAwDMMwDGM42MBhGIZhGMZwsIHDMAzDMIzhYAOnnTz33HMYNWoUvL29ERgY2KbXUB73k08+KZQXvby8cOGFF2L//v0NtqmoqMDdd9+N0NBQ+Pj44JJLLsGZM2dgb/Ly8jBnzhwhpEQPWs/Pz2/1NaQG3dzjv//9r2Ub+s6N//+aa66BGnTkO950001Nxj9ixAhDHMOqqio89NBD6Nu3rxg3nac33nijUPuuj5rH8O2330bnzp3h6emJwYMHY+3ata1uv3r1arEdbd+lSxe8++67TbZZtGgRevXqJao8aLl48WKoRXu+3/fff4/JkycjLCxMCKiNHDkSS5cubbDNggULmv1NlpeXQ+vfb9WqVc2O/dChQ5o9fu39js1dT+jRu3dvTR7DNWvWYObMmeLaQGP44YcfzvkaTfwGqYqKaTuPP/646eWXXzbdf//9poCAgDa95t///rfJz8/PtGjRItPevXtNs2fPNkVFRZkKCwst29xxxx2mmJgY0/Lly007duwwjR8/3tS/f39TdXW1XQ/PtGnTTH369DFt2LBBPGh9xowZrb4mPT29weOjjz4yOTk5mY4fP27ZZty4cabbbrutwXb5+fkmNejId/zLX/4iXld//Lm5uQ220esxpOMwadIk08KFC02HDh0ybdy40TR8+HDT4MGDG2yn1jH8+uuvTW5ubqYPPvjAdODAAdO9995r8vHxMSUnJze7/YkTJ0ze3t5iO9qeXkev/+677yzb0H5xcXExPf/886aDBw+Kpaurq2nTpk0me9Pe70f//+KLL5q2bNliOnLkiOnhhx8Wr6dzTuHjjz82+fv7N/ltqkF7v9/KlSupstd0+PDhBmOv/zvS0vHryHek303975aSkmIKDg42PfHEE5o8hkuWLDE98sgj4h5Gx2bx4sWtbq+V3yAbOB2ETr62GDi1tbWmyMhIYeQolJeXi9e+++67lpOdDj79SBRSU1NNzs7Opt9//91kL+hEpJO3/glGNzt6jm58beXSSy81TZgwocnNkU52tenodyQDh75XSxjtGNLNk15T/wKt1jEcNmyYMB7rk5iYaPrnP//Z7Pb/+Mc/xP/X5/bbbzeNGDHC8vfVV18tDMH6TJ061XTNNdeYtP79mqNXr16mp556qt3XJy1+P8XAycvLa/E9tXT8rHEMyWCgSeGpU6c0eQzr0xYDRyu/QQ5R2ZiTJ08iIyMDU6ZMsTxH7rhx48Zhw4YN4u/t27eLMEH9bcgV2KdPH8s29mDjxo0ipDF8+HDLcxSGoefaOo7MzEz8+uuvuPXWW5v83xdffCHCN+SGffDBBy1d3e3J+XxHcp2Hh4ejR48euO2225CVlWX5PyMdQ6KgoEC4ohuHYe19DCsrK8W+rb9fCfq7pe9D+6Dx9lOnTsW2bdvEMWptG3seq45+v+aaB9NxCA4ObvB8cXExEhISRN+9GTNmYOfOnbA35/P9Bg4ciKioKEycOBErV65s8H9aOX7WOobz58/HpEmTxPHS2jHsCFr5DTpks017QsYNERER0eB5+js5Odmyjbu7O4KCgppso7zeXmOlG3hj6Lm2juOTTz4RXdovv/zyBs9ff/31Ij4dGRmJffv24eGHH8bu3buxfPly2JOOfsfp06fjqquuEhcbMlofe+wxTJgwQVzYyGA10jGkGP8///lPXHfddQ2a5KlxDHNyclBTU9Ps76el70PPN7d9dXW1eD+6aba0jT2PVUe/X2P+97//oaSkBFdffbXlucTERJHDQXlV1NX5tddew+jRo8Xx6t69O7T8/ej4vP/++yJ/g/LaPvvsM2Hk0ARj7NixYhutHD9rHMP09HT89ttv+PLLLxs8r5Vj2BG08htkAwcQCcBPPfVUqztq69atGDJkSId3NM2G60OevsbPNaYt21jz+zU3zvaO46OPPhI3Qkosqw95PBTIq0E/UNqfO3bswKBBg6D17zh79uwG46exk7FD3qrGxlx73ldrx5BmV5Q4TF4BSpq05zG05u+nue0bP9+R36St6OhYvvrqK3Fu/Pjjjw0MW/La1U+CpxsjHaM33ngDr7/+OrT8/Xr27CkeCpREnZKSgpdeesli4LT3Pe1BR8dDRgx5SmfNmtXgea0dw/aihd8gGzgA7rrrrnNWg3Tq1KlDO5hmuwRZpWS1KlB4Q7FeaRtyc1L1S30PAG1DFVv2+n579uwRIabGZGdnN7G0m4OqBg4fPoyFCxeec1v6obq5ueHo0aNWuTna6zsq0LEkA4fGb5RjSMYNeQHIQ/Xnn3828N7Y4xg2B4XDXFxcmszq6v9+GkPHorntXV1dERIS0uo27TkH1Pp+CvQ7o1Dwt99+K8IbreHs7IyhQ4dazlc9fL/60I3+888/t/ytleN3vt+Rbug0KaRKR/IAa/EYdgTN/Aatls3jYLQ3yZiqHhQqKiqaTTKmKhaFtLQ01RJUN2/ebHmOklXbmqBKibiNK29agqrJ6H1Xr15tsifn+x0VcnJyTB4eHqZPPvnEEMewsrLSNGvWLFPv3r1NWVlZmjqGlMD5t7/9rcFzSUlJrSYZ0//XhxJAGyc4Tp8+vcE2lPCoVpJxe74f8eWXX5o8PT3PmexZ/zo0ZMgQ080332zSw/drzBVXXCGqErV4/M7nOyoJ1fRb0vIx7EiSsRZ+g2zgtBOqKtm5c6eoWPD19RXr9CgqKrJs07NnT9P3339v+ZsqqMigoefoRL722mubLROPjY01rVixQpR7UhWSWiXG/fr1E5U39Ojbt2+TEuPG348oKCgQZYHvvPNOk/c8duyY2F9bt241nTx50vTrr7+KDPuBAwfa/ft15DvSsX3ggQdEWSONny5KI0eOFCXhRjiGVVVVpksuuUSMfdeuXQ1KUskYV/sYKiW48+fPFwbcvHnzRAmuUnFCN5E5c+Y0KVG97777xPb0usYlquvXrxclqvTbpBJVWqpdJt7W70fGDY31rbfearFk/8knnxSGNUk10PWJbor0mvqGr1a/3yuvvCJuoFQCv2/fPvH/dFOlEmUtHr+OfEeFG264QUgyNIeWjmFRUZHlXkfHgqRSaF2pstTqb5ANnHZCXgo6wI0fdNOz7FRAeHjqW96kb0CeHJr1jx07tonFXlZWZrrrrruEFoKXl5e4IZ0+fdpkb0jb5frrrxe6PfSg9cblmo2/H/Hee++JcTeni0Lfg74zfTd3d3dT165dTffcc08THRmtfsfS0lLTlClTTGFhYeJHGh8fL86DxsdHr8eQDJbmzun657Xax5Bu5gkJCeKzBw0a1MBrRMeCStjrs2rVKmF80fadOnVq1vD+9ttvhaFHx5SMtfo3UHvTnu9H680dK9pOgW6wdJ7S+9F5S+cvGeh6+H7k7abzizxUQUFBpgsuuEAY1Fo+fh05R+laSdeJ999/v9n309IxXGn2NLV0zmn1N+hE/1gv4MUwDMMwDKM+rIPDMAzDMIzhYAOHYRiGYRjDwQYOwzAMwzCGgw0chmEYhmEMBxs4DMMwDMMYDjZwGIZhGIYxHGzgMAzDMAxjONjAYRiGYRjGcLCBwzAMwzCM4WADh2EYQzFv3jzMmjVL7WEwDKMybOAwDGMotm7dimHDhqk9DIZhVIZ7UTEMYwiqqqrg4+Mjlgpk6GzevFnVcTEMow6uKn0uwzCMVXFxccG6deswfPhw7Nq1CxEREfD09OS9zDAOChs4DMMYAmdnZ6SlpSEkJAT9+/dXezgMw6gM5+AwDGMYdu7cycYNwzACNnAYhjEMFJpi7w3DMAQbOAzDGIa9e/eiX79+ag+DYRgNwAYOwzCGoba2Fnv27BG5OAUFBWoPh2EYFWEDh2EYw/Dss89i4cKFiImJwdNPP632cBiGURHWwWEYhmEYxnCwB4dhGIZhGMPBBg7DMAzDMIaDDRyGYRiGYQwHGzgMwzAMwxgONnAYhmEYhjEcbOAwDMMwDGM42MBhGIZhGMZwsIHDMAzDMIzhYAOHYRiGYRjDwQYOwzAMwzCGgw0chmEYhmFgNP4f07Uw+7rlXlgAAAAASUVORK5CYII=", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "\n", "degree = 4\n", "\n", "t = np.linspace(-1, 1, num=10000)\n", "\n", "gen = NodesGenerator(\"CHEBY-1\")\n", "alpha, beta = gen.getOrthogPolyCoefficients(degree+1)\n", "\n", "pi1, pi2 = gen.evalOrthogPoly(t, alpha, beta)\n", "\n", "plt.plot(t, pi1, label=r\"$\\pi_{M-1}(t)$\")\n", "plt.plot(t, pi2, label=r\"$\\pi_{M}(t)$\")\n", "plt.legend(); plt.xlabel(\"$t$\");" ] }, { "cell_type": "markdown", "id": "8d07a8e5", "metadata": {}, "source": [ "> 📣 Note that the `quadType` argument does not matter when generating orthogonal polynomials, and can simply be left to its default value." ] } ], "metadata": { "language_info": { "name": "python" } }, "nbformat": 4, "nbformat_minor": 5 }