{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# The trapezoidal rule: vectorization and perormance\n", "\n", "Here we discuss in detail the performance profile of various solutions to the exercise from our Numpy introduction, the trapezoidal rule for the numerical approximation of definite integrals.\n", "\n", "If we denote by $x_{i}$ ($i=0,\\ldots,n,$ with $x_{0}=a$ and $x_{n}=b$) the abscissas\n", "where the function is sampled, then\n", "\n", "$$\\int_{a}^{b}f(x)dx\\approx\\frac{1}{2}\\sum_{i=1}^{n}\\left(x_{i}-x_{i-1}\\right)\\left(f(x_{i})+f(x_{i-1})\\right).$$\n", "\n", "The common case of using equally spaced abscissas with spacing $h=(b-a)/n$ reads:\n", "\n", "$$\\int_{a}^{b}f(x)dx\\approx\\frac{h}{2}\\sum_{i=1}^{n}\\left(f(x_{i})+f(x_{i-1})\\right).$$\n", "\n", "One frequently receives the function values already precomputed, $y_{i}=f(x_{i}),$ so the formula becomes\n", "\n", "$$\\int_{a}^{b}f(x)dx\\approx\\frac{1}{2}\\sum_{i=1}^{n}\\left(x_{i}-x_{i-1}\\right)\\left(y_{i}+y_{i-1}\\right).$$\n", "\n", "In this exercise, you'll need to write two functions, `trapz` and `trapzf`. `trapz` applies the trapezoid formula to pre-computed values, implementing equation trapz, while `trapzf` takes a function $f$ as input, as well as the total number of samples to evaluate, and computes the equation above.\n", "\n", "Test it and show that it produces correct values for some simple integrals you can compute analytically or compare your answers against `scipy.integrate.trapz`, using as test function $f(x)$:\n", "\n", "$$\n", "f(x) = (x-3)(x-5)(x-7)+85\n", "$$\n", "\n", "integrated between $a=1$ and $b=9$." ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": true }, "outputs": [], "source": [ "%matplotlib inline\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "import scipy.integrate as sint" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Reference value (trapezoid): 680.0\n" ] } ], "source": [ "def f(x):\n", " return (x-3)*(x-5)*(x-7)+85\n", "\n", "a, b, n = (1, 9, 200)\n", "x = np.linspace(a, b, n)\n", "y = f(x)\n", "\n", "print(\"Reference value (trapezoid):\", sint.trapz(y, x))" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Our functions\n", "680.0\n", "680.0\n" ] } ], "source": [ "def trapz(y, x):\n", " return 0.5*((x[1:]-x[:-1])*(y[1:]+y[:-1])).sum()\n", "\n", "def trapzf(f, a, b, n):\n", " x = np.linspace(a, b, n)\n", " # Sample the input function at all values of x\n", " y = f(x)\n", " # Note that in the equations above, the list had n+1 points. Since here we're using\n", " # n as the total number of points, the actual size of the interval h should \n", " # use n-1 in the formula:\n", " h = (b-a)/(n-1)\n", " # Compute the trapezoid rule sum for the final result\n", " return (h/2)*(y[1:]+y[:-1]).sum()\n", "\n", "print(\"Our functions\")\n", "print(trapz(y, x))\n", "print(trapzf(f, a, b, n))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Some performance notes\n", "\n", "Let's compare the time it takes to compute the fully vectorized form in our `trapz` with a more manual implementation that does the computation with a pure Python loop." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "680.00000000000011" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "def trapzl(y, x):\n", " \"Pure python version of trapezoid rule.\"\n", " s = 0\n", " for i in range(1, len(x)):\n", " s += (x[i]-x[i-1])*(y[i]+y[i-1])\n", " return s/2\n", "\n", "trapzl(y, x)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's compute this over a much finer grid to really see the performance difference" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "35.2 µs ± 2.61 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)\n", "10.8 ms ± 444 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)\n" ] } ], "source": [ "a, b, n = (1, 9, 20_000)\n", "x = np.linspace(a, b, n)\n", "y = f(x)\n", "\n", "tv = %timeit -o trapz(y, x)\n", "tl = %timeit -o trapzl(y, x)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The output of `%timeit` is a special `TimeitResult` object that contains all the information about the timing run in a set of attributes. Try typing `tv.` to see some of its attributes." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "tv" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In particular, we can look at the average run time, which as we can see, is in seconds:" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "3.5238153685661796e-05" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "tv.average" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "With this, we can now see how much faster the Numpy version was at 20,000 sample points:" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "306.7110725043031" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "tl.average/tv.average" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's see what the performance difference is across a range of sizes (note this will take a while to run). We're actually going to need the same computation of $x$ and $f(x)$ over and over, so let's do it once and retrieve them later on:" ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "collapsed": true }, "outputs": [], "source": [ "a, b = 1, 9\n", "\n", "xy = []\n", "sizes = np.logspace(1, 7, 7)\n", "for n in sizes:\n", " x = np.linspace(a, b, int(n))\n", " y = f(x)\n", " xy.append((x, y))" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "3.43 µs ± 147 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)\n", "5.46 µs ± 136 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)\n", "3.72 µs ± 219 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)\n", "52.4 µs ± 1.63 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)\n", "5.49 µs ± 89 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)\n", "514 µs ± 11.6 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)\n", "17.3 µs ± 338 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)\n", "5.29 ms ± 143 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)\n", "144 µs ± 3.32 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)\n", "51.7 ms ± 1.03 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)\n", "3.56 ms ± 225 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)\n", "548 ms ± 10.5 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)\n", "71.7 ms ± 672 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)\n", "5.18 s ± 50.6 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)\n" ] } ], "source": [ "tv, tl = [], []\n", "tve, tle = [], []\n", "for x, y in xy:\n", " t1 = %timeit -o trapz(y, x)\n", " tv.append(t1.average)\n", " tve.append(t1.stdev)\n", " t2 = %timeit -o trapzl(y, x)\n", " tl.append(t2.average)\n", " tle.append(t2.stdev)\n", "\n", "tva = np.array(tv); tvea = np.array(tve)\n", "tla = np.array(tl); tlea = np.array(tle)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now let's plot these, along with proper error bars (remember your error propagation formulas, assuming here that the measurements are normal IID for simplicity):" ] }, { "cell_type": "code", "execution_count": 11, "metadata": { "collapsed": true }, "outputs": [], "source": [ "def rat_err(p, dp, q, dq):\n", " ratio = p/q\n", " err = ratio*np.sqrt((dp/p)**2 + (dq/q)**2)\n", " return ratio, err" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEaCAYAAAAL7cBuAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3Xd8VfX5wPHPkxAgIUAII2zCVnAA\nYQsKtopbap3UjYJWa1utq621y5+rrbWt1sF0IlUE66KIBAQBAUGQjSSMsEcCgQQynt8f5wQv4Sa5\nGfeee2+e9+t1Xjn3zOebm3ufnHO+Q1QVY4wxprQYrwMwxhgTnixBGGOM8csShDHGGL8sQRhjjPHL\nEoQxxhi/LEEYY4zxyxKEqTYRyRSRH4bgPL8XkTeCfZ7qEpF0Ebmjivu2F5FcEYmt6bjKOaeKSJdQ\nnc9EDksQtZT7pZ7nfhntFpGJIpIYwH6TROTPoYixNiidXFV1q6omqmqRl3EZA5YgarvLVTUR6AP0\nA37rcTxhT0TqeB2DMaFiCcKgqlnAJ8AZInKNiCzzXS8iD4jIdBEZA/wEeMi98vivz2a9RGSliOSI\nyDsiUt9n/ztFZJOIHBCRD0Sktc86FZG7RGSjiBwUkRdERAKJW0SuEJHVIpLt3tY53Wfd6e6ybHeb\nK3zWTRKRl0RklogcFpG5ItKhjHOkujGOFpGtwOfu8oEi8qV7/G9EZFgZ+3cWkc9FZL+I7BORN0Uk\nyV33OtAe+K/7+3zI53x1ROR6EVla6ni/FJEP3Pl6IvIXEdnqXgW+JCLxZcTRxS1njhvHO2Vs11hE\nXhORvSKyRUR+KyIx7rotIpLmzt/oxtnDfX2HiEx3538vIlPd4xx2f/99fc7RWkTec8+RISL3+azr\nLyJLReSQW6a/ucvri8gb7u8xW0SWiEiKvzKYGqSqNtXCCcgEfujOtwNWA38C6gEHgNN9tl0O/Nid\nnwT82c+xvgJaA8nAWuAud935wD6cq5R6wD+BeT77KvAhkITzZbkXuKiMmH8PvOHOdwOOABcAccBD\nwCagrvt6E/Br9/X5wGGgu08ZDgPnujE9D8wv45ypboyvAQ2AeKANsB+4BOefrAvc183dfdKBO9z5\nLu76ekBzYB7wd3/vQ6nz1QES3Di7+qxfAlzvzv8d+MD9nTcE/gs8WUY53gZ+48ZbHxhS6j3o4s6/\nBsxwj5cKbABG+6x7wJ1/BfgOuNtn3S993qd89/cTCzwJLHLXxQDLgN+5700nYDMwwl2/ELjJnU8E\nBrrzY93yJbjHTAMaef05ivbJ8wBs8uiNd76YcoFsYAvwIhDvrvs38IQ73xM4CNRzX0/Cf4K40ef1\nM8BL7vx44BmfdYlAAZDqvtZSX1ZTgUfKiPn3fJ8gHgOm+qyLAbKAYcBQYBcQ47P+beD3PmWYUiqm\nIqCdn3OmujF28ln2MPB6qe1mAre48+m4CcLP8UYCy0v97vwmCPf1G8Dv3PmuOAkjARCcBNnZZ99B\nQEYZ530N50u9rZ91ipPIYoFjQA+fdWOBdHd+NPCBO78WuKPk9+j+DfXxeZ8+8zlGDyDPnR8AbC11\n/keBie78POAPQLNS29wOfAmc5fVnpzZNdoupdhupqkmq2kFVf6qqee7yycAo91bPTThfxMcqONYu\nn/mjOF+64FxVbClZoaq5OP9ttwlg3/KUPm4xsM09bmtgm7usxJZS59xWKqYD7n5l2eYz3wG4xr3V\nkS0i2cAQoFXpnUSkhYhMEZEsETmE84XfLIDylXgLuMGdHwVMV9WjOFcjCcAynxg+dZf78xBOUvnK\nveVzu59tmuH8V7/FZ5nv720uMFREWuIkk3eAc0QkFWgMrPDZr/R7Wl+c5zcdgNalfne/BkpuF43G\nuTpc595Gusxd/jpOEp4iIjtE5BkRiSujrKaG2AM3cwpVXSQix3H+Ex/lTidWV/JwO3C+FAAQkQZA\nU5z/9qtjB3Cmz3EF51ZZFu7VgIjE+CSJ9ji3S0q089k3Eec2zY5yzudb7m04VxB3BhDnk+6+Z6nq\nfhEZCfyrjOP68z+gmYj0wkkUv3SX7wPygJ7qPEMql6ruAu4EEJEhwGciMk9VN/lstg/n6q4DsMZd\n1h73vVLVTSJyFLgP5zbhYRHZBYzBuUXnm5DLsg3nKqdrGXFuBG5wn3tcBbwrIk1V9QjOlcUf3IT0\nMbAe5wrVBIldQZiyvIbzRVaoqvN9lu/GuW8cqLeA20Skl4jUA/4PWKyqmdWMbypwqYj8wP1P8gGc\n2yNfAotxbr88JCJx7gPky4EpPvtfIiJDRKQuzrOXxaq6jcC8AVwuIiNEJNZ9gDpMRNr62bYh7q08\nEWkDPFhqfbm/T1UtBN4FnsVJYrPc5cXAq8BzItICQETaiMgIf8cRp/JBSXwHcRLTSVVp1alaOxV4\nQkQaivPg/n63vCXmAve6P8G5neb7uiJfAYdE5GERiXd/f2eISD83zhtFpLlbvmx3nyIRGS4iZ4rT\nPuQQTiKzqsBBZgnClOV14Az3p6/xQA/39sD0ig6iqrNxnhe8B+wEOgPXVzc4VV0P3Ijz0HsfTgK4\nXFWPq+px4ArgYnfdi8DNqrrO5xBvAY/j3FpKw6mdFei5twFX4twa2YvzX/GD+P88/QHnAX0O8BEw\nrdT6J4Hfur/PX5VxyreAHwL/cRNGiYdxHsYvcm9ffQZ0L+MY/YDFIpKL82D756qa4We7n+Ek183A\nfPfcE3zWz8VJevPKeF0uNwldDvQCMnDen3E4t6gALgJWu3E+j/NAPh9oiZMoD+E8/5jLyYnLBIGo\n2oBB5lRudck9OA8eN3odT00SkUnAdlW1dh/GlMOuIExZ7gaWRFtyMMYEzh5Sm1OISCZOjZeRHodi\njPGQ3WIyxhjjl91iMsYY45clCGOMMX5F9DOIZs2aaWpqapX2PXLkCA0aNKjZgDxiZQlP0VKWaCkH\nWFlKLFu2bJ+qltXq/oSIThCpqaksXbq04g39SE9PZ9iwYTUbkEesLOEpWsoSLeUAK0sJEdlS8VZ2\ni8kYY0wZLEEYY4zxyxKEMcYYvyxBGGOM8csShDHGGL8sQRhjjPHLEoQxxhi/IrodhDEmuK57eSHZ\n2XlESdMBU0l2BWGMMcYvSxDGGGP8sgRhjPFr+vIslm/NZv3BYs556nOmL8/yOiQTYpYgjDGnmL48\ni0enreJ4UTEAWdl5PDptlSWJWsYShDHmFM/OXE9eQdFJy/IKinh25nqPIjJesARhjDnFjuy8Si03\n0cmquRpjTigsKuadpdsQAX+jEbdOig99UMYzQbuCEJH6IvKViHwjIqtF5A/u8kkikiEiK9ypl7tc\nROQfIrJJRFaKSJ9gxWaMOVX6+j1c/PwX/Ob9b0lt2oB6dU7+eqhXJ4YHR3T3KDrjhWBeQRwDzlfV\nXBGJA+aLyCfuugdV9d1S218MdHWnAcC/3Z/GmCBav+swT3y8lnkb9tKhaQIv3ZjGiJ4pzFixg4fe\nXcnxomIEOLttY0b2buN1uCaEgpYgVFWBXPdlnDv5uWg94UrgNXe/RSKSJCKtVHVnsGI0pjbbe/gY\nf5u1gXeWbCWxXh0eu6wHNw3sQF33ymFk7za8/dVWsrOzGXR6e95cvIWdOXm0amy3mWoLUX83Gmvq\n4CKxwDKgC/CCqj4sIpOAQThXGLOBR1T1mIh8CDylqvPdfWcDD6vq0lLHHAOMAUhJSUmbMmVKlWLL\nzc0lMTGxagULM1aW8BSuZTlepMzMLOCjzQUUFMP57etwZee6JNaVU7Z9cnEeRUVFjO2VwMNf5HFh\nhziuP62uB1HXjHB9T6qiOmUZPnz4MlXtW+GGqhr0CUgC5gBnAK0AAeoBk4Hfudt8BAzx2Wc2kFbe\ncdPS0rSq5syZU+V9w42VJTyFW1mKior1/a+366D/+0w7PPyh3jl5iX6353CF+5WU4763v9aev/tU\ns48eD3KkwRNu70l1VKcswFIN4Ls7JNVcVTUbSAcuUtWdbozHgIlAf3ez7UA7n93aAjtCEZ8x0W5J\n5gF+9OICfvHOCpIT6zJlzEBeubkvnZoH/h/omHM7kXuskDcXBzTevYkCwazF1FxEktz5eOCHwDoR\naeUuE2Ak8K27ywfAzW5tpoFAjtrzB2OqZcv+I9z9xjKueWkhuw8d42/Xns0H9wxhYKemlT5Wz9aN\nObdbcyYuyCS/VCM6E52CWYupFTDZfQ4RA0xV1Q9F5HMRaY5zm2kFcJe7/cfAJcAm4ChwWxBjMyaq\n5Rwt4J+fb2TywkziYmO4/4Ju3Dm0E/F1Y6t13LvO7cSocYuZvjyL6/u3r5lgTdgKZi2mlUBvP8vP\nL2N7Be4JVjzG1AYFRcW8sWgLz8/eSE5eAdemteOBC7vRolH9Gjn+oM5NObNNY16Zt5lr+7YjJubU\nB9smelhLamOigKoya81unvpkHZv3HWFIl2b8+pLT6dG6UY2eR0QYe14n7n1rObPW7mZEz5Y1enwT\nXixBGBPhvs3K4c8frWHR5gN0bt6Aibf2Y1j35jiP+WreRT1b0j45gZfmfseFPVKCdh7jPUsQxkSo\nXTn5PDtzPdOWb6dJQl3+dGVPru/fnrjY4FZOrBMbw51DO/LYjNUsyTxI/47JQT2f8Y4lCGMizJFj\nhbw8bzOvzPuO4mKn+uk9w7vQqH5cyGK4pm87/v7ZRl6e+50liChmCcKYCFFUrLy3bDt/+d969hw+\nxmVnteLhi06jXXJCyGOpHxfLLYNT+dusDWzYfZhuKQ1DHoMJPhsPwpgIsGDTPi7753weem8lbZrE\n897dg/nXqD6eJIcSNw3sQHxcLK/M2+xZDCa47ArCmDC2aU8uT368ltnr9tC2STz/vKE3l53VKiwe\nDDdpUJfr+rXjzcVbeODCbtaJXxSyKwhjwtCBI8f53YxvGfH3eXyVcYBHLj6Nz+4/j8vPbh0WyaHE\n6CEdKVaYuCDT61BMENgVhDFh5FhhEZMWZPKvOZs4eryIUf3b84sfdqVpYj2vQ/OrXXICl53VircW\nb+We4V1oHB+6B+Um+CxBGBMGVJWPV+3iqU/Xsu1AHuef1oJfX3IaXVqE/8PfMed2YsaKHby5eAs/\nHdbF63BMDbIEYYzHvt56kCc+WsuyLQc5rWVDXh/dn6Fdm3sdVsB8O/G7/ZyO1I+rXn9PJnzYMwhj\nPLLtwFF+9vZyrnrxS7YeOMrTPz6Tj+4bGlHJocRd53Zi7+FjTF+e5XUopgbZFYQxQXDdywvJzs5j\n2LBT1x3KL+DFOd8xYUEGMQL3nd+Fsed1pkG9yP04Wid+0Sly/yKNiTCFRcVMWbKN52ZtYP+R41zV\npw0PjugeFdVDrRO/6GQJwpggU1XS1+/l/z5ey8Y9uQzomMykS3twZtvGXodWo6wTv+hTYYIQkXuB\nN1X1YAjiMSbiTV+exfKt2RwvKqb/E5/RJCGO9btzSW2awMs3pUXtl6d14hd9AnlI3RJYIiJTReQi\nica/bGNqyPTlWTw6bRXHi4oB2HP4GOt35/KjXq353y/PY0TPllGZHEpcndaO5AZ1eXnud16HYmpA\nhQlCVX8LdAXGA7cCG0Xk/0Skc3n7iUh9EflKRL4RkdUi8gd3eUcRWSwiG0XkHRGp6y6v577e5K5P\nrWbZjAm5Z2euJ8/PeM1fZR6kbp3orzQYXzeWWwenMnvdHjbsPux1OKaaAvqLdYcD3eVOhUAT4F0R\neaac3Y4B56vq2UAv4CIRGQg8DTynql2Bg8Bod/vRwEFV7QI8525nTETZkZ1XqeXRyDrxix4VJggR\nuU9ElgHPAAuAM1X1biAN+HFZ+6kj130Z504KnA+86y6fDIx05690X+Ou/4HdzjKRpmVj/2M/t06K\n/JpKgSrpxG/Giix25tSexBiNArmCaAZcpaojVPU/qloAoKrFwGXl7SgisSKyAtgDzAK+A7JVtdDd\nZDvQxp1vA2xzj10I5ABNK1keYzw1wM+D2fi4WB4c0d2DaLxjnfhFB3HuHlWwkUgfYAjOFcACVf26\nUicRSQLeB34HTHRvIyEi7YCPVfVMEVkNjFDV7e6674D+qrq/1LHGAGMAUlJS0qZMmVKZUE7Izc0l\nMTGxSvuGGytLeCgqVh6al0dcjLIvHwqLlab1Y/hxtzgGt47cTuyq+p689E0+K/YU8ddhCTSIC4+b\nAZH891VadcoyfPjwZarat6LtAqnm+hhwLTDNXTRRRP6jqn8ONBhVzRaRdGAgkCQiddyrhLbADnez\n7UA7YLuI1AEaAwf8HOsV4BWAvn376jB/TVUDkJ6eTlX3DTdWlvDw0cqd7M//mlduSmP8/Ayys7OZ\n+fDFXodVbVV9T5p3y+HSf8xnS1y7sOnEL5L/vkoLRVkCucU0Cuinqo+r6uM4X/I/qWgnEWnuXjkg\nIvHAD4G1wBzganezW4AZ7vwH7mvc9Z9rIJc3xoSJcfM3k9o0gR+cnsI7Ywfx6IDa89zBn56tGzO0\nazMmLsgk30/NLhP+AkkQmYDvk7d6OM8SKtIKmCMiK4ElwCxV/RB4GLhfRDbhPGMY724/HmjqLr8f\neCSgEhgTBpZtOcjyrdncPqQjsdYP0Ql3n9fZOvGLYIF0tXEMWC0is3CeQVwAzBeRfwCo6n3+dlLV\nlUBvP8s3A/39LM8Hrgk8dGPCx/j5m2kcH8fVaW29DiWsWCd+kS2QBPG+O5VID04oxkSmbQeO8um3\nuxh7XmcS6lr3Zr6sE7/IVuFfs6pOrmgbY2qziQsyiRHhlkGpXocSlqwTv8gVSEO5y0RkuYgcEJFD\nInJYRA6FIjhjwt2h/ALeWbKVy89uXWYjudqupBO/5VuzWZJpfX5GkkAeUv8dp3ZRU1VtpKoNVbVR\nkOMyJiK889U2jhwvYvSQjl6HEtasE7/IFEiC2AZ8a1VOjTlZYVExExdkMLBTMme0ia6xHWqadeIX\nmQJJEA8BH4vIoyJyf8kU7MCMCXeffLuLHTn53DGkk9ehRATrxC/yBJIgngCO4rSFaOgzGVNrqSrj\nvthMx2YNOP+0Fl6HExGsE7/IE0idvGRVvTDokRgTQZZtOcg323P408gzrG5/JYwe0pHXF21hwvwM\nfnNpD6/DMRUI5AriMxGxBGGMj3FfZJCUEMeP+7SpeGNzQrvkBC47qxVvLd5KTl6B1+GYCgSSIO4B\nPhWRPKvmagxs2X+EmWt28ZMB7a1hXBWMObcTR44X8ebiLV6HYioQyJCjDVU1RlXjrZqrMU7DuDox\nws3WMK5KrBO/yBFIQ7lz/U2hCM6YcJOTV8DUpdu4/OzWpDSyhnFVZZ34RYZAro8f9Jmvj9PR3jKc\noUONqVWmfLWVo9YwrtqsE7/IEMgtpst9pguAM4DdwQ/NmPBSUFTMpC8zGdy5KT1bW8O46ijpxG/z\nviPMWmtfJ+EqkIfUpW3HSRLG1Cofr9rJzpx87hhqVw81wbcTP+uoITwFMuToP3HGgQAnofQCvglm\nUMaEG1Vl/PwMOjVvwLBu1jCuJpR04vfYjNUsyTxI/47JXodkSgnkCmIpzjOHZcBC4GFVvTGoURkT\nZpZkHmTl9hxGD+lo98trkHXiF95sPAhjAjDui800SYjjqt42YlxNiq8byy2DUnnusw1s2H2YbinW\ni084KfMKQkTmiMjnZUyzKzqwiLRzj7FWRFaLyM/d5b8XkSwRWeFOl/js86iIbBKR9SIyomaKaEz1\nZLoPUm8c2IH4urFehxN1bh5knfiFq/KuIH7lZ9lAnN5d9wRw7ELgAVX9WkQaAsvcca0BnlPVv/hu\nLCI9gOuBnkBrnC4+uqmqtaQxnpq4IIO4mBhuGtTB61CiUkknfm8u3sIDF3ajVeN4r0MyrjKvIFR1\nWckEJAJP43yB36Wq/So6sKruVNWv3fnDwFqgvI5rrgSmqOoxVc0ANuG0uTDGMzlHC5i6dDtX9GpN\ni4bWMC5YRg/pSLHChPkZXodifJT7kFpERojIfOAx4AlVHaqqn1T2JCKSCvQGFruL7hWRlSIyQUSa\nuMva4AxOVGI75ScUY4Lura+2kldgDeOCzTrxC09SVv1jEVkCNAeexam9dJKSq4MKTyCSCMzFSTDT\nRCQF2IdTdfZPQCtVvV1EXgAWquob7n7jgY9V9b1SxxsDjAFISUlJmzJlSkAFLS03N5fExMQq7Rtu\nrCzBUVis/GpuHm0ShQf7Vf62RziVpTpCVY4th4p4/Mt8ru4Wx2Wd6gblHNHynkD1yjJ8+PBlqtq3\nwg1V1e8EpANzypg+L2u/UseIA2YC95exPhVnOFOAR4FHfdbNBAaVd/y0tDStqjlz5lR533BjZQmO\n97/erh0e/lA/X7e7SvuHU1mqI5TluHHcIu3751mad7wwKMePlvdEtXplAZZqAN/h5d1iukFVh5cx\nVdgPk4gIMB5Yq6p/81neymezHwHfuvMfANeLSD0R6Qh0Bb6q6DzGBIOqMm7+Zrq0SOS8rs29DqfW\nsE78wkt5tZhKng+kA58C81W1sBLHPge4CVglIivcZb8GbhCRXji3mDKBsQCqulpEpgJrcGpA3aNW\ng8l4ZHHGAb7NOsSTV51pDeNCyDrxCy9lJghVvVhE6gPDcP7T/4uIbMVJFp+q6tbyDqyq8wF/7+7H\n5ezzBM4Y2MZ4atwXGSQ3qMuPels9iVAq6cTv3reWM2vtbkb0bOl1SLVaubWYVDVfVT9V1Z+r80Dj\nAZyk8i8Rsds/Jipt3pvL7HVOw7j6cdYwLtSsE7/wUaneXFU1Q1VfVNUrgCFBiskYT01ckOk0jBto\nDeO8UNKJ3/Kt2SzJPOh1OLVaICPKXSUiG0Ukx3dMalU9HooAjQml7KPH+c+ybYzs3ZrmDet5HU6t\nZZ34hYdAriCeAa5Q1cZqY1KbKPfm4q3kFxQzekgnr0Op1Uo68Zu9bg8bdh/2OpxaK5AEsVtV1wY9\nEmM8drywmMlfZjK0azO6t7ReRb1mnfh5L6DxIETkHRG5wb3ddJWIXBX0yIwJsQ9X7mDP4WPcMdSu\nHsJBSSd+M1ZksTMnz+twaqVAEkQj4ChwIXC5O10WzKCMCTVVZdwXGXRtkci5XZt5HY5xWSd+3gpk\nwKDbQhGIMV5auHk/a3Ye4ukfn4nTCYAJB76d+N17flcax8d5HVKtEkgtprYi8r6I7BGR3SLynojY\nsFomqoz/IoOmDepyZS9rGBduxpzbiSPHi3hz8RavQ6l1ArnFNBGnn6TWON1v/9ddZkxU2LQnl9nr\n9nDTIGsYF456tm7M0K7NmLggk/wC630nlAJJEM1VdaKqFrrTJJxuwI2JChMWZFC3Tgw3WsO4sHWX\ndeLniUASxD4RuVFEYt3pRmB/sAMzJhQOHDnOe8u2c1XvNjRLtIZx4WqwTyd+xcXW/UaoBJIgbgeu\nBXYBO4Gr3WXGRLw3F23hWGExt9uIcWGtpBO/zfuO8L81u70Op9YIpBbTVuCKEMRiTEgdKyxi8sIt\nnNetOd1SrGFcuPPtxG9EzxSrbRYCFSYIEWkO3Ikz+tuJ7VXVriJMRPtgxQ725R7jjqF29RAJSjrx\ne2zGapZkHqR/x2SvQ4p6gdximgE0Bj4DPvKZjIlYqsr4+Rl0T2nIkC7WMC5SWCd+oVXhFQSQoKoP\nBz0SY0Jowab9rNt1mGeuPstuVUSQkk78nvtsAxt2H7Zbg0EWyBXEhyJySdAjMSaExs3fTLPEelzZ\nq7XXoZhKsk78QqfMBFEy7gPwc5wkkec7HkRFBxaRdiIyR0TWishqEfm5uzxZRGa5Y0zMcse9Rhz/\nEJFNIrJSRPrUVCGN8bVx92HS1+/l5kEdqFfHGsZFGuvEL3TKTBAl4z64P2NUNb6S40EUAg+o6unA\nQOAeEekBPALMVtWuwGz3NcDFQFd3GgP8uxrlMqZMExZkUK9ODD8Z0N7rUEwVWSd+oRFIX0yzA1lW\nmqruVNWv3fnDwFqcrjquBCa7m00GRrrzVwKvqWMRkCQirQIqhTEB2p97jPe+zuKqPm1pag3jIpZv\nJ345eQVehxO1yrvFVF9EmgLNRKSJe2soWURScfplCpi7T29gMZCiqjvBSSJAC3ezNsA2n922u8uM\nqTFvLNrK8cJiRg9J9ToUU03WiV/wlVeLaSzwC5xk8LXP8kPAC4GeQEQSgfeAX6jqoXJqjPhbcUqb\nehEZg3MLipSUFNLT0wMN5SS5ublV3jfcWFkCc7xIGT/vKGc1j2X7mmVsXxOU05wQLe9LOJfjjKax\nvDxnA12KtlE3tuLaaOFclsoKSVlUtdwJ+FlF25SzbxwwE7jfZ9l6oJU73wpY786/DNzgb7uyprS0\nNK2qOXPmVHnfcGNlCcw7X23VDg9/qPM37g3aOXxFy/sSzuWYv3Gvdnj4Q3178ZaAtg/nslRWdcoC\nLNUAvsMDqeZ6i4j8VESSKpN4xLlUGA+sVdW/+az6ALil5Ng4DfFKlt/s1mYaCOSoeyvKmOpSVcbN\n38xpLRsyuHNTr8MxNcQ68QuuQBLE9TjPApaKyBQRGSGBtSw6B7gJOF9EVrjTJcBTwAUishG4wH0N\n8DGwGdgEvAr8tJJlMaZMX2zcx4bdudwxtJM1jIsi1olfcAXSWd8m4Dci8hjOWNQTgGIRmQA8r6oH\nythvPv6fKwD8wM/2CtwTaODGVMa4+Rk0b1iPy8+2inHRxjrxC55AriAQkbOAvwLP4jxwvhrnYfXn\nwQvNmJqxftdh5m3Yyy3WMC4qlXTit2JbNksyD3odTlQJpB3EMuA5YAlwlqrep6qLVfWvOLeEjAlr\nE+ZnUD8uhlEDbMS4aGWd+AVHIFcQV6vqD1T1LVU95rtCVa8KUlzG1Ii9h4/x/oosftynLckN6nod\njgmSkk78Zq/bw4bdh70OJ2qU11BugIh8A6wSkYVuNxnGRJQ3Fm3huI0YVytYJ341r7wriBeAXwFN\ngb/h3GYyJmLkFxTxxqIt/OC0FnRunuh1OCbIrBO/mldegohR1VmqekxV/wM0D1VQxtSE6cuz2H/k\nOKNtxLhawzrxq1nlVXNNEpGrynqtqtOCF5Yx1eM0jMugR6tGDOpkDeNqC99O/O49vyuN4+O8Dimi\nlXcFMRe43GfyfX1Z8EMzpurmbtjLpj253DG0o9WLr2WsE7+aU+YVhKreFspAjKlJ4+dn0KJhPS47\ny0aMq216tm7M0K7NmLggk9ux6/pzAAAcNElEQVTP6Uj9OGv7UlUBNZQzJpKs23WILzbu45bBqdSt\nY3/itdFd53Vm7+FjTF+e5XUoEc0+PSbqjP8ig/i4WBsxrhYb3LkpZ7RpZJ34VZMlCBNV9hzOZ8aK\nHVyd1pakBGsYV1uJCHed19k68aumQLraiBOR+0TkXXf6mYhY1QATlt5YuIWC4mJuOyfV61CMx3w7\n8XP6AjWVFcgVxL+BNOBFd+rjLjMmrOQXFPH6oi384LQUOlnDuFrPOvGrvgq7+wb6qerZPq8/d7vg\nMCasTPs6i4NHC7jDGsYZ19Vp7Xjus428PPc7/vq/9WRn5zFsmNdRRY5AriCKRKRzyQsR6QQUBS8k\nYyqvuFgZP38zZ7RpxICOyV6HY8KEbyd+R48Xeh1OxAkkQTwIzBGRdBGZizMGxAPBDcuYypm7YS/f\n7T3CHUNsxDhzspsHdSAuRli94xDrDxZzzlOfW/XXAAUyotxsEekKdMcZIW5d6W6/jfHauPmbadmo\nPpecaSPGmZPN3bCXYqCktmtWdh6PTlsFwMjebbwLLAKU1933+e7Pq4BLgS5AZ+DSUn00lbX/BBHZ\nIyLf+iz7vYhklRqjumTdoyKySUTWi8iI6hTK1C5rdhxiwab91jDO+PXszPUUlWoLkVdQxLMz13sU\nUeQo7wriPJzbSZf7WadARZ31TQL+BbxWavlzqvoX3wXuWBPXAz2B1sBnItJNVe1Zh6nQ+PlOw7hR\n/a1hnDnVjmz/XX+Xtdx8r7y+mB53Z/+oqif1nSsiFVYTUdV5IpIaYBxXAlPcW1cZIrIJ6A8sDHB/\nU0vtOZTPB99kMap/exonWPMcc6rWSfFk+UkGrRrX9yCayBJINdf3cNo++HoXp21EVdwrIjcDS4EH\nVPUg0AZY5LPNdnfZKURkDDAGICUlhfT09CoFkZubW+V9w01tLst7G45TWKT0qLMn7H4H0fK+RHo5\nLm1fxKRDcLz45OXFBfl8+L85JNaNzEoNoXhfykwQInIazi2fxqWeOTQCqpp6/w38CecW1Z+AvwK3\n4zz8Ls1v00dVfQV4BaBv3746rIqVmtPT06nqvuGmtpYl73gRv5g3mwt6pHDdpX2DG1gVRMv7Eunl\nGAb0WJ7FQ++u5HhRMW2S4jn/9Ba889U2/rpSmHRbPzo0beB1mJUWivelvCd63XHGfUji5HEh+gB3\nVuVkqrpbVYtUtRh4Fec2EjhXDO18Nm0L7KjKOUzt8d7X28k+WsAdQzt5HYoJcyN7t6F3+yS6N4lh\nwSPn86crz+DNOwdw8OhxfvTilyzbYi2t/SkzQajqDHdMiMtU9Taf6T5V/bIqJxMR3zqIPwJKajh9\nAFwvIvXc5xtdga+qcg5TOxQXKxPmZ3BW28b0S23idTgmAvVLTWba3YNpWL8ON7y6iI9W7vQ6pLAT\nyDOI5SJyD87tphO3llT19vJ2EpG3ca7umonIduBxYJiI9MK5fZQJjHWPtVpEpgJrgELgHqvBZMoz\nZ/0eNu87wvPX97KGcSYg74wddMo9+07NE5l292DGvL6Me976mu0HT2PMudbYskQgCeJ1YB0wAvgj\n8BNgbUU7qeoNfhaPL2f7J4AnAojHGMZ9kUGrxtYwzlRf08R6vHnHAB6Y+g1PfrKOrQeO8ocrelIn\n1trUBPIb6KKqjwFHVHUyTqO5M4MbljFl+zYrh4Wb93Pr4FTi7ENsakD9uFj+eUNv7jqvM28u3sod\nry0l95j13RTIp6vA/ZktImcAjYHUoEVkTAUmzM8goW4s11vDOFODYmKERy4+jf/70Zl8sXEf17y0\nkF05+V6H5alAEsQrItIE+C3Ow+Q1wNNBjcqYMuzKyeeDb3Zwbd92NI63hnGm5o0a0J7xt/Rl6/4j\njHxhAWt2HPI6JM+UmyBEJAY4pKoHVXWeqnZS1Raq+nKI4jPmJK8tzKRIldvPsTEfTPAM696C/9w1\nGIBrXvqSuRv2ehyRN8pNEG57hXtDFIsx5Tp6vJA3F29lRI+WtG+a4HU4Jsr1aN2I9+8ZTPumDbh9\n0hLeWrzV65BCLpBbTLNE5Fci0k5EkkumoEdmTCnvLdtOTp6NGGdCp1XjeP5z1yCGdGnGr99fxVOf\nrKO4uPaMbx1INdeS9g73+CxTwJqvmpBxRozL4Ox2SaR1sIZxJnQS69Vh/C19+d0Hq3lp7ndsO3iU\nv15zNvXjYr0OLegCGTDI/l0znpu9bg+Z+4/yzwu7WyMmE3J1YmN4YuQZdEhO4MlP1rE7J59Xbu5L\ncoO6XocWVFaJ3ESEcV9spk1SPBef0dLrUEwtJSKMPa8zL4zqw8qsHK56cQEZ+454HVZQWYIwYW/V\n9hwWZxzg1sGp1rrVeO7Ss1rx9p0DyMkr4KoXF7A084DXIQWNfdpM2Bs/fzMN6sZyXf92FW9sTAik\ndUjm/Z+eQ1JCXUaNW8x/v4nOzqcrTBAi0sfP1FlEAnnAbUy17MzJ48OVO7muX3sa1beGcSZ8pDZr\nwLS7B3N228b87O3l/Dv9O1Sjq4ZTIFcQL+KM9vYKzhgOC4EpwAYRuTCIsRnD5C+3UKzKbeekeh2K\nMado0qAur48ewOVnt+bpT9fx6/dXUVBUXPGOESKQBJEJ9FbVvqqaBvTGGcfhh8AzQYzN1HJHjhXy\n1uItXHxGK9olW8M4E57qx8Xy/HW9uGd4Z97+ahujJy/lcH5BxTtGgEASxGmqurrkhaquwUkYm4MX\nljHw7rLtHMovZLQ1jDNhLiZGeHDEaTz94zNZsMnp6G9nTp7XYVVbIAlivYj8W0TOc6cXcW4v1eP7\nnl6NqVFFxcqEBRn0aZ9En/bWMM5Ehuv6tWfirf3YfjCPkS8sYPWOHK9DqpZAEsStwCbgF8Avgc3u\nsgJgeLACM7XbZ2t3s2X/URtv2kScc7s15927BxErwrUvLWTOuj1eh1RlFSYIVc1T1b+q6o9UdaSq\n/kVVj6pqsarmhiJIU/uM/yKDtk3iubBHitehGFNpp7VsxPv3nENqswaMnryENxZt8TqkKgmkmus5\nIjJLRDaIyOaSKYD9JojIHhH51mdZsnusje7PJu5yEZF/iMgmEVkpIn2qVywTyb7Zls1XmQe47ZyO\n1jDORKyURvWZOnYQw7q34LfTv+XJj9dGXEd/gXz6xgN/A4YA/XymikwCLiq17BFgtqp2BWa7rwEu\nBrq60xjg3wEc30Sp8fMzaFivDtf2bet1KMZUS4N6dXjlpjRuGtiBl+dt5t63vya/oMjrsAIWSILI\nUdVPVHWPqu4vmSraSVXnAaXboF8JTHbnJwMjfZa/po5FQJKI2Gj0tdD+vGI+WrWT6/u3o6E1jDNR\noE5sDH+8sie/vfR0Pvl2Fze8uoj9uce8DisgUlHLPxF5CogFpgEnSqWqX1d4cJFU4ENVPcN9na2q\nST7rD6pqExH5EHhKVee7y2cDD6vqUj/HHINzlUFKSkralClTKgrDr9zcXBITE6u0b7iJprK8sSqX\nz3cIz5wbT7P4yL69FC3vS7SUA7wvy9Jdhby88hhJ9YT70+rTKrHqf+PVKcvw4cOXqWrfirYLpLuM\nAe5P34MpcH5VAiuDv/6b/WYuVX0Fp1U3ffv21WHDhlXphOnp6VR133ATDWWZvjyLpz9dx84cIT4u\nhjotuzOsdxuvw6qWaHhfIHrKAd6XZRhw/uCD3Dl5KU8tK+TVm/vSv2PVxl8LRVkCqcU03M9U1eSw\nu+TWkfuzpP7XdsC3J7a2QHT2fmVOMX15Fo9OW8XOnHwA8gqKeXTaKqYvz/I4MmNqXp/2TXj/p+fQ\nNLEuN45bzIwV4ft3XmaCEJEb3Z/3+5uqeL4PgFvc+VuAGT7Lb3ZrMw3Eee6xs4rnMBHm2ZnrySv1\n4C6voIhnZ673KCJjgqt90wSm3T2Y3u2T+PmUFfzr841h2dFfeVcQDdyfDcuYyiUib+N07NddRLaL\nyGjgKeACEdkIXOC+BvgYpwHeJpwOAX9a+aKYSJWV7b9Lgh1lLDcmGiQl1OW10f0Z2as1f/nfBh5+\nb2XYdfRX5jMIVX3Z/fmHqhxYVW8oY9UP/GyrnDzmtakFVJUX5mwqc33rpPgQRmNM6NWrE8tz1/Wi\nfXIC//h8Eztz8nnhJ33Cpmv7Ch9Si0hH4GdAqu/2qnpF8MIy0e5YYRGPvLeK95dn0ad9Emt2HiK/\n4Pv/nuLjYnlwRHcPIzQmNESE+y/sTrvkBB6dtopr/r2QCbf1o00Y/IMUSC2m6TiN5f4LhNf1j4lI\n+3KPMfb1ZSzbcpAHLujGved3YcaKHTw7cz1Z2Xm0SYrnwRHdGRnhtZiMqYxr+rajdVI8d72+jJEv\nLGDirf04o01jT2MKJEHkq+o/gh6JqRU27D7M7ZOWsPfwMf41qjeXndUagJG92zCydxvPqyEa46Vz\nujTjvZ8O5raJS7j25YX884be/OB07/ojC6SVxvMi8riIDPIddjTokZmok75+D1e9+CXHCot5Z+yg\nE8nBGPO9bikNef+ewXRunsidry3l9YWZnsUSyBXEmcBNOA3jSm4x1XRDORPFVJXJX2byxw/X0L1l\nI8bf0tceQBtTjhYN6/PO2IHc9/ZyHpuxmi37j/LrS04nJsZfm+LgCSRB/AjopKrHgx2MiT4FRcX8\n4b+reWPRVn54egrPX9+LBvUC+bMzpnZLqFuHl2/qy58+XMO4+RlsP5jHc9f1Ir5ubMhiCOST+g2Q\nxPetno0JSE5eAfe+9TVfbNzH2HM78dBFpxEb4v+AjIlksTHC76/oSfvkBP700Rquf3URV6e14aX0\nzU6FjkWfB7VCRyAJIgVYJyJLOLmzPqvmasqUue8IoycvYeuBozxz9Vlc27ddxTsZY/y6fUhH2jSJ\n5543l7FyW/aJjuqysvN4dNoqgKAkiUASxOM1flYT1RZt3s9dbywD4PXRAxjYqanHERkT+Ub0bEmT\nhHrsLdVVeEm3NJ4kCFWdW+NnNVFr6tJt/Ob9VbRLTmDCLf1Ibdag4p2MMQHZV8Y4EsHqlqbMBCEi\n81V1iIgc5uSutwWnd4xGQYnIRKTiYuXpmet4ee5mhnRpxguj+tA4ITy6CzAmWrROivfbd1mwagVW\n2FmfqjZU1UY+U0NLDsbXkWOFjH1jGS/P3cxPBrRn4m39LDkYEwQPjuhOfNzJtZiC2S1NebeYwq/v\nWRN2dubkMXrSUtbtOsTjl/fg1sGpiFhNJWOCoeQ5Q6i6pSkvQbQob9wHVf1bEOIxEeSbbdnc+dpS\njh4vYvyt/RjevYXXIRkT9ULZLU15CSIWSMT/cKCmlvto5U7un7qCZon1eO/uAXRvWeEQIcaYCFNe\ngtipqn8MWSQmIqgq//p8E3+dtYG0Dk14+aY0miXW8zosY0wQlJcg7MrBnCS/oIhH3lvJ9BU7GNmr\nNU/9+Czqx4Wu2b8xJrTKSxCnjPxmai9/YzjYw2hjolt5Q44eCNZJRSQTOAwUAYWq2ldEkoF3cEau\nywSuVdWDwYrBBG79rsOMnuyM4fDCqD5celYrr0MyxoRAIONBBMtwVe2lqn3d148As1W1KzDbfW08\nNmf9Hn78b2cMh6ljB1lyMKYW8TJBlHYlMNmdnwyM9DCWWk9Vmbggg9GTltA+OYEP7j2Hs9sleR2W\nMSaERDX07eFEJAM4iNMY72VVfUVEslU1yWebg6raxM++Y4AxACkpKWlTpkypUgy5ubkkJiZWad9w\nU9NlKSxW3lx7nDnbCunTIpYxZ9Wjfp3QPG+w9yX8REs5wMpSYvjw4ct87t6UTVVDPgGt3Z8tcMab\nOBfILrXNwYqOk5aWplU1Z86cKu8bbmqyLNlHjuuoVxdqh4c/1Cc/XqtFRcU1duxA2PsSfqKlHKpW\nlhLAUg3gu9qTob1UdYf7c4+IvA/0B3aLSCtV3SkirbABikIuc98Rbp+8hG0HjvLs1WdxjY3hYEyt\nFvJnECLSQEQalswDFwLfAh8At7ib3QLMCHVstdmizfsZ+eICDhw5zhujB1hyMMZ4cgWRArzv1qGv\nA7ylqp+6I9ZNFZHRwFbgGg9iq5WmLtnGb6avon1yAhNu7UeHpjaGgzHGgwShqpuBs/0s3481zgup\nomLlmU/X8fK8zQzt2ox/jepD43jrptsY4/DkGYTx3pFjhfzinRXMWrObGwe25/HLexIXG061no0x\nXrMEUQvtyM5j9OSlrN91iN9f3oNbbAwHY4wfliBqmW+2ZXPHa0vJszEcjDEVsARRi3y4cgcPTP2G\n5g3r8cZoG8PBGFM+SxC1gNoYDsaYKrAEEeV8x3D4Ue82PHnVmTaGgzEmIJYgoti+3GOMeW0pX2/N\n5lcXduOe4TaGgzEmcJYgotT6XYe5fdIS9h85xos/6cMlZ1o33caYyrEEEYXmrNvDz95eTkLdWKaO\nHcRZba2bbmNM5VmCiCKqyoQFmTzx0RpOb9WIcbf0pVXjeK/DMsZEKEsQEWz68iyenbmerOw8Wi+a\nTYemCSz87gAX9kjh79f3IqGuvb3GmKqzb5AINX15Fo9OW0VeQREAO7Lz2ZGdzw9Ob8FLN6YRE2MP\no40x1WOd70SoZz5ddyI5+Fq387AlB2NMjbAriDBVXKzsPpzP1v1H2XrgKNsOHGXbwTy2HnBe7z18\nzO9+O7LzQhypMSZaWYLw0OH8ArYdyDuRAEq+/LcdOMr2g3kcLyo+sW2MQKvG8bRPTuD87i345Nud\nHMovPOWYrZPsobQxpmZYggiiwqJidubkn/TF7/vz4NGCk7ZvHB9H++QETmvVkAt6ptA+OeHE1Dop\n/qTuuAd1bnrSMwiA+LhYHhzRPWTlM8ZEN0sQ1aCqZB8tYNvBo36SQB5Z2XkUFeuJ7evECG2bxNMu\nOYFLzmxF++QE2rkJoF2TBBonBD5Yz8jebQBO1GJqkxTPgyO6n1hujDHVVesShG/V0DaLPq/wS/VY\nYRFZ7r3/bQfznATg81zg8LGTb/M0S6xLu+QEerVL4oqzW3+fBJom0LJRfWJr8AHyyN5tGNm7Denp\n6QwbNqzGjmuMMRCGCUJELgKeB2KBcar6VE0du3TV0KzsPB6dtpJD+QX0bN34lOcA2w4cZeehfPT7\niwDq1Yk58V9//47JtEtOoF2TeNo3da4CGtQLu1+pMcZUSVh9m4lILPACcAGwHVgiIh+o6pqaOP6z\nM9efUjU0r6CY381YfdKylo3q0z45gYGdm570HKB9cgLNEutZNVJjTK0QVgkC6A9sUtXNACIyBbgS\nqJEEUV4V0Im39qNdcgJtm8Rbd9jGGEP4JYg2wDaf19uBAb4biMgYYAxASkoK6enpAR88ub6wP19P\nWd60viC71rB9l3PCSJObm1up30M4s7KEn2gpB1hZKivcEoS/ezcnfaOr6ivAKwB9+/bVyjycfaxx\nlt+qoY9deSbDIrj2TzQ9pLayhJ9oKQdYWSor3BLEdqCdz+u2wI6aOrhVDTXGmMCFW4JYAnQVkY5A\nFnA9MKomT2BVQ40xJjBhlSBUtVBE7gVm4lRznaCqqyvYzRhjTBCEVYIAUNWPgY+9jsMYY2o76+7b\nGGOMX5YgjDHG+GUJwhhjjF+iemrDsUghIjnARp9FjYGcMl6XzJf8bAbsq8bpS5+rMtv4Wx5I7GXN\nV6cs1SlHWesisSyVLUfp16X/viByyhLM96S8OAPZJpzKEg6flZr6++qqqo0r3EpVI3YCXgn0dcm8\nz8+lNXnuymzjb3kgsZdTpiqXpTrliKayVLYcFf19RVJZgvmeRFNZwuGzEsq/L1WN+FtM/63E6/+W\nsU1Nnbsy2/hbHkjs5c1XVXXKUda6SCxLZctR+rX9fZUtWsoSDp+VUL4nkX2LqTpEZKmq9vU6jppg\nZQlP0VKWaCkHWFkqK9KvIKrjFa8DqEFWlvAULWWJlnKAlaVSau0VhDHGmPLV5isIY4wx5bAEYYwx\nxi9LEMYYY/yyBOESkU4iMl5E3vU6luoSkZEi8qqIzBCRC72Op6pE5HQReUlE3hWRu72Op7pEpIGI\nLBORy7yOpTpEZJiIfOG+N8O8jqc6RCRGRJ4QkX+KyC1ex1NVIjLUfT/GiciXNXXcqE4QIjJBRPaI\nyLelll8kIutFZJOIPAKgqptVdbQ3kVaskmWZrqp3ArcC13kQbpkqWY61qnoXcC0QdlUTK1MW18PA\n1NBGGZhKlkWBXKA+YThKbyXLciXOUMcFhFlZKvlZ+cL9rHwITK6xIKrTQjLcJ+BcoA/wrc+yWOA7\noBNQF/gG6OGz/l2v467BsvwV6ON17NUpB3AF8CUwyuvYq1MW4Ic4A2DdClzmdezVLEuMuz4FeNPr\n2KtZlkeAse42YfXZr+JnfirQqKZiiOorCFWdBxwotbg/sEmdK4bjwBSc/yLCWmXKIo6ngU9U9etQ\nx1qeyr4nqvqBqg4GfhLaSCtWybIMBwbijJB4p4iE1WevMmVR1WJ3/UGgXgjDDEgl35ftOOUAKCKM\nVPazIiLtgRxVPVRTMYTdgEEh0AbY5vN6OzBARJoCTwC9ReRRVX3Sk+gqx29ZgJ/h/MfaWES6qOpL\nXgRXCWW9J8OAq3C+hCJlECm/ZVHVewFE5FZgn8+XbDgr6325ChgBJAH/8iKwKijrs/I88E8RGQrM\n8yKwSiqrHACjgYk1ebLamCDEzzJV1f3AXaEOpprKKss/gH+EOphqKKsc6UB6aEOpNr9lOTGjOil0\noVRbWe/LNGBaqIOpprLKchTnizVSlPn3paqP1/TJwuoyN0S2A+18XrcFdngUS3VFS1mipRxgZQlX\n0VKWkJajNiaIJUBXEekoInVxHhx+4HFMVRUtZYmWcoCVJVxFS1lCWw6vn9QHuRbA28BOvq/CNtpd\nfgmwAac2wG+8jrM2lSVaymFlCd8pWsoSDuWwzvqMMcb4VRtvMRljjAmAJQhjjDF+WYIwxhjjlyUI\nY4wxflmCMMYY45clCGOMMX5ZgjC1lojkBuGYqSIyqpL7fCwiSTUdizHVZQnCmJqVitNja8BU9RJV\nzQ5OOMZUnSUIU+u5I6SluyPXrRORN0VE3HWZIvK0iHzlTl3c5ZNE5GqfY5RcjTwFDBWRFSLyy1Ln\naSUi89x137o9iJaco5mI3OWuWyEiGSIyx11/oYgsFJGvReQ/IpIYit+LMZYgjHH0Bn6BM4hMJ+Ac\nn3WHVLU/TtfWf6/gOI8AX6hqL1V9rtS6UcBMVe0FnA2s8F2pqi+56/rhdK3wNxFpBvwW+KGq9gGW\nAvdXpYDGVFZt7O7bGH++UtXtACKyAudW0Xx33ds+P0t/6VfGEmCCiMQB01V1RRnbPQ98rqr/FWf8\n6h7AAveipi6wsBoxGBMwSxDGOI75zBdx8mdD/cwX4l6Bu7ej6lZ0AlWdJyLnApcCr4vIs6r6mu82\n7oBCHYB7SxYBs1T1hsCLYkzNsFtMxlTsOp+fJf+9ZwJp7vyVQJw7fxho6O8gItIB2KOqrwLjccYb\n9l2fBvwKuFG/H3FuEXCOz7OPBBHpVt0CGRMIu4IwpmL1RGQxzj9UJf/JvwrMEJGvgNnAEXf5SqBQ\nRL4BJpV6DjEMeFBECoBc4OZS57kXSAbmuLeTlqrqHe5VxdsiUjL+829xuns2Jqisu29jyiEimUBf\nVd3ndSzGhJrdYjLGGOOXXUEYY4zxy64gjDHG+GUJwhhjjF+WIIwxxvhlCcIYY4xfliCMMcb4ZQnC\nGGOMX/8POywGG5WOrgwAAAAASUVORK5CYII=\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "ratio, err = rat_err(tla, tlea, tva, tvea)\n", "\n", "plt.errorbar(sizes, ratio, yerr=err, fmt = 'o-')\n", "plt.xscale('log')\n", "\n", "plt.grid(True)\n", "plt.xlabel(\"Input size\")\n", "plt.ylabel(\"Timing ratio Python/Numpy\")\n", "plt.title(\"Python loop relative slowness\");" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "A few observations:\n", "\n", "- At very small sizes the difference isn't very dramatic. The overhead of the Python function call and basic arithmetic dominates.\n", "\n", "- For input sizes ranging from the thousands to the million elements, the vectorized version can be ~ 100x to ~300x fastr than the pure Python loop. This is not uncommon in real-world applications.\n", "\n", "- At very large sizes (at a million) the performance profile changes.\n", "\n", "- There's significant uncertainty for one of the values but little for the others, and even accounting for that uncertainty, the conclusion doesn't fundamentally change.\n", "\n", "In general it's extremely hard to predict exact performance ratios, which often vary in hard to understand ways with problem size, as hardware, cache and algorithm specifics intersect. But the results of this example aren't suprising: whether exactly 100x or 350x slower, the fact remains that the vectorized solution vastly outperforms a pure Python loop. This is a key aspect of Numpy as a key piece of high-performance Python.\n", "\n", "For an in-depth look at the \"vectorized mindset\" for Numpy, Nicolas Rougier has written the excellent (and free!) book [\"From Python to Numpy\"](http://www.labri.fr/perso/nrougier/from-python-to-numpy)." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Numba\n", "\n", "Let's have a look at how [Numba](https://numba.pydata.org) can help in this type of loop-intensive code." ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "680.0000000000001" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "from numba import jit\n", "\n", "trapzn = jit(trapzl)\n", "\n", "x = np.linspace(a, b, 200)\n", "y = f(x)\n", "\n", "trapzn(y, x)" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "34.6 µs ± 1.39 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)\n", "22.5 µs ± 732 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)\n" ] } ], "source": [ "x = np.linspace(a, b, 20_000)\n", "y = f(x)\n", "\n", "%timeit trapz(y, x)\n", "%timeit trapzn(y, x)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "That's pretty impressive! By simply *compiling* our slow, Python code, we get it to beat the Numpy version at a size of 200.\n", "\n", "Let's collect Numba performance numbers similarly to the code above, and make a plot comparing Numba to Numpy at the same sizes:" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "413 ns ± 4.92 ns per loop (mean ± std. dev. of 7 runs, 1000000 loops each)\n", "512 ns ± 11.7 ns per loop (mean ± std. dev. of 7 runs, 1000000 loops each)\n", "1.56 µs ± 25.3 ns per loop (mean ± std. dev. of 7 runs, 1000000 loops each)\n", "11.3 µs ± 150 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)\n", "110 µs ± 1.3 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)\n", "1.16 ms ± 14.6 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)\n", "13 ms ± 948 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)\n" ] } ], "source": [ "tn, tne = [], []\n", "for x, y in xy:\n", " t = %timeit -o trapzn(y, x)\n", " tn.append(t.average)\n", " tne.append(t.stdev)\n", "\n", "tna = np.array(tn)\n", "tnea = np.array(tne)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's look at the Numba data alongside the Python one:" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY4AAAEaCAYAAAAG87ApAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3Xd4VFX6wPHvmx5aAqEnQOi9JASC\nuggIKhYE7IoiUpTd1bXsT1dXXV3L6q7uurruWkFEBQREFBuiiF16lQ5SEkAETGgBUs7vj3MDQ0yZ\nSWZyJ5P38zzzZObeO3feMzO579xzzj1HjDEopZRS3gpzOwCllFJViyYOpZRSPtHEoZRSyieaOJRS\nSvlEE4dSSimfaOJQSinlE00cIUBEHhWRfSKyx+1YgpWIbBORQeV8bl8R2eDvmEKViLwgIg/4eZ+j\nRORrf+5TlZ8mDhc4B7EcETksIj+JyKsiUquc+2oG/BHoZIxp7N9IqycRMSLSpvCxMeYrY0z7ALzO\nQyLyhr/3WxEiMklEHq3IPowx440xj/grprKISLLzmX1QZPkbIvJQZcVRnWjicM8QY0wtIBXoBdzv\n6w5EJAJoAew3xuwt5/OrvFApRygQkXAXX76PiJzl4utXG5o4XGaMyQQ+AroAiEiciEwQkd0ikulU\nQ4U760aJyDci8rSIHAAWAPOAps7ZyyRnu0tE5AcRyRKRBSLSsfD1nLOdP4nIKuCIiEQ4y+4SkVUi\ncsR5/UYi8pGIHBKRT0Wkrsc+ZojIHhHJFpEvRaSzx7pJIvJfEfnAee5CEWntsb6ziMwTkQPO2daf\nneVhInKPiGwRkf0iMl1E6hX3nolIfxHJcMqxB3jVWX6xiKxwyv2tiHQr4fm9ReQ7Z7vdIvKciEQ5\n6750NlvpvKdXFb6es/4eEZlZZH/PiMizZX1+RZ4zGPgzcJXzOiud5U1F5D3n/dksIuNKKEMf5zMI\n91g23Plcy3w/ReQ3znuUJSI7ne/WTcAI4G4npjnOth2d71GW8726pMjn/byIfCgiR4AB4nHWIiJz\nnH0V3gpEZJSzroPHd2GDiFzpsd8E5304KCKLgJPfoVL8Ayj2bEmKqeoSjzNLJ+b/Od/5w2L/zxqL\nyL9F5BcRWS8iKR7P3SYi94rIWmf9qyIS46xbIyJDPLaNFFuV3MOLMlQNxhi9VfIN2AYMcu43A34A\nHnEezwZeBGoCDYFFwM3OulFAHnArEAHEAv2BDI99twOOAOcCkcDdwGYgyuO1VzivG+ux7HugEZAI\n7AWWASlANDAfeNDjNUYDtZ11/wZWeKybBBwAejsxvglMc9bVBnZjq9ZinMfpzrrbnRiSnP2+CEwt\n4f3r77wPf3e2jcWeue0F0oFw4AanXNHFvOc9gT5OfMnAOuB2j/0boE2R18tw7rcAjgJ1nMfhTpn6\nlPX5FVOOh4A3iiz7Avif8/70AH4GBpbw/C3AuR6PZwD3lPV+As2BQ8A1znckAejh8fk96rHPSOz3\n589AFHCO89z2HttnA2dhf4jGFN2Hx74GA7uw372awE7gRudzSAX2AZ2dbacB053tugCZwNclvA/J\nzmdWy9mu8HN+A3jI43/n6yLPO/k5OzHvc74bMdjv/I/ASOczfhT4vMj/8BqnLPWAbwrLjP2fe8tj\n26HAarePO349hrkdQHW8OV+6w0AWsN05UMRiD9zHcQ7ozrbXFH5hnS//jiL76s/pieMBYLrH4zDn\nn6m/x2uPLiaeER6P3wae93h8KzC7hLLEO/+Acc7jScArHusvBNZ7lGV5CftZh8cBEmgC5AIRxWzb\nHzgBxHgsex4n+Xos2wD08yjjoBJe+3bgHY/HJSYO5/HXwEjn/rnAFud+qZ9fMa/7EB6JwzkI5QO1\nPZY9Dkwq4fmPAhOd+7WxPxhalPV+Avd6lrfIPidxeuLoC+wBwjyWTeXUAXkSMLm0fTjL2mETe1/n\n8VXAV0W2eRF4EHugzgU6eKz7G2Unjgjgd8D3znJfE8fLRb7z6zwedwWyivzPjC/yPS/8HjTFJtfC\nHxczgbuLi72q3rRu2D3DjDGfei4Qka7YX3i7RaRwcRj2l1khz/vFaYpNRgAYYwpEZCf2TKK0ffzk\ncT+nmMe1nBjDgceAK4AGQIGzTX3sL0+wB5pCRwufiz0wbikh7hbAOyJS4LEsH3swzixm+5+NMceK\nPP8GEbnVY1kU9v04jYi0A/4FpAE1sAecpSXEVZwp2IQwGbjWeVwYQ1mfX2maAgeMMYc8lm134iwp\njm9F5LfApcAyY0zhZ1/a+1na51BcTDuNMZ772U7Z36eTRCQOeBd4wBjzlUd86SKS5bFpBPA69nsV\nUWS/2/HOy8BdnlVFPvDqf8BD0fiaAhhjdonIN8BlIvIOcAFwWzniCVqaOILLTuwv1vrGmLwStilr\nOONd2F9HAIg9gjXj9INvRYZEvhZ76j0I+6srDvgFkFKeU2gn9oBb0rrRxphvvIyjaBl2Ao8ZYx7z\n4rnPA8uBa4wxh0TkduByL18XbJXQP0UkCRgOnOERQ1mfn6eiZdgF1BOR2h7JoznFJ06MMWtFZDv2\nwOSZwApjKfb9dH5I9PYhpmYiEuaRPJoDG0t5judrhTlxfW6MebFIfF8YY84t5jnh2KrIZsB6j9cs\nkzEmV0T+CjyCrQIudAT7I6HwNfzRA7GZx/3m2Peq0GvAWOwx9jtj2zJDhjaOBxFjzG7gE+xBqY7T\nwNlaRPr5sJvpwEUiMlBEIrHtCceBb/0UZm1nf/ux/4h/8+G57wONReR2EYkWkdoiku6sewF4TERa\nAIhIAxEZ6sO+XwbGi0i6WDVF5CIRqV1CGQ4Ch0WkA/DbIut/AlqV9ELGmJ+xHRNeBX40xqxzlvv6\n+f0EJDsHV4wxO7Gf0+MiEiO2cX8Mtp2oJFOAPwBnYxNaodLezzeBQSJypdjOEQkeDbdFy74Qe9C9\n22nk7Q8MwbZBeOMxbDtF0V/c7wPtROR6Z7+RItJLRDoaY/KBWcBDIlJDRDph26y89Tq2XWewx7KV\nQGcR6eE0Yj/kw/5K8nsRSRLb6eDPwFse62Zj221uw56ZhhRNHMFnJLaKZS32l/xMbP20V4wxG4Dr\ngP9gG/uGYLv+nvBTfJOxp+WZTozf+xDbIWybwBBsddYmYICz+hngPeATETnk7De9uP2UsO8lwDjg\nOez7thlbr12c/8P+Qj+ETThvFVn/EPCa04voSoo3BXvWNaXIcl8+v8ID/X4RWebcvwZbZ78LeAfb\nKWFeCc8H297QH5hvjNnnsbzE99MYswNbJ/9HbEeGFUB353kTgE5O2Wc735tLsGc1+7DtcSONMevx\nzjXYjgi/ePSsGuF8F84DrnbKuodTnR0AbsFWDe3Btj+86uXr4SSeB7GN1oXLNgIPA59iv3f+uJhw\nCvaHwlbndrJHlzEmB9tW2BKbBEOKOI03SimlvCQi24CxRdspi2zzF6CdMea6Sguskmgbh1JK+ZlT\nfTUGuN7tWAJBq6qUUsqPxF60uRP4yBjzZVnbV0VaVaWUUsonesahlFLKJ5o4lFJK+SQkG8fr169v\nkpOTy/XcI0eOULNmTf8G5JJQKUuolAO0LMEqVMpSkXIsXbp0nzGmgTfbhmTiSE5OZsmSJeV67oIF\nC+jfv79/A3JJqJQlVMoBWpZgFSplqUg5nFEIvBJSVVUiMkREXsrOzi57Y6WUUuUSUonDGDPHGHNT\nXFyc26EopVTICqnEoZRSKvBCso2jOLm5uWRkZHDs2LFSt4uLi2PdunWVFFVgeVuWmJgYkpKSiIyM\nrISolFJVXbVJHBkZGdSuXZvk5GQ85kr4lUOHDlG7dnEDqlY93pTFGMP+/fvJyMigZcuWlRSZUqoq\nqzZVVceOHSMhIaHUpFEdiQgJCQllnokpVWj28kzOemI+oz4+wllPzGf28pCaakJ5IaQSR1m9qnxN\nGle9+B1XvfidP0ILappMlbdmL8/k3lmryczKASAzK4d7Z63W5FHNhFTiCPZeVeHh4fTo0YMuXbpw\nxRVXcPTo0RK33bZtG1OmnJrqYdKkSdxyyy2VEaZSAOTmF/DzoeNs/OkQC7fu5+M1u3novR/Iyc0/\nbbuc3HyenLvBpSiVG6pNG4evZi/PZPmOLE7kF3DWE/O56/z2DEtJLPuJpYiNjWXFihUAjBgxghde\neIE777yz2G0LE8e1115boddUCiAvv4CsnFx+OXKCA0dO8MvRXH456twv+vioXXbwmDez31q7nDMQ\nVT1o4ihG4en4iXw7xXLh6ThQ4eRRqG/fvqxatYoHHniA+vXrc9ttdmbN++67j0aNGjFlyhTWrVtH\njx49uOGGG6hbty67du1i8ODBbNmyheHDh/OPf/wDgKlTp/K3v/0NYwwXXXQRf//73wFo0qQJt912\nG++//z6xsbG8++67NGrUyC/xK9/MXp7Jk3M3kJmVQ+L3FfshkpdfQHZO4YE+lwNHTpB19AQHnAP+\ngSO5RR6XngRiI8OpVzOK+BqR1KsZRfN6NU57XLdG1MnHYyYtZs/B47/ahwGun7CQEektGNSxIRHh\nIVWZoYqolonjr3N+YO2ug8Wuy8/PZ1XmoZNJo1BObj53z1zF1EU7in1ep6Z1eHBIZ69ePy8vj48+\n+ojBgwdzwQUXcOmll3LbbbdRUFDAtGnTWLRoEd26deOpp57i/fffB2xV1YoVK1i+fDnR0dG0b9+e\nW2+9lfDwcP70pz+xdOlS6taty3nnncfs2bMZNmwYR44coU+fPjz22GPcfffdvPzyy9x///0+vFPK\nHwp/iBRW8Xj+EBnSvSlZR8s4Ayg8Czhqk0R2Tm6JrxUbGU7dGpHUrWkP9kl1a1DP43F8jSjq1Yii\nbs1TSSEmMtzrstxzQcfTygIQExFG//YNWZmRxfg3ltKoTjRX92rO1b2b0SQutpzvmgpm1TJxlKVo\n0ihrubdycnLo0aMHYM84xowZQ1RUFAkJCSxfvpyffvqJlJQUEhISin3+wIEDKWy/6dSpE9u3b2f/\n/v3079+fBg3s2GQjRozgyy+/ZNiwYURFRXHxxRcD0LNnT+bNK23qahUoT87dUGy7wB3TV3DH9BWU\nNCVOdEQYCTWjqOsc4BOdJBDvnAHUrWmTgOeZQWyU90mgPArPkk6ePcXHnjx7yssvYP76vby5cAfP\nzt/Ec59vZmCHhozo04K+beoTFqadMEJFtUwcpZ0ZHDp0iMH/XXyy14inxPhY3rr5jHK/rmcbh6ex\nY8cyadIk9uzZw+jRo0t8fnR09Mn74eHh5OXlUdpEXJGRkSd7TBVuryqXMabY75JdB7cNbFtstVBl\nJIHyGpaSyLCUxF8NqBcRHsZ5nRtzXufG7Nh/lCmLdjBjyU4+WfsTzevV4Nr05lzRM4mEWtEl71xV\nCSFVEemvQQ7vOr89sUVO32Mjw7nr/PYV2m9Jhg8fzscff8zixYs5//zzAahduzaHDh0q87np6el8\n8cUX7Nu3j/z8fKZOnUq/fv0CEqfyzYY9hxg5cVGJ6xPjY7nj3HbccGYyQ3sk0rdtA7okxtE0PjZo\nk4a3mifU4J4LOvDtvefwzNU9aBwXwxMfreeMx+dz27TlLPrxQKk/elRwC6kzDmPMHGBOWlrauIrs\np/B0/O6ZqziRX3Da6XggREVFMWDAAOLj4wkPtweMbt26ERERQffu3Rk1ahR169Yt9rlNmjTh8ccf\nZ8CAARhjuPDCCxk6dGhA4lTeOXDkBP+at4EpC3dQOyaSS1Oa8uGaPRzLPVXVGcgfIsEkOiKcoT0S\nGdojkU0/HeLNhTt4e2kG767YRbtGtRiR3oLhqYnUidHhbqqSkJxzPC0tzRSdj2PdunV07NixzOd6\nDtNRePFfRaqnvFFQUEBqaiozZsygbdu2ftuvL8OnePv+uKGqzJVwIq+Ayd9t45nPNnH0RD7XpTfn\n9kHtqFsz6vReVQH+IVJZyvu5HD2Rx5yVu3hz4Q5WZWQTGxnO0B5Nua5PC7okunMNVlX5jpWlgvNx\nLDXGpHmzbUidcfhboBMGwNq1a7n44osZPny4X5OGqjzGGD5bt5fHPlzHj/uOcHa7BjxwUUfaNjqV\ntEtqF6iOakRFcFWv5lzVqzmrMrJ48/sdzF6RybTFO+meFMeI9BYM6d60ylfXhTJNHC7r1KkTW7du\ndTsMVU7r9xzk0ffX8fXmfbRuUJNXb+zFgPYN3Q6ryuiWFE+3y+P580UdeWdZBm8u3MHdb6/ikQ/W\ncllqEtf1aU6bhqEx6Ggo0cShVDnsP3ycpz/deLId48EhnbiuTwsi9cK3comLjWTUWS254cxkFv14\ngDcX7uDNhduZ9O020lvW47o+LTi/c2OiIvT9DQaaOJTyQdF2jJFnJHPbwLbUrRnldmghQURIb5VA\neqsE9h3uxIwlGUxZtJ1bpy6nfq0orkxrxjW9m9OsXg23Q63WNHEo5YWi7Rj92jXggYs7ajVKANWv\nFc1v+7fm5rNb8eWmn3lz4Q5e+GILz3+xhX7tGnBdegsGdGhIuF5YWOk0cZTm1Yvs3xs/cDcO5Spt\nx3BXWJjQv31D+rdvyK6sHKYt3sm0RTsYO3kJTeNiuKZ3c67q1YyGdWLcDrXaCKkKQ39dABgoIsIf\n//jHk4+feuopHnroIb/se9SoUcycOdMv+1LW/sPHue+d1Vz4zFeszszmoSGd+Pj2szVpuKhpfCx3\nntuOb+45hxeuS6V1w1r8c95GznxiPr97cynfbN5HQUHoXWIQbELqjMNfFwACsGo6ZCyG/OPwdBcY\n+BfodmWFdhkdHc2sWbO49957qV+/foVDVIFRXDvG7YPaEl9D2zGCRWR4GIO7NGFwlyZs23fk5PAm\nH67eQ8v6NRmR3pzLUpO07SlAQuqMw29WTYc5f7BJAyB7p328anqFdhsREcFNN93E008//at1Rc8Y\natWqBdgLevr168eVV15Ju3btuOeee3jzzTfp3bs3Xbt2ZcuWLSef8+mnn9K3b1/atWt3clTdbdu2\n0bdvX1JTU0lNTeXbb7+tUBlCmTGGeWt/4rynv+DRD9bRs0Vd5t7el4cu6axJI4gl16/Jny/syHf3\nDuTpq7qTUDOKRz9YR/rjn3Hn9BUs3f6LDm/iZyF1xuG1j+6BPauLXRWbnwe7l59KGoVyc+DdW2Dp\na8Xvs3FXuOCJMl/697//Pd26dePuu+/2OtyVK1eybt066tWrR6tWrRg7diyLFi3imWee4T//+Q//\n/ve/AZskvvjiC7Zs2cKAAQNYvnw5DRs2ZN68ecTExLBp0yauueYail5Vr2w7xiPvr+Wbzfu1HaOK\niokMZ3hKEsNTkli3+yBTFu7gneWZzFqWSccmdRiR3pxhKYnUiq6ehz1/0newOEWTRlnLfVCnTh1G\njhzJs88+S2ysd3MV9OrViyZNmgDQunVrzjvvPAC6du3K559/fnK7K6+8krCwMNq2bUurVq3YuHEj\nXbp04ZZbbmHFihWEh4ezcePGCpchlOw/fJx/zdvI1EX2eoyHhnRihF6PUeV1bFKHR4Z14Z4LOvDu\nil288f127p+9hsc/XMewlERGpLegU9M6fp1gqzqpnomjlDODnEOHqP3KGbZ6qqi4Zn7pYXX77beT\nmprKjTfeeHJZREQEBQV2EDxjDCdOnDi5znM49bCwsJOPw8LCThsqvXAIdc/HTz/9NI0aNWLlypUU\nFBQQE6M9T8C2Y7z27Tae/WwTR3O1HSNU1YyO4Nr05lzTuxkrdmbx5sIdzFxqr1BvUS+WXdnHyM23\n1ViBmOkzVOnPquIM/AtEFjkbiIy1y/2gXr16XHnllUyYMOHksuTkZJYuXQrAu+++S25uybO8lWTG\njBkUFBSwZcsWtm7dStu2bcnOzqZJkyaEhYXx+uuvk5+fX/aOQphnO8ZjH64jLVnbMaoDESGleV2e\nuqI7i/48iAcu7kRm1qmkUSgnN58n525wKcqqQxNHcbpdCUOehXDnl35cM/u4gr2qPP3xj39k3759\nJx+PGzeOL774gt69e7Nw4UJq1qzp8z7bt29Pv379uOCCC3jhhReIiYnhd7/7Ha+99hp9+vRh48aN\n5dpvqFi3+yDXTVjIuMlLiAgPY9KNvXj1xt56EV81E1cjkjG/aUl+Cd12M7Ny2HngaCVHVbXosOpF\nnDYUeRW/AFCHVbf2Oe0Y0xbtoE5sJHcMase16c1daccIpdFxq3pZznpifomzMwL0blmPy1ITuaBr\nkyozX4gOqx4MqmjCUJa2Y6jS3HV+e+6dtfq0+eDtBFvtOHoin1nLM/nT26v5y7s/cG6nRlyaamdp\n1I4TmjhUCCpsx/jbh+vYtv8oA9o34L6LdFwpdbrCBvCSJtj6/YA2rMzI5p1lGby3chfvr9pN/VpR\nXNI9kUtTE+nctM6vOqRUF5o4VEhZt/sgj35gr8do07AWk27sRX+9HkOVoLQJtkSEHs3i6dEsnvsu\n6sSCDXt5Z3kmb3y/nYnf/Ei7RrUYnpLEsJSmNInzrmt9qKhWicMYU21/IZQmFNq5irZj/PWSzq61\nY6jQExURxnmdG3Ne58ZkH83l/dW7mLUsk79/vJ5/zF3Pma0TuDQlicFdGlOzGlxgGPoldMTExLB/\n/34SEhI0eXgwxrB///4qe33H8bx8Xvt2G//5bDM5ufnccKadH0PbMVSgxNWIZER6C0akt2DbviP2\n6vTlGfxxxkrun72GwV0ac2lqIme2rh+yQ74HfeIQkWHARUBD4L/GmE/Ks5+kpCQyMjL4+eefS93u\n2LFjVfYgWpS3ZYmJiSEpKakSIvKfwnaMxz5cx/aT7RidaNOwltuhqWokuX5N7ji3HbcPasvS7b/w\n9rJMPli1i3eWZ9KoTjTDeiQyPDWRDo3ruB2qX7mSOERkInAxsNcY08Vj+WDgGSAceMUY84QxZjYw\nW0TqAk8B5UockZGRtGzZssztFixYQEpKSnleIuiEUlk8rdttx5X6dou2Y6jgICKkJdcjLbkeDw7p\nxPz1e5m1LJMJX//Ii19upVOTOlyamsglPZrSsHbV/2Hq1hnHJOA5YHLhAhEJB/4LnAtkAItF5D1j\nzFpnk/ud9aqaKDqO0Pj+rVi76xBvLbbtGA8P7cy1vZsToe0YKojERIZzYdcmXNi1CfsPH2fOSnsG\n8ugH63j8o/X0bVuf4SmJnNepMbFR4W6HWy6uXQAoIsnA+4VnHCJyBvCQMeZ85/G9zqZPOLd5xphP\nS9nfTcBNAI0aNeo5bdq0csV1+PDhk0OaV3VVuSzf7spl0poTnCg4fbkAg1pEMKxNFDUjq179cVX+\nTIrSsvhm1+ECvt2Vx7e78jhwzBATDr0aR3Bm0wja1wsjzA9trxUpx4ABA6rkBYCJgOfIghlAOnAr\nMAiIE5E2xpgXinuyMeYl4CWwV46X9+rJqn41rKeqXJb7npj/q6QB0KB2NC//dlDlB+QnVfkzKUrL\n4rtrgYICw/c/7uedZZl8uHo3X2UeIzE+lmEpTRieklShdrrKKkcwJY7i0q0xxjwLPFvZwSh37Sph\nKIifD1V8aHul3BQWJpzZuj5ntq7Pw0O78MnaPcxalsnzC7bw38+30D0pjktTkxjSvSn1gnQGw2BK\nHBlAM4/HScAuX3YgIkOAIW3atPFnXMoF8TUi+eXor0cIbhpfvS60UqEtNiqcoT0SGdojkb0Hj/He\nyl28vSyTB9/7gUfeX0v/9g25LDWRczo2JDoieNpDgilxLAbaikhLIBO4Gntm5zW/zjmuXDPpmx/5\n5WguYQKeA5jacYTauxeYUgHUsE4MY/u2YmzfVqzbfZB3lmcye3kmn677iToxEVzcvSmXpiTSs0Vd\n169Fc6s77lSgP1BfRDKAB40xE0TkFmAutjvuRGPMD27Ep9xhjOHpeRt5dv5mzuvUiHM7NeLfn24q\ndhwhpUJZxyZ16NikDn8a3IFvNu9j1rIM3lmWyZSFO2iRUINhPex4WS0S7DQJlT2ToSuJwxhzTQnL\nPwQ+LO9+taqq6sovMDzw7hqmLNzBVWnNeGx4FyLCw7girVlINcIq5YvwMOHsdg04u10DDh/P4+M1\ne3hneQbPzt/EM59tIq1FXVo2qMmclbs4lmt7k1TGTIYh1QHeGDPHGHNTXFyc26EoHxzPy+eWKcuY\nsnAHv+vfmicu66rXZihVRK3oCC7vmcSbY/vwzZ/O4e7B7cnOyWXGkoyTSaNQoGcy1P9O5apDx3K5\n8dXFfLRmD/df1JG7B3dwvf5WqWDXND6W3/Vvwyd3nF1sd1QouWeiP4RU4hCRISLyUnZ2ttuhKC/s\nO3yca17+nkU/HuDpq7oztm8rt0NSqkoRkRJ7GgayB2JIJQ6tqqo6dh44yuXPf8vmvYd5eWQaw1Oq\n1iCLSgWLu85vT2zk6V11A90DMZi646pqYv2eg4ycsIjjeQW8ObYPPVvUdTskpaqssmYyDARNHKpS\nLd52gDGTFlMjKoIZ48+gXSOdzlWpiiptJsNACKmqKm3jCG6frv2J615ZSP1a0cz8rSYNpaqqkEoc\n2sYRvGYuzeDmN5bSvnFtZow/g6S6NdwOSSlVTlpVpQLupS+38LcP1/ObNvV54fqe1KoGczIrFcr0\nP1gFjDGGJz5ez4tfbOWibk3415Xdg2qgNqVU+YRU4tAhR4JHXn4B985azYylGVzfpwUPXdKZ8DC9\nsE+pUKBtHMrvjuXmM/6NZcxYmsHtg9ry8FBNGkqFkpA641Duy87JZdxrS1i8/QAPD+3MyDOS3Q5J\nKeVnmjiU3+w9eIyRExex5efDPHt1CkO6N3U7JKVUAGjiUH6xbd8Rrp+4kP2HTzBxVC/6tm3gdkhK\nqQDRxKEqbE1mNqNeXUR+gWHquD50bxbvdkhKqQAKqcZxvXK88n23ZT9Xv/Q90RHhzBh/piYNpaqB\nkEoc2quqcn28Zjc3TFxEk7gYZv72DNo0rOV2SEqpSqBVVapcpi7awX3vrKZHs3gmjupFfI0ot0NS\nSlUSTRzKJ8YY/rdgC0/O3UD/9g3434hUakTp10ip6kT/45XXCgoMj3ywlle/2cbwlET+cXk3InVu\ncKWqHU0cyisn8gq4a+ZK3l2xi9FnteT+izoSpleDK1UtaeJQZTp6Io/fvrGMLzb+zN2D2/Pbfq0R\n0aShVHUVUolDBzn0v1+OnOCVRsY1AAAgAElEQVTGSYtZlZHFE5d25erezd0OSSnlspCqoNbuuP61\nKyuHK178jrW7D/L8dT01aSilgBA741D+s3nvYUZOWMihY3lMHt2bPq0S3A5JKRUkNHGoX1mxM4sb\nX11EeFgY027uQ+emeganlDpFE4c6zVebfubm15dSv1Y0r4/pTYuEmm6HpJQKMpo41ElzVu7izukr\naNOwNq/d2IuGdWLcDkkpFYQ0cSgAJn+3jQff+4FeyfV4eWQacbGRboeklApSZfaqEpF/iEgdEYkU\nkc9EZJ+IXFcZwanAM8bwr3kb+cu7PzCwQyMmj+6tSUMpVSpvuuOeZ4w5CFwMZADtgLsCGpWqFPkF\nhgfeXcOzn23iip5JvHBdKjGR4W6HpZQKct5UVRX+/LwQmGqMOaBXDVd9x/PyufOtlXywejfj+7Xm\nT4Pb69XgSimveJM45ojIeiAH+J2INACOBTas8tErx71z+HgeN7++hG827+e+Czsy7uxWboeklKpC\nyqyqMsbcA5wBpBljcoEjwNBAB1YeeuV42fYfPs41L33P91sP8M8rumvSUEr5rMwzDhGJBK4Hznaq\nMr4AXghwXCoAdh44yg0TF7ErO4eXR/bknA6N3A5JKVUFeVNV9Ty2neN/zuPrnWVjAxWU8r8New4x\ncuJCck7k88aYdNKS67kdklKqivImcfQyxnT3eDxfRFYGKiDlf0u2HWD0pMXERoUzY/yZtG9c2+2Q\nlFJVmDeJI19EWhtjtgCISCsgP7BhqYqYvTyTJ+duIDMrh4Sv5pGdc4Jm9WoyeXRvmtWr4XZ4Sqkq\nzpvEcRfwuYhsBQRoAdwY0KhUuc1ensm9s1aTk2tz+/4jJxCB0b9J1qShlPKLMhOHMeYzEWkLtMcm\njvXGmOMBj0yVy5NzN5xMGoWMgRcWbOX6PsnuBKWUCiklJg4ROccYM19ELi2yqrWIYIyZFeDYVDns\nysrxablSSvmqtDOOfsB8YEgx6wygiSMIxcVGkpWT+6vlTeNjXYhGKRWKSkwcxpgHnbsPG2N+9Fwn\nIi0DGpUql3eWZ5CVk0uYQIE5tTw2Mpy7zm/vXmBKqZDizSCHbxezbKa/A1EV8/GaPfzfjFWc2TqB\nv1/WjUTnDCMxPpbHL+3KsJRElyNUSoWK0to4OgCdgbgi7Rx1AJ3hJ4gs2LCXW6cuo3tSHC+PTKNm\ndARXpDVjwYIF9O/f3+3wlFIhprQ2jvbYodTjOb2d4xAwLpBBeXKuG7kPiDPGXF5Zr1tVLNy6n5tf\nX0rbhrV59cbe1IzWubmUUoFVWhvHu8C7InKGMeY7f76oiEzEJqW9xpguHssHA88A4cArxpgnjDFb\ngTEiotVjRazYmcWY15bQrF4NXh+jEzAppSqHNz9Pl4vI77HVVierqIwxoyvwupOA54DJhQtEJBz4\nL3AudsKoxSLynjFmbQVeJ2St232QGyYuom7NSN4Yk05CrWi3Q1JKVRPeNI6/DjQGzseOjJuEra4q\nN2PMl8CBIot7A5uNMVuNMSeAaQTp8O1u2/rzYa6fsJDYyHCmjO1D4zhtclJKVR4xxpS+gchyY0yK\niKwyxnRzhlmfa4w5p0IvLJIMvF9YVSUilwODjTFjncfXA+nAg8Bj2DORV4wxj5ewv5uAmwAaNWrU\nc9q0aeWK6/Dhw9SqVatcz60MPx8t4PFFx8gtMNzbO5amtUrO/cFeFm+FSjlAyxKsQqUsFSnHgAED\nlhpj0rzZ1puqqsKrybJEpAuwB0guV2SlK27eUmOM2Q+ML+vJxpiXgJcA0tLSTHl7EwVzT6S9B49x\nxYvfkUc4b40/g05N65S6fTCXxRehUg7QsgSrUClLZZXDm6qql0SkLnA/8B6wFvh7AGLJAJp5PE4C\ndvmyAxEZIiIvZWdn+zWwYHDgyAlGvLKQfYeOM2l07zKThlJKBUqpiUNEwoCDxphfjDFfGmNaGWMa\nGmNeDEAsi4G2ItJSRKKAq7GJymuhOnXswWO5jJy4kB0HjvLKDb1IbV7X7ZCUUtVYqYnDGFMA3OLv\nFxWRqcB3QHsRyRCRMcaYPOe15gLrgOnGmB/8/dpVzdETedz46mI27DnEC9f15IzWCW6HpJSq5rxp\n45gnIv8HvAUcKVxojCnaK8prxphrSlj+IfBhefcrIkOAIW3atCnvLoLKsdx8xk1ewvIdv/Dfa1MZ\n0KGh2yEppZRXbRyjgd8DXwJLnduSQAZVXqFUVZWbX8AtU5bxzeb9PHl5dy7o2sTtkJRSCvBuIicd\nCbeS5RcY7py+kk/X7eWRoZ25rGeS2yEppdRJ3pxxqEpUUGD486zVzFm5i3su6MD1ZyS7HZJSSp0m\npBJHVe+Oa4zh4ffX8taSnfzhnDaM79fa7ZCUUupXQipxVPU2jn9+spFJ325j9FktuePcdm6Ho5RS\nxSqzjUNEUotZnA1sd7rQKj/434LNPPf5Zq7u1YwHLu6ISHEX0iullPu86Y77PyAVWIUdFqSLcz9B\nRMYbYz4JYHw+qardcV/7dhv/+HgDQ3s05bHhXTVpKKWCmjdVVduAFGNMmjGmJ5ACrAEGAf8IYGw+\nq4pVVTOW7OTB937g3E6NeOqK7oSHadJQSgU3bxJHB88ruJ35MVKcCZZUBXywajd/ensVfdvW5z/X\npBAZHlJNTkqpEOVNVdUGEXkeOz8GwFXARhGJ5tTIucpH89f/xG3TltOzRV1evL4nMZHhboeklFJe\n8eYn7ihgM3A7cAew1VmWCwwIVGCh7NvN+xj/xjI6NqnDhFG9qBGl84QrpaoOb64czwH+6dyKOuz3\niCqgKjSOL93+C2MnLyE5oQaTR/emTozOE66UqlrKPOMQkbNEZJ6IbBSRrYW3ygjOV8HeOL4mM5tR\nry6iYe1o3hiTTt2aUW6HpJRSPvOmjmQCtopqKZAf2HBC1+a9hxg5cRG1oyN4Y2w6DevoPOFKqarJ\nm8SRbYz5KOCRhLAd+48y4pWFhInw5rg+JNWt4XZISilVbt4kjs9F5ElgFnC8cKExZlnAogohu7Nz\nuPaV7zmeV8C0m/rQsn5Nt0NSSqkK8SZxpDt/0zyWGeAc/4dTMcHWOL7v8HFGvLKQrKO5TBmXTofG\nOk+4Uqrq86ZXVZXpcmuMmQPMSUtLG+d2LNlHc7l+wiJ2ZeUweXQ63ZLi3Q5JKaX8osTEISLXGWPe\nEJE7i1tvjPlX4MKq2g4fz+OGVxexZe9hXrkhjd4t67kdklJK+U1pZxyFlfG1KyOQUHEsN5+xry1m\ndWY2/xuRytntGrgdklJK+VWJicMY86Lz96+VF07VdiKvgPFvLGXhjwd4+soenN+5sdshKaWU33kz\nH0dL4FYg2XN7Y8wlgQur6snLL+D2t5azYMPP/G14V4alJLodklJKBYQ3vapmYy8CnAMUBDacqqmg\nwHD326v4cPUe7r+oI9emN3c7JKWUChhvEscxY8yzAY/ED9zojmuM4cH3fmDWskzuGNSOsX1bVdpr\nK6WUG7wZHfcZEXlQRM4QkdTCW8AjK4fKHqvKGMMTH6/n9e+3c/PZrfjDwOC4fkQppQLJmzOOrsD1\n2Av+CquqgvICwMr23PzNvPjFVq7r05x7LuigU74qpaoFbxLHcKCVMeZEoIOpSiZ8/SP/nLeRS1MS\nefiSLpo0lFLVhjdVVSsBvezZw7RFO3jk/bUM7tyYf1zejTCdJ1wpVY14c8bRCFgvIos5fZDDatkd\n990Vmdz7zmr6tWvAs9ekEKHzhCulqhlvEseDAY+iivjkhz3cOX0lvZPr8cJ1PYmK0KShlKp+vBnk\n8IvKCCTYfbXpZ26ZspwuiXFMGNWL2Khwt0NSSilXlDbI4dfGmN+IyCFsL6qTqwBjjKk2Y4Qv3naA\ncZOX0KpBTV67sRe1or05UVNKqdBU5iCHxphqPcjhqowsRr+6mKZxsbw+Jp34GjpPuFKqeiutkt6U\nsi4oicgQEXkpOzvbL/vbsMfOE14nNpI3xqbToHa0X/arlFJVWWlnHA1LmosDgnM+Dn9O5PTjviNc\nN2EhUeFhTBmXTtP4WD9EqJRSVV9piSMcqIVt06hWMrNyGPHy9+QXGN66qQ8tEnSecKWUKlRa4tht\njHm40iJx2ezlmTw5dwOZWTmEf/I5EQJv/+4s2jaq1k08Sin1K6W1cVSbM43ZyzO5d9ZqMrNyAMgv\nMCDC5r2HXY5MKaWCT2mJY2ClReGyJ+duICc3/7Rlx/MKeHLuBpciUkqp4FVi4jDGHKjMQNy0yznT\n8Ha5UkpVZzpmBpTYY0p7Uiml1K9p4gDuOr89sZGnDyESGxnOXee3dykipZQKXjp2BjAsJRHgZK+q\nxPhY7jq//cnlSimlTtHE4RiWksiwlEQWLFhA//793Q5HKaWCllZVKaWU8okmDqWUUj7RxKGUUson\nQd/GISI1gf8BJ4AFxpg3XQ5JKaWqNVfOOERkoojsFZE1RZYPFpENIrJZRO5xFl8KzDTGjAOq5Tzn\nSikVTNw645gEPAdMLlwgIuHAf4FzgQxgsYi8ByQBq53NTh8XpCQbNkA5e0b1yMqC+PhyPTfYhEpZ\nQqUcoGUJVqFSlsoqhyuJwxjzpYgkF1ncG9hsjNkKICLTgKHYJJIErKCUMyQRuQm4CaBLZCRZWVnl\nii0/P7/czw02oVKWUCkHaFmCVaiUpbLKEUxtHInATo/HGUA68CzwnIhcBMwp6cnGmJeAlwDS0tJM\n/JIl5QoilK7jCJWyhEo5QMsSrEKlLBUqh3g/IHowJY7iojbGmCPAjZUdjFJKqeIFU3fcDKCZx+Mk\nYJcvO/D3nONKKaV+LZgSx2KgrYi0FJEo4GrgPV92YIyZY4y5KS4uLiABKqWAVdPh6S70WzAMnu5i\nH6tqxa3uuFOB74D2IpIhImOMMXnALcBcYB0w3RjzgxvxKaVKsGo6zPkDZO9EMJC90z7W5FGtuNWr\n6poSln8IfFje/YrIEGBImzZtyrsLpaqn/Dw4lgU5v8DRA/Zvzi+Q43H/6AHY8AHkHT/9ubk58NnD\n0O1Kd2JXlS6YGscrzBgzB5iTlpY2zu1YlDrNqunw2cP0y86A5Ukw8C+BOdB6JgDPA35JiSDnF8jJ\nguOltAtKGMTEQ416v04ahbJ32v3VqOf/MqmgE1KJQ6mgVFi9k5tjuw4WVu9AycmjIN8e0Es94Bdd\n/ot3CSC2rr3VaggN2juP651aXqPuqfux9SC6DoQ5tdpPd7HxF+dfHaHzpdBrDCT29Kl7p6paQipx\naFWVCkqfPWyrczzl5sAHf4SMxcUngmOl9QwUiPVIADXqQ0Jb+2vf84BfNBFEx51KAOU18C8nk+BJ\nkbFw9t02oayaDiunQONuNoF0vQKialbsNVXQCanEoVVVKugUFEB2RvHrjh+0B9rYuvagX6MeJLQu\nctD3TAbOLSYOwsKL32egFZ4hffYwJjsDiStS7Xbuw7ZMSybCnNvgkweg21U2iTTs6E7Myu9CKnEo\nFTSO7IPlb8DSSYApfpu4JLijCnYc7HYldLuSL4q7Sjm6tk0SaaNh5yJYMgGWvQaLX4bmZ9p1HYdA\nRLQroSv/0MShlL8YA9u/sb+2182B/BPQ4ixoPRBWvvnr6p2BD7oXa6CJQPN0ezv/cVjxBix5Fd4e\nY6vWUq+HnqOgbrLbkapyCKnEoW0cyhVHD8DKabD0Vdi30VYlpY2xB8aGHew2zdNLrt4JdTUT4Kzb\n4IxbYevnNrF+8wx8/W9oM8iehbQ9z73qN+WzkEoc2sahKo0xtmF7yUT44R3IOwZJvWDY89BpGETV\nOH370qp3qouwMGgz0N6yM20V1tLXYOrVENcMet4AKSOhdiO3I1VlCKnEoVTAHct2Gn9fhb0/QFRt\n6DEC0m6Exl3djq7qiEuEAX+Gs++CDR/C4gkw/1FY8AR0uNiehST31S69QUoTh1LeyFxmzy7WvA25\nR6FJdxjyDHS5HKJruR1d1RUeCZ2G2tu+zba6b/kbsHY21G9nG9m7X2O7H6ugEVKJQ9s4lF8dPwxr\nZtqEsXslRNaArpdDzxshMdXt6EJP/TZw/mNwzv22+m/xBPj4Hvj0r9D1MptEEnu6HaUixBKHtnEo\nv9iz2lZFrZoOJw5Bw85w4VO2nSJGR14OuMhY6HGtve1eaRPI6hn2TKRJD1uN1eXyX7cjqUoTUolD\nqXLLzbG/cpdMtI3e4dHQ5VJ7dtGst9a1u6VJd7jkWTjvEZvIF0+A926FufdDj2vsWUiD9m5HWe1o\n4lDV288b7NnFyim24Tuhrb3uoPvVOmBfMImJg97joNdY2PGdTSCLJ8DCF2wjetqN0GEIRES5HWm1\noIlDVT95x2Hte7Yhdvs3EBYJnS6xv15bnKVnF8FMBFqcaW+Hn4Dlr9vPceZoqNnw1IWF8c3djjSk\naeJQ1cf+LfYgs2IKHN0PdVvCoL/a7rS1GrgdnfJVrQbQ904463bY8pk9A/n6aXtre569CLPNQL2w\nMABCKnForyr1K/m5sP4D23bx4xcg4dDhIlu10bJ/xUeLVe4LC4O259pb1k47PtiyybDxY3vm0XOU\nvbBQfxz4TUglDu1VpU76Zbu9MnnZ63Bkr70y+Zz7IeV6qN3Y7ehUoMQ3g4EPQL8/wfr37Q+Gzx6G\nzx93qiPH2GourY6skJBKHKqay8+DTXNtY/fmT+3Boe359uyizSCtsqhOIqJsr7gul8LPG20CWTHF\nXsDZoINzYeHVsHFu5czMGGI0caiqLzvTVk0smwyHdkHtJtDvbkgdaYcuV9Vbg3ZwwRM2Kax52w71\n/tHdMPc+MAVg8r2fmVEBmjhUVVWQD1vm21+SGz+2gw62PgcufBLaDYZw/WqrIqJq2F5XqdfbIWQm\nXQy5R07fJjcHPn1IE0cZ9L9LBa9V039djdCyn9MF8zXI3gE1G9ghu1NvgHot3Y5YVRWJqXbMseIc\nzLRJpc1AW8XZqIu2iRShiUMFp1XTT85tfbIa4Z3xtmoBAy3PhnP/akdS1Yu+VHnEJdnvVVHRtSEn\ny555fPoQ1GpsE0ibgdCqv14YSoglDu2OG0I+e/j0GfMATL4dxvymBXZAPKUqYuBfTv44OSkyFi76\nl62qOrjbXh+y+VNYP8fOYihhkJhmu/62GQhNUqpll+6QShzaHTcEHNwNP8wq/pcgwInDmjSUfxS2\nY5Q0M2OdJpBynb3l58GuZTaJbJoHn/8NPn8MaiTYtrU259q/1eRakZBKHKqKOnoA1r0Hq2fCtq8B\nY4cBKcj99bbaS0r5k7czM4ZH2MEum/W2E1Ad2QdbPreJZPOndvResKP3thlkb0m9QraTRmiWSgW/\nE0dgw0f2H27zZzZJJLSxF251vRx2LS++GmHgX9yLWalCNetDtyvsraAA9qx0kshndsiTr56C6Dho\n3d8mkdYD7ayHIUITh6o8eSdsnfHqGTZp5B6F2k0h/WboeoUdQruw90r9tvZvSdUISgWLsDBommJv\nZ99lG9a3LjiVSNa+a7dr2OnU2UjzPhAR7WrYFaGJQwVWQb4dgXb1TPsPdCwLYutCt6vsmUXzM0tu\nXPS2GkGpYBIbD52H2ZsxsHftqSqt75+Hb5+FyJq2Z2BbJ5HUTXY7ap9o4lD+Z4xtSFz9tr1S9/Ae\n+4/S4SKbLFoN0C60qnoQgUad7e2s2+x0xNu+sg3sm+fBxo/sdgltTp2NJP/GVssGMU0cyn9+3mDP\nLNbMhANbITzK9jbpehm0u0Cn+lQquha0v8DejLFD/W/+1CaRpZPsxFQRMXZemMJEUr9t0F2AqIlD\nVUzWTntWsXom/LTa9nNP7gu/uQM6DrHVUkqpXxOxXcvrt4E+421HkO3f2HaRTfNg7r32Ft/8VBJp\neba9QNFlmjiU747ss/Nzr54JO7+3yxLTYPAT0Hm4DluuVHlExp5KEIMfh1+22SSy+TM7ksKSibab\nevM+p7Zr1NkmoOKG5wlgR5KQShx65XgAHTtoJ0RaPcP2GDH5dnjqc+6HLpdBvVZuR6hUaKmbDL3G\n2FveCfsjrbCn1qcP2lvtJlC3FWQuhvwTlTbKb0glDr1y3M9yj8GmT2yy2PQJ5B2DuOZw1h+gy+Wn\nfu0opQIrIspWU7U8G859GA7usqNDb5rndPc1p2+fm2OH7dHEoSpFfp6dYnX1TDuD2vGDdgTa1JH2\nWoukXposlHJbnaanhkN5KL74bbIzAvbymjiU7d2xc5HtDfXDO3DkZ4iuYxu3u1xmhzIP0aETlKry\nShrlN4DD8+jRoLoyBn76wSaL1W/buS0iYqDd+bYaqu15EBnjdpRKqbKUNMpvAIfn0cRR3Rz40UkW\nM+Hn9SDh0HqAHbitw0UQU8ftCJVSvihrlN8A0MQRiop2zTvrNijIs43cmUvtNs3PgAufst1na9Z3\nN16lVMVU8vA8mjhCTXEz5334f3Zd464w6K+23SK+mZtRKqWqME0cwaqgwPZoOpZtBwbMyTp1/1i2\n87iY+/s3O9OrFlGrMYz/uvLLoZQKOZo4CgXiysu8EyUc+LNKPvAXbnf8YPEJoJCEQUycc4u3I3LW\naQr7Nha//eGfKlYWpZRyaOKA4qt3Cq+8bH+hb7/4PbfLPVr660bEnH7gr9UI6re3y2Lj7fJi78fZ\nubeLG4786S6V3jVPKVW9aOIAe4WlZ1c2sI9neXEBenSd0w/qCa1PJYKiB/uTy537geju6kLXPKVU\n9aKJA0q/wnLQX0s+8MfEQVh45cXpDRe65imlqpegTxwi0gq4D4gzxlwekBcp8crLZvCb2wPykgGl\nM+cppQKohDk7/UNEJorIXhFZU2T5YBHZICKbReSe0vZhjNlqjBkTyDgZ+Jdfz7il1TtKKVWsQJ9x\nTAKeAyYXLhCRcOC/wLlABrBYRN4DwoHHizx/tDFmb4Bj1OodpZTyQUAThzHmSxFJLrK4N7DZGLMV\nQESmAUONMY8DFwcynlJp9Y5SSnlFjDFlb1WRF7CJ431jTBfn8eXAYGPMWOfx9UC6MeaWEp6fADyG\nPUN5xUkwxW13E3ATQKNGjXpOmzatXPEePnyYWrVqleu5wSZUyhIq5QAtS7AKlbJUpBwDBgxYaoxJ\n82ZbNxrHi5vMocTsZYzZD4wva6fGmJeAlwDS0tJMec8aFoTQGUeolCVUygFalmAVKmWprHIEtHG8\nBBmA50BJScAuf+xYRIaIyEvZ2dn+2J1SSqliuJE4FgNtRaSliEQBVwPv+WPHxpg5xpib4uLi/LE7\npZRSxQh0d9ypwHdAexHJEJExxpg84BZgLrAOmG6M+SGQcSillPKfQPequqaE5R8CH/r79URkCDCk\nTZs2/t61UkopR8B7VblBRLKBTR6L4oBsL+/XB/aV86U99+frNsUtL7qstMeF9z2XVcWy+PszKS1O\nb7bxtSzB+v0qaZ2vZanq3y/P+1WxLIH8frXAXjs3p8wtjTEhdwNeKulxWfeBJf56XV+2KW55aeUo\nJX7PZVWuLP7+TCq7LMH6/fJXWar696uqlyWQ3y9vy2KMcaVxvDIUzZhzfLzvr9f1ZZvilpdWjqKP\n55SwTXm5VRZ/fybe7sdfZQnW71dJ63wtS1X8TIo+rsplCeT3y+v9hGRVVUWIyBLj5UUwwS5UyhIq\n5QAtS7AKlbJUVjlC9YyjIl5yOwA/CpWyhEo5QMsSrEKlLJVSDj3jUEop5RM941BKKeUTTRxKKaV8\noolDKaWUTzRxlEFEWonIBBGZ6XYsFSEiw0TkZRF5V0TOczueihCRjiLygojMFJHfuh1PRYlITRFZ\nKiLuzUfjByLSX0S+cj6b/m7HU14iEiYij4nIf0TkBrfjqQgR6et8Hq+IyLf+2m+1TBy+TGlrKmPq\n2nLysRyzjTHjgFHAVS6EWyofy7LOGDMeuBIIui6U5Zgy+U/A9MqN0js+lsUAh4EY7CjYQcPHcgwF\nEoFcgqwc4PP/ylfO/8r7wGt+C6IiVxlW1RtwNpAKrPFYFg5sAVoBUcBKoJPH+plux+2ncvwTSHU7\n9oqWBbgE+Ba41u3YK1IWYBB2hOhRwMVux17BsoQ56xsBb7odewXKcQ9ws7NNqPzfTwfq+CuGannG\nYYz5EjhQZPHJKW2NMSeAadhfHkHLl3KI9XfgI2PMssqOtSy+fibGmPeMMWcCIyo30rL5WJYBQB/g\nWmCciATV/6QvZTHGFDjrfwGiKzHMMvn4mWRgywCQX3lResfX/xURaQ5kG2MO+isGN2YADFaJwE6P\nxxlAusfUtSkicq8pYeraIFJsOYBbsb9u40SkjTHmBTeC81FJn0l/4FLswcnvoywHSLFlMc6UySIy\nCtjncfANZiV9LpcC5wPxwHNuBOajkv5XngH+IyJ9gS/dCKwcSioLwBjgVX++mCaOU4qd0tZ4OXVt\nECmpHM8Cz1Z2MBVUUlkWAAsqN5QKK3XKZGPMpMoLpcJK+lxmAbMqO5gKKKkcR7EH26qkxO+XMeZB\nf79YUJ0WuyxgU9pWslApB2hZglWolCVUygGVXBZNHKcEbErbShYq5QAtS7AKlbKESjmgssvidg8B\nl3olTAV2c6q73Rhn+YXARmzvhPvcjrO6lEPLEry3UClLqJQjWMqigxwqpZTyiVZVKaWU8okmDqWU\nUj7RxKGUUsonmjiUUkr5RBOHUkopn2jiUEop5RNNHEoVISKHA7DPZBG51sfnfCgi8f6ORamK0sSh\nVOVIxo6A6zVjzIXGmKzAhKNU+WniUKoEzox2C5yZBteLyJsiIs66bSLydxFZ5NzaOMsnicjlHvso\nPHt5AugrIitE5I4ir9NERL501q1xRmUtfI36IjLeWbdCRH4Ukc+d9eeJyHciskxEZohIrcp4X5TS\nxKFU6VKA27ET/LQCzvJYd9AY0xs7hPi/y9jPPcBXxpgexpini6y7FphrjOkBdAdWeK40xrzgrOuF\nHWLiXyJSH7gfGGSMSQWWAHeWp4BK+UqHVVeqdIuMMRkAIrICW+X0tbNuqsffosnAF4uBiSISCcw2\nxqwoYbtngPnGmDli59Y1QrsAAAEsSURBVCfvBHzjnARFAd9VIAalvKaJQ6nSHfe4n8/p/zOmmPt5\nOGfyTrVWVFkvYIz5UkTOBi4CXheRJ40xkz23cSZ6agHcUrgImGeMucb7oijlH1pVpVT5XeXxt/DX\n/jagp3N/KBDp3D8E1C5uJyLSAthrjHkZmICdT9pzfU/g/4DrzKkZAr8HzvJoW6khIu0qWiClvKFn\nHEqVX7SILMT+ACv85f8y8K6ILAI+A444y1cBeSKyEphUpJ2jP3CXiOQCh4GRRV7nFqAe8LlTLbXE\nGDPWOQuZKiKF83vfjx1WW6mA0mHVlSoHEdkGpBlj9rkdi1KVTauqlFJK+UTPOJRSSvlEzziUUkr5\nRBOHUkopn2jiUEop5RNNHEoppXyiiUMppZRPNHEopZTyyf8DCcCCNXkwqcMAAAAASUVORK5CYII=\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "rp, ep = rat_err(tla, tlea, tva, tvea)\n", "rn, en = rat_err(tna, tnea, tva, tvea)\n", "\n", "plt.errorbar(sizes, rp, yerr=ep, fmt = 'o-', label=\"Python\")\n", "plt.errorbar(sizes, rn, yerr=en, fmt = 'o-', label=\"Numba\")\n", "plt.xscale('log')\n", "plt.yscale('log')\n", "\n", "plt.grid(True)\n", "plt.axhline(1, color='red')\n", "plt.legend()\n", "plt.xlabel(\"Input size\")\n", "plt.ylabel(\"Timing ratios\")\n", "plt.title(\"Performance relative to vectorized Numpy\");" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's quickly see the Numba/Numpy ratios (the raw numbers):" ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[ 0.12029788 0.13762593 0.28367041 0.65578777 0.76074653 0.32554678\n", " 0.18141704]\n", "[ 8.31269874 7.26607237 3.52521785 1.52488357 1.31449828 3.07175515\n", " 5.51216148]\n" ] } ], "source": [ "print(tna/tva)\n", "print(tva/tna) # Let's look at the reverse relation, which may be a bit easier to interpret" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Some quick thoughts:\n", "\n", "- The Numba performance is *always* better than Numpy. Note that this is a pretty ideal case for Numpy, and we did exactly *nothing* to the Python code other than compile it with `numba.jit`. We're basically getting better-than-numpy performance for free. Nubma can be a bit like magic ([this blog post](https://jakevdp.github.io/blog/2015/02/24/optimizing-python-with-numpy-and-numba/) by Jake VanderPlas contains much more detail).\n", "\n", "- In some cases, the difference is quite remarkable. On the other hand, where Numpy performs best comapred to raw Python, Numba still beats it but not by too much. This probably indicates Numpy is already nearly optimal in that range.\n", "\n", "\n", "Let's look at the actual times instead, which in this case can be informative:" ] }, { "cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY4AAAEaCAYAAAAG87ApAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3Xd4lFX2wPHvSQgh9CYBEpSOIkSB\nACqLBJUFCyuiYgEsCOiuZS2Lyuqqy66rrr39VhEQK4iACIiiIsVGLwFEkCoJvSQQSJ/z++OdQAgp\nM8lMZjI5n+eZJ/P2c1Pm5L73vveKqmKMMcZ4KizQARhjjKlYLHEYY4zxiiUOY4wxXrHEYYwxxiuW\nOIwxxnjFEocxxhivWOIwppREpJqIqIjEBjCG70TkhkBd31ROljhMSBGRtHwvl4ik51seXMKx/URk\nc3nF6i0ReVZExuVfp6qXqOongYrJVE5VAh2AMb6kqjXz3ovIdmC4qn4buIhOxBIGoKquQMdiTFlZ\njcNUKiISJSJvishuEUkSkedFJEJEGgCfAS3z1VAaiEgPEVkiIqkisktEXhYRj/7hEpHFIjJGRJYA\nx4GmIrJHRP6Qb58TtQgROVtEckTkdnds+0VklHvbAOBB4FZ3bEvzXWOI+/1d7ltXb7jj/U1E4kVk\npIgki8heEbmxwPfiFRHZ6Y7rdRGJ9M132oQySxymsvknEAd0BLoACcDDqnoQuAbYqqo13a+DQDZw\nD1Af6An0B4Z7cb0hwC1ALWCPB/uHA/FAa+AK4GkRaamqM4CXgPfcsXUr4viewE9AA2AGMA04B2gB\njAD+JyLV3Pu+DMTifC/aAW2BR70om6mkLHGYymYw8KSqHlDVvcC/gaFF7ayqS1V1marmquoWYBzQ\ny4vrjVPVjaqarao5Hh7zpKpmqOoy4FecROepX1X1Y/e1pgBnAk+papaqzgSqAs3dtaZhwF9VNUVV\nU4FngRuLPLMxbtbGYSoNERGgMbAj3+odQEwxx7QHXgQ6A1E4fzM/enHZnV6GmauqB/ItHwdqFrVz\nIfbme58OZLqTQv51NYGmQASw3vm2ACCAp8nNVGJW4zCVhjpDQe8Bzsq3+kwgOW+XQg57B1gJtFLV\n2sAYnA9Yjy9bYPkYUD3fcuMynKssduMkiVaqWtf9qqOqDXx4DROiLHGYymYS8KS74bsR8BjwoXvb\nXqCRiOT/D78WkKqqaSJyLk47QVmsBm4SkSoicgFwtRfH7gVaSL4qQmmpajYwAXhVRBqKo5mI9Cnr\nuU3os8RhKpsngF+A9Tgf4j8C/3VvWwPMBHaISIqI1AceAIaLSBrwJlDWZyb+jtMYnQKMBiZ7cexk\nnNrKIRH5qYxxANwP7AKWA6nAVziN8sYUS2wiJ2OMMd6wGocxxhivWOIwxhjjFUscxhhjvGKJwxhj\njFcscRhjjPFKSD453rBhQ23evHmpjj127Bg1atTwbUABYmUJPqFSDrCyBKOylmPFihUHVPWMkvYL\n+sQhIi1xHtKqo6rXeXJM8+bNWb58eamut2DBAhISEkp1bLCxsgSfUCkHWFmCUVnLISI7St7Lz7eq\nRGSCiOwTkXUF1vcTkY0isllEih2NU1W3quod/ozTGGOM5/xd45gIvAG8n7dCRMJxnsDtAyQBy0Rk\nJs5w0s8UOH6Yqu7zc4zGGGO84Pcnx0WkOTBbVTu4ly/EGea5r3t5NICqFkwaBc8ztbhbVSIyEhgJ\nEB0d3WXyZG9GcjgpLS2NmjW9GYw0eFlZgk+olAOsLMGorOXo3bv3ClWNL2m/QLRxxHDqUNNJQPei\ndnbPzPY00ElERheVYFR1LDAWID4+Xgve58vOziYpKYmMjIxig6tTpw7VqlUrdp+KwtOyVKtWjdjY\nWCIiIsohqtKxe9DBx8oSfMqrHIFIHIWN7Flktcc9C9tdHp1YpD/Qv3Xr08dpS0pKolatWjRv3pzi\nBhc9evQotWrV8uRyQc+TsqgqBw8eJCkpiRYtWpRTZMaYiiwQz3EkAc3yLcfijNBZZqo6S1VH1qlT\n57RtGRkZNGjQoNikURmJCA0aNCixJmaMCV4zViXT49nvuO2rY/R49jtmrEou+aAyCETiWAa0EZEW\nIlIVZ6rKmb44sYj0F5GxqampRW336nw3vP0zN7z9sy9CC2qWTI2puGasSmb09LUkp6QDkJySzujp\na/2aPPzdHXcS8DPQTkSSROQO91zI9wBzgQ3AFFVd74vrFVfjCAYiwkMPPXRi+YUXXuCpp54KXEDG\nmArv+bkbSc/OPWVdenYuz8/d6Ldr+jVxqOpNqtpEVSNUNVZVx7vXz1HVtqraSlWf9tX1SqpxeGPG\nqmRW/Z7Ckm2HfFb1i4yMZPr06Rw4cKDknY0xxgO73DUNT9f7QkiNVeWrGkde1S8r1wX4rupXpUoV\nRo4cycsvv3zatttuu42pU6eeWM7rUrdgwQJ69erFoEGDaNu2LY8++igfffQR3bp1o2PHjmzZsuXE\n8XfddRc9e/akbdu2zJ49G4CePXuyevXqE+ft0aMHiYmJZSqHMSY4qCrVI8ML3da0bpTfrhv0Q474\nwz9nreeXXUcK3Zabm0ti8tETSSNPenYuD09NZNLS3ws9rn3T2jzZ/9wSr3333XcTFxfHww8/7HG8\na9asYcOGDdSvX5+WLVsyfPhwli5dyquvvsrrr7/OK6+8AsD27dtZuHAhW7ZsoXfv3qxatYrhw4cz\nceJEXnnlFTZt2kRmZiZxcXEeX9sYE5xcLuWJmes4lplLeJiQ6zrZOTUqIpxRfdv57dohVePw1a2q\ngkmjpPXeqF27Nrfccguvvfaax8d07dqVJk2aEBkZSatWrfjjH/8IQMeOHdm+ffuJ/QYNGkRYWBht\n2rShZcuWbNq0ieuvv57Zs2eTnZ3NhAkTuO2228pcBmNMYOW6lEenJ/Lh4t+58+KWvHBdHDHuGkZM\n3SieGdiRAZ1i/Hb9kKpxqOosYFZ8fPyI4vYrrmZw9OhR+r257EQPhfxi6kbxyZ0XljnO+++/n86d\nO3P77befWFelShVcLicxqSpZWVkntkVGRp54HxYWdmI5LCyMnJycE9sK9o4SEapXr06fPn34/PPP\nmTJlSqkHfzTGBIecXBcPfbqGz1fv4r5L2/DAZW0QEa7pHFtuDwBajaMQo/q2Iyri1PuGvqz61a9f\nn0GDBjF+/PgT65o3b86KFSsA+Pzzz8nOzvb6vJ9++ikul4stW7awdetW2rRpA8Dw4cO577776Nq1\nK/Xr1/dJGYwx5S8rx8W9k1bx+epdjOrbjgf7tA1Id/qQShy+ahwf0CmGZwZ2pGq48+3xR9XvoYce\nOqV31YgRI1i4cCHdunVjyZIlpRpTv127dvTq1YvLL7+ct95668RwI126dKF27dqn1HCMMRVLZk4u\nf/loBV+u28PjV57D3b1PHyGjvITUrSpfGtAp5kRDuC9uT4EzAFme6Ohojh8/fsry4sWLTyw/84wz\nJFdCQsIpVc8FCxaceF9wW48ePU7psXX06FEAdu3ahcvlOtE2YoypWNKzcrnzwxUs2rSff119LkMv\nbB7QeCxxFMNXCSOQ3n//fR577DFeeuklwsJCqoJpTKVwLDOH4e8tZ/G2gzx3bUdu6HpmoEMKrcRR\n3CCHoW7ixImFrr/lllu45ZZbyjcYY4xPHM3I5vZ3l7Hy98O8NOg8rukUG+iQAGvjMMaYoJR6PJsh\n45eyemcKr9/UOWiSBoRYjcMYY0LBoWNZDBm3hM370vjfkC70aR8d6JBOYYnDGGOCyP6jmQwet5gd\nB48z9pYuJLRrFOiQTmOJozjvXul8vf2LwMZhjKkU9qRmcPO4xexOyeDd27pyUeuGgQ6pUCHVxuHL\n0XH9wZ/DqhccJNEYU7EkHT7OoLd/Zt+RTN6/o1vQJg0IscTh08bxxCmQtAx2/AAvd3CWy8iGVTfG\nFGbHwWPc8PZiUo5n8cEd3ejaPLhHeAipxOEziVNg1n2Qm+ksp+50lsuYPPw5rDrAt99+e9qw6tu3\nb6dnz5507tyZzp0789NPP5WpDMYY39q8L41Bb//M8awcPh5xAZ3OrBfokEpUOds4vnwU9qwtdFNU\nbg7sXnUyaeTJTofP74EV7xV+zsYd4fJnS7x0eQ+r3qhRI7755huqVavGb7/9xk033WQDHRoTJDbu\nOcrgcc6IEZNHXki7xrUCHJFnKmfiKEnBpFHSei/kH1Y9KsqziVbyhlUHThtWff78+Sf2K2xY9Q4d\nOnDPPfewevVqwsPD2bRpU5nLYIwpu3XJqQwdv4SqVcL4aPgFtG5UM9AheaxyJo5iagbpR49Sa9yF\nzu2pguo080kPq/IcVv3ll18mOjqaNWvW4HK5Tgx8aIwJnFW/H+bWCUupVS2Cj0d056wG3g9qGkgh\n1cbhs15Vlz4BEQVqAxFRznofKM9h1VNTU2nSpAlhYWF88MEH5ObmlnwiY4zfLNt+iKHjl1K3elU+\nufOCCpc0IMQSh896VcUNgv6vQbj7P/06zZzluEFlD9KtvIZV/8tf/sJ7773HBRdcwKZNm0p1XmOM\nb/y0+QC3jF9Ko1qRTLnzQmLrVQ90SKVSOW9VeSJu0MmGcB89AOjPYdULG+Tw6NGjtGnThsTExNPO\na4wpXws27uPOD1ZwVoPqfDi8O41qVdzbxpY4imNPjBtjfOCbX/Zy90crad2oJh8O7079GlUDHVKZ\nhNStKmOMCTZfJO7mzx+u4JwmtZg04oIKnzTAEocxxvjNZ6uSuHfSSs5vVpcPh3enTvUI/1wocQq8\n3IFeCwb4bKSL4lSqW1WqGpCJ3YOdqgY6BGNCzpRlO3lkeiIXtGjAuFvjqRHpp4/bvJEustMRODnS\nBfi0Q09+labGUa1aNQ4ePGgfkgWoKgcPHrTnO4zxoQ9+3s7D0xLp2eYM3r29q/+SBsC8Mc7IFvll\npzvr/aRC1DhEZABwJdAIeFNVv/b2HLGxsSQlJbF///5i98vIyAiZD1FPy1KtWjViY4NndjFjKrJx\n32/l319s4LJzGvHm4M5EVgn37wVTk7xb7wN+TxwiMgG4Ctinqh3yre8HvAqEA+NUtcjHuVV1BjBD\nROoBLwBeJ46IiAhatGhR4n4LFiygU6dO3p4+KIVSWYypCN6cv5nn527kio6NeeWGTlSt4uebOus/\nK3pbHf/9M1geNY6JwBvA+3krRCQceBPoAyQBy0RkJk4SKfigwTBV3ed+/7j7OGOMCRqqysvf/sZr\n835jwPlNeeH686gS7sekkZsN3z4FP78B9VrA0d2Qk3Fyuw9HuiiM3xOHqi4SkeYFVncDNqvqVgAR\nmQxcrarP4NROTiFOi/azwJequtK/ERtjjOdUlWe/+pW3F25lUHwszwyMIzzMj51wju6FqbfDjh+h\n6wjo+x/4ZQbMG4OmJiF1Yp2k4aeGcQApj8Zid+KYnXerSkSuA/qp6nD38lCgu6reU8Tx9wG3AsuA\n1ar6ViH7jARGAkRHR3eZPHlyqWJNS0s7MRdGRWdlCT6hUg6wsoCTND7+NYtvduRwSbMqDGlflTA/\n9tysnbqBc9c/R5WcY2xqezd7Gyecsr2sP5PevXuvUNX4EndUVb+/gObAunzL1+O0a+QtDwVe99X1\nunTpoqU1f/78Uh8bbKwswSdUyqFqZcnNdemj0xL1rEdm65hZ69Xlcvk+sDwul+rP/1P9Z33VV85T\n3b220N3K+jMBlqsHn7GB6lWVBDTLtxwL7CrrSUWkP9C/devWZT2VMcYUKdelPDw1kWkrk/hLQitG\n9W3nv2fEMtNg1l9h3VRoezlc8xZE1fXPtTwUqMSxDGgjIi2AZOBG4OYAxWKMMR7LznXx4JQ1zFqz\niwcua8t9l7b2X9I4sBk+GQL7f4VL/gF/eBDCAv/4nd8jEJFJwM9AOxFJEpE7VDUHuAeYC2wApqjq\n+rJeS301rLoxxhQiK8fFPR+vZNaaXTx6+dn89bI2/ksaG2bB2ARI2wtDp8PFfwuKpAHl06vqpiLW\nzwHm+Pv6xhjjCxnZufzlo5V89+s+nuzfntt7lPxcWKnk5sB3/4IfX4GmnWDQB1C3WcnHlaMK8eS4\np6yNwxjjD+lZuYz8YDnf/3aAp6/pwODuZ/nnQmn7Ydow2LYIutwG/Z6DiOAbySI46j0+YreqjDG+\nlpaZw63vLuXHzQd44frz/Jc0di6Dty+GnUvh6jeh/6tBmTQgxBKHz+YcN8YY4EhGNreMX8KKHYd5\n+Ybzua6LH4bxUIVl4+DdyyE8Au74GjoN8f11fCikEofVOIwxvpJyPIsh45awNjmVN2/uxNXnx/j+\nIlnH4bO74IuHoGUCjFwATc7z/XV8LKTaOIwxxhcOpmUyeNwSth44xttDu3DJ2dG+v8ihrfDJUNi7\nHhJGw8UPB02vqZKEVOKwxnFjTFntO5LBzeOWkHT4OONvjadnmzN8f5GNX8L0O0EEBn8Kbfr4/hp+\nVDHSm4fsVpUxpix2paQz6O2f2ZWSzsTbu/k+abhy4bt/w6Qbod5ZcOfCCpc0IMRqHMYYU1o7Dx3n\npncWk3o8mw/u6EaXs+r79gLHDsL04bDlOzh/CFz5gjP8eQVkicMYU2nNWJXM83M3kpySTtjc+URW\nCeOTOy8kLtbHY0Elr4QptzhPgfd/FTrf6tymqqBC6laVdcc1xnhqxqpkRk9fS3KKM1+3S53X1v3H\nfHuhFe/BhL7O+2FfOQ/2VeCkASGWOKyNwxjjqefnbiQ9O/eUdZk5Lp6fu9E3F8hOh8/vhln3wVk9\nYORCiOnim3MHmN2qMsZUSrvcNQ1P13vl8Hbn1tTuNdDzb9D77xAWXvbzBglLHMaYSmfnoeOECeQW\nMgFq07plbLD+7RuYNtx5IvymydDu8rKdLwhZ4jDGVCq7U9O5edxiIsKFKgiZOa4T26IiwhnVt13p\nTuxywaLnYcEzEH0uDHofGrTyUdTBJaQShz0AaIwpzr6jGQx+ZwmHj2XzyZ0Xse3AsRO9qmLqRjGq\nbzsGdCrF0CLHD8Fnd8JvX0PcDXDVK1C1uu8LECRCKnGo6ixgVnx8/IhAx2KMCS4H0zIZ/M4S9hzJ\n4P1h3TivWV3Oa1aXAZ1iWLBgAQkJCaU78e41ztAhR3bBFS9A1+EVvtdUSUIqcRhjTGFSjmcxdPxS\nfj90nIm3dyO+uY8e7lv1EXzxIETVh9u/hGZdfXPeIGeJwxgT0o5kZHPrhKVs3pfGO7fGc2GrBmU/\naU4mfPkIrHgXmveE696Fmn4Y0ypIWeIwxoSsY5k5DHt3Get3HeGtIV3o1dYHH+4pO52utrtWQo/7\n4ZJ/QHjl+iitXKU1xlQaGdm5DH9vOSt/P8wbN3fmsvY+GBp9y3yYOgxys+GGD+Gc/mU/ZwUUUk+O\n25AjxhiAzJxcRn6wgsXbDvLSoPO5omOTsp3Q5YJFL8CHA6FmtDPhUiVNGhBiicOGHDHGZOW4uPuj\nlSzatJ/nBsaVrnttfukp8Mlg+O5fcO41MPxbaFi5u/zbrSpjTMjIyXVx/yer+HbDPv519bkM6tqs\nbCfcsw6mDIWU36Hfc9D9zpDvausJSxzGmJCQ61L+9uka5qzdw+NXnsPQC5uX7YRrPoFZf4VqdeC2\nL+DMC3wSZyiwxGGMqfBcLuXv09cyY/UuRvVtx/CeLUt/spws+PoxWDrWGdX2unehlh/mHK/ALHEY\nYyo0VeXJmev5ZPlO7rukNXf3LkP7Q2oyfHobJC2FC++By56C8AgfRRo6LHEYYyosVeXpLzbwweId\n3HlxSx7o09a7EyROgXlj6JWaBMsaQnYGqMupZXQY6J+gQ4AlDmNMhfXi15sY98M2bruoOY9efjbi\nTcN14hRnkqXsdATg2H5AoM8/LWmUIOi744rIOSLylohMFZE/BzoeY0xweH3eb7wxfzM3dWvGE1e1\n9y5pAMwb48zSdwqFpe/4LMZQ5dfEISITRGSfiKwrsL6fiGwUkc0i8mhx51DVDap6FzAIiPdnvMaY\nimHsoi28+M0mBnaK4ekBHQkLK0UX2dQk79abE/xd45gI9Mu/QkTCgTeBy4H2wE0i0l5EOorI7AKv\nRu5j/gT8AMzzc7zGmCD33k/b+c+cX7kyrgn/vS7O+6ShCismAoVM/wdQJ7asIYY8US3im+erC4g0\nB2aragf38oXAU6ra1708GkBVn/HgXF+o6pVFbBsJjASIjo7uMnny5FLFm5aWRs2aNUt1bLCxsgSf\nUCkHBKYsC3ZmM3F9Fp0ahXP3+ZFU8TJphOVm0HbT/2i8dwFp1ZsRlbGXcFfWie25YZFsbHc3+6J7\n+Tr0clHWn0nv3r1XqGqJd3YC0TgeA+zMt5wEdC9qZxFJAAYCkcCcovZT1bHAWID4+Hgt7aQsZZrQ\nJchYWYJPqJQDyr8s01cm8d4va+jV9gzG3tKFyCrh3p1g/0ZnVNv9GyFhNDUvHgXrpsG8MWhqElIn\nlvBLn6B93CDa+6cIfldeP5NAJI7C/kUostqjqguABR6d2KaONSYkzU7cxd8+XcOFLRvw9tBSJI3E\nT52nwCOiYOhn0Kq3sz5uEMQNYmEIJfTyEIheVUlA/gFkYoFdvjixDXJoTOj5ev0e7p+8mi5n1WPc\nrfFUi/AiaWRnwKz7YfpwaBIHd31/MmmYUgtEjWMZ0EZEWgDJwI3Azb44sdU4jAkt8zfu4+6PV9Ih\npg4TbutK9apefGQd2gaf3urMCd7jr+4Jl+wpcF/wd3fcScDPQDsRSRKRO1Q1B7gHmAtsAKao6npf\nXM9qHMaEjp82H+CuD1bQNroW793ejVrVvPjQ3zAL3u4Fh7fDjZOgzxhLGj7k1xqHqt5UxPo5FNPQ\nXVpW4zAmNCzbfog73ltO8wY1+OCO7tSp7uGHfm42fPsU/PwGNO0E178H9c7ya6yVUdA/Oe4Nq3EY\nU/Gt+v0wt7+7jCZ1q/Hh8O7Ur1HVswNTk+DdK5yk0W0kDJtrScNPbKwqY0zQWJecyq0TllK/RlU+\nHn4BZ9SK9OzAzd/CtBGQmwXXTYAO1/o30EoupGocNue4MRXXxj1HGTp+CbWqRfDxiO40rlOt5INc\nufDdv+HD66BWExi50JJGOQipxGG3qoypmDbvS2PwuMVUrRLGR8O7E1uveskHHd0L718Ni56H8wfb\nXODlqMRbVe4hQoYAPYEmQDqwDvgC+FBV7d97Y0yp7Th4jMHjFgPw0fALaN6wRskHbf8Bpg6DjCNw\n9ZvQaYifozT5FVvjEJEvgeE4XWf74SSO9sDjQDXgc/cAhEHBblUZU7EkHT7Oze8sITPHxYfDu9O6\nUQnjLLlc8P2L8F5/iKwFI+ZZ0giAkmocQ1X1QIF1acBK9+tFEWnol8hKQVVnAbPi4+NHBDoWY0zx\n9qRmMHjcEo5kZDNpxAWc3bh28QccPwSf3QW/zYVzB8KfXnOShyl3xSaOvKQhIjWAdFV1iUhb4Gzg\nS1XNLiSxGGNMsfYfzWTwuMUcOJrJh8O70yGmhHbJpOXOXOBpe+GKF6DrcPB24ibjM542ji8CqolI\nDM6cGLfjzLURVOxWlTHB79CxLIaMW8KulAzevb0bnc6sV/TOqrD4LZjQz0kUw+ZCtxGWNALM08Qh\nqnocZ3jz11X1Ggi+kYetV5UxwS01PZuh45ew7eAxxt0aT7cW9YveOeOIM9bUV49A68vgzkUQ07n8\ngjVF8vQBQHH3rhoM3OHlscYYQ1pmDrdOWMqmvUcZe0s8PVoX0zy6Z60zd8bhHc44UxfdZ7WMIOLp\nh/9fgdHAZ6q6XkRaAvP9F5YxJpQcz8ph2LvLWJucyv8N7kzvdo0K31EVVn0Ac0ZBVD247Qs468Ly\nDdaUyKPEoaqLcNo58pa3Avf5KyhjTOjIyM5lxPvLWb7jEK/e2Im+5zYufMesY/DFQ7BmErRMgIHj\noOYZ5Rmq8VBJz3GMFZGORWyrISLDRGSwf0LznjWOGxNcMnNy+fOHK/hpy0Gev+48+p/XtPAd92+C\ndy6FNZMhYTQMmW5JI4iVVOP4P+Af7uSxDtiP8+BfG6A2MAH4yK8ResGe4zAmeGTnurj341XM37if\n/1zTkWu7xBa+49qpMPM+97Su06HVJeUbqPFaSc9xrAYGiUhNIJ6TQ45sUNWN5RCfMaYCynUpD3yy\nmq9/2ctT/dtzc/czT98pOwPm/h2Wj4czL3RGta1dRI3EBBVP2zjSgAX+DcUYEwpcLmXU1DXMTtzN\n6MvP5rYeLU7fyaZ1rdCsS60xxmdUlcdmrGX6ymQe7NOWO3u1On2nDbNhxl9AcKZ1PfuKco/TlI0l\nDmOMT6gq/5z1C5OW7uQvCa2495ICQ5zbtK4hw6vEISI1VPWYv4IxxlRMqsqzX/7KxJ+2c8cfWjCq\nbzsk/wN7qckw9XbYucSZ1vWP/4YqHs7uZ4KOR0OOiMhFIvILsMG9fJ6I/J9fIysF645rTGC8/O1v\nvL1oK0MuOJPHrzzn1KSxeR683RP2rncawK943pJGBefpWFUvA32BgwCquga42F9BlZaNVWVM+Xtz\n/mZem/cbg+JjGfOnDieThisXvnsaPrwWajaGkQtsWtcQ4fGtKlXdKaeOFZPr+3CMMcFuxqpknp+7\nkeSUdOosnEtqeg5Xn9+UZwbGERbm/oxI2wfThsO2hXD+EKeWUdWD6WBNheBp4tgpIhcBKiJVcYYb\n2eC/sIwxwWjGqmRGT19Lerbzf2Nqeg5hAr3aNCQ8L2ls/9E9rWuqTesaojy9VXUXcDcQAyQB57uX\njTGVyPNzN55IGnlcCi9+85szresPL7unda1p07qGME8fADyAM6S6MaYS25WSXuj6Yyn7YfJNsOkr\nm9a1EvAocYhIC+BeoHn+Y1T1T/4JyxgTbGYn7nIe2tNT158nm3kr8nXYnGLTulYSnrZxzADGA7MA\nl//CMcYEm8PHsnhi5npmrdlFs/pRdDs6jwdkMk3lAEeoQU2Ok1m9KQyeCzFdAh2uKQeeJo4MVX3N\nr5EUQ0Rq4MwH8qSqzg5UHMZUNvN/3cfD0xI5fCyLv/2xLX+uvxKdNY4quRkA1OUYLsKocckoSxqV\niKeN46+KyJMicqGIdM57lXSQiEwQkX0isq7A+n4islFENovIox5c/xFgioexGmPKKC0zh9HTE7l9\n4jLqV6/KjLt7cM8lbQif/686VarzAAAdDElEQVQTSSNPGC74/sUARWoCwdMaR0dgKHAJJ29VqXu5\nOBOBN4D381aISDjwJtAHp4fWMhGZCYQDzxQ4fhgQB/yCMw+IMcbPFm89yN8+XcOulHTu6tWKB/q0\nIbJKuNNrKnVn4QelJpVvkCagRFVL3knkVyBOVbO8voBIc2C2qnZwL18IPKWqfd3LowFUtWDSyDv+\naaAG0B5nLpBrVPW0dhYRGQmMBIiOju4yefJkb0MFIC0tjZo1a5bq2GBjZQk+wVyOrFxl2qYsvt6R\nwxnVhREdI2lTLxyAiKwUzv71NRocWlHosRmRZ7D4wnHlGa5PBfPPxRtlLUfv3r1XqGp8Sft5WuNY\nA9QF9pU6opNigPz/tiQB3YvaWVUfAxCR24ADhSUN935jgbEA8fHxmpCQUKrgFixYQGmPDTZWluAT\nrOVITErhwSlr2Lwvh6EXnMWjl59NjUj3x8PWBTD9YUhPgfNvhvWfQXa+brkRUVS78j8kxCUEInSf\nCNafi7fKqxyeJo5o4FcRWQZk5q0sZXfcwvrplVjtUdWJJZ5YpD/Qv3Xr1iXtaozBmd719e828+b8\nzZxRM5L3h3Xj4rbuub5zs2HBM/D9S9CwrTMPeOMO0LI3zBuDpiYhdWLh0icgblBgC2LKlaeJ40kf\nXjMJaJZvORbY5YsT25zjxnhu096jPDhlNeuSjzCwUwxP/ulc6kS5Z+E7vMMZayppKXS+Bfo9C1Vr\nONviBkHcIBaGyH/pxnuePjm+0IfXXAa0cT9UmAzcCNzsixNbjcOYkuW6lPE/bOWFrzdRK7IKbw3p\nTL8OTU7usH4GzLwPUGcYdBvR1hRQbOIQkR9U9Q8icpRTbycJoKpau4TjJwEJQEMRScJ5DmO8iNwD\nzMXpSTVBVdeXpRB5rMZhTPF+P3ichz5dzbLth/lj+2j+M7AjDWu658bIToevRsOKdyEmHq4bD/Wa\nBzReE5xKqnHUAFDVUg06o6o3FbF+DjCnNOcsjtU4jCmcqvLx0t95+osNhIvw4vXnMbBzzMm5M/b+\n4oxou38D9LgfLnkcwiMCG7QJWiUljpL76gYRq3EYc7o9qRk8Mi2RhZv284fWDfnvdXE0rRvlbFR1\nahhfjYbI2jD0M2hV0uNZprIrKXE0EpEHi9qoqi/5OJ4ysRqHMSepKjPX7OIfM9aRletizNXnMqT7\nWScnW0o/7LRlbJjpJItr3oaajQIbtKkQSkoc4UBNCu9CG3SsxmGM42BaJo/PWMeX6/bQ+cy6vDjo\nfFo0rHFyh9+XwLQ74Ohu6DMGLrwXwjwdgchUdiUljt2qOqZcIjHG+MQ3v+xl9PREUtOzebhfO+68\nuNXJ2flcuc5kS/P/A3WbwbCvIdYGJzTeKSlxVIiaRh67VWUqsyMZ2YyZ9QtTVyRxTpPafHBHd85p\nkq/j45Hd8NlI2LYIOlwHV70M1YrtGGlMoUpKHJeWSxQ+YreqTGX10+YDjJqayO7UdO7u3Yq/XtqW\nqlXy3Xra9DXMuMvpcnv1m3D+YJtsyZRasYlDVQ+VVyDGGO+lZ+Xy3Fe/MvGn7bRsWIOpf76IzmfW\nO7lDThbM+yf8/AZEd4Dr3oUz2gYuYBMSPB1yxBgTZFb+fpi/TVnD1gPHuO2i5jzS72yiqoaf3OHg\nFufZjN2rodtI6PMviLDZCUzZhVTisDYOUxlk5bh4dd4m/rdgC03qRPHx8O5c1LrhqTut+QS+eNB5\niO/Gj+HsKwMTrAlJIZU4rI3DhLoNu4/w4JQ1bNh9hOu7xPKP/u2pXS3fE96ZaTDnb7BmEpx5EVz7\nDtSJDVzAJiSFVOIwJlTlupS3F23h5W82UScqgnduiadP++hTd9q12rk1dXgb9HoULh4F4fYnbnzP\nfquMCXLbDhzjoSmrWfl7Cld0bMy/B3Skfo2qJ3dQhSVvwTdPQPWGcOssaP6HwAVsQp4lDmOClMul\nfLhkB8/M+ZWIcOHVG8/nT+c1PTkwIcCxg/D5X2DTV9D2chjwf1C9fuCCNpVCSCUOaxw3oWJXSjoP\nT03kh80HuLjtGfz32jga1ynQI2rb9zB9BBw/CJf/1+k5Zc9mmHIQUonDGsdNRaeqTF+ZzFOz1pPr\nUp6+pgM3dzvz1FpGbg4sfA4WPQ8NWsPNU6BJXOCCNpVOSCUOYyqyA2mZ/H36Wr7+ZS9dm9fjhevP\n46wGNU7dKWWnU8v4/Wfn6e/L/wuRNQMTsKm0LHEYEwS+WreHxz5by9GMHP5+xdnc8YeWJwcmzLNh\nFnx+jzNQ4cBxEHd9YII1lZ4lDmMCKDU9m6dmruezVcl0iKnNpEHn0za6wISb2enw9eOwbBw07eTM\nA16/ZWACNgZLHMYEzKJN+3l4aiL70zK579I23HtJayLCC8yJse9X59mMfevhonvhkiegStXCT2hM\nOQmpxGG9qkywmrEqmefnbiQ5JZ0mP8+jZcPq/LjlEK3OqMHbQy/ivGZ1Tz1AFVa+D18+AlVrwOBp\n0OaywARvTAEhlTisV5UJRjNWJTN6+lrSs3MB2J2awe7UDHq1bcjbQ+OpFhF+6gEZqTDrr7D+M2iZ\n4EzpWqtxucdtTFFCKnEYE4z++9WvJ5JGfpv3HTs9aexcBtOGQWoyXPok9LjfpnQ1QccShzF+oKqs\nSz7CtJVJ7ErNKHSfXSnpJxdcLvjpVfju31C7KQybC826llO0xnjHEocxPrT3SAafrUpm+sokNu1N\no2p4GFERYaRnu07bt2ndKOfN0b3w2Z2wdT60HwD9X4Wouqftb0ywsMRhTBmlZ+Xy9S97mLYymR9+\n249LofOZdXn6mg5c1bEp8zfuO6WNAyAqIpxRfdvB5m/hs7uc4dD7vwqdb7VhQ0zQs8RhTCmoKku3\nHWL6ymS+WLubtMwcYupGcXfv1gzsHEuLhief+B7QKQbgRK+qmLpRPNynBVfvfws+fw0atYdbZ0Oj\nswNVHGO8YonDGC/8fvA401YmMX1VEjsPpVO9ajhXdGzCwM4xXNCiAWEFn/Z2GxD+IwMix6DVkpCI\nxvB9VUjZAfHDoO9/ICKqnEtiTOkFfeIQkQTgX8B6YLKqLghoQKbSOZKRzZzE3UxbmcSy7YcRgR6t\nGvLAZW3p16Ex1auW8GeUOAVm3QfZ6QjA0d3O+u53weXP+Tt8Y3zOr4lDRCYAVwH7VLVDvvX9gFeB\ncGCcqj5bzGkUSAOqAUl+DNeYE3Jdyve/7WfaymS+Xr+HzBwXLc+owai+7bimU8zJhm1PzBvjDBtS\n0K9fWOIwFZK/axwTgTeA9/NWiEg48CbQBycRLBORmThJ5JkCxw8DvlfVhSISDbwEDPZzzKYS27jn\nKNNXJvHZqmT2Hc2kTlQEg+KbcW2XWM6LrXPq8OaeUIXUnYVvS7X/g0zF5NfEoaqLRKR5gdXdgM2q\nuhVARCYDV6vqMzi1k6IcBiL9Eaep3A6mZTJzzS6mrUxiXfIRqoQJCe0acW3nGC45pxGRVcJLPklh\nDm2D2fcXvb1ObOnOa0yAiar69wJO4pidd6tKRK4D+qnqcPfyUKC7qt5TxPEDgb5AXeB/RbVxiMhI\nYCRAdHR0l8mTJ5cq3rS0NGrWDI35DawsRct2KWv25fLjrhwS9+eSq3BW7TB6NK3CBU2qUDuy9F1i\nxZVLTPIsWmz7CJVw9je8iEb7fyDclXlin9ywSDa2u5t90b18UZyAsN+v4FPWcvTu3XuFqsaXtF8g\nGscL+4ssMnup6nRgekknVdWxwFiA+Ph4TUhIKFVwCxYsoLTHBhsry6lUlTVJqUxfmcTMNbtIOZ7N\nGbUiuaPnWQzsHMPZjWuXPdDdiTDzXti92pkD/MoXaVInxmkgnzcGTU1C6sQSfukTtI8bRPuyXzFg\n7Pcr+JRXOQKROJKAZvmWY4FdvjixjY5rCrM7NZ3PViUzbUUSW/YfI7JKGH88tzHXdo7hD60bUqXg\nUOalkZ0OC56Fn16H6vXh+onOU+B5bSJxgyBuEAtD5APKVG6BSBzLgDYi0gJIBm4Ebg5AHCaEHc/K\nYe76PUxbkcyPWw6gCl2b12NEz5ZcEdeE2tUifHexbYuc0WwPbYVOQ6DPv5zkYUyI8nd33ElAAtBQ\nRJKAJ1V1vIjcA8zF6Uk1QVXX++J6Nqx65eZyKUu2HWLayiS+XLubY1m5xNaL4t5L2nBt55jT5+8u\nq/TD8PU/YNUHUK8F3PK5Mwy6MSHO372qbipi/Rxgjj+vbSqPbQeOMX1lEtNXJpOckk7NyCpcGdeE\nazvH0rV5/SKf5i41Vfjlc5gzCo4fhB5/hV6PQtXqvr2OMUEq6J8c94a1cVQeqenZzE7cxfSVyazY\ncZgwgR6tG/Jwv3b8sX1joqqWsgttiRdOhjl/g41zoMl5MGSq89WYSiSkEofdqgo9+adcbbp4Hld1\nbEJyagbf/LKXrBwXbRrV5NHLz2bA+TE0rlPNf4G4XLBiAnzzFLhynHaMC/4C4SH1J2SMR0Lqt95q\nHKFlxqpkHp2eSIZ7LotdKRmM/X4b1SPCuKmr8zR3x5hSPM3trf2bnLGmfv8ZWvSC/q9A/Zb+vaYx\nQSykEofVOCqejOxcklPSSTqcTtLh4yQdTmfnIedrYlIKrkKe8KlbvSr/vLrD6Rt8LScLfnwFFj0P\nEdXh6v+D82+2+TJMpRdSicMEn8ISQ/73+49mnrJ/RLgQUzeK2HrVC00aALuLmIrVp3Yucx7k278B\nzh3oDEZYs5H/r2tMBRBSicNuVZW/jOxcdqWks9OLxNC0bhSx9aK4pF0jYutFEVvfSRSx9aJoVKsa\n4e5eUD2e/Y7klNNHlfVqZFpvZR515v1e8rYz9/dNk6Hd5f67njEVUEglDrtV5Xt5iaFgQsj7uq9A\nYqgS5nliKMmovu2KnnLVHzZ9DV886Ixa220EXPIPqOaDYUiMCTEhlTiMI39PpJjF3zGqb7sT05cW\nVJbEkNDujBMJIe9rdG3PE0NJCptytbiylFrafvjqUVg3FRq2g2Fz4czuvr2GMSHEEoebNx+2wWzG\nquRT/ktPTknnkWmJ/LrnCGc1qHHa7aS9R05NDOFhQtO61YitW51ebU8mhmb1fZ8YPDGgUwwDOsX4\nZ/A2VVgzGeaOhsw0SBgNf3gAqtjo/cYUJ6QSR2nbOAr7sB09fS1AmZOHy6Vk5bqcV87JV3aui8yc\nU9dn5713b8su7Jgi9s/KcY5Zsu0QWTmuU2LIzHHx1sKtwKmJ4eI2+WsMUcTWr050rUjfDPoX7A5v\nh1n3w9b5ENsN/vQ6NDo70FEZUyGEVOIobRvH83M3nnIfHSA9O5fHZqxl+Y5DZOfoKR/OBT+wC0sE\nefvkFNU1qJSqVgkjMjyMqlXCiHB/rVoljKrhYURUCTstaeQR4PtHetO4drXKkRiKkpsDS96C+U+D\nhMEVL0D8HRBWib8nxngppBJHae0qpOcOwLHMXL5cu6fQD+jI8DCiIsKpXa2Ke1s4VcPDqFpF3F8L\nfLiHhxFZxId9pIf7VwmTEh92K64nUmy9Sj6W0ilzZfSDK1+0WfiMKQVLHDgfqoV92MbUjeLHRy8J\nQESlV+49kSqC7HRY+Bz8+Joz3Pl1E5xnM+xBPmNKJaTq5yLSX0TGpqamenXcqL7tiIo4dVC8ivph\nO6BTDM8M7EiM+1mHmLpRPDOwY4Vs6PeJbd/D/y6CH16G826Cu5dCh2staRhTBiFV4yhtG0e5dfss\nJ37tiVRRpB+Gb56Ale9DveYwdAa06h3oqIwJCSGVOMrCPmxDhCpsmOnMlXFsP1x0n9PN1ubKMMZn\nLHGY0HFkl5Mwfp0NjePg5inQ9PxAR2VMyLHEYSo+lwtWvAvfPgW5WdBnDFxwt82VYYyf2F+WqdhO\nmSvjYrjqFWjQKtBRGRPSLHGYiiknC358FRb9FyKi4E9vQKch1lvKmHIQUonDhlWvJJKWOw/y7fsF\nzr0G+j0HtaIDHZUxlUZIJQ4bVj0EJU6BeWPolZoEK5vCGWfDlu+gVhO4cRKcfUWgIzSm0gmpxGFC\nTOIUp/0iOx0BOJLsvFokwA0f2FwZxgRISD05bkLMvDHOcCEFHdpiScOYALIahwk+B7c4tY3UnYVv\nT00q33iMMaewxGGCw7EDsG46JH4CycsBgfBIyM08fV8b0daYgLLEYQIn6zhsnOMki83zQHMhuoPz\nAF+H62DHjyfaOE6IiIJLnwhczMYYSxymnLlyYdtC51bUhlmQlQa1Y+Cie6DjIGjc4eS+cYOcr/PG\noKlJSJ1YJ2nkrTfGBETQJw4RCQP+BdQGlqvqewEOyXhLFXavcZLFuqmQthci6zjPYMTdAGf1KHoG\nvrhBEDeIhTb4pDFBw6+JQ0QmAFcB+1S1Q771/YBXgXBgnKo+W8xprgZigEOAtYpWJId3wNopkPgp\nHNgIYRHQtq+TDNr0hYhqgY7QGFMK/q5xTATeAN7PWyEi4cCbQB+cRLBMRGbiJJFnChw/DGgH/Kyq\nb4vIVGCen2M2ZXH8EPwyw6ld/P6zs+7Mi+Cql6H9AGcGPmNMhSaq6t8LiDQHZufVOETkQuApVe3r\nXh4NoKoFk0be8UOALFWdIiKfqOoNRew3EhgJEB0d3WXy5MmlijctLY2aNWuW6thgU15lCcvNosHB\nZUTvXUj9QysI0xyOVY9lb3QC+xpdTEZU2YcDCZWfS6iUA6wswais5ejdu/cKVY0vab9AtHHEAPk7\n6CcB3YvZfzrwuoj0BBYVtZOqjgXGAsTHx2tp74eH0kROfi2LywU7fnB6RP0yEzKPQM1ouOAuiBtE\njcZxtBShpY8uFyo/l1ApB1hZglF5lSMQiaOw4UuLrPao6nHgDo9ObIMc+t/e9U6yWDvVGf6jak04\n509Ou0WLiyEsvORzGGMqtEAkjiSgWb7lWGCXL05sgxz6SWoyrP3UabfYtx7CqkDry5znLdpdYdOy\nGlPJBCJxLAPaiEgLIBm4EbjZFye2GocPZaQ6t6ASP4HtPwAKsV3hihecbrQ1GgY6QmNMgPi7O+4k\nIAFoKCJJwJOqOl5E7gHm4vSkmqCq631xPatxlFFOFmz+xkkWG79yhvuo3woSRkPH62xmPWMM4OfE\noao3FbF+DjDH19ezGkcpuFywc4m7kXsGpB+G6g0h/nbnSe6YzjarnjHmFEH/5Lg3rMbhhf0b3Y3c\nn0LK7xBRHc6+0nmSu2UChEcEOkJjTJAKqcRh3PLPmrcq3/hOR/fAumlOwti9BiQMWvaG3o87SSOy\n4vdjN8b4X0glDrtVxemz5qXuhM/vhu9fcob9UBc07QT9noVzB9pc3cYYr4VU4qi0t6pcuU4vqIwU\n+Pqx02fNy82CA5ug50NOu8UZbQMTpzEmJIRU4qjwcrKcxun0w04SyHt/2qvAtoxUinmG0qEuuOTx\ncimGMSa0hVTiKNOtqqLaBbylClnHCnywF5UEUk5NAtnHiilcGFSrC1H1nFf1+k732LzlvNfcx+D4\ngdOPt1nzjDE+ElKJo9S3qgprF5h1nzPJUIteJz/cPa0FuLKLvlZ4VYiqf/KDvm4zaBLnXq57ahLI\nnygiaxc9Z0V+Emaz5hlj/CqkEkepzRtzertAdjrMfqDoY6rWOvXDvtHZp//3X1gCiIjy73MRNmue\nMcbPQipxlPpWVWox80Nd83YhCaBucD/nYLPmGWP8KKQSR6lvVdWJdW5Pnba+GZx3o2+CM8aYEOHB\nTfNK4NInnFtI+Vm7gDHGFMoSBzi3dvq/BnWaoYhT0+j/mrULGGNMIULqVlWZWLuAMcZ4JKRqHCLS\nX0TGpqamBjoUY4wJWSGVOFR1lqqOrFOnTqBDMcaYkBVSicMYY4z/WeIwxhjjFUscxhhjvCKqJYyq\nWgGJSCrwW75VdYBUD983BAoZJdAj+c9Xmn0K21ZwXUUoi7flKLic9z7/uopSFn/+TIqL05N9gqks\nvv5bKalsofL7VXC5YFnK+vvVRlVLbiRW1ZB7AWOLWi7pPbDcV9f1dp/CtlXEsnhbjmLiz7+uQpTF\nnz+TUCqLr/9WSipbqPx+lVSW8vj9UtWQvVU1q5hlT9776rre7lPYtopYFm/LUXB5VhH7lFZ5lsWf\nPxNPz1MRyuLrv5WSyhYqv18FlwNRltC8VVUWIrJcVeMDHYcvWFmCT6iUA6wswai8yhGqNY6yGBvo\nAHzIyhJ8QqUcYGUJRuVSDqtxGGOM8YrVOIwxxnjFEocxxhivWOIwxhjjFUscJRCRliIyXkSmBjqW\nshKRASLyjoh8LiJ/DHQ8pSUi54jIWyIyVUT+HOh4ykpEaojIChG5KtCxlIWIJIjI9+6fTUKg4ykt\nEQkTkadF5HURuTXQ8ZSFiPR0/zzGichPvjpvpUwcIjJBRPaJyLoC6/uJyEYR2SwijwKo6lZVvSMw\nkZbMy7LMUNURwG3ADQEIt0helmODqt4FDAKCrgulN2VxewSYUr5ResbLsiiQBlQDkso71uJ4WY6r\ngRggmyArB3j9t/K9+29lNvCez4Ioy1OGFfUFXAx0BtblWxcObAFaAlWBNUD7fNunBjpuH5blRaBz\noGMvSzmAPwE/ATcHOvaylAW4DLgRJ5lfFejYy1iWMPf2aOCjQMdehnI8Ctzp3ifo/u5L+Tc/Bajt\nqxgqZY1DVRcBhwqs7gZsVqeGkQVMxvnPI6h5UxZxPAd8qaoryzvW4nj7M1HVmap6ETC4fCMtmZdl\n6Q1cANwMjBCRoPqb9KYsqupybz8MRJZjmCXy8meShFMGgNzyi9Iz3v6tiMiZQKqqHvFVDDZ17Ekx\nwM58y0lAdxFpADwNdBKR0ar6TECi806hZQHuxfkPt46ItFbVtwIRnBeK+pkkAANxPpzmBCCu0ii0\nLKp6D4CI3AYcyPfhG8yK+rkMBPoCdYE3AhGYl4r6O3kVeF1EegKLAhFYKRRVFoA7gHd9eTFLHCdJ\nIetUVQ8Cd5V3MGVUVFleA14r72DKoKhyLAAWlG8oZVZoWU68UZ1YfqGUWVE/l+nA9PIOpgyKKsdx\nnA/biqTI3y9VfdLXFwuqanGAJQHN8i3HArsCFEtZhUpZQqUcYGUJRqFSDijnsljiOGkZ0EZEWohI\nVZwGy5kBjqm0QqUsoVIOsLIEo1ApB5R3WQLdQyBAvRImAbs52d3uDvf6K4BNOL0THgt0nJWpLKFS\nDitLcL5CpRzBUhYb5NAYY4xX7FaVMcYYr1jiMMYY4xVLHMYYY7xiicMYY4xXLHEYY4zxiiUOY4wx\nXrHEYUwBIpLmh3M2F5GbvTxmjojU9XUsxpSVJQ5jykdznBFwPaaqV6hqin/CMab0LHEYUwT3jHYL\n3DMN/ioiH4mIuLdtF5HnRGSp+9XavX6iiFyX7xx5tZdngZ4islpEHihwnSYissi9bZ17VNa8azQU\nkbvc21aLyDYRme/e/kcR+VlEVorIpyJSszy+L8ZY4jCmeJ2A+3Em+GkJ9Mi37YiqdsMZQvyVEs7z\nKPC9qp6vqi8X2HYzMFdVzwfOA1bn36iqb7m3dcUZYuIlEWkIPA5cpqqdgeXAg6UpoDHesmHVjSne\nUlVNAhCR1Ti3nH5wb5uU72vBZOCNZcAEEYkAZqjq6iL2exX4TlVniTM/eXvgR3clqCrwcxliMMZj\nljiMKV5mvve5nPo3o4W8z8Fdk3ff1qpa0gVUdZGIXAxcCXwgIs+r6vv593FP9HQWcE/eKuAbVb3J\n86IY4xt2q8qY0rsh39e8//a3A13c768GItzvjwK1CjuJiJwF7FPVd4DxOPNJ59/eBfgbMERPzhC4\nGOiRr22luoi0LWuBjPGE1TiMKb1IEVmC8w9Y3n/+7wCfi8hSYB5wzL0+EcgRkTXAxALtHAnAKBHJ\nBtKAWwpc5x6gPjDffVtquaoOd9dCJolI3vzej+MMq22MX9mw6saUgohsB+JV9UCgYzGmvNmtKmOM\nMV6xGocxxhivWI3DGGOMVyxxGGOM8YolDmOMMV6xxGGMMcYrljiMMcZ4xRKHMcYYr/w/8JaeP6Cp\nQNcAAAAASUVORK5CYII=\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.errorbar(sizes, tva, yerr=tvea, fmt='o-', label=\"Numpy\")\n", "plt.errorbar(sizes, tna, yerr=tnea, fmt='o-', label=\"Numba\")\n", "plt.xscale('log')\n", "plt.yscale('log')\n", "\n", "plt.legend()\n", "plt.grid(True)\n", "plt.xlabel(\"Input size\")\n", "plt.ylabel(\"Time (s)\")\n", "plt.title(\"Total runtime\");" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We see that the runtime profile is more smoothly linear (as it should be) for the Numba code." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Exercise: Numba atop Numpy?\n", "\n", "Let's make `trapznv` which is a Numba-fied version of `trapz`, the Numpy version. The question is: can Numba be of any benefit on a function that is already a single Numpy expression? Let's see..." ] }, { "cell_type": "code", "execution_count": 19, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "680.0000000000001" ] }, "execution_count": 19, "metadata": {}, "output_type": "execute_result" } ], "source": [ "trapznv = jit(trapz)\n", "\n", "x = np.linspace(a, b, 200)\n", "y = f(x)\n", "\n", "trapznv(y, x)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now that we see it's correct, let's look at performance in just one case:" ] }, { "cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "4.19 µs ± 97.7 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)\n", "1.04 µs ± 22.1 ns per loop (mean ± std. dev. of 7 runs, 1000000 loops each)\n" ] } ], "source": [ "%timeit trapz(y, x)\n", "%timeit trapznv(y, x)" ] }, { "cell_type": "code", "execution_count": 21, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "631 ns ± 3.56 ns per loop (mean ± std. dev. of 7 runs, 1000000 loops each)\n", "839 ns ± 13.1 ns per loop (mean ± std. dev. of 7 runs, 1000000 loops each)\n", "2.84 µs ± 30.6 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)\n", "23.5 µs ± 287 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)\n", "228 µs ± 3.59 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)\n", "3.49 ms ± 125 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)\n", "43.7 ms ± 842 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)\n" ] } ], "source": [ "tnv, tnve = [], []\n", "for x, y in xy:\n", " t = %timeit -o trapznv(y, x)\n", " tnv.append(t.average)\n", " tnve.append(t.stdev)\n", "\n", "tnva = np.array(tnv)\n", "tnvea = np.array(tnve)" ] }, { "cell_type": "code", "execution_count": 22, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEaCAYAAAAL7cBuAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3Xd4VFX6wPHvm5AGhAQSaugkgPQA\nSnXFCipIYLEssi5r766rWNaG3f25dkUUdW0rIqwbwF4woogimNB7T0JLIAGSACnn98ed4BAmySSZ\nyZ3yfp4nT2buvTP3PdPee8659xwxxqCUUkpVFGJ3AEoppXyTJgillFIuaYJQSinlkiYIpZRSLmmC\nUEop5ZImCKWUUi5pglDHicjbIvJYFesPi0hnx+0oEZkvIvkiMrs+41D1R0Smisj7dseh7KEJwoeJ\nyDYR2SMijZyWXS0iaXbEY4xpbIzZ4rg7AWgJxBljLq6vGERksogYEZlSYXmmiIyorzhU7YhImogc\nEZF2TsvOEZFtNoalKqEJwvc1AG6zOwgXOgAbjDElNux7P3C3iDSxYd8+RUQa2B1DZUQktJJVBcAD\n9RmLqh1NEL7vaeBOEYmtuEJEOjqOphs4LUsTkasdtyeLyCIReU5E8kRki4gMdSzfKSJ7ReQvFZ42\nXkS+FpFDIvK9iHRwem4jIoki8jDwIHCpo9npOhHZLyK9nbZtISJFItLccX+0iGQ44vhJRPo4bZss\nIr859jkLiKzmNVkLLAZud7WyYhOViIwQkUyn+9tEZIqIrBCRAhF5U0Raisjnjhi+EZGmFV7ja0Uk\nW0R2icgdjnWtRKRQROKcnnuAiOwTkbAKMbVxvB7NKpQ7R0TCHK/r944muxzH6+CqbOXxXCUiO4AF\njuWDHa9rnogsd65NiUgnEVnoVLZXypuNKr42Tq/POZXsf7aI7HbEuVBEelZ43V8Vkc9EpAA409Vz\nAC8CfxKRxEr2YZzXOb+f5fGKyF2Oz+8uEUkRkQtEZIPjc/gPp8dOFZE5IjLLUf7fRKSvY90UEflv\nhX2/JCLPVxJ30NEE4fuWAmnAnbV8/CBgBRAHfAB8CJwKJAKTgJdFpLHT9pcDjwLxQAbwn4pPaIx5\nCHgCmOVodnrN8byTnDb7E/CNMWafiPQH3gKuc8TxGjBPRCJEJBxIBd4DmgGzgT+6Ua4HgNudf3Br\n6I/AuUBXYAzwOfAPrHKHALdW2P5MIAk4D7hHRM4xxuzGem8ucdpuEvChMabY+cHGmGyspOZctonA\nHMe2jwJfAU2BtsBL1cR/BnAKMFJEEoBPgcewXsM7gf+WJ2es930J1ms/FfhzNc9dlc+xXocWwG+c\n/PmYCDwORAM/VvIcWcAMRyy10QrrICIB60BlBtbrPgA4HXhQHH1lDmOxPlfNsF6LVEcCfx8YJY6D\nL7EOtC7F+iwqNEH4iweBW5y+8DWx1Rjzb2NMKTALaAc8Yow5aoz5CjiGlSzKfWqMWWiMOQrcBwwR\np/biKrwDTBSR8s/Un/n9i3YN8Jox5hdjTKkx5h3gKDDY8RcGPG+MKTbGzAF+rW5nxpgMrB/Uu92I\nzZWXjDF7jDFZwA/AL8aYdEe5/wckV9j+YWNMgTFmJfBvrARYXu5JcLxJ5U9U/gPzQfnjRESAyxzL\nAIqxmu3aGGOOGGMq+3EtN9URT5Fj/58ZYz4zxpQZY77GOrC4QETaYx0QPGiMOeZ43nnVvTiVMca8\nZYw55HidpgJ9RSTGaZO5xphFjjiOVPFUTwJjnGsgNVAMPO5IrB9iJfUXHHGtBlYDfZy2X2aMKU/E\nz2Ill8HGmF3AQqC8D20UkGOMWVaLmAKSJgg/YIxZBXwC3FOLh+9xul3keL6Ky5xrEDud9nsYq72/\njRsx/oLVtnyGiHTHSjrlP0QdgDsczR95IpKHlajaOP6yzImjRm53s2wPAjeISCs3t3dW8TWo6jUB\np9cFK77y12Qu0MNxxHoukG+MWVLJPudgJdw2wB8Ag5WcAO4CBFgiIqtF5Mpq4neOpwNwcYXXdzjQ\n2hHnfmNMYSWPdZuIhIrIUyKyWUQOAtscq+Jr+tzGmH3Ay8AjtQgl13HAA47PNFW/f86f6TIgk9/f\nv+MJ3vFfaw9OfLaDS53kIawq/TNOywoc/xsCBx23a/Nj6cz57JLGWNXybDcfW/5l243VdFJ+BLkT\n64jv8YoPEJEzgAQREack0R7YXN3OjDHrRORjrKYhZwVYr0m5ur4mYL0u65ziy3bEcEREPsJqmutO\nFT8wxpg8EfkKq0nqFGBmeZkdzVXXAIjIcOAbEVlojNlU2dM53d4JvGeMuabiRmL1ITUTkYZOScK5\nRnjCa+WoBVVWU52I1VxzDlZyiAEOYCU2V3FV52lgC1bzl7NCTn7/Mqk95890CFYTXvlnOhV4VUR6\nAaOxErVy0BqEn3D8UMzCqW3ccRSWBUxyHN1dCXSp464uEJHhjr6BR7GaXtw94nwPGIeVJN51Wj4D\nuF5EBomlkYhcKCLRWO3yJcCtItJARMYDp9Ug3oeBvwLOnfgZjnI0c9Qu/laD56vMAyLS0NEk8les\n96Lcu8Bk4CKsdu2qfABcgdUXUd68hIhcLCJtHXcPYP3Qlp78cJfex2quGen4HEQ6OnPbGmO2YzU3\nTRWRcBEZgtXnUm4DEOl4P8KA+4GISvYTjdU0mIv1A/6Em/G5ZIzJwzrgqfijnIHVXBkqIqOw+lvq\nYoCIjHf0MfwNqww/O2I4glWz+wBYYozZUcd9BRRNEP7lEaBRhWXXAFOwvrQ9gZ/quI8PsGor+7E6\n/S5394HGmEysWo5z0wnGmKWOOF/G+vHbhPWDijHmGDDecf8AVifhxzXY51asxOT8urwHLMc6yv2K\nE3/Ma+t7rLi/Bf7l6L8pj2ERUAb8ZozZVs3zzMPq5N1jjFnutPxU4BcROezY5jZH2arlSOBjsWpS\n+7BqFFP4/ft9OTAE6zPyGNbrcdTx2HzgRuANrIONAio/Wn8Xq3ktC1iD40e2jl7g5ER4G1YSy3PE\nnlrHfczF+lwdwOobG1/hJIJ3gN5o89JJRCcMUp4kIm8B2caY++2OxRNEpCOwFQir6poPEVkAfGCM\neaOeQqs1sU6hXec4Gy2gichUINEYM6mKbdpjNR+2MsYcrGy7YKR9EMpjHD+m4zn5DKCAJiKnAv2x\njuJ9jiO+/ViJ7jysOJ+yNSgf4eiT+DvWqcmaHCrQBKE8QkQexbpw7Ul3m0YCgYi8A6RgNQkdsjue\nSrTCaraLw2o+usEYk25vSPYTawibPVjNZqNsDscnaROTUkopl7STWimllEuaIJRSSrnk130Q8fHx\npmPHjrV6bEFBAY0aVTxj1D9pWXxPoJQDtCy+qi5lWbZsWY4xptqhe/w6QXTs2JGlS5fW6rFpaWmM\nGDHCswHZRMviewKlHKBl8VV1KYuIuDWcjTYxKaWUckkThFJKKZc0QSillHLJr/sglFIKoLi4mMzM\nTI4cqWoKCoiJiWHt2rX1FJV3uVOWyMhI2rZtS1hYWJXbVUYThFLK72VmZhIdHU3Hjh2x5mJy7dCh\nQ0RHR9djZN5TXVmMMeTm5pKZmUmnTp1qtQ9tYlJK+b0jR44QFxdXZXIINiJCXFxctbWqqmiCUEoF\nhJomh0tfW8ylry32UjS+oa4JUxOEUsql1PQsut73OZO/KGDYUwtITc+yOySfFhoaSr9+/ejVqxcX\nX3wxhYWFlW67bds2Pvjg+HxRvP3229x88831EWaNaIJQSp0kNT2Lez9eybHSMgCy8oq49+OVAZMk\nUtOzSN+Rxy9b93ss+UVFRZGRkcGqVasIDw9n+vTplW5bMUH4Kk0QSqkTGGN44rO1FBWfONFbUXEp\nT3+53qaoPKc+kt/pp5/Opk2beOCBB3jhhReOL7/vvvt48cUXueeee/jhhx/o168fzz33HADZ2dmM\nGjWKpKQk7rrr91lYZ86cSe/evenVqxd333338eWtW7fmvvvuo2/fvgwePJg9e/Z4LP5yehaTUord\n+Uf4YeM+ftyUw6JNueQcPupyu+y8onqOrOYenr+aNdmu5/4pLS1lRdah48mhXFFxKXfNWcHMJa6n\npO7RpgkPjenp1v5LSkr4/PPPGTVqFOeffz7jx4/ntttuo6ysjA8//JAlS5bQp08f/vWvf/HJJ58A\nVhNTRkYG6enpRERE0K1bN2655RZCQ0O5++67WbZsGU2bNuW8884jNTWVlJQUCgoKGDx4MI8//jh3\n3XUXM2bM4P77PTuRoyYIpYLQoSPF/LJlPz9uyuHHTTls2nsYgPjG4QxLjOf7DfvIKyw+6XHNGoXX\nd6geVzE5VLfcXUVFRfTr1w+wahBXXXUV4eHhxMXFkZ6ezp49e0hOTiYuLs7l488++2xiYmIA6NGj\nB9u3byc3N5cRI0bQvLk1rt7ll1/OwoULSUlJITw8nNGjRwMwYMAAvv766zrF74omCKWCQHFpGct3\n5lkJYWMOGTvzKCkzRIaFcFqnOC4d2I7hSfF0axlNSIgcb4ZxbmYSgdyCY7zxwxauGt7JZ08prepI\n/9ChQ4x65VeyXNSEEmKjmHXdkFrvt7wPoqKrr76at99+m927d3PllVdW+viIiIjjt0NDQykpKaGq\nCd3CwsKOvwfl23uaJgilApAxhs37Cli0KYcfNubw85ZcDh8tQQT6JMRw7R86MzwpngEdmhLRIPSk\nx6ckJwDw9JfrycorIiE2itvOTuK79Xt57NO1bNp7mEfG9iK8gf91Y04Z2e2k5BcVFsqUkd28sr9x\n48bx4IMPUlxcfLxjOjo6mkOHqp+hdtCgQdx2223k5OTQtGlTZs6cyS233OKVOF3RBKFUgMg5fJRF\njhrCj5ty2JVvXSDVvllDxvRtw+lJ8QztEkdsQ/eaiVKSE0hJTjhhWOkJA9ry7NcbePm7TWzLLeDV\nywfQ1M+ancqT311zVnCstIyE2CimjOx2fLmnhYeHc+aZZxIbG0toqJWM+/TpQ4MGDejbty+TJ0+m\nadOmLh/bunVrnnzySc4880yMMVxwwQWMHTvWK3G6oglCKT9VdKyUJdv28+PGffy4KZe1u6yO2Zio\nMIYlxnFLYnOGJ8bTPq6hx/YZEiLcObIbXVo04u45Kxk3bRFv/OVUEls09tg+6kNKcsLxDum6NCs5\nO3z4sMvlZWVl/Pzzz8yePfv4srCwML799tsTtps8efLx2+Wd1wATJ05k4sSJJz3vrl27jt+eMGEC\nEyZMqG3oldIEoZSfKC0zrM7O54eNVi1h2fYDHCstIzw0hAEdmjJlZDdOT4qnZ5sYQkO82z8wLrkt\n7Zs14rr3ljJu2iJevXwAw5PivbpPT/NUYqjKmjVrGD16NOPGjSMpKcnr+/M0TRBK+bAduYWOM432\n8dPm3ONnFp3Sugl/GdqB4UnNOa1jM6LCT+5H8LYBHZqSetMwrnp7KX/59xKmXtSTPw/uUO9x+LIe\nPXqwZcsWu8OoNU0QSvmQvMJjLN6cyw+OvoQd+63hGlrHRHLuKS0ZnhTP0C7xNI+OqOaZ6kfbpg35\n741DuXVmOg+krmLz3sPcf+EpNAj1v85rdTJNEErZ6GhJKcu2HzjeubwiKx9joHFEAwZ3juOq4Z0Y\nlhhPl+aNfPa00sYRDZhxxUCe/Gwtb/y4lS05Bbw8MZkmkbWbg6C2jDE++xrZparTZN2hCUIpD0pN\nz/r91NCfF5x0dowxhnW7Dx0/02jJ1v0UFZcSGiIkt4vltrOTOD0pnj5tYwnzo6Pw0BDh/tE9SGzR\nmPtTVzF+2k+8+ZeBdIhrVC/7j4yMJDc3V4f8dlI+H0RkZGStn0MThFIeUvHisvIxfvKKjtEovMFJ\nw1h0ad6IS09tx/DEeAZ1bkZ0PR9xe8Nlp7WnfVxDbnj/N1JeWcT0SQMY1Nn1lcOe1LZtWzIzM9m3\nb1+V2x05cqROP5i+xJ2ylM8oV1uaIJTykKe/XO9ygLup89YAEN84guGJcQxLjGd4UjytY6LsCNPr\nhnaJd3Re/8qkN3/h8XG9uWRgO6/uMywszK1Z09LS0khOTvZqLPWlPsqiCUIpD6lqILvPbzud7q2i\ng6b5o1N8I/534zBu+uA37pqzgs17D3PXqO5eP/1WeZb/NHIq5ePaxLquESTERnFK6yZBkxzKxTQM\n499/PZVJg9vz2sItXPfeMgqOen68IOU9miCU8pApI7sRWiEJeHOMH38QFhrCo2N7MXVMDxas28OE\n6YtdDpSnfJMmCKU8pF2zKEqNoXGE1XKbEBvFk+N7e22MH38hIkwe1om3Jp9K5v5Cxr68iPQdB+wO\nS7lBE4RSHlBaZnho3mpaNYnkl3+czdujGrHonrOCPjk4G9GtBR/fOJSo8BAuff1n5i3PtjskVQ1N\nEEp5wKxfd7Iq6yD/uPAUGkXouR+VSWoZTeqNw+jbNoZbZ6bz7Ncb6nwxl/IeTRBK1VFe4TGe/nId\ngzo1Y0yf1naH4/PiGkfw/tWDmDCgLS9+u5GbZ6ZzpMLpwco36KGOUnX0zFcbyC8qZupFPYPuTKXa\nimgQytMT+pDUojFPfbGOzP2FzLhiIC2aBMZFbIFCaxBK1cGa7IP855ftXDGkI6e0bmJ3OH5FRLju\njC68NmkAG/ceZuwri1iVlW93WMqJJgilaskYw9R5q4ltGM7t53S1Oxy/dV7PVsy+fggCXDx9MV+u\n3m13SMpBE4RStTRveTZLtu3nrpHdiGno/+Mo2alnmxhSbx5G11bRXP/+Ml5N26yd1z5AE4RStXD4\naAlPfLaWPm1jvD7OULBoER3JrGsHM7pPG/75xTrunL2CoyXaeW0n7aRWqhZeWrCRPQePMn3SAEJ0\nfCGPiQwL5cXL+tGleSOe/2YjO/YXMH3SAOIa+8YEScFGaxBK1dDmfYd568etTBjQluT2Te0OJ+CI\nCH87pysv/SmZFZn5pExbxIY9h+wOKyhpglCqBowxPDJ/DZENQrl7VHe7wwloY/q2YdZ1QzhSXMYf\np/1E2vq9docUdDRBKFUD36zdy/cb9vG3c7v6zLzQgaxfu1jm3jSMds0acuXbv/LvRVu187oeaYJQ\nyk1Hikt59JM1JLVozBVDOtgdTtBoExvF7OuHcM4pLXl4/hruT11FcWmZ3WEFBU0QSrlpxsIt7Nhf\nyMMX9fSr+aIDQaOIBkyfNIDrz+jCf37ZweR/LyG/sNjusAKefsqVckPmgUJeSdvEBb1bMTQx3u5w\nglJIiHDP+d15ekIflmzdz7hpi9iaU2B3WAFNE4RSbnjis7UA3HdhD5sjURcPbMd/rh7MgcJjpLyy\niJ825dgdUsDSBKFUNRZtyuGzlbu5aUQiCZVMK6rq12mdmjH3puG0iI7gireW8MEvO+wOKSBpglCq\nCsWlZUydt5r2zRpyzR862x2OctI+riH/vXEowxLj+cf/VvLI/DWUlukZTp6kCUKpKry7eDsb9x7m\ngdE9iAwLtTscVUGTyDDe/MtAJg/tyFuLtnL1O79y6Ih2XnuKJgilKrHv0FGe/3oDZ3RtzjmntLA7\nHFWJBqEhTL2oJ4+l9GLhxhwmvLqYnfsLT9ru0tcW8+QvRTZE6L98JkGISCMReUdEZojI5XbHo9Q/\nv1jHkZJSHhrTQycC8gOTBnfgnb+exq78IlJeWcTSbfvtDsnveTVBiMhbIrJXRFZVWD5KRNaLyCYR\nucexeDwwxxhzDXCRN+NSqjq/7TjAnGWZXDm8E52bN7Y7HOWm4Unx/O+mYTSJCmPijF/4+LdMu0Py\na96uQbwNjHJeICKhwCvA+UAP4E8i0gNoC+x0bKZj/CrblJVZEwG1bBLBLWcl2R2OqqEuzRvzvxuH\nMqBDU/7+0XL+74t1/G9ZJuk78lh/oIxhTy0gNT3L7jD9gleH+zbGLBSRjhUWnwZsMsZsARCRD4Gx\nQCZWksjA3cS1fj2MGFGr2Prl5UFsbK0e62u0LJ71UfPerOgyihc2fkLjkY/V6jl8oRye4o9liQXe\nlRAe7HgO09IgxJRRJtbPSlZeEfd+8Cs89igpuetsjbMu6uN9sWM+iAR+rymAlRgGAS8CL4vIhcD8\nyh4sItcC1wL0CgsjLy+vVkGUlpbW+rG+RsviOQcbRPLP/qfTf/9W/rB5MbWNxO5yeJI/l+WuAx/x\n6VlJHAxveMLyotAw/tl2OCM2/2xTZHVXH++LHQnCVW+fMcYUAH+t7sHGmNeB1wEGDhxoYpcurVUQ\naWlpjKhl7cPXaFk857m5q8j/eTuP3XkFTdvcXOvnsbscnuTvZTl0z6cul++OakpsRkY9R+M5dXpf\n3Dzpwo6zmDIB5zka2wLZNsSh1AnW7jrIez9v5/JBHejRpond4SgPaVPJ1e+VLVe/syNB/AokiUgn\nEQkHLgPm2RCHUscZY3ho3mpiosK447yudoejPGjKyG5EubjIsVebJjq3RDW8fZrrTGAx0E1EMkXk\nKmNMCXAz8CWwFvjIGLPam3EoVZ35K3axZOt+pozsTmzDcLvDUR6UkpzAk+N7E+4Yor1NTCSndWzK\nl2v2cPPMdIqO6UmTlfH2WUx/qmT5Z8Bn3ty3Uu4qOFrCE5+upVdCEy49tV31D1B+JyU5gZlLdpCX\nl8eXd5+NMYbXF27hqS/WsXN/ITOuGEjLJpF2h+lzfOZKaqXs8sp3m9h98AgPX9SL0BC9YjoYiAjX\nndGF1/88kM17D3PRyz+yMjPf7rB8jiYIFdS25hQw44ctjO+fwIAOTe0OR9Wzc3u0ZM4NQ2kQEsLF\nr/3Epyt22R2ST9EEoYLaI/NXE9EglHvO7253KMrLZl03hHsHnXzm0imtm5B60zB6tonhpg9+48Vv\nN2rntYNfJggRGSMir+fna5VQ1d63a/fw3fp93HZ2Ei2itf05mDWPjuA/Vw9ifHICz369gVs/zOBI\nsXZe+2WCMMbMN8ZcGxMTY3coyk8dKS7lkU/WkNiiMZOHdbQ7HOUDIsNCeeaSvtw9qjufrMjm0tcW\ns/fgEbvDspVfJgil6urNH7eyPbeQqWN6EhaqXwNlERFuGNGF6ZMGsHHvYca+sohVWcHbUqHfDBV0\nsvOKeHnBJkb1bMXwpHi7w1E+aGTPVsy+fggCXDx9MV+sCs7Oa00QKug8/tlayozhvgtPsTsU5cN6\ntokh9eZhdGsVzfXv/8bLC4Kv81oThAoqP23O4dMVu7hhRBfaNWtY/QNUUGsRHcmH1w4mpV8b/vXV\nBm6fFVyd13aM5qqULUpKy3h43hraNo3i+jO62B2O8hORYaE8d2k/klpG8/SX69mWW8jrVwwIijPf\ntAahgsZ7P29n/Z5DPDC6B5EuBm9TqjIiwk1nJjJ9Un/W7z5EysuLWJ0d+J3XmiBUUMg5fJRnv97A\n6UnxnNejpd3hKD81qldrZl8/BANMeHUxX67ebXdIXqUJQgWF//tiHUXHSnloTE/EzclSlHKlV0IM\nc28aRtdW0Vz33jKmpW0K2M5rTRAq4GXszOOjpZlcObwTiS0a2x2OCgAtmkQy69rBjOnbhv/7Yj13\nfLQ8IDuv/TJB6FAbyl1lZYaH5q6ieXQEt5yVaHc4KoBEhoXy4mX9uOPcrnycnsXEGT+z79BRu8Py\nKL9MEDrUhnLXnGWZLM/M597zuxMdGWZ3OCrAiAi3nJ3EtMv7s2bXQVJeWcTaXQftDstj/DJBKOWO\n/KJi/vnFOgZ2aMq45AS7w1EB7ILerZl93VBKywx/fPUnvl6zx+6QPEIThApYz3+zgf2Fx5h6kXZM\nK+/r3TaGuTcPI7FFY659bynTv9/s953XmiBUQFq/+xDvLt7OxNPa0ytBmyJV/WjZJJJZ1w7hgt6t\neerzddw5ewVHS/y381qvpFYBxxjDQ/NWER3ZgDvP62Z3OCrIRIWH8vKfkklq0Zjnv9nI9twCpv95\nAPGNI+wOrca0BqECzqcrd/Hzlv3ccV43mjYKtzscFYREhL+d05WXJyazMiufsS8vYt1u/+u8rjZB\niEiYiNwqInMcf7eIiJ4OonxS4bESHv90LT1aN2Hiae3tDkcFudF92vDRdUMoLi3jj9N+4tu1/tV5\n7U4N4lVgADDN8dffsUwpnzPtu83syj/CI2N7EhqiHdPKfn3bxTLv5uF0at6Iq99dyoyFW/ym89qd\nPohTjTF9ne4vEJHl3gpIqdrallPA6wu3MC45gYEdm9kdjlLHtYqJZPZ1Q7ljdgaPf7aWDXsO8fi4\n3oQ38O1WfneiKxWR42Mji0hnwH+75VXAevSTNYSFCvee393uUJQ6idV53Z9bz05i9rJMJr3xC7mH\nffvKa3cSxBTgOxFJE5HvgQXAHd4Nq2o61Iaq6Lt1e/l23V5uPTuJFk0Cf5x+5Z9CQoS/n9uVFy7r\nR0ZmHinTFrFhzyG7w6pUtQnCGPMtkATc6vjrZoz5ztuBVROTDrWhjjtaUsrD81fTuXkj/jqsk93h\nKFWtsf0S+Oi6IRwpLmP8tJ/4bt1eu0NyqdIEISJnOf6PBy4EEoEuwIWOZUr5hDd/3Mq23EIeGtPT\n59t0lSrXr10s824eRoe4hlz1zq+88YPvdV5X1Ul9BlZz0hgX6wzwsVciUqoGduUX8fKCTZzXoyVn\ndG1udzhK1UjrmChmXz+Ev89azmOfrmXT3sM8MraXzxzoVJogjDEPOW4+YozZ6rxORLQer3zCE5+t\no7TM8MDoHnaHolStNAxvwLTL+/Ps1xt4+btNbM0p4NVJA2jmAxd5upOm/uti2RxPB6JUTf28JZf5\ny7O57owutGvW0O5wlKq1kBDhzpHdeOGyfqTvzCPllUVs9IHO60prECLSHegJxFToc2gC6GkiylYl\npWVMnbeahNgobjijS/UPUMoPjO2XQLtmDbn23WWMn/YTL01MZkS3FrbFU1UNohswGojF6oco/+sP\nXOP90JSq3H9+2cG63Ye4/8JTiAoPtTscpTymf/umzL15GG2bNeTKt3/lrR+32tZ5XVUfxFxgrogM\nMcYsrseYlKpS7uGjPPPVeoYlxjGqVyu7w1HK4xJio5hz/RBun5XBI5+sYePewzwytidhofXbee3O\nUBvpInITVnPT8aYlY8yVXotKqSr866v1FB4rZeoYnQhIBa5GEQ2YPmkA//pqPdPSNrMtp4Bpl/ev\n1xGK3UlH7wGtgJHA90BbwP7GpZpBAAAUy0lEQVTeExWUVmTm8eGvO5k8tCNJLaPtDkcprwoJEe4a\n1Z1nL+nLsu0HSJm2iNe+38ywpxYw+YsChj21gNT0LO/t341tEo0xDwAFxph3sC6a6+21iJSqRFmZ\n4cG5q4lrFMFt5yTZHY5S9WZ8/7bMvHYQuYeP8uTn68jKKwIgK6+Iez9e6bUk4U6CKHb8zxORXkAM\n0NEr0ShVhf/+lknGzjzuOb870ZE6JYkKLgM6NKNRxMm9AkXFpTz95Xqv7NOdPojXRaQpcD8wD2gM\nPOCVaJSqxMEjxfzzi3Ukt49lfHKC3eEoZYu9B12P/prtqFF4WpUJQkRCgIPGmAPAQqCzV6KoIREZ\nA4xJTEy0OxRVT174ZiO5Bcf49+TTCNGJgFSQahMbdbx5qeJyb6iyickYUwbc7JU914GO5hpcNuw5\nxNs/beOyU9vTu62+5yp4TRnZjaiwE6/7iQoLZcrIbl7ZnztNTF+LyJ3ALKCgfKExZr9XIlLKiTGG\nqfNW0ziigde+BEr5ixRH8+rTX64nK6+IhNgopozsdny5p7mTIMqvd7jJaZnBR5qbVGD7fNVuftqc\nyyNje/rE4GVK2S0lOYGU5ATS0tIYMWKEV/dVbYIwxujIrcoWRcdKefzTtXRvFc3E09rbHY5SQced\nGoRStng1bRNZeUXMunYwDep5iAGllHvXQShV73bkFjJ94RYu6tuGQZ3j7A5HqaCkCUL5pEc+WUOD\nEOEfF5xidyhKBa1qm5hEpL+LxfnAdmNMiedDUsEubf1evlm7h7tHdadVjE49opRd3OmDmIY1B8QK\nQIBejttxInK9MeYrL8angkxJmeGJ+WvoFN+IK4d3tDscpYKaO01M24BkY8xAY8wAIBlYBZwD/J8X\nY1NBJDU9i2FPLeDqrwrZklPAuT1aENFAJwJSyk7uJIjuxpjV5XeMMWuwEsYW74Wlgklqehb3frzy\nhCEE3lu8w6vDGCulqudOglgvIq+KyBmOv2nABhGJ4PeRXpWqtae/XE9RcekJy7w5QqVSyj3uJIjJ\nwCbgb8DtwBbHsmLgTG8FpoJHZSNRemuESqWUe9y5kroIeMbxV9Fhj0ekgk6rmEh25R85abm3RqhU\nSrmn2hqEiAwTka9FZIOIbCn/q4/gVHBo2STipGXeHKFSKeUed05zfROraWkZUFrNtvVC54MIHPOW\nZ5OxM58LerVieWZ+vYxQqZRyjzsJIt8Y87nXI6kBY8x8YP7AgQOvsTsWVXu784/wQOoq+rWL5cU/\nJdMgNKReRqhUSrnHnQTxnYg8DXwMHJ/vzhjzm9eiUgHPGMPd/13B0ZJSnr2krw7Gp5QPcidBDHL8\nH+i0zABneT4cFSw+WLKD7zfs45GxPencvLHd4SilXHDnLCY9lVV51LacAh77ZC2nJ8UzaVAHu8NR\nSlWi0gQhIpOMMe+LyN9drTfGPOu9sFSgKi0z3DF7OQ1Chf+b0IeQELE7JKVUJaqqQTRy/I+uj0BU\ncHh94RaWbT/A85f2o3WMXueglC+rNEEYY15z/H+4/sJRgWztroM8+/V6LujdirH92tgdjlKqGu7M\nB9EJuAXo6Ly9MeYi74WlAs3RklJun5VBTFQ4j6X0RkSblpTyde6cxZSKdbHcfKDMu+GoQPX8NxtZ\nt/sQb/5lIM0ahdsdjlLKDe4kiCPGmBe9HokKWEu37ee17zdz2antOPuUlnaHo5RykzsJ4gUReQj4\nCr1QTtVQwdES7pi9nDaxUdw/uofd4SilasCdBNEb+DPWhXHlTUx6oZxyyxOfrWXH/kI+vGYwjSPc\n+bgppXyFO9/YcUBnY8wxbwejAst36/fyn192cM3pnRjUOc7ucJRSNeTOADjLgVhvB6ICS17hMe6e\ns4KkFo254zwdtlspf+RODaIlsE5EfuXEPgg9zVVV6oG5q9lfcIy3Jp9KZFio3eEopWrBnQTxkNej\nUAFl/vJs5i/P5o5zu9IrIcbucJRSteTOYH3f10cgKjDsOXiEB+auom+7WG4Y0cXucJRSdVDVYH0/\nGmOGi8ghrLOWjq8CjDGmidejU36lfI6HI8U6x4NSgaDawfqMMTpYn3LLzCU7SVu/j4cv6kkXneNB\nKb9X1SGeqWKdrURkjIi8np+fb3coymF7bgGPfbqG4Ynx/HmwzvGgVCCoqgbRorK5IMDe+SB0Tmrf\nUlpmuOOj5YSG6BwPSgWSqhJEKNAYq89BqUrN+GELS7cf4LlL+9ImVud4UCpQVJUgdhljHqm3SJRf\nWrvrIM9+tYHze7UipV+C3eEopTyoqj4IrTmoKpXP8dAkKozHUnrpHA9KBZiqahBn11sUyi+94Jjj\n4Y0rBhLXOMLucJRSHlZpDcIYs78+A1H+Zdn2/Uz/fjOXDmzHOT10jgelApFeyaRqrOBoCX//aDmt\nY6K4f/QpdoejlPISHaBf1diTn1tzPMy8ZjDRkWF2h6OU8hKtQaga+X7DPt7/eQdXDevEYJ3jQamA\npglCuS2/sJi75iwnqUVj7hypczwoFei0iUm57cF5q8g9fIw3rtA5HpQKBlqDUG75ZEU2czOyufXs\nJHq31TkelAoGmiBUtfYePML9qdYcDzfqHA9KBQ1NEKpK5XM8FB3TOR6UCjb6bVdV+vDXnXy3fh/3\nnt9d53hQKshoglCV2pFbyKOfrGFYYhxXDOlodzhKqXqmCUK5VFpmuGN2BqEhwtMT+uocD0oFIT3N\nVbn0xg9b+HXbAZ69ROd4UCpYaQ1CnWTd7oM889UGRvVsxbhkneNBqWClCUKd4FhJGbfPWk6TqAY8\nPk7neFAqmGkTkzrBC99uYO2ug8zQOR6UCnpag1DHLdt+gFfTNnPxgLacq3M8KBX0NEEoAAqPlXDn\nbGuOhwfH9LA7HKWUD9AmJgXAU5+vY2tOgc7xoJQ6TmsQioUb9vHu4u1cNbwTQ7roHA9KKYtfJggR\nGSMir+fn59sdit+z5nhYQWKLxkzROR6UUk78MkEYY+YbY66NidFhp+vqoXmryDl8lOcu6adzPCil\nTuCXCUJ5xqcrdpGakc0tZ+kcD0qpk2mCCFLWHA8r6ds2hhvP1DkelFIn0wQRhIwx3PPxSgqPlfLM\nJf0I0zkelFIu6C9DEJr1604WrNvLPed3J7GFzvGglHJNE0SQKZ/jYWiXOP6iczwopaqgCSKIlJYZ\n7py9nBARnr5Y53hQSlVNr6QOIm/+uIUl2/bzzMV9SdA5HpRS1dAaRJBYv/sQ//pyAyN7tmR8f53j\nQSlVPU0QQeBYSRl//yiD6MgGPDGut87xoJRyizYxBYGXFmxkdfZBXv/zAJ3jQSnlNq1BBLjfdhzg\nle82MWFAW87r2crucJRSfkQTRAArOlbKHR/pHA9KqdrRJqYA9tTna9maU8AH1wyiic7xoJSqIa1B\nBKgfNu7jncXbuXJYJ4Z2ibc7HKWUH9IEEYDyC4uZMnsFXZo34q5ROseDUqp2NEEEoKnzV7Pv8FGe\nu1TneFBK1Z4miADz2cpd/C89i1vOSqRP21i7w1FK+TFNEAFk76Ej3Pe/lfRpG8NNZybaHY5Sys/p\nWUx+7tLXFpOXV8QZZxju/a81x8Ozl/TVOR6UUnWmvyIB4qOlO/l23V7uHtWdxBbRdoejlAoAWoPw\nY6npWaTvyONYaRn3/Hclic0bMXloR7vDUkoFCK1B+KnU9Czu/Xglx0rLADDAzgNFzFuebW9gSqmA\noQnCTz395XqKiktPWHa0pIynv1xvU0RKqUCjTUx+puBoCV+v2UNWXpHL9dmVLFdKqZrSBOEHikvL\n+GHjPuZmZPPV6j0UFZcSKlBqTt62jc4Up5TyEE0QPsoYw287DpCans2nK3exv+AYsQ3DGNc/gZR+\nCWTtL+QfqatOaGaKCgtlykgdWkMp5RmaIHzMxj2HSM3IYm5GNpkHiohoEMK5PVoytl8CZ3RtTngD\nR7dRp2ZIiHDXnBUcKy0jITaKKSO7kZKs04kqpTxDE4QP2JVfxPzl2aSmZ7Nm10FCBIYlxnP7OV0Z\n2asVjSNcv00pyQnMXLKDvLw8vrz7rHqOWikV6DRB2CS/sJjPV+0iNSOLX7buxxjo2y6Wh8b0YHSf\nNjSP1qlBlVL20gRRj44Ul7Jg3V5S07NIW7+PY6VldIpvxG1nJzG2XwKd4hvV+DlnXTeEtLQ0zwer\nlAp6miC8rLTM8POWXFLTs/hi1W4OHS2heXQEkwZ3ICW5Db0TYhARu8NUSqmTaILwAmMMq7IOkpqR\nxfzl2ew9dJTGEQ0Y1asVY/u1YWiXeEJDNCkopXybJggP2p5bwNyMbFIzstiyr4CwUGFEtxak9Evg\n7FNa6OQ9Sim/ogmijnIOH+WT5dmkZmSTsTMPgEGdmnHN6Z05v1crYhuG2xyhUkrVjl8mCBEZA4xJ\nTLRnUpzDR0v4avVu5mZk8+OmHErLDN1bRXPP+d25qG8bvZpZKRUQ/DJBGGPmA/MHDhx4TX3t81iJ\nNdxFakY2X6/ZzZFi6+K06/7QmbH9EujWSudgUEoFFr9MEPWlrMywbMcBUtOz+GzlLg4UFhPbMIw/\n9m9LSnICA9o3JUQ7m5VSASooE0T5NJ0jRrhev2HPIVLTreEusvKKiAwL4dwerUjp14bTk5yGu1BK\nqQAWdAnCeRa2YU8tOD5+UXaeNdlOanoW63YfIjREGJ4Yz50ju3Juj8qHu1BKqUAVVL96FWdhy8or\nYsqc5bz07Ua25BZgDPRrF8vUMT24UIe7UEoFuaBKEK5mYSsuNWzfX8jt53Tlor5t6FiL4S6UUioQ\nBVWCqGy2tdIyw61nJ9VzNEop5duCqre1susT9LoFpZQ6WVAliCkjuxFVYbgLnYVNKaVcC6ompvLZ\n1nQWNqWUql5QJQiwkkRKcgJpaWmMqOxCCKWUUsHVxKSUUsp9miCUUkq5pAlCKaWUS5oglFJKuaQJ\nQimllEuaIJRSSrmkCUIppZRLmiCUUkq5pAlCKaWUS2KMsTuGWhORfGCj06IYIN/N2/FATh127/yc\nNV3val3FZfVVlurKUd02VcVd3f3y287L7CpLTd+TivcrlsXbn6+qtgnkz5erZf5QFk9/vqBuZUky\nxsRUu5Uxxm//gNcru1/dbWCpJ/ddk/Wu1tlVlurKUdOy1OS+U/zOy2wpS03fk+rK4u3PlyfL4k+f\nL38ti6c/X/VRFmOM3zcxza/ivju3Pbnvmqx3tc6usrjzHDUpS03uz69km9qqS1lq+p5UvO/PZfGn\nz5erZf5QFn/8fPl3E1NdiMhSY8xAu+PwBC2L7wmUcoCWxVfVR1n8vQZRF6/bHYAHaVl8T6CUA7Qs\nvsrrZQnaGoRSSqmqBXMNQimlVBU0QSillHJJE4RSSimXNEE4iEhnEXlTRObYHUtdiUiKiMwQkbki\ncp7d8dSWiJwiItNFZI6I3GB3PHUlIo1EZJmIjLY7lroQkREi8oPjvRlhdzy1JSIhIvK4iLwkIn+x\nO566EJHTHe/HGyLyk6eeN6AThIi8JSJ7RWRVheWjRGS9iGwSkXsAjDFbjDFX2RNp9WpYllRjzDXA\nZOBSG8KtVA3LsdYYcz1wCeBzpybWpCwOdwMf1W+U7qlhWQxwGIgEMus71qrUsBxjgQSgGB8rB9T4\nu/KD47vyCfCOx4Ko7ZV4/vAH/AHoD6xyWhYKbAY6A+HAcqCH0/o5dsftwbI8A/S3O/a6lAO4CPgJ\nmGh37HUpC3AOcBlW0h5td+x1LEuIY31L4D92x16HctwDXOfYxue+97X8zn8ENPFUDAFdgzDGLAT2\nV1h8GrDJWDWGY8CHWEcSPq0mZRHLP4HPjTG/1XesVanpe2KMmWeMGQpcXr+RVq+GZTkTGAxMBK4R\nEZ/67tWkLMaYMsf6A0BEPYZZrRq+J5lYZQAorb8o3VPT74qItAfyjTEHPRVDA089kR9JAHY63c8E\nBolIHPA4kCwi9xpjnrQluppxWRbgFqwj1hgRSTTGTLcjuBqo7D0ZAYzH+hH6zIa4asNlWYwxNwOI\nyGQgx+lH1pdV9r6MB0YCscDLdgRWQ5V9T14AXhKR04GFdgRWC5WVBeAq4N+e3FkwJghxscwYY3KB\n6+s7mDqqrCwvAi/WdzB1UFk50oC0+g2lzlyW5fgNY96uv1DqrLL35WPg4/oOpg4qK0ch1o+qP6n0\n82WMecjTO/Opam49yQTaOd1vC2TbFEtdBUpZAqUcoGXxRYFSDqjnsgRjgvgVSBKRTiISjtVxOM/m\nmGorUMoSKOUALYsvCpRyQH2Xxe6eei+fBTAT2MXvp7Fd5Vh+AbAB62yA++yOM5jKEijl0LL45l+g\nlMNXyqKD9SmllHIpGJuYlFJKuUEThFJKKZc0QSillHJJE4RSSimXNEEopZRySROEUkoplzRBqKAl\nIoe98JwdRWRiDR/zmYjEejoWpepKE4RSntURa8RWtxljLjDG5HknHKVqTxOECnqOGdLSHDPXrROR\n/4iIONZtE5F/isgSx1+iY/nbIjLB6TnKayNPAaeLSIaI3F5hP61FZKFj3SrHKKLl+4gXkesd6zJE\nZKuIfOdYf56ILBaR30Rktog0ro/XRSlNEEpZkoG/YU0k0xkY5rTuoDHmNKyhrZ+v5nnuAX4wxvQz\nxjxXYd1E4EtjTD+gL5DhvNIYM92x7lSsoRWeFZF44H7gHGNMf2Ap8PfaFFCpmgrG4b6VcmWJMSYT\nQEQysJqKfnSsm+n0v+KPfk38CrwlImFAqjEmo5LtXgAWGGPmizV/dQ9gkaNSEw4srkMMSrlNE4RS\nlqNOt0s58bthXNwuwVEDdzRHhVe3A2PMQhH5A3Ah8J6IPG2Medd5G8eEQh2Am8sXAV8bY/7kflGU\n8gxtYlKqepc6/S8/et8GDHDcHguEOW4fAqJdPYmIdAD2GmNmAG9izTfsvH4AcCcwyfw+49zPwDCn\nvo+GItK1rgVSyh1ag1CqehEi8gvWAVX5kfwMYK6ILAG+BQocy1cAJSKyHHi7Qj/ECGCKiBQDh4Er\nKuznZqAZ8J2jOWmpMeZqR61ipoiUz/98P9Zwz0p5lQ73rVQVRGQbMNAYk2N3LErVN21iUkop5ZLW\nIJRSSrmkNQillFIuaYJQSinlkiYIpZRSLmmCUEop5ZImCKWUUi5pglBKKeXS/wMKgdMrOShGSQAA\nAABJRU5ErkJggg==\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "r, e = rat_err(tnva, tnvea, tva, tvea)\n", "\n", "plt.errorbar(sizes, r, yerr=e, fmt='o-', label=\"Python\")\n", "plt.xscale('log')\n", "plt.yscale('log')\n", "\n", "plt.grid(True)\n", "plt.axhline(1, color='red')\n", "plt.legend()\n", "plt.xlabel(\"Input size\")\n", "plt.ylabel(\"Timing ratio\")\n", "plt.title(\"Numbifyed Numpy vs regular Numpy\");" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Numerical accuracy\n", "\n", "We should finally take a quick look at accuracy considerations. In this case we know the exact answer (680), so we can see how each method compares. As $n \\rightarrow \\infty$ we expect them all to converge to 680, but let's see if any of them accumulates error worse than the others. Remember, since the Numba version is algorithmically nothing but the Python code, any differences between Python and Numba would be due to Numba's compilation/optimizations:" ] }, { "cell_type": "code", "execution_count": 23, "metadata": { "collapsed": true }, "outputs": [], "source": [ "def err(val, exact=680):\n", " return abs((val-exact)/exact)\n", "\n", "# This time, let's make empty arrays and fill them in\n", "epy = np.empty_like(sizes)\n", "enpy = np.empty_like(sizes)\n", "enba = np.empty_like(sizes)\n", "\n", "for i, (x, y) in enumerate(xy):\n", " epy[i] = err(trapzl(y, x))\n", " enpy[i] = err(trapz(y, x))\n", " enba[i] = err(trapzn(y, x))" ] }, { "cell_type": "code", "execution_count": 24, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAAEOCAYAAACjJpHCAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3Xl8VOXZ//HPNdmBEPYdgUBCBcIO\n4hIJBUFEFKhSi4ooS6miaKmKj9WmpT72+WlR0LauyCKiQkEUVEQlLKKiJGyKRAhBw6bsW0K2+/fH\nmYQQs0wyy5mZXO/Xa16ZOTNzzvdOMrlyn3Pu+4gxBqWUUqq6HHYHUEopFdi0kCillHKLFhKllFJu\n0UKilFLKLVpIlFJKuUULiVJKKbdoIVFKKeUWLSRKKaXcooVEKaWUW7SQKKWUckuo3QF8oVGjRqZt\n27bVfv/Zs2epXbu25wLZJFjaAdoWfxQs7QBtS5HNmzcfMcY0rux1QV1IRGQ4MLxDhw58/fXX1V5P\nSkoKSUlJHstll2BpB2hb/FGwtAO0LUVEZJ8rrwvqXVvGmPeMMZNiYmLsjqKUUkErqAuJUkop79NC\nopRSyi1BfYykInl5eWRlZZGTk1Ppa2NiYti5c6cPUnmXq+2IjIykVatWhIWF+SCVUirQ1dhCkpWV\nRXR0NG3btkVEKnzt6dOniY6O9lEy73GlHcYYjh49SlZWFu3atfNRMhXIVqY8xqyMZRxyQLM5MDV2\nJMOSZtgdS/lQjd21lZOTQ8OGDSstIjWNiNCwYUOXempKrUx5jOS9yzgYIhgRDoYIyXuXsTLlMbuj\nKR+qsYUE0CJSDv2+KFfNylhGjuPi35cchzArY5lNiZQdanQhsZuIMG3atOLHTz/9NMnJyfYFUqqK\nDpXzF+SQA3jvftj0Muz7HHJO+jSX8q0ae4ykqt5J289Tq3Zx4EQ2LepF8eCQjozo0dKtdUZERLB0\n6VIeeeQRGjVq5KGkSvlOs0LDwZBf9mCbFRTCN8tg82sXFta7BJp2gaadnbcEaNAOHCE+TKy8QXsk\nLli54zCPLN3O/hPZGGD/iWweWbqdd9L2u7Xe0NBQJk2axDPPPPOL58aNG8eSJUuKH9epUwewRqn2\n79+f0aNHEx8fz/Tp01m4cCF9+/YlISGBPXv2FL9/8uTJJCYmEh8fz4oVKwBITExky5Ytxeu98sor\n2bZtm1vtUDVQQT58+nemHj1GZGHhRU9FFhqmtv8NPJwJD3wLYxbDwL9Aqz5wLAPWz4TF4+D5XvC/\nLeGlAbB8Cnz5ImRugHPHbGmSqj7tkQB/fe8bvj1wqtzn0344Tm6BuWhZdl4BDy3ZxqJNP5T5nk4t\n6vKX4Z0r3fY999xD165deeihh1zOu3XrVnbu3EmDBg2IjY1lwoQJbNq0iVmzZvHcc8/x7LPPApCZ\nmcnatWvZs2cPAwYMIC0tjQkTJjB37lyeffZZ0tPTOX/+PF27dnV520px6iD8dzzs+4xhPW6DOnWZ\nte8966ytwlJnbcW0tG7xgy+8Py8HjuyCQzvg8DdweAfseh/SFlx4Td2WpXovXaBhBwjRP1n+SH8q\nLihdRC4sLyxzeVXUrVuXsWPHMnv2bKKiolx6T58+fWjevDkA7du3Z/Bg60OakJDAmjVril83evRo\nHA4HcXFxxMbGkp6ezs0338yMGTN46qmnmDNnDuPGjXO7DaoG2f0xLP095GXDyBeh2y0MA4bxpOtz\nOoVFQvNu1q2IMXDmsFVUDn9j3Q7tgD2fQGG+9ZqQCGjyqxIFpot1q93QGy1VVaCFBCrtOVz+vx9z\n8NT5XyxvWS+Kt35/udvbv//+++nZsyd33nln8bLQ0FAKnbsMjDHk5uYWPxcREVF83+FwFD92OBzk\n5+cXP1f67CsRoVatWlxzzTUsX76ct99+263JLFUNUpAPa56ADTOhSSe4eR40jvfc+kUgupl16zDo\nwvL8XDiS7iwu262vuz+GLQsvvKZOM6uwNOtyocg0jIPQcM/lUxXSQuKCqQPa8tf3d5OdV1C8LCos\nhAeHdPTI+hs0aMDo0aN59dVXueuuuwBo27YtmzdvZvTo0Sxfvpy8vLwqr3fx4sXccccd7N27l4yM\nDOLi4gCYMGECw4cPJzExkQYNGnikDSqIndxv7cr64XPoORau/T8Ir+WbbYeGWwWiWRfgtxeWn/n5\n4t7L4e3wxX+gwPkPlyMMGv+qxK6xztAsAeo08U3uGkYLiQuGdWlKZGSUx8/aKmnatGk8//zzxY8n\nTpzIjTfeSN++fRk4cGC1rifQsWNH+vfvz+HDh3nhhReIjIwEoFevXtStW/eiHpBSZfp+NSydZP2B\nHvUKdL3Z7kSWOo2hzgBoP+DCsoI8OLrbuVvM2XvZuw62vXnhNbUbX7xbrGlnaNwRQq1evY7Srx4t\nJC4a0aOlRwsHwJkzZ4rvN23alHPnzl30+Isvvih+/OSTTwKQlJR00X7olJSU4vuln7vyyisvOiPs\n9OnTABw4cIDCwsLiYytK/UJBHnz6d/jsWesP7s1zoVGc3akqFhIGTS61bgk3XVh+7liJ3ovz61ev\nQL5z9gYJgUbxrIwMIdlxgpwQ62TWgyGQvNcaWKnFpGJaSGqY+fPn8+ijjzJz5kwcDj37W5XhZBYs\nuQt+/BJ63QnXPglhrp0I4pdqNYB2V1u3IgX51qnIh3cUF5dZ2dvIcVz8J7FolP6w/n+zjuOoMmkh\nCVJz584tc/nYsWMZO3asb8OowJG+Cpb93uqR/ObVi/+zDyYhodbJAo3jocsoAA7N7VLmSw85gGcT\nIP5a6DgU2l5VvCtMWbSQKKWswvHJ32DjbGvE+eh50LC93al8qlmhtTvrF8sLDDTrBmmvw1cvQ3gd\n6DAQ4odC3GA9/RgtJEqpEz9au7KyNkHv8TDkf62xHjXM1NiRJO+9eBJKa5T+KEiaYY2d2bsOdn0A\n6R/Ct8tBHND6Mmdv5TrrOFIN3AWmhUSpmmzXB7BsMhQWwE2vFe/mqYmKDqgXn7VVepR+WBTED7Fu\nxsDBLdb3b9cH8PFfrFuDWKugxF8Ll1xeY0bi14xWKqUulp8Ln/wVPn/eGmF+02s1bldWWYYlzWBY\n0ozKR+mLQIse1m3A/1gnKKR/aBWVTS9Z39fIGGvXV/y11iDLqHo+a4evaSGxUUhICAkJCeTn53Pp\npZcyb948atUqe6BXZmYmGzduZMyYMYB1MP3rr7++aOyJUi45vg+W3An7N0PfSTD473rw2F0xraDP\nBOt2/gxkrHHuAlsF2xeDIxTaXGEdV+k41Jr1OIjo+Z+u2vY2PNMFkutZX7e97fYqo6Ki2LJlCzt2\n7CA8PJwXXnih3NdmZmbyxhtvuL1NVcPtXAEvJsKR761pTq57SouIp0XUgUuHw4h/w5/S4a6P4Ip7\n4cxPsOoRmN0d/nUZfJwMP3xp7VYMcFpIXBC6cxm8dx+c/BEw1tf37vNIMSmSmJjI7t27eeyxx5g1\na1bx8kcffZTZs2czffp01q9fT/fu3YsHGR44cIBrr72WuLi4i2YPXrRoEQkJCXTp0oWHH364eHnz\n5s159NFH6datG/369ePw4cMey6/8XH4ufDAd3roV6reD36+FziPsThX8HCFwyWUwKBnu+RLuS4Mh\nT1pTtWx8DuYMhqfj4Z274dt3rd5MAPL7XVsiEgs8CsQYY24qsbw2sA74izFmhVsb+WC6NaVCOSKz\nNl2Yw6dIXrZ1DYXN88p+U7MEGPoPlzafn5/PBx98wLXXXsvQoUMZNWoUU6dOpbCwkDfffJNNmzbR\ntWtXnn766eLrisydO5ctW7aQlpZGREQEHTt25N577yUkJISHH36YzZs3U79+fQYPHsw777zDiBEj\nOHv2LP369eOJJ57goYce4uWXX+bPf/6zSxlVADueCYvvhAOpcNlkuOZv2guxS4NYuPxu65Z9wpqA\nMv1D+G6FNRFlSLg1cLLjUOvYSkwruxO7xKuFRETmANcDPxljupRYfi0wCwgBXjHGlPsX1xiTAYwX\nkSWlnnoY8FyXoCKli0jx8l/OCFwV2dnZdO/eHbB6JOPHjyc8PJyGDRuSlpbG4cOH6dGjBw0bln2e\n+sCBA4mJiQGgU6dO7Nu3j6NHj5KUlETjxo0BuPXWW1m3bh0jRowgPDyc66+/HrDm21q9erVb+VUA\n+PZd6x8egNELoNMN9uZRF0TVswZ8JtxkjeP54QvnWWDvw8pp1q1ZwoWzwJp3Bz+djcLbPZK5wPPA\n/KIFIhIC/Au4BsgCvhKRd7GKypOl3n+XMean0isVkUHAt4BnTnavpOdg/tkJOV3G1RBjWsOdK6u9\n2aJjJKUVXXzq0KFDxbMBl6XkdPIhISHk5+djTNnXTgEICwsrnlq+6PUqSOWfh48eg00vQouecPNr\nUL+t3alUeULCoF2idRvyhDV1ftGpxeuegrX/B9HNrVOPO15n9Vr8aNoarxYSY8w6EWlbanFfYLez\np4GIvAncaIx5Eqv34ooBQG2gE5AtIu8bY9y/ylQ5zidOJ2r1Q9burCJhUTDwca9sb+TIkTz++OPk\n5eUVH2CPjo4unnSxIpdddhlTp07lyJEj1K9fn0WLFnHvvfd6JafyU8cyrF1ZB7dAv7th0F/12hyB\nRMSakbhxR7jqfjh7FL7/yOqpbF8Cm+dCaJQ183HHoRA3BKKb2hrZjmMkLYEfSzzOAi4r78Ui0hB4\nAughIo8YY540xjzqfG4ccKSsIiIik4BJYM2kW3KWXICYmBiX/jADFMRbuwMi1v8DOX0AE92C84nT\nyW83FFxcR3nKy3DVVVcRExNTPCNwu3btEBESEhIYM2YM9evXJzc3t/j9+fn5nDt3jjp16vD444/T\nv39/jDEMHjyYX//618WvK/qanZ1NXl5eudvPycn5xffMn5w5c8av81WFJ9vS+KfP6LjreYwI33X5\nH45GXgYbNnpk3ZXRn4k3NYem45HGY6l3YjuNjnxFw8xNRO56H4BT0fEcadSHow37crZ2GxBh1/4F\nLMn5gsMhQtNXDTdF9qNjy9u9kk4q2hXikQ1YPZIVRcdIRORmYIgxZoLz8e1AX2OM1/5t7t27tyl9\nJcCdO3dy6aWXuvT+06dPEx0d7Y1oZSosLKRnz54sXry4+GJUnlCVdlTl+2MHly/rGgA80pa8HPjo\nUWt69Ja94aY5UL+NR/K5Sn8mPmaM83r3zl1gB1Kt5TGXsLJBE5ILD5FT4phKZKEhuV3Vrq8iIpuN\nMb0re50dR26ygNYlHrcCDtiQwy99++23dOjQgYEDB3q0iKggdnQPvHqNVUQunwJ3fuDzIqJsIGId\njO//EExaA9N2wfBZ0LQzs/IOXFRE4MKU+N5gx66tr4A4EWkH7AduAcbYkMMvderUiYyMDLtjqECx\n47/w7lRrvMLv3rT2mauaKboZ9BoHvcZVPCW+F3i1RyIii4DPgY4ikiUi440x+cAUYBWwE3jbGPON\nN3MoFXTycmDFA9asvU1+BZM3aBFRxZqVc+pRecvd5e2ztn5XzvL3gfe9uW0AERkODO/QoYO3N6WU\n7xzZDYvHweHtcMV91tmDIWF2p1J+pNwp8WNHemV7/jm6xUOMMe8ZYyYVDdpTKuBtXwIv9YdTWTDm\nbRg8Q4uI+oVhSTNIbjeS5gUGMYbmBVU/0F4Vfj9FilIKawzTh9OtMQSt+8FNrwbM9BnKHi5Pie8B\nQd0j8XciwrRp04ofP/300yQnJ3tk3ePGjWPJktKzyqiAdOR7eGWQVUSuegDGrdAiovyKFhIXrcxY\nyeAlg+k6ryuDlwxmZUb1p0YpEhERwdKlSzly5IgHEqqgtO1teLE/nDoAty6xZpHVXVnKzwR1IRGR\n4SLy0smTJ91az6ofVpG8MZmDZw9iMBw8e5DkjcluF5PQ0FAmTZpUPC18SaV7FHXq1AGsgVL9+/dn\n9OjRxMfHM336dBYuXEjfvn1JSEhgz549xe/5+OOPSUxMJD4+vnjW4MzMTBITE+nZsyc9e/Zk40bf\njHpWVZR7zppscelEaN7VOisr7hq7UylVpqAuJJ462P7CNy+QU5Bz0bKcghxmpc4q5x2uu+eee1i4\ncCFVKXZbt25l1qxZbN++nQULFpCens6mTZuYMGECzz33XPHrMjMzWbt2LStXrmTy5Mnk5OTQpEkT\nVq9eTWpqKm+99Rb33Xef221QHvbzLnhlIKQtgMRpcMcKiGlpdyqlyqUH213wU/YvJiAG4NDZQ26v\nu27duowdO5bZs2cTFeXabJ59+vShefPmALRv357BgwcDkJCQwJo1a4pfN3r0aBwOB3FxccTGxpKe\nnk6XLl2YMmUKW7ZsISQkhPT0dLfboDxoyyJY+UcIqwW3/de61rdSfk4LiQuaRDXhcPYvrybYrHYz\nj6z//vvvp2fPntx5553Fy0JDQykstEYPGWPIzb1wTZSS08c7HI7ixw6H46Kp4YumjC/5+JlnnqFp\n06Zs3bqVwsJCIiM9MxO/qrqVKY8xK2MZhxzQbA5MlYYM27cF2lwFv3kF6ja3O6JSLgnqXVueMrnz\nZCJDLv6DGxkSydSeUz2y/gYNGjB69GheffXV4mVt27Zl8+bNACxfvpy8vLwqr3fx4sUUFhayZ88e\nMjIyiIuL4+TJkzRv3hyHw8GCBQsoKAj860UHopUpj5G8dxkHQwQjwsEQIZkjrOxwBYxdrkVEBZSg\nLiSeOtg+5JIhJF+RTPPazRGE5rWbk3xFMsNih3koKUybNu2is7cmTpzI2rVr6du3L19++SW1a9eu\n8jo7duxI//79GTp0KC+88AKRkZHcfffdzJs3j379+pGenl6t9Sr3zcq4eNQxQI7DwazcHyFEdxSo\nwOL1aeT9QaBNI+8tOo28/+g6twum1K5HADGGbeN22JDIfYH+MylJ22Lx52nklarxfD2pnlLepIVE\nKRtMbTeSyMKLq4Y3J9VTypu0kChlg2GX/JrkI8donl/ok0n1lPKmGn1Uzxjzi1NklfV9UV6WOp9h\nBeEMG7uZlM++CJr98apmqrE9ksjISI4ePap/NEsxxnD06FEdX+JNZ4/Cdyug2y0Qpt9nFfiCukdS\n0YWtWrVqRVZWFj///HOl68nJyQmKP6yutiMyMpJWrXR2Wa/Z9hYU5EKP2+1OopRHBHUhMca8B7zX\nu3fviaWfCwsLo127di6tJyUlhR49eng6ns8FSzsCmjGQOh9a9oJmZV9XW6lAU2N3bSlli6yv4eed\n0HOs3UmU8hgtJEr5Uuo8CKsNXX5jdxKlPEYLiVK+cv407FgKXUZCRODPlKBUES0kSvnKjqWQdxZ6\n3mF3EqU8SguJUr6StgAa/wpa9bE7iVIepYVEKV84/C1kfWUdZNdBsCrIBHUh8dQ08kq5LW0BOMKg\n6y12J1HK44K6kHjqmu1KuSX/PGxdBJdeD7Ub2p1GKY8L6kKilF/4bgVkH9exIypoaSFRyttS50PM\nJdAuye4kSnmFFhKlvOl4JmSkQM/bwaEfNxWc9DdbKW9Kex3EAd3H2J1EKa/RQqKUtxTkQ9pC6DAI\nYnQ2ZRW8tJAo5S17PoHTB/Qguwp6WkiU8pbU+VC7McRfa3cSpbxKC4lS3nD6MOz6ALr9DkLC7E6j\nlFcFdSHRke3KNlvfAFOgu7VUjRDUhURHtitbFF0F8ZIroFGc3WmU8rqgLiRK2WLfZ3AsQ3sjqsbQ\nQqKUp6XOh4i60OlGu5Mo5RMVFhIRcYjIaF+FUSrgZZ+Ab5dDws0QXsvuNEr5RIWFxBhTCEzxURal\nAt/2xZCfo7u1VI3iyq6t1SLyJxFpLSINim5eT6ZUIEqdD826QovudidRymdCXXjNXc6v95RYZoBY\nz8dRKoAd2AKHtsF1T9udRCmfqrSQGGPa+SKIUgEvdT6ERlrHR5SqQSotJCISBvwBuNq5KAV40RiT\n58VcSgWW3HPW8ZFOIyCqnt1plPIpV3Zt/QcIA/7tfHy7c9kEb4VSKuB8uxzOn9KD7KpGcqWQ9DHG\ndCvx+FMR2eqtQEoFpNT50KA9tLnC7iRK+ZwrZ20ViEj7ogciEgsUeC+S5+hcW8onjnwPP2y0eiMi\ndqdRyudcKSQPAmtEJEVE1gKfAtO8G8szdK4t5ROp88ERas30q1QNVOGuLRFxANlAHNAREOA7Y8x5\nH2RTyv/l58LWRdY1R6Kb2p1GKVtUWEiMMYUi8k9jzOXANh9lUipwpH8IZ3+GnnfYnUQp27iya+sj\nEfmNiO78VeoXUudDdAvoMNDuJErZxpWztv4I1AbyRSQHa/eWMcbU9WoypfzdySzY/TFc/SA4QuxO\no5RtKjtGIkBnY8wPPsqjVOBIW2h97XGbvTmUsllls/8aYJmPsigVOAoLIG0BxCZB/TZ2p1HKVq4c\nI/lCRPp4PYlSgSQjBU7+CD1vtzuJUrZz5RjJAGCyiGQCZ7lwjKSrN4Mp5ddS50NUffjV9XYnUcp2\nrhSSoV5PoVQgOXsEvlsJfSdCaITdaZSyXaW7towx+4DWwK+d98+58j6lgta2t6AwD3robi2lwIWC\nICJ/AR4GHnEuCgNe92YopfyWMdZurVZ9oGknu9Mo5Rdc6VmMBG7AOj6CMeYAEO3NUEr5rayv4Ofv\ndLp4pUpwpZDkOk8DNgAiUtu7kZTyY6nzILwOdB5ldxKl/IYrheRtEXkRqCciE4GPgZe9G0spP5Rz\nCnYshS6jIKKO3WmU8huuXLP9aRG5BjiFNQPw48aY1V5PppS/+WYp5J3TCRqVKsWV039xFo6AKx4i\nMhwY3qFDB7ujqGCQOh+adIKWvexOopRfCerTePXCVspjDu2A/Zv1KohKlSGoC4lSHpO2AELCoetv\n7U6ilN9xqZCISJSIdPR2GKX8Ul4ObH0TLh0OtRrYnUYpv+PKgMThwBbgQ+fj7iLyrreDKeU3vlsB\nOSd07IhS5XClR5IM9AVOABhjtgBtvRdJKT+TOg/qtYG2V9udRCmXvZO2nyv/8SnjPjzLlf/4lHfS\n9nttW64UknxjzEmvJVDKnx3LgL3rrOniHXpIUQWGd9L288jS7ew/kQ3A/hPZPLJ0u9eKiSufjB0i\nMgYIEZE4EXkO2OiVNEr5m7TXQRzQ/Va7kyjlsqdW7SI7r+CiZdl5BTy1apdXtudKIbkX6AycB94A\nTgL3eyWNUv6kIN+6nG7cYKjbwu40SrnsgLMn4upyd7kyILGjMeZR4FGvJFDKX+1eDWcO6UF2FXBa\n1Isq3q1Verk3uNIjmSki34nIDBHp7JUUSvmj1PlQp6nVI1EqgIy9/JJfLIsKC+HBId4ZxeHKha0G\nAEnAz8BLIrJdRP7slTRK+YvThyB9FXQfAyFhdqdRqkq27T9FWIjQLCYSgJb1onhyVAIjerT0yvZc\nnWvrEDBbRNYADwGPA3/3SiKl/MGWN8AU6FUQVcDZefAUK7cd5J4B7XlwyK9ISUkhKSnJq9t0ZUDi\npSKSLCI7gOexzthq5dVUStmp6CqIba6Chu3tTqNUlTyzOp3oiFAmJsb6bJuu9EheAxYBg51XR1Qq\nuGVugON7IemRyl+rlB/ZnnWSj749zAOD4qlXK9xn23XleiT9fBFEKb+ROh8iYqDTDXYnUapKZq7e\nRb1aYdx1VVufbrfcQiIibxtjRovIdpyX2S16CjDGmK5eT6eUr2Ufh2+XW6f8hnnnVEmlvGHzvuOs\n2fUzD13bkehI354gUlGPZKrz6/W+CKKUX9i2GArO69gRFXCeWZ1Ow9rh3HF5W59vu9yD7caYg867\ndxtj9pW8AXf7Jp5SPmSMNUFj8+7QXDvcKnB8kXGUDbuP8Iek9tSOcOlkXI9yZUDiNWUsG+rpIErZ\n7kAaHN6hvREVUIwxzPwonSbREdzWr40tGSo6RvIHrJ5HrIhsK/FUNPCZt4Mp5XOp8yE0ChJusjuJ\nUi7bsPsImzKP8dcbOhMZFmJLhor6QG8AHwBPAtNLLD9tjDnm1VRK+VruWdi+BDqPhMgYu9Mo5RJj\nDP/8KJ0WMZHc0re1bTkqOkZy0hiTaYz5nfO4SDbW2Vt1ROSXE7koFci+eQdyT+tuLRVQ1uz6iS0/\nnuDegXFEhNrTGwEXL7UrIt8De4G1QCZWT0Wp4JE6HxrGwSU6bEoFBmMMM1enc0mDWtzUy97JRlw5\n2P53oB+QboxpBwxEj5GoYPLzLvjxC6s3ImJ3GqVcsuqbw+zYf4r7BsYRFmLv1Ttd2XqeMeYo4BAR\nhzFmDdDdy7mU8p3U+eAIhW6/szuJUi4pLDQ8szqd2Ea1GdHd/ouuuVJITohIHWAdsFBEZgH53o11\ngYjEisirIrKkxLIkEVkvIi+ISJKvsqgglJ8LWxdBx+ugTmO70yjlkpXbD7Lr8GmmDooj1ObeCLhW\nSG7EOtD+APAhsAcY7srKRWSOiPzknDm45PJrRWSXiOwWkenlvR/AGJNhjBlfejFwBogEslzJolSZ\ndr0P545CzzvsTqKUSwoKDc9+nE580zoM72p/bwRcm7TxbImH86q4/rlYU8/PL1ogIiHAv7AGOmYB\nX4nIu0AI1qnGJd1ljPmpjPWuN8asFZGmwEzg1irmUsqSOh/qtoL2A+xOopRLlm/Zz56fz/KfW3vi\ncPjHMb2KBiSepozJGrkwaWPdylZujFknIm1LLe4L7DbGZDi38yZwozHmSVyc18sYU+i8exyIcOU9\nSv3CiR9hz6fQ/2Fw2HfqpFKuyiso5NmPv6dT87oM6dzM7jjFyi0kxphoL22zJfBjicdZwGXlvVhE\nGgJPAD1E5BFjzJMiMgoYAtTD6vGU9b5JwCSApk2bkpKSUu3AZ86ccev9/iJY2gGeaUubzDdpC3xx\nvj3nbfy+BMvPJVjaAf7blrU/5vHDsVym9oxg3bq1Lr3HJ20xxlR6A64C7nTebwS0c+V9zte3BXaU\neHwz8EqJx7cDz7m6vurcevXqZdyxZs0at97vL4KlHcZ4oC0F+cbM7GzM/JEeyeOOYPm5BEs7jPHP\ntuTk5ZsrnvzE3PD8BlNYWOjy+9xpC/C1ceFvrCsDEv8CPAwUXS4uHHjdjdqVBZQcy98K0CsvKt/K\nWAMnf9SR7CpgvP3Vj+w/kc20a+IRPxvv5MpZWyOBG4CzAMa63K47u72+AuJEpJ2IhAO3AO+6sT6l\nqi51PtRqaJ32q5Sfy8kr4Pkjt1qtAAATdklEQVQ1u+nTtj6JcY3sjvMLrhSSXGcXxwCISG1XVy4i\ni4DPgY4ikiUi440x+cAUYBWwE3jbGPNN1aMrVU1nfobv3rcGIIb67rrWSlXXwi9/4PCp8/zxmo5+\n1xsBF07/Bd4WkReBeiIyEbgLeMWVlRtjyhwqbIx5H3jf5ZTVJCLDgeEdOnTw9qZUINn2JhTmQY/b\n7U6iVKXO5ebzn5TdXNG+IZe3b2h3nDJV2iMxxjwNLAH+C3QEHjfGzPZ2ME8wxrxnjJkUE6PTgisn\nY6zdWq0vgya/sjuNUpWa//k+jpzJZdrgeLujlMulazIaY1YDq8EaUCgitxpjFno1mVLe8OOXcCQd\nbvyX3UmUqtTpnDxeXLuH/vGN6dWmgd1xylVuj0RE6orIIyLyvIgMFssUIAMY7buISnlQ6nwIj4ZO\nI+xOolSlXvssk+Pn8vjjNf7bG4GKeyQLsEaOfw5MAB7EOvX3RmPMFh9kU8qzck7CN8ug62iIqGN3\nGqUqdPJcHi+vz2DQpU3p1rqe3XEqVFEhiTXGJACIyCvAEeASY8xpnyTzAD3Yri6y47+Qd07HjqiA\n8MqGDE7n5Pt9bwQqPtieV3THGFMA7A2kIgJ6sF2VkjofmnaBFj3tTqJUhY6dzWXOhr0MS2hOpxaV\nTmtou4p6JN1E5JTzvgBRzscuT9qolN84uA0OpMHQ/6dXQVR+78V1eziXV8D9g+LsjuKSiiZt1OlQ\nVfBIWwAhEZBws91JlKrQz6fPM3/jPm7s1oK4pt6aO9ezXDr9VwW2lSmPMStjGYcc0GwOTI0dybCk\nGXbH8p28bNj2FnS6AWr57ymUSgH8J2UPuQWFTB3k/8dGimghCXIrUx4jee8yckKs3TkHQyB57zKA\nmlNMdr5nnbGlB9mVnzt0MofXv9zHqB4tadfI5dmobGf/xX69SESGi8hLJ0+etDuKbWZlLCOn1FXU\nchzCrIxlNiWyQep8qN8O2lxldxKlKvSvNbspLDTcNzAwjo0UCepComdtwaFyfsLlLQ86R/dA5nro\neTs4akqjVSDKOn6ON7/6gdF9WtO6QS2741SJfrKCXLPCcpYXFELuWd+GsUPa6yAh0G2M3UmUqtBz\nn+xGEKYMCLxxb1pIgtzU2JFEFpqLlkUWFjL12HF4vo81SM+Yct4d4AryYctCiB8CdZvbnUapcmUe\nOcuS1CzGXHYJLepF2R2nyrSQBLlhSTNIbjeS5gUGMYbmBYbkdqMY9ttl1oWdltwFc6+HQzvsjup5\n338EZw7rQXbl92Z/8j1hIcLdA9rbHaVa9KytGmBY0gyGJc0gJSWFpKSkC09MSoHUefDJDHgxEXqP\nhwH/EzynyKbOhzrNoMM1didRqly7fzrDO1v2MyExlibRkXbHqZag7pHoWVuVcIRA77vg3s1WEfn6\nVXiuF3w9BwoL7E7nnlMH4PtV0ONWCNH/l5T/evbjdCLDQvj91bF2R6m2oC4ketaWi2o1gGFPw+/X\nQ5NLYcUD8FIS/PCF3cmqb8sbYAqhx212J1GqXN8dOsWKbQe588q2NKwTYXecagvqQqKqqFkXGLcS\nbpoD547CnCGwdBKcOmh3sqopLLSmRGl3NTQI3P/yVPB7ZnU60RGhTEwM7N9TLSTqYiLQ5Tcw5StI\n/JN1/Y7nesGGZyD/vN3pXJO5Ho5nQg89yK781479J1n1zWHGJ7ajXq1wu+O4RQuJKlt4bRj4GNzz\nJcT2h4+T4d+Xw/er7U5WudT5EFkPLh1udxKlyjVzdToxUWHcdVU7u6O4TQuJqliDWPjdIrj1v1Zv\nZeFN8MZvrRHj/ujcMdj5LnT9LYQF5hkwKvil/nCcT7/7iUlXx1I3MszuOG7TQqJcEzcI/vA5XPM3\nyNwA/+4HH/8Vzp+xO9nFtr0NBbk6dkT5tZkfpdOwdjjjrmhrdxSP0EKiXBcaDldOtU4X7jwKNsy0\nRsdvX+Ifo+ONscbFtOhpnTiglB/6MuMoG3Yf4Q9J7akdERynpgd1IdFxJF4S3QxGvQh3fQR1GsN/\nx8Nr18Gh7fbm2p8KP32rvRHlt4wx/HN1Ok2iI7itXxu743hMUBcSHUfiZZdcBhPXwPBZcGQXvHg1\nrPijdZzCDqnzIKyWddaZUn7os91H2bT3GPcM6EBkWPBchDaoC4nyAUcI9Bpn7e7qMxE2z4XnesJX\nr/h2dPz5M9YElJ1HQWRd321XKRdZvZFdtIiJ5Ja+re2O41FaSJRnRNWH6/4fTF4PTbvAymnwYn/Y\nt9E32/9mGeSe0d1aym+l7PqZtB9OMOXXcUSEBk9vBLSQKE9r2hnueA9ungvZx+G1ofDfCdbcV96U\ntgAadYTWfb27HaWqwRjDzNXptG4Qxc29W9kdx+O0kCjPE4HOI2HKJrj6Ifj2XXiuN6yf6Z3R8T99\nBz9+afVGRCp/vVI+9tG3h9m+/yT3/TqOsJDg+7MbfC1S/iO8Nvz6UWt0fPsB8MlfrfEn6as8u520\nBeAIg263eHa9SnlAYaHhmdXpxDaqzcgeLe2O4xVaSJT3NWgHtyyE25Zal719YzQsvNkzo+Pzz8PW\nRfCrYVC7kfvrU8rDVm4/yHeHTjN1UByhQdgbAS0kypc6DIQ/bITBf4d9n8O/LoPVf3FvdPyu962Z\nivUgu/JDBYWGZz9OJ75pHYZ3bWF3HK/RQqJ8KzQcrrjXOl044Wb47Fl4vrc1tUl1RsenzoeYSyB2\ngOezKuWm5Vv2s+fnszwwKB6HI3iP3wV1IdGR7X4suimM/A+M/9gaKb90Isy5Fg5udX0dx/fBnjXW\nxascQf2rrAJQXkEhsz75nk7N6zKkczO743hVUH/6dGR7AGjdByZ8Cjc8B0d3W2NPVjzg2uj4LQut\nrz1u9W5GpaphaWoW+46e44/XBHdvBIK8kKgA4XBYxzju3QyXTYbN82B2D9j0MhTkl/0eUwBpr0OH\nQRATfOflq8CWm1/I7E920611PQZe2sTuOF6nhUT5j6h6MPQfMHkDNEuA9/8EL/WHzM9+8dIGx7bA\nqf16kF35pbe+/pH9J7L54zXxSA0Y26SFRPmfpp2co+PnQc5JmHsdLLkLTu5nZcpjDJ7ThdtOv8bg\n1i1YeehLu9MqdZGcvAL+9eluerepz9VxNeOU9OCYDF8FHxHoPALiBltndm14lpWZH5HcMIYc57n4\nB0NDSd63HFIcDEuaYXNgpSxvfPkDh07lMPO33WpEbwS0R6L8XXgtGPA/MGUTs+rXJafU2Vk5DmFW\nxjKbwil1sXO5+fw7ZQ9XtG/IFe1rRm8EtJCoQFG/LYfKGRV8SH+LlZ+Y//k+jpw5z7TB8XZH8Sn9\nCKqA0aywasuV8qUz5/N5ce0e+sc3plebBnbH8SktJCpgTI0dSWThxaPfIwsNU2NH2pRIqQte27CX\n4+fy+OM1Nas3AnqwXQWQogPqszKWcchh9USmxo7UA+3Kdiez83h5fQaDLm1Kt9b17I7jc1pIVEAZ\nljSDYUkzSElJISkpye44SgHw6voMTuXk18jeCAT5ri2da0sp5W3Hz+Yy57NMrktoRqcWde2OY4ug\nLiQ615ZSytteXJfB2dx87h9UM3sjEOSFRCmlvOnn0+eZtzGTG7u1IL5ptN1xbKOFRCmlqumFtXvI\nLShkag3ujYAWEqWUqpZDJ3N4/Yt9jOrRknaNatsdx1ZaSJRSqhr+tWY3BYWG+wbG2R3FdlpIlFKq\nirKOn+PNr35gdJ/WtG5Qy+44ttNCopRSVfT8p7sRhCkDOtgdxS9oIVFKqSrYd/QsizdnMeayS2hR\nL8ruOH5BC4lSSlXBrE++J9Qh3J3U3u4ofkMLiVJKuWj3T2d4J20/Yy9vQ5O6kXbH8RtaSJRSykWz\nPvmeyLAQJvfX3khJWkiUUsoFuw6dZsW2A9x5ZVsa1omwO45f0UKilFIueGZ1OnXCQ5mYGGt3FL+j\nhUQppSqxY/9JPvzmEOMT21GvVrjdcfyOFhKllKrEzNXpxESFcddV7eyO4pe0kCilVAVSfzjOp9/9\nxKSrY6kbGWZ3HL8U1IVEL2yllHLXM6vTaVg7nHFXtLU7it8K6kKiF7ZSSrlj095jrP/+CJP7t6d2\nhF6ZvDxBXUiUUqq6jDH886NdNI6O4LZ+beyO49e0kCilVBk27jnKl3uPMWVAB6LCQ+yO49e0kCil\nVClFvZEWMZHc0re13XH8nu70qwHeSdvPU6t2sf9ENi2/+JQHh3RkRI+WdsdSyu+U/KzAOW7u3YqI\nUO2NVEZ7JEHunbT9PLJ0u/ODAftPZPPI0u28k7bf5mRK+ZfSnxWA97Ye0M+KC7SQBLmnVu0iO6/g\nomXZeQU8tWqXTYmU8k9lfVZy8gr1s+ICLSRB7kCJ/65cWa5UTaWflerTQhLkyruCm17ZTamL6Wel\n+rSQBLkHh3QkKuzig4VRYSE8OKSjTYmU8k/6Wak+PWsryBWdnVV81la9KD1rS6ky6Gel+rSQ1AAj\nerRkRI+WpKSkkJSUZHccpfyWflaqR3dtKaWUcosWEqWUUm7RQqKUUsotWkiUUkq5RQuJUkopt4gx\nxu4MXiciJ4HvSyyKAU6W87jofslljYAj1dx86W1V5TVlLXcle3n33WlHRTldeY0/tcWddpT1XGVt\nC5bfr9KPS7dFf78qzujqa/ypLXHGmMqvDGiMCfob8JKrj4vul1r2tae2XZXXlLXclewVtKna7Qim\ntrjTDld+n1zMH3C/X5W1RX+/vPP75e9tMcbUmF1b71Xh8XvlvMZT267Ka8pa7kr2iu67I1ja4k47\nynqusrYFy+9X6ceB3JZA+v0qa5k/taVm7Npyl4h8bYzpbXcOdwVLO0Db4o+CpR2gbamqmtIjcddL\ndgfwkGBpB2hb/FGwtAO0LVWiPRKllFJu0R6JUkopt2ghUUop5RYtJEoppdyihaSKRCRWRF4VkSV2\nZ3GXiIwQkZdFZLmIDLY7jztE5FIReUFElojIH+zO4w4RqS0im0XkeruzuENEkkRkvfPnkmR3HneI\niENEnhCR50TkDrvzuENEEp0/k1dEZKMn1qmFBBCROSLyk4jsKLX8WhHZJSK7RWQ6gDEmwxgz3p6k\nlatiW94xxkwExgG/tSFuharYlp3GmMnAaMCvTtusSjucHgbe9m1K11SxLQY4A0QCWb7OWpkqtuVG\noCWQR4C3xRiz3vlZWQHM80gAd0aiBssNuBroCewosSwE2APEAuHAVqBTieeX2J3bg235J9DT7uzu\ntgW4AdgIjLE7e3XbAQwCbsEq7tfbnd3NtjiczzcFFtqd3c22TAd+73yN3332q/m5fxuo64nta48E\nMMasA46VWtwX2G2sHkgu8CbWfyV+rSptEcv/AR8YY1J9nbUyVf25GGPeNcZcAdzq26QVq2I7BgD9\ngDHARBHxq89oVdpijCl0Pn8ciPBhTJdU8eeShdUOgALfpXRNVT8rInIJcNIYc8oT29dL7ZavJfBj\nicdZwGUi0hB4AughIo8YY560JV3VlNkW4F6s/4BjRKSDMeYFO8JVUXk/lyRgFNYfrPdtyFVVZbbD\nGDMFQETGAUdK/DH2Z+X9TEYBQ4B6wPN2BKuG8j4rs4DnRCQRWGdHsGoory0A44HXPLUhLSTlkzKW\nGWPMUWCyr8O4qby2zAZm+zqMm8prSwqQ4tsobimzHcV3jJnruyhuK+9nshRY6uswbiqvLeew/vgG\nknJ/x4wxf/Hkhvyq2+xnsoDWJR63Ag7YlMVd2hb/EyztAG2Lv/JZW7SQlO8rIE5E2olIONYB0Hdt\nzlRd2hb/EyztAG2Lv/JdW+w+28AfbsAi4CAXTu0b71x+HZCOdebDo3bn1LYEZluCpR3aFv+92d0W\nnbRRKaWUW3TXllJKKbdoIVFKKeUWLSRKKaXcooVEKaWUW7SQKKWUcosWEqWUUm7RQqKUC0TkjBfW\n2VZExlTxPe+LSD1PZ1HKHVpIlLJPW6xZfl1mjLnOGHPCO3GUqh4tJEpVgfOqfynOKzF+JyILRUSc\nz2WKyP+JyCbnrYNz+VwRuanEOop6N/8AEkVki4g8UGo7zUVknfO5Hc5ZZ4u20UhEJjuf2yIie0Vk\njfP5wSLyuYikishiEanji++Lqtm0kChVdT2A+7EueBQLXFniuVPGmL5Y06Y/W8l6pgPrjTHdjTHP\nlHpuDLDKGNMd6AZsKfmkMeYF53N9sKbEmCkijYA/A4OMMT2Br4E/VqeBSlWFTiOvVNVtMsZkAYjI\nFqxdVBuczy0q8bV0caiKr4A5IhIGvGOM2VLO62YBnxpj3hPrGu+dgM+cnaRw4HM3MijlEi0kSlXd\n+RL3C7j4c2TKuJ+Ps/fv3A0WXtkGjDHrRORqYBiwQESeMsbML/ka58Wv2gBTihYBq40xv3O9KUq5\nT3dtKeVZvy3xtag3kAn0ct6/EQhz3j8NRJe1EhFpA/xkjHkZeBXretwln+8F/Am4zVy4iuIXwJUl\njs3UEpF4dxukVGW0R6KUZ0WIyJdY/6QV9QxeBpaLyCbgE+Csc/k2IF9EtgJzSx0nSQIeFJE84Aww\nttR2pgANgDXO3VhfG2MmOHspi0Sk6Brpf8aaRlwpr9Fp5JXyEBHJBHobY47YnUUpX9JdW0oppdyi\nPRKllFJu0R6JUkopt2ghUUop5RYtJEoppdyihUQppZRbtJAopZRyixYSpZRSbvn/BJQNw8aHxi8A\nAAAASUVORK5CYII=\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.loglog(sizes, enpy, 'o-', label='Numpy')\n", "plt.loglog(sizes, epy, 'o-', label='Python')\n", "plt.loglog(sizes, enba, 'o', label='Numba')\n", "plt.legend()\n", "plt.grid()\n", "plt.xlabel(\"Input size\")\n", "plt.ylabel(\"Relative error\");" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Fortunately, we see that Numba produces identical results to Pure Python (good news: its optimizations don't change the results). Furthermore, while we see a small loss of accuracy relative to Numpy's [high accuracy summation algorithm](https://github.com/numpy/numpy/pull/3685), the error is $\\cal{O}(10^{-14})$, which is still close to floating-point precision." ] } ], "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.6.3" } }, "nbformat": 4, "nbformat_minor": 2 }