{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Introduction to Numpy\n", "[Numpy](http://numpy.scipy.org/) is the fundamental library for scientific computing in Python. It contains list like objects that work like arrays, matrices, and data tables. This is how scientists typically expect data to behave. Numpy also provides linear algebra, Fourier transforms, random number generation, and tools for integrating C/C++ and Fortran code.\n", "\n", "[Matplotlib](http://matplotlib.org/) is the reigning library for 2D (with budding support for 3D) scientific plotting in Python. It produces publication quality figures in a variety of hardcopy formats and interactive environments across platforms. Between Numpy and Matplotlib, much of MATLAB's functionality can be replaced with Python.\n", "\n", "If you primarily want to work with tables of data, [Pandas](pandas), which depends on Numpy, is probably the module that you want to start with." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Why Numpy?\n", "\n", "Python was never designed originally for scientific computing, and contains many high-level abstractions necessary to enable its enormously flexible object-oriented interface. In Python, storing most integers requires more than just 4-8 bytes. It also requires at least a couple pointers per-integer. Performing a calculation on two numbers requires one or two bytecode operations, each of which can take dozens of CPU instructions for each pass through the Python eval loop. And when it comes to looping and index operations of Python lists the situation is even more dire.\n", "\n", "### A basic example\n", "\n", "`Z = A + B * C`\n", "\n", "In pure Python:" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "# Create 3 lists of a million ints\n", "A = range(1000000)\n", "B = range(1000000)\n", "C = range(1000000)" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1.36 s ± 500 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)\n" ] } ], "source": [ "%%timeit\n", "# Time doing the operation with a for loop\n", "Z = []\n", "for idx in range(len(A)):\n", " Z.append(A[idx] + B[idx] * C[idx])\n", " \n", "# print(Z) # DON'T DO THIS! It will print all 10000 array items" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "# Create 3 Numpy arrays with a million ints\n", "\n", "import numpy as np\n", "\n", "A = np.arange(1000000)\n", "B = np.arange(1000000)\n", "C = np.arange(1000000)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Using Numpy:" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "18.9 ms ± 475 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)\n" ] } ], "source": [ "%%timeit\n", "# Time the operation with Numpy\n", "Z = A + B * C" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[ 0 2 6 ... 999995000006 999997000002\n", " 999999000000]\n" ] } ], "source": [ "# Print the result\n", "Z = A + B * C\n", "print(Z)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In addition to just *looking* simpler, the Numpy version is significantly faster. The for loop disappears completely and is replaced by vectorized array operations. When printing very large arrays the output is also truncated by default, and for machine integers much less memory is used." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Memory Usage\n", "3 x 1000000 lists of Python ints: ~96 MB\n", "\n", "3 x 1000000 Numpy arrays of 64-bit ints: ~32 MB" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### More human-friendly interfaces to numerical libraries\n", "\n", "" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Numpy and Matplotlib" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAD8CAYAAAB0IB+mAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nO3dd3xU55U38N/RqPeOhApCBYRERzRjwKa4YnCysWP7jTfJGwdn146deMsnebPxZpN3P9l9dzebTZxk7diJneK+NrhgG2PANIORQPQm0STUJdTrzJz3j7m6c0eoT7kz957v56MPM9Jo7iPOzJnnPvd5zkPMDCGEEMYXpHcDhBBC+IYkfCGEMAlJ+EIIYRKS8IUQwiQk4QshhElIwhdCCJPwSMInot8RUQMRnRzh57cQURsRlStfT3viuEIIIcYv2EPP8yKAZwD8YZTH7GXmDR46nhBCiAnySA+fmfcAaPHEcwkhhPAOT/Xwx2M5ER0DUAPgb5n51HAPIqLNADYDQFRU1KLCwkIfNlEIIQJbWVlZEzOnDPczXyX8IwCmMXMnEd0FYAuAguEeyMzPAXgOAEpKSri0tNRHTRRCiMBHRFdG+plPZukwczszdyq3twEIIaJkXxxbCCGEg08SPhGlEREpt5cox232xbGFEEI4eGRIh4heAXALgGQiqgbwjwBCAICZ/xvAlwD8FRFZAfQAeIClTKcQQviURxI+Mz84xs+fgWPaphBCCJ3ISlshhDAJSfhCCGESkvCFEMIkJOELIYRJSMIXQgiTkIQvhBAmIQlfCCFMQhK+EEKYhCR8IYQwCUn4QghhEpLwhRDCJCThCyGESUjCF0IIk5CEL4QQJiEJXwghTEISvhBCmIQkfCGEMAlJ+EIIYRKS8IUQwiQk4QshhElIwhdCCJOQhC+EECYhCV8IIUxCEr4QQpiEJHwhhDAJSfhCCGESkvCFEMIkPJLwieh3RNRARCdH+DkR0S+IqIKIjhPRQk8cVwghxPh5qof/IoA7Rvn5nQAKlK/NAH7joePeoN9qBzN76+mFECJgBXviSZh5DxHljPKQTQD+wI5MfJCI4okonZlrPXF8rX969xS2HL2GjIQIFEyJwT1z07Fu1hQEW2T0yggaOnrx6udV2F/RhNq2XoRYCMVT43DXnDSsL0qDJYj0bqLwgIb2XrzyeRUOVDrjPDsjDnfNScf6WVMQJHGeFI8k/HHIAFCluV+tfM/jCb+mtQdd/Tacr+/E+fpOvH+8FoVpMfjpF+dgQXaCpw8nfMRqs+PZPRfxX59cQL/V7vKzysYuvHOsBjOnxOA/7p+H2RlxOrVSuMtqs+M3uyvxy10Vw8Z5a3kNCtNi8O/3SZwnw1fd3uE+jocddyGizURUSkSljY2NEz5Qc1f/Dd87W9eB+5/9DK+XVg3zG8LfdfZZ8bXfH8a/fXTuhiSgda6+A1/49X68daTah60TntLRO4C//N3n+I+Pz48a57N1jjhvOXrNh60zBl/18KsBZGnuZwKoGe6BzPwcgOcAoKSkZMKD8VsfW4HW7gFcbenGR6fq8Pv9l9EzYMOAjfH3bx6Hzc54cEn2ZP4GoYOuPiu+8vwhlFe1qt8rnhqLR1ZOx9zMeHT2WvHhqTq8qInzU68fQ++AHQ8tlTgHik4lzseq29Tvzc6IxSM352JOZhw6eq348GQdXjrgjPN3XitHn9WGLy+WOI8XeeoCpzKG/x4zzx7mZ3cDeBzAXQCWAvgFMy8Z6zlLSkq4tLTUrXZVtXRj8x/LcKa2XWkL8NuHS7CuaIpbzyu8z2qzY/Mfy7DzbIP6vW+vyceTawtuuCZzuakL3/pTGc7WdQAAggh4/qslWFMocfZ3Vpsdj/yhFLvPOc/on1xbgCfWFtxwTeZSUxe+9ccynKt3xNkSRHj+qyW4dWaqT9vsz4iojJlLhvuZp6ZlvgLgMwAziaiaiL5BRN8iom8pD9kG4CKACgC/BfDXnjjueGQlRuK1R5dhbqZjvI8Z+O5r5bjU1OWrJohJ+tWuSpdk/6N7ivA3t80c9gJ8TnIUXtu8HHOUcV07A0++Uo6qlm6ftVdMzi92Vrgk+59sKsZ3188Y9gL89OQovPboMhRPjQUA2OyMJ145iurrEufx8EjCZ+YHmTmdmUOYOZOZX2Dm/2bm/1Z+zsz8GDPnMfMcZnav2z5BseEhePHrS5ARHwEA6Oiz4ruvlcNqG3mcUOirvKoVv9h5Qb3/rdV5+NqK6aP+TlxkCF74WolLnL/9ylGJsx87cvU6frWrQr3/2K15eHh5zqi/Ex8Zit9/bTGmxoUDADp6rXjy1XLY7DIdeyymmauYGBWKZx9ehBCLo9dQXtWK3+69pHOrxHB6B2z47mvON/DinAT83e0zx/W7qTHheOahBQgOcsb5jweveK2tYvK6+60ucV4yPRFPrR9nnGPD8cuHFqpnAWVXruPPhyTOYzFNwgeA2Rlx+M66Ger9X3xyAXVtvTq2SAznuT0X1SG36LBg/Oz++ROaX78gOwFPrC1Q7//H9vNoaJc4+5tnP72IK82OoZiYsGD87P55E4rzomkJ+PaafPX+v310Dg0dEufRmCrhA8Cjq3JRmBYDAOgZsOH/fXhW5xYJrbq2Xvxmd6V6//t3FSIrMXLCz/Po6lzkpkQBcMwA+b/vn/FYG4X7alp78OweZ5x/cPcsZCZMPM5/dUsepic74tzRa8VPt8n7eTSmS/jBliA8fU+Rev+to9dw9Op1HVsktP71w7PoGbABAGalx+KBSU65Cwu24CebnBPG3jlWgxOaKX9CX//64Vn0DjiurRRPjcV9JVlj/MbwwoIt+PGmYvX+20ev4VSNxHkkpkv4AHBTXjJuL3ZO1/v37ed0bI0YdK6uA29rFtM8vaHIrVIJK/KTcefsNPX+zz6WOPuDM7Xt2FruXIbjbpxXFqS4vJ//8+PzbrXPyEyZ8AHg+3fOUl9k+yuaUXZFevl6e0YzW2PdrFQsz0ty+zm/u34GSMklu841Spz9wDM7nXFeXzQFS3Pdj/NT62eqcd5xpsFloZ5wMm3Cz0mOwqb5U9X7v9RMARS+V9nYifeOO3t9T66dMcqjx2/GlBhsmueM8893SO9PTxUNHdh20llC60nNxXV3zEyLwYa5zjhLL394pk34APDYrflqr2D3uUYcr5ZegV5+tasCg4u+b5mZgjmZniuM9eS6GRgcMdh7oUlddS1875mdzjivKUz1aAG076wrUOP86flGnFNWXQsnUyf8vJRol17B8zIvXxf17b14RzOmq51q5wnTk6Nw5+x09f4L+yTOeqht68G7x529e0/HOS8lGrcVOa/ZvLDvokef3whMnfAB4Furc9Xb207Uol7ma/vcnw9dhVVZfFMyLQGLpiV6/BiPrHSu0t1afk3m5evgzwevuiyy8ka58m+ucsZ5y9EaNHb0efwYgcz0Cb94ahyW5DgSjNXO+POhqzq3yFz6rDa8rFkh+fUxyidM1oLsBCya5kgwAzbGHz6TVZm+1Dtgw8ufO99b/3tFjleOszA7AfOz4gEA/Ta7rLIewvQJHwC+elOOevvlQ1fQZ7Xp1xiT2XaiFk2djj0M0uPCcVux96pbPnKz88PktdIqDEiNHZ9573gtWpS9KjLiI7BulnfiTEQuZ3Ovfn5VailpSMIHcFvxFKQrhZiaOvvx4ck6nVtkHn866Oz1fWXZNIR4cSvK9UVTkBoTBgBo7OjDJ2caxvgN4Sl/0vS0v7Jsmle3HL29OA3J0Y44N3T0uVRcNTtJ+ABCLEF4SLMpyhulsmOSL1Q0dKrz4oODCF9ePLnVluMVbAlyOcYrn8vwnS9cqO9Q58WHDomBN4RYgnB/SaZ6X+LsJAlf8aWSTHWK5v7KJqmj7gNvljk/WNfOSlV7Zd50f0mWGuc9Fxolzj7whibO64pSkRgV6vVjakty7D7fiGutPV4/ZiCQhK9Ij4vAqoIUAI5NUrTJSHie1WZ32Xv2vkXe7fUNykqMdInzGxJnrxqw2fHWEWe5DF/FOTspEisLkgEocZb9rAFIwnehPdV8s6xaNlTwor0XmtCgTJlLjg7DLTNTfHZsbZy3HL0GT23zKW706blGNHU64pwaE6YmYV+4X1OQbWt5jcQZkvBdrJ2VioTIEADAtdYefH6pRecWGdcbZc4e1xcXZnj1It5QawpTERMeDAC42tKNo1J3xWtc45zp0zivL5qC6DBHnC81dblskG5WkvA1woItLitv3zlWM8qjxWRd7+rHjtPOmRP3Lcoc5dGeFx5iwV2albdbNBU6hee0dPW7zIS6r8T3cb5DUy1V4iwJ/wYbNQXVPjhZi36rzOH1tA9P1aFfmRs9LzMOBVNifN6GTQuccX7veK3MyfeCD07WqiuoF2THIy8l2udtuHd+hnr7veM1pp+TLwl/iEXZCeom2K3dA9h7oVHnFhmPtirmPZpKlr60bHqSuvaipatf4uwF7x1z1s3ZqFOcl+clqWsvmjr7sa+iSZd2+AtJ+EMEBRE2zHOe7suwjmc1dfbhs8pm9f7dc9NHebT3BAWRSxLaclTi7EkNHb04dMkRZyLgrjn6xNkyJM7ajVfMSBL+MLQvkI9P16OnX0oteMoHJ+swOPlpcU4C0uMidGvLvQucp/vbT9dJnD3ogxPOOC/JScSU2HDd2qKN80en6tA7YN44S8IfRlF6LPKUDbC7+23YcaZe5xYZx3uaMybtBXI9zEqPRX6qY1y5d8COT8/LEnxPeV9TBnmDTsM5g4qnxiI32fl+3nPevMN3kvCHQUTYOM/ZK9h2onaUR4vxqm/vxeeXHVNdgwi4c07aGL/hfdo9bz+QGkoeUdfWi8NXNHGerW+cichlto6Za2VJwh/B3XOdL5Dd5xpNfRroKdtO1Kq7HS2dnoTUGP1O8wdpE8HOMw1SKdUD3tfE+aa8ZJ+UzBiLdgOcj8/Um3b2nST8EeSnxiBXGdbpGTD3aaCnaHtWel2sHaooPRbZiZEAgI4+K/abfBaHJ3zkh3GenRGLzATH9aKOXisOVJozzpLwR3F7sbP399EpGcd3x/WufhxWhnOI4NW69xNBRK7DOifMe7rvCS1d/Si94ozz+iL/ifMdxTKs45GET0R3ENE5Iqogou8N8/OvEVEjEZUrX4944rjepn2BfHK23vSLNtyx82yDOmtjfla8XwznDNIO63x8pl4WYbnhkzP1apwXZSf4xXDOIO01o+2nzfl+djvhE5EFwK8A3AmgCMCDRFQ0zENfY+b5ytfz7h7XF+ZmxqmLc1q7B6S2jhs+Pu08Q/KXXt+geZnxLnE+dFHiPFn+HOcFWQmYEuv4AGrp6lcnEJiJJ3r4SwBUMPNFZu4H8CqATR54Xt0REW7TvGg/PGXO00B39Q7YsEezkvU2P0sEQUHkMnwn03Anp3fAhr0XnGPj/pbwb4jzafNNw/VEws8AoC02Xa18b6i/IKLjRPQmEY1YFJuINhNRKRGVNjbqf6H0ds3p/vZT9VJidRIOVDahW1nUND05SpeaKmPRJqdPzkqcJ2PfhSb0KLPZ8lKikOuHcdbupWvGOHsi4dMw3xv6v/gugBxmngtgB4CXRnoyZn6OmUuYuSQlxXc10keyJCcR8UrJ5Lr2Xpyqade5RYFn6Gk+0XAvGX0tzklEjFJKt6qlBxcaOnVuUeBxjbP+ayyGszQ3EVGhFgDAleZuVDZ26dwi3/JEwq8GoO2xZwJwKVjBzM3M3Kfc/S2ARR44rk8EW4JwywznB49siDwxdjtjh6ZErr+d5g8KDQ7CKk2cZYPzibHZGZ+c9d/x+0FhwRasLNC+n801fOeJhH8YQAERTSeiUAAPAHhH+wAi0k7G3QjgjAeO6zO3FqaqtyXhT0x5dSsalZ2tkqJCsTA7QecWjWyNJs6fyDj+hJRXXUdTZz8Axw5mC7LidW7RyNbMcsZ5h8k+2N1O+MxsBfA4gI/gSOSvM/MpIvoxEW1UHvYEEZ0iomMAngDwNXeP60urZ6QgSBmFOFbdqm7ZJsamPc1fU5gKS5D/DecMumVmirrB+ZGr19HS1a9vgwLIdk2c181KRZAfx/nWmalqnMuuXEdrt3ni7JF5+My8jZlnMHMeM/+z8r2nmfkd5fb3mbmYmecx863MfNYTx/WV+MhQLJrm6JkyO0otiPHR/l+t89PT/EFJ0WHqGYidgd3nzNX7c8en2jjP8u84p8SEYV6m4wzEZmd8aqJV9LLSdpzWFDpfxLtkWGdcGjp6cabWcZE7OIiwIt93G1hPlsuwjsR5XBrae3G2rgMAEGIh3JSfpHOLxrbWZfjOPHGWhD9O2kSw53yjrMYcB21dmoXTEtQNpf3ZWs347p5zjaYtsjUR2rn3JdMSERnq/3HWjuPvPtdgmvezJPxxmjElWt36sKPPqtaFESPbc96ZCFYV+H/vHgBmTomROE+QdlHdyhmBEeei9Fh1dXV7rxWll6/r3CLfkIQ/TkSEWwud07lkWGd0dju79Py0Ux79GRG59PJlHH90djtjnzbOBYETZ+1Z+26TbH4jCX8C1mrG8WV65ujO1nWos5kSIkNQPDVO5xaN3y0znUlLe5YibnS6th3NymympKhQFKXH6tyi8btlpnaY1hxxloQ/AcvzkhAe4vgvq2zsQlVLt84t8l97Naf5K/KT/Xo65lDLcpMQanHE+Vx9B2rbenRukf/SDufcXJDs19Mxh1qel4QQi6O9Z2rb0dDeq3OLvE8S/gSEh1iwLNc5A8FM07kmSpsIAuU0f1BkaDAWT3cuENtrkt7fZGj/b1YGWJyjw4LV6dYAsOeC8eMsCX+CtC9qbS9WOPX023D4kvMiWKBcyNNarbnmIB/sw+vut6qbnQCBc2Fea/UM57COGeIsCX+CVmuS14GKZtNM55qIQ5ea0a/8vxSkRiM9LkLnFk2c9iLzvoomU26WMZZDF1swYHPUSSxMi0FqrP9sajNeqzTv530XGmGzG7t6piT8CcpLicZUZTpXR58V5VWtOrfI/2hn5wTaaf6gmVNi1M0y2noGcKy6TecW+R9tj3hlAPbuAcf0zJQYR5yvdw/gxDVjx1kS/gQRkeuwjglOAydKu+H7qgAczgEccdZee5BN7G+kHdIM1A92x/vZ+Ro1epwl4U+C9nT/UxNc6JmI2jZnLflQSxCWTvf/ZfYjWT1TxvFHcq21R60lHxYchCXTE3Vu0eSZ6XqNJPxJWJGfpFbPPF7diutSVVGlHc5ZPD0BEcpmE4Ho5vxkifMItGe2S6YnIjwkcOO8ssBZJfXo1eto6x7Qt0FeJAl/EuIjQzFXqbbHDOyvlF7+ICOM3w+KjwzFPKWuu50dF2+Fw94AXF07ksSoUMzNcCwMtBv8/SwJf5K0wzpGH/cbL5udsS+A598PZ7XE+QY2O7t8+AXitNuhzBJnSfiTpJ2eued8k+k2Qx7OqZo2XFdOh5Ojw1CYFqNzi9y3asj4rsQZOHGtDW09jjinxoRh5hSJc6CQhD9J8zLjERPuKANb196LCtn0eshwTmAtsx/JvMx4xEU4NrFv6OhT676b2Z7zrrNz/HFT+oman+V8P9e29Rp2E3tJ+JMUbAnCijxnL9/oV/fH41MDTMccyhJEuLlAuzjHuOO746WdjmmUOAdbgnCzZoOevQaNsyR8N7iM4xv0BTJenX1WHLniLKcQCLtbjZe2ZMAek5fT6OgdwJGrzsWGRoqzGcqmSMJ3g3bBxqGLzegdsOnYGn0drGyGVVmWPis9FqkxgbfMfiQ3axLB55daTB3nA5XNavmB2RmxSI4O07lFnqN9Px806PtZEr4bshIjkZscBQDos9pNvTuSy2l+gC6zH0lGfATyUpxx/vySxBkI/Gm3Q2UlRmK68n7uHbCj7IrxdsGShO8m7bCOUcf9xmNPAO5uNRGucTbm6f54DL0wbzRGH76ThO8m7YverAm/qqUbl5ocy+zDQ4JcaowbxaoC+WC/0tyFK82OTX8iQiyGjPPKAmPvdiYJ301Lc5MQHOTcNaexo0/nFvmeNgEunZ4U0MvsR7I0N1HdHelsXQfqTbA70lDas7jleUkICzZenJflub6fGzqMFWdJ+G6KDgvGQk1PZ78Jl9+7jusa7zQfcOyCVTLNWSDMjL38vQYohzwWo7+fJeF7wMp8Y4/7jcZqs7sss19twPH7QWYexx+w2fFZZbN632gXbLW0r2GjbW8pCd8DVmp3R7pgrjILx6rb0NFrBQCkxYYjPzVa5xZ5z8ohC7DsBt8dSau8qhUdfY44T40LV2ctGZFLfXyDxdkjCZ+I7iCic0RUQUTfG+bnYUT0mvLzQ0SU44nj+os5GXEuy+/P1xtzWfZwhg7nGGGZ/UiK0mORFBUKAGju6sfp2nadW+Q7e11WURujnMJIiqfGISHS8X5u6jRWOQ23Ez4RWQD8CsCdAIoAPEhERUMe9g0A15k5H8B/AvhXd4/rTyxBhBX5zo0+zHS6v9fg0zG1goKG7I5kojjvMVDZ67E4ymkYc/jOEz38JQAqmPkiM/cDeBXApiGP2QTgJeX2mwDWksG6CCtNOG2vrWdA3dOXyFjL7Efiur2lOeLc2t2P49XaOAfuLmbjZdQPdk8k/AwAVZr71cr3hn0MM1sBtAEY9lVDRJuJqJSIShsbA+c/Wlt46dAlYy7LHuqzyiZ1mf2cjDgkKsMdRqZNBKVXWtCljGsb2f6KZgwOY8/NjEd8pLnifPjSdfT0G+P97ImEP1xPfehVjvE8xvFN5ueYuYSZS1JSAufU0QzLsofaY/BVl8NJjQ1X6/wP2BiHLjWP8RuBTzuksdokcU6Pi0CBMgGh32Y3TJw9kfCrAWRp7mcCqBnpMUQUDCAOgOEKkphp1S0zu9RFN8LuVuPlutuZ8ePsUk7B4NdptIxYNsUTCf8wgAIimk5EoQAeAPDOkMe8A+Cryu0vAdjJBpy7aIbyqoMuN3ej+noPACAq1IIF2cZbZj8S1w92Y8e5srEL11odcY4OC8Z8ZY9fM3AZxzfIfhduJ3xlTP5xAB8BOAPgdWY+RUQ/JqKNysNeAJBERBUAngJww9RNI1iWmwiLsiz7VE07mjqNW2ZBm+iW5yUhNNg8SzoW5yQiTPl7tQnRiLRxvikvCSEW88R56fQkhCp/74WGTtS2BX6cPRI9Zt7GzDOYOY+Z/1n53tPM/I5yu5eZ72PmfGZewswXPXFcfxMTHoKF2c4ekNGWZWtphzKMPh1zqPAQC5ZM15RZMEjvbzhmHc4BgIhQCxZPd565GmFYxzwf1z5yc77xxv2G6rfa8VmleeZlD2e1Acd3h+qz2lzKKRhtn4PxWOVSPTPwP9gl4XvYyhmuy+8NeKkCR69eR5cyTS0zIQI5SZE6t8j3tB9y+yqc01ONpOzKdfQo04uzEyMxLcm45RRGYrQ4S8L3sLkZcYgNDwYA1LX3oqLBeGUW9g5ZdWmwNXTjMmNKNFJjHNv7tfUMqAuTjMR1FbX5evcAUJgWo27j2No9gFM1bTq3yD2S8D0s2BKEm/Jciy8Zjcu8bJMmAiIy/OpqI29nOF5Dy2kEepwl4XuB67BO4I/7abV09eP4NUcvJ4iA5XnmTPiAa6/XaNMzmzr7cPKaozicJYiwPM/45RRGoo3zpwE+ji8J3wu0F3oOXmxBn9UYy7IBx8yjwcsS87Pi1SqhZqQtp3Hkais6egd0bI1naWeYLciKR2y4eeOsrRF15Mp1dAZwOQ1J+F6QlRiJacqFzJ4Bm6HKLGh7smabjjlUUnQYZmfEAgBsdsaBSmMsvwfMPe12qNSYcMxKd8TZamccDOA4S8L3kqGbZRjBDcvsTTquq2XE1dWOOBt/O8OJWGWQ1dWS8L3EiPPxHasNHZs6x4YHY15mnM4t0t8qA164PVvXgYYOxyrxuIgQzM00TzmFkbjUTwrgOEvC95LleUlqmYWTNW1o6erXuUXu+/Scs2dzc0Eygk20zH4kC6fFIzLUAgC40tyNK81dOrfIfdoLkzcXJKuvYzNbNC0B4SGO1/ulpi5UtXTr3KLJkXesl8RFhKiFppiNUWZhj8t0TBnOAYCwYAuW5TpnsARy72+QdkWpxNkhPMSCpdO1u9oFZpwl4XuRdhZHII/7AUB3vxWHLjorWpv9Qp6Wy/hugE/b6+qz4vBlTZzlOo3KtSx2YMZZEr4XrTJQmYVDF1vQb7MDcKwyTY+L0LlF/kNbVOyzymYMKP9PgejgxWYM2Byv08K0GKTFhevcIv+h/WDfX9kEawDGWRK+F83LjEdMmKPMQk1bLyobA3d891OTbnYyHrnJUciId3wAdvRZcawqcMssfCrDOSPKT41GWqzjA7Cj14pj1YFXZkESvhcFW4JcVigG8rCOy7juTEkEWo4yC8bYLMNlFzNJ+C6GxjkQ38+S8L1Me7ofqPPxq1q6cbHJcXYSHhKExTmJY/yG+Rhh2t6V5i5cbnbMPokIsaAkxzy7mI1XoI/jS8L3Mu2432cXm9FvDbxxP+1p/rLcJISHWHRsjX+6KS8Jg7MXj1e3orU78KbhahPY8rwkhAVLnIdakZ+MweKw5VWtaOsJrHIakvC9bFpSFLISHeO73f02HL0aeGUWzLpZ+UTER4aqC5TsjIAssyDj92NLjArFnAzHgkM7w2UjoEAgCd8HtMvvA63a3oDN7pK8ZPx+ZKsCeBzfsYuZZncrSfgj0o7jf3peEr4YQttb2nm2QceWTJy2OmBGfARyk82369F4rRqy7WEgTcMtu+LcxSw7MdKUu5iN19BtDwMpzpLwfeDm/GSEKmUIztZ1oKa1R+cWjd/uIbNzzLi71XjNy3JOw73W2qNe6A4Eu887OyKrZ0icR7MgOwFRSjmNa6096oXuQCAJ3weiwoKxNNc5syWQevk7z7gmAjGykKHTcANoWEfiPH6hwYE73VoSvo+sKUxVb+8KkIRf1dKNc/UdABwvcimTO7aVATg982pzNy4oey+HBQe5bPghhreyIDCnZ0rC9xFtwt9f2YTeAf/fBUt7JnJTXhIiQ4N1bE1gWF3gWmYhEHY7++RsvXp7RX4yIkJlOuZYVg0ppxEo060l4fvItKQo5KY4Lnj2DrjOiPBXn2gS/lrNB5YYWXbSkN3OLvv/NNxPNJHlUUYAABMnSURBVMM5a2dJnMcjJykSmQmO6dZd/YGzq50kfB/SJk1/H8fv6rO6bOW2ZtYUHVsTWLRj4DvO+HecO3oHcOiSJs7ywT4uROQS50/O1I/yaP8hCd+Hbh2S8P15OtfeC01qdczCtBi1OJgY2zrNh+P203V+Hed9F5rU6pjFU2OlCuoErCtyxvnjM/V+HedBbiV8Ikokoo+J6ILy77DFN4jIRkTlytc77hwzkC3OSXSZtne+vlPnFo1sp2ZcV07zJ2ZZbpIa5+rrPThb16Fzi0amPQORYbuJuSkvSZ2eeUVz4dufudvD/x6AT5i5AMAnyv3h9DDzfOVro5vHDFghliCXiz3ai2X+xG5n7DzrnHmwVoZzJiQ0OMhlRfLHp/0zzjY7Y/c5Z8KXYbuJCQu2BESctdxN+JsAvKTcfgnAvW4+n+Fpx0i3n/LPF8jxa21o6nRsYp0UFYp5son1hK3Xnu77aSIor2pFs7LXcnJ0GOZmyKb0E6WN83Y/jbOWuwl/CjPXAoDy70jnhOFEVEpEB4lo1A8FItqsPLa0sTFw5reO19pZqeqm0OVVrahr69W5RTfSXoC6tTBVNrGehFtmpiJY+X87ca0NtW3+t7paG+c1hSkIkjhP2K0zne+PY1WtqG/3v/ez1pgJn4h2ENHJYb42TeA42cxcAuAhAD8noryRHsjMzzFzCTOXpKQYb8VffGQolmlW3W4/Xadja4b34Ulnm9bJ+P2kxEWEuGxuvsPPen/MjA9POeMsw3aTEx8ZiiWa/SF2+PlsnTETPjOvY+bZw3xtBVBPROkAoPw77Bw0Zq5R/r0IYDeABR77CwLQ7cVp6u2PTvlXwq9o6FAvPoWHBEnVRDf48+n+hYZOXFS23IwIsUjZazcEwvDdIHeHdN4B8FXl9lcBbB36ACJKIKIw5XYygBUATrt53IB2W5Ez4R+82ILrXf6zWcYHJ5wfQLfMSJXVtW7QTts7eLEZ7b3+s1mGNs63FqbI6lo3aBP+gYpmtbqsP3I34f8LgPVEdAHAeuU+iKiEiJ5XHjMLQCkRHQOwC8C/MLOpE35aXDjmZzkuhNrs7LKiVW8faIZz7pyTNsojxVgy4iNQPDUWADBgY7+qofTByVr19h2z03VsSeDLSoxEYVoMAKDfZneZ+eRv3Er4zNzMzGuZuUD5t0X5fikzP6LcPsDMc5h5nvLvC55oeKDzx2Gdq83dOF3bDgAItQTJqksPuEMT5/eP147ySN+53NSlrg0IDZY4e8Ids/0vzsORlbY6ub3YeRq453wjOvzgdF/b67u5IBkx4SE6tsYY7p7r7D3v9ps4OzsYqwqSER0mw3bu2qCJ886zDX47rCMJXye5KdHqaWCf1e4XF3u0iUDbYxGTl5sSjaJ0x7BOv5/E+UMZzvG4/NQYl/ezv9bWkYSvo43zp6q33zlWo2NLHLXvy6taAQCWIMJ6mabnMRvmOZPqezqf7l9t7sax6jYAQLDE2aO0vfx3j/nnsI4kfB3dM9eZ8PddaEKLjrN1tB84KwuSkRAVqltbjGbDHGec915oRFu3fsM6W8uvqbdXzUhBXKQM23nKBs37ec/5RrT16D98N5QkfB1lJUZiYbZjto7Vzth2Qp9eATPj7aPORHDv/Axd2mFU2UmRmJfpKFswYGPdLtIzM97WJPx7F0icPSknOQqzM5ThO5t/DN8NJQlfZxvn6T+sc6qmHRXKYqvIUAtuK5bTfE/TXrx997g+cT55rV1dbBUVapHhHC+4W3M2967Ow7TDkYSvs7vnTsVgCZPDl1t0qbmiPc2/rWiKLLbygg1zp4KUOO+raEJNq+/jrD2Lu312miy28gLtOP7eC41+VytLEr7OUmLCcFOeY9NoZmDLUd/2Cmx2xtZy5zHlNN87psZH4KY8R20dZuCtI9U+Pb7VZnc5g5RhO+/ISoxUa2XZGXjrqG/jPBZJ+H7gC5ok+3pplU93zjl4sRkNHY5SyMnRobg5P9lnxzab+xZlqbffLKv2aZwPVDarJa8dnYykMX5DTJZLnEt9G+exSML3A3fNSVd3SLrU1IVDl1p8duzXDleptzfMnYpgi7wkvOX24jQ1zpebu3HYhxucv1bqjPM9EmevunNOmrqY7WJTF45c9Z8NziXqfiAi1IJNC5wXe7RJ2JtauvpdSiHfX5I1yqOFuyJCLdiguUj/Rqlv4tzU2YftmplB9y/O9MlxzSoyNNhlLP+NUv8Z1pGE7yceWJyt3t52otYnc7XfOlKtblQ+LyseRUqhL+E995U4k+37J2rR5YMl+P9TVq1uVL4wOx6FaRJnb9PG+b3jteju949SC5Lw/cTsjDi1smKf1Y4tmpkz3sDMePnzq+r9BxdL794XFmTFIy8lCgDQ3W/zSZxf1ZwxPrAke5RHC09ZmJ2AXCXOnX1WvFPuH1M0JeH7Ee2b8eVDV716sefgxRaXOdn3aIYahPcQER5aOk29/9KBy16N84HKZlxqcsQ5Jsx1qEF4DxHhIc37+UUvx3m8JOH7kY3zpiIixDE3+lx9B/ZXNHvtWC/su6jevndBBqKkYqLP3FeSiUhlDvz5+k58Vum9OD+/1xnnLyzMkDUWPnRfSZb6fj5b14GDF303GWMkkvD9SFxECO7XjP09r0nKnnSxsRM7zjg3afj6iuleOY4YXmx4CL60yBnn3x+47JXjVDR0YNe5RgAAkcTZ1+IiQvDFhc4p1y95Kc4TIQnfz3x9xXR1Rebuc424UN/h8WP8bv8l9fbawlTkp0Z7/BhidH+5PEe9veNMPSobOz1+jBf2XVZvr5s1BdOTozx+DDG6r92Uo97efroOF70Q54mQhO9ncpKjXGqcPLfHs738ps4+vFnmnCb2jZXS69NDfmo0bp3p2DicGfjVrgqPPn9DR6/Lat5HbpY466FgSgxWzXDE2c7Ab3ZX6toeSfh+6JurctXbbx29pl5084Tf7K5E74BjKmbx1Fgsz5UVl3p5fE2BentreQ2uNnd77Ll/vasSfVZHnOdmxmHJ9ESPPbeYmG+vyVdvv330GqpaPBfniZKE74cW5yRiRb4jEdvsjP/acd4jz1vX1os/Hbyi3v/OuhmgwfEj4XOLpiW4xNlTvfya1h68fMg55fbJtQUSZx0tzklU6+tY7Yxf69jLl4Tvp55aP1O9vfVYDc57YCz/V7sq1F7fvMw4rJslm1fr7duaXv4bZVU4V+d+nJ/ZVaEuqJufFS+blPuBJzRxfr20yivX5sZDEr6fWjQtAbdoxnh//O5pt+bxnqvrcFlo9d310rv3B8tyk9SCdXYGfvKee3E+XdPuUprjKYmzX1iel6QWrLPZGT95/4wu8/Il4fuxv7t9plorf19FE96f5I5YzIx/2HICNrvjBbYsNxGrlQtJQn8/3FDkEmftlNmJsNsZP9x6Uo3zTXlJWFkg1U/9ARG5xHnP+UbsOje5OLtDEr4fK54a5zJ978fvnkZ778Rr7LxZVq1WZgwOIvxk02zp9fmRmWkxeGipc1XmD7ecRGv3xPc3fqOsCmVXHHEOsRB+vKlY4uxHZqXH4suamln/8PZJn+97Kwnfzz112wwkR4cBABo6+vDDLScndCpY1dKNH797Wr3/yMpcFEyJ8Xg7hXueWj8TicrG8XXtvfjBBON8pbkLP3nvjHr/mytzkZ8qcfY3f3PbDCQoG8fXtPXi6a0nfXp8Sfh+LjY8BD/aWKTe31pe4zIWP5qefhsee/kIOpSKjNmJkXhibf4YvyX0kBgVip9+cY56//3jtS5Fz0bT3W/F4y8fRacS52lJkXh8jcTZHyVHh7nEeWt5DV73UTl0QBJ+QNgwdyru0yzFf3rrKew8Wz/q7wzY7Hjy1aM4Xt0GwDGU84sHF0gtFT92e3EavqzZk+CHW06OOc7bb7XjiVfKceKaI84hFsIvJc5+7Y7Z6S6lNf7P2yew53yjT47tVsInovuI6BQR2YmoZJTH3UFE54iogoi+584xzeqfNhWr5ZNtdsajfyzDlqPDl9Zt7x3Ao38sw/bTzg+Fp+8pwvyseJ+0VUze0/cUoSjdEWernbH5D6Uum8xrtfcOYPMfS7HjjDPO/3hPMeZmSpz93Y82FqMwzTHkZrUzHnmpFO8e834JZXJnahARzQJgB/AsgL9l5tJhHmMBcB7AegDVAA4DeJCZTw997FAlJSVcWnrDU5pWfXsvvvjrA7jW2qN+7+456Xh0dS5mT41DV78V20/V42cfn3d5zDdXTscP7i4a7imFH6pr68Vf/MY1zhvmpuPRVXkonhqLzn4rPjpZh5/vuODymEdX5+L7d87So8liEmrbevAXvz6AmrZe9Xsb503F5lW5KEqPRVDQ5C64E1EZMw/bAXcr4WsOsBsjJ/zlAH7EzLcr978PAMz807GeVxL+jWpae/DwC4dQ2ehabsESROp0PK3Hbs3D3942U2ZrBJhrSpwvjjPO316TL3PuA1D19W48/MLnN5RPiQ0PxprCVPz8gQUTfs7REr4vxvAzAGivSlQr3xOTMDU+Am8/tgJfXOD6Xzg0CcRFhOCZhxbg724vlCQQgDLiI/D2X6/AvfNdN6YZGuf4yBD8+n8txN/Ih3pAykyIxJa/XoGNQzYgau+1or3X89sijnllh4h2AEgb5kc/YOat4zjGcK/CEU8riGgzgM0AkJ0t27ENJzY8BD/78nw8vHwaXjxwGfsrmtHU2YcgAgrTYnHXnDQ8vCwHccr0LxGY4iJC8PMHFuAvb8rBi/sv40ClM86z0mNx15x0fGXpNIlzgIuLDMEvHlyAr96UgxcPXMZnlU1o6uxHSU6Cx481ZsJn5nVuHqMagHbD1EwAI16dYObnADwHOIZ03Dy2oS3ITsCCbMeLot9qhyWIYJnkuJ/wXwuzE7AwOwHMjAEbS5wNatG0BCya5ojz5eZuRIVZPH4MX8zdOgyggIimA7gG4AEAD/nguKYSGiwzbI2OiBAaLIne6IjIa5vVuDst8wtEVA1gOYD3iegj5ftTiWgbADCzFcDjAD4CcAbA68x8yr1mCyGEmCi3evjM/DaAt4f5fg2AuzT3twHY5s6xhBBCuEfGAYQQwiQk4QshhElIwhdCCJOQhC+EECYhCV8IIUxCEr4QQpiEJHwhhDAJSfhCCGESkvCFEMIkJOELIYRJSMIXQgiTkIQvhBAmIQlfCCFMQhK+EEKYhCR8IYQwCUn4QghhEpLwhRDCJCThCyGESUjCF0IIk5CEL4QQJiEJXwghTEISvhBCmIQkfCGEMAlJ+EIIYRKS8IUQwiQk4QshhElIwhdCCJOQhC+EECbhVsInovuI6BQR2YmoZJTHXSaiE0RUTkSl7hxTCCHE5AS7+fsnAXwRwLPjeOytzNzk5vGEEEJMklsJn5nPAAAReaY1QgghvMbdHv54MYDtRMQAnmXm50Z6IBFtBrBZudtJROcmecxkAGY7o5C/2RzkbzaHyf7N00b6wZgJn4h2AEgb5kc/YOat42zACmauIaJUAB8T0Vlm3jPcA5UPgxE/EMaLiEqZecTrCkYkf7M5yN9sDt74m8dM+My8zt2DMHON8m8DEb0NYAmAYRO+EEII7/D6tEwiiiKimMHbAG6D42KvEEIIH3J3WuYXiKgawHIA7xPRR8r3pxLRNuVhUwDsI6JjAD4H8D4zf+jOccfJ7WGhACR/sznI32wOHv+biZk9/ZxCCCH8kKy0FUIIk5CEL4QQJmG4hE9EdxDROSKqIKLv6d0eXyCi3xFRAxGZ4mI4EWUR0S4iOqOU9nhS7zZ5GxGFE9HnRHRM+Zv/Se82+QoRWYjoKBG9p3dbfMGbpWgMNYZPRBYA5wGsB1AN4DCAB5n5tK4N8zIiWgWgE8AfmHm23u3xNiJKB5DOzEeUGWBlAO41cpzJsZw9ipk7iSgEwD4ATzLzQZ2b5nVE9BSAEgCxzLxB7/Z4GxFdBlDijVI0RuvhLwFQwcwXmbkfwKsANuncJq9TFrG16N0OX2HmWmY+otzuAHAGQIa+rfIuduhU7oYoX8bprY2AiDIB3A3geb3bYgRGS/gZAKo096th8ERgdkSUA2ABgEP6tsT7lKGNcgANAD5mZsP/zQB+DuDvAdj1bogPDZaiKVNKzXiM0RL+cFXcDN8LMisiigbwPwC+w8zterfH25jZxszzAWQCWEJEhh6+I6INABqYuUzvtvjYCmZeCOBOAI8pQ7YeYbSEXw0gS3M/E0CNTm0RXqSMY/8PgD8z81t6t8eXmLkVwG4Ad+jcFG9bAWCjMqb9KoA1RPQnfZvkfdpSNAAGS9F4hNES/mEABUQ0nYhCATwA4B2d2yQ8TLmA+QKAM8z8M73b4wtElEJE8crtCADrAJzVt1XexczfZ+ZMZs6B4728k5m/onOzvMrbpWgMlfCZ2QrgcQAfwXEh73VmPqVvq7yPiF4B8BmAmURUTUTf0LtNXrYCwMNw9PjKla+79G6Ul6UD2EVEx+Ho2HzMzKaYpmgyXi1FY6hpmUIIIUZmqB6+EEKIkUnCF0IIk5CEL4QQJiEJXwghTEISvhBCmIQkfCGEMAlJ+EIIYRL/H7EmypiyXEucAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "%matplotlib inline\n", "import numpy as np\n", "from matplotlib import pyplot as plt\n", "\n", "ax = plt.subplot(111)\n", "x = np.arange(0.0, 5.0, 0.01)\n", "y = np.cos(x * np.pi)\n", "plt.ylim(-1.5, 1.5)\n", "lines, = plt.plot(x, y, lw=3)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Numpy Array Basics\n", "#### Creating a Numpy array" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([1, 2, 3])" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Create an array from a list of ints and show the array\n", "import numpy\n", "\n", "vals = [1, 2, 3]\n", "arr = numpy.array(vals)\n", "arr" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[1 2 3]\n" ] } ], "source": [ "# Print the array--notice any difference?\n", "print(arr)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Unlike Python lists, NumPy arrays are homogeneous: all values must have exactly the same type. This allows values to be packed together as shown here, which saves memory and is much faster to process.\n", "\n", "\n", "\n", "If we give NumPy initial values of different types, it finds the most general type and stores all the values in the array using that type. For example, if we construct an array from an integer and a float, the array's values are both floats: " ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([1. , 2.3])" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Create an array from a heterogeneous list\n", "arr = numpy.array([1, 2.3])\n", "arr" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "If we want a specific type, we can pass an optional argument to array called dtype (for \"data type\"). For example, we can tell NumPy to create an array of 32-bit floats even though all the initial values are integers: " ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([1., 2., 3., 4.], dtype=float32)" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Create an array of floats from a list of ints\n", "numpy.array([1, 2, 3, 4], dtype=numpy.float32)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "NumPy provides many basic numerical data types, each of which is identified by a name like float32.\n", "The three called `int`, `float`, and `complex` are whatever the underlying hardware+OS platform uses as its native C types: this will usually be 32- or 64-bit.\n", "\n", "Note: Changing the dtype of an array is usually not going to yield anything useful. Instead use `arr.astype()` or similar:" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[1 2 3 4]\n", "[1. 2. 3. 4.]\n", "[1.+0.j 2.+0.j 3.+0.j 4.+0.j]\n" ] } ], "source": [ "# Create arrays of ints, floats, and complex from lists of ints\n", "print(numpy.array([1, 2, 3, 4], dtype=numpy.int))\n", "print(numpy.array([1, 2, 3, 4], dtype=numpy.float))\n", "print(numpy.array([1, 2, 3, 4], dtype=numpy.complex))" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[1 2 3 4]\n", "[1.e-45 3.e-45 4.e-45 6.e-45]\n", "[1. 2. 3. 4.]\n" ] } ], "source": [ "# Create an array of ints and show the dtype\n", "arr = numpy.array([1, 2, 3, 4], dtype=numpy.int32)\n", "print(arr)\n", "\n", "# Try reassigning the dtype to float\n", "arr.dtype = numpy.float32\n", "print(arr)\n", "\n", "# Restore the correct dtype and use .astype() instead\n", "arr.dtype = numpy.int32\n", "print(arr.astype(numpy.float32))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "There are many other ways to create arrays besides calling the basic `array()` constructor. For example, the `zeros()` function takes a tuple specifying array dimensions as an argument and returns an array of zeros of that size: " ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[0., 0., 0.],\n", " [0., 0., 0.]])" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Create a 2x3 array of zeros\n", "z = numpy.zeros((2, 3))\n", "z" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The array's dtype defaults `float` unless something else is specified with the `dtype=` argument. This is typical in most functions in Numpy that create/return arrays." ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[0, 0, 0],\n", " [0, 0, 0]], dtype=int32)" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Create a 2x3 array of integer zeros\n", "z = numpy.zeros((2, 3), dtype=numpy.int32)\n", "z" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[1. 1. 1.]\n", " [1. 1. 1.]]\n", "\n", "[[1. 0.]\n", " [0. 1.]]\n" ] } ], "source": [ "# The ones and identity functions work much the same way:\n", "print(numpy.ones((2, 3)))\n", "print()\n", "print(numpy.identity(2))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "It's also possible to create NumPy arrays without filling them with data using the empty function. This function does not initialize the values, so the array contains whatever bits were lying around in memory when it was called:" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[1., 0.],\n", " [0., 1.]])" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Create a 2x2 empty array\n", "arr = numpy.empty((2, 2))\n", "arr" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This might not seem particularly useful, but if a program is going to overwrite an array immediately after creating it, perhaps by filling it with the result of some computation, there's no point taking the time to fill it with zeroes or ones." ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[10 11 12 13 14 15 16 17 18 19]\n" ] } ], "source": [ "# Another frequently useful array creation function is arange:\n", "\n", "arr = np.arange(10, 20)\n", "print(arr)" ] }, { "cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[10. 10.1 10.2 10.3 10.4 10.5 10.6 10.7 10.8 10.9 11. 11.1 11.2 11.3\n", " 11.4 11.5 11.6 11.7 11.8 11.9 12. 12.1 12.2 12.3 12.4 12.5 12.6 12.7\n", " 12.8 12.9 13. 13.1 13.2 13.3 13.4 13.5 13.6 13.7 13.8 13.9 14. 14.1\n", " 14.2 14.3 14.4 14.5 14.6 14.7 14.8 14.9 15. 15.1 15.2 15.3 15.4 15.5\n", " 15.6 15.7 15.8 15.9 16. 16.1 16.2 16.3 16.4 16.5 16.6 16.7 16.8 16.9\n", " 17. 17.1 17.2 17.3 17.4 17.5 17.6 17.7 17.8 17.9 18. 18.1 18.2 18.3\n", " 18.4 18.5 18.6 18.7 18.8 18.9 19. 19.1 19.2 19.3 19.4 19.5 19.6 19.7\n", " 19.8 19.9]\n" ] } ], "source": [ "arr = np.arange(10, 20, 0.1)\n", "print(arr)" ] }, { "cell_type": "code", "execution_count": 19, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[ 0. 5.26315789 10.52631579 15.78947368 21.05263158\n", " 26.31578947 31.57894737 36.84210526 42.10526316 47.36842105\n", " 52.63157895 57.89473684 63.15789474 68.42105263 73.68421053\n", " 78.94736842 84.21052632 89.47368421 94.73684211 100. ]\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX0AAAEICAYAAACzliQjAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAR4klEQVR4nO3dfYxddZ3H8fcXasHism1hINjSDiaND2tUyETxIcal/iFoLH9oojtZG0PSbOKu+JAI2j+Mf5DoxviUNW4moFs2DatbiTTGuCEFY/YPuw5CeLDslkVbCpUOqTysjQHS7/5xzoRhuHemc++5c+89v/crae7cM+fe+5szZz7z6znnfiYyE0lSGc4a9gAkSavH0Jekghj6klQQQ1+SCmLoS1JBDH1JKsiyoR8R34+IExHx4IJlGyPizog4XN9uqJdHRHwnIh6JiPsj4opBDl6StDKx3HX6EfFe4P+AWzPzzfWyfwROZuZXI+JGYENm3hAR1wD/AFwDvAP4dma+Y7lBXHjhhTk5OdnfVyJJbXXyJBw5AqdPv7TsrLN48PTp5/6cef5KnmrNcitk5i8jYnLR4h3A++qP9wC/AG6ol9+a1W+SX0XE+oi4JDOPL/Uak5OTzM7OrmTcklSOycmXBz7A6dO8Fc5Z6VP1ekz/4vkgr28vqpdvAh5bsN6xetkrRMSuiJiNiNm5ubkehyFJBTh6tOPiV8HalT5V0ydyo8OyjsePMnMmM6cyc2piYqLhYUhSi2zZ0nHxC/D8Sp+q19B/MiIuAahvT9TLjwGXLlhvM/BEj68hSQK46SZYt+7ly9at4wl4fKVP1Wvo7wd21h/vBO5YsPwT9VU8VwLPLHc8X5K0jOlpmJmBrVshorqdmeEpOLnSpzqTq3duozppeyHwJPBl4CfAj4AtwFHgo5l5MiIC+CfgA8Ap4JOZuewZ2qmpqfREriStTETck5lTK3nMmVy98/Eun9reYd0EPrWSAUiSVo/vyJWkghj6klQQQ1+SCmLoS9Ig7d1bvaP2rLOq2717hzqcZU/kSpJ6tHcv7NoFp05V948cqe5DdRnmEDjTl6RB2b37pcCfd+pUtXxIDH1JGpQunTldl68CQ1+SBqVLZ07X5avA0JekQenSmcNNNw1nPBj6kjQ4XTpzhnUSF7x6R5IGa3p6qCG/mDN9SSqIoS9JBTH0Jakghr4kdTNiFQpN8ESuJHUyghUKTXCmL0mdjGCFQhMMfUnqZAQrFJpg6EtSJyNYodAEQ1+SOhnBCoUmGPqS1MkIVig0wat3JKmbEatQaIIzfUkqiKEvSQUx9CWpIIa+pHZqYYVCEzyRK6l9Wlqh0ARn+pLap6UVCk0w9CW1T0srFJpg6Etqn5ZWKDTB0JfUPi2tUGiCoS+pfVpaodAEr96R1E4trFBogjN9SSpIX6EfEZ+NiIci4sGIuC0izo2IyyLiYEQcjogfRsTapgYrSepPz6EfEZuATwNTmflm4GzgY8DXgG9m5jbgj8B1TQxUktS/fg/vrAFeHRFrgHXAceAqYF/9+T3AtX2+hqSSWJ8wUD2HfmY+DnwdOEoV9s8A9wBPZ+aL9WrHgE2dHh8RuyJiNiJm5+bmeh2GpDaZr084cgQyX6pPMPgb08/hnQ3ADuAy4LXAecDVHVbNTo/PzJnMnMrMqYmJiV6HIalNrE8YuH4O77wf+F1mzmXmC8DtwLuA9fXhHoDNwBN9jlFSKaxPGLh+Qv8ocGVErIuIALYDvwXuBj5Sr7MTuKO/IUoqhvUJA9fPMf2DVCdsfwM8UD/XDHAD8LmIeAS4ALilgXFKKoH1CQPX1ztyM/PLwJcXLX4UeHs/zyupUPPvoN29uzqks2VLFfi+s7Yx1jBIGi3WJwyUNQySVBBDX5IKYuhLUkEMfUkqiKEvqTn25ow8r96R1Iz53pz5GoX53hzwapwR4kxfUjPszRkLhr6kZtibMxYMfUnNsDdnLBj6kpphb85YMPQlNWN6GmZmYOtWiKhuZ2Y8iTtivHpHUnPszRl5zvQlqSCGviQVxNCXpIIY+pIqVigUwRO5kqxQKIgzfUlWKBTE0JdkhUJBDH1JVigUxNCXZIVCQQx9SVYoFMSrdyRVrFAogjN9SSqIoS9JBTH0Jakghr7UBlYo6Ax5Ilcad1YoaAWc6UvjzgoFrYChL407KxS0Aoa+NO6sUNAKGPrSuLNCQSvQV+hHxPqI2BcRD0fEoYh4Z0RsjIg7I+JwfbuhqcFK6sAKBa1AvzP9bwM/z8w3AG8FDgE3AgcycxtwoL4vaZCmp+H3v4fTp6tbA19d9Bz6EXE+8F7gFoDMfD4znwZ2AHvq1fYA1/Y7SElSM/qZ6b8OmAN+EBH3RsTNEXEecHFmHgeoby9qYJySpAb0E/prgCuA72Xm5cCfWMGhnIjYFRGzETE7NzfXxzAkSWeqn9A/BhzLzIP1/X1UvwSejIhLAOrbE50enJkzmTmVmVMTExN9DEMaY9YnaJX1HPqZ+QfgsYh4fb1oO/BbYD+ws162E7ijrxFKbTVfn3DkCGS+VJ9g8GuAIjN7f3DE24CbgbXAo8AnqX6R/AjYAhwFPpqZJ5d6nqmpqZydne15HNJYmpysgn6xrVurK3CkZUTEPZk5tZLH9FW4lpn3AZ1ecHs/zysVwfoEDYHvyJWGxfoEDYGhLw2L9QkaAkNfGhbrEzQE/hEVaZimpw15rSpn+pJUEENfkgpi6EtSQQx9qVdWKGgMeSJX6sV8hcL8HySfr1AAT8xqpDnTl3qxe/dLgT/v1KlquTTCDH2pF1YoaEwZ+lIvrFDQmDL0pV5YoaAxZehLvbBCQWPKq3ekXlmhoDHkTF+SCmLoS1JBDH1JKoihL0kFMfRVJntzVCiv3lF57M1RwZzpqzz25qhghr7KY2+OCmboqzz25qhghr7KY2+OCmboqzz25qhgXr2jMtmbo0I505ekghj6klQQQ1+SCmLoa/xYoSD1zBO5Gi9WKEh9caav8WKFgtQXQ1/jxQoFqS99h35EnB0R90bET+v7l0XEwYg4HBE/jIi1/Q9TqlmhIPWliZn+9cChBfe/BnwzM7cBfwSua+A1pIoVClJf+gr9iNgMfBC4ub4fwFXAvnqVPcC1/byG9DJWKEh96ffqnW8BXwD+or5/AfB0Zr5Y3z8GbOr0wIjYBewC2OJ/zbUSVihIPet5ph8RHwJOZOY9Cxd3WDU7PT4zZzJzKjOnJiYmeh2GJGkF+pnpvxv4cERcA5wLnE81818fEWvq2f5m4In+hylJakLPM/3M/GJmbs7MSeBjwF2ZOQ3cDXykXm0ncEffo5QkNWIQ1+nfAHwuIh6hOsZ/ywBeQ+PI+gRp6BqpYcjMXwC/qD9+FHh7E8+rFrE+QRoJviNXq8P6BGkkGPpaHdYnSCPB0NfqsD5BGgmGvlaH9QnSSDD0tTqsT5BGgn9ERavH+gRp6JzpS1JBDH1JKoihL0kFMfR1ZqxQkFrBE7lanhUKUms409fyrFCQWsPQ1/KsUJBaw9DX8qxQkFrD0NfyrFCQWsPQ1/KsUJBaw6t3dGasUJBawZm+JBXE0Jekghj6klQQQ1+SCmLol8DeHEk1r95pO3tzJC3gTL/t7M2RtICh33b25khawNBvO3tzJC1g6LedvTmSFjD0287eHEkLePVOCezNkVRzpi9JBTH0Jakghr4kFcTQH2XWJ0hqmCdyR5X1CZIGoOeZfkRcGhF3R8ShiHgoIq6vl2+MiDsj4nB9u6G54RbE+gRJA9DP4Z0Xgc9n5huBK4FPRcSbgBuBA5m5DThQ39dKWZ8gaQB6Dv3MPJ6Zv6k/fg44BGwCdgB76tX2ANf2O8giWZ8gaQAaOZEbEZPA5cBB4OLMPA7VLwbgoi6P2RURsxExOzc318Qw2sX6BEkD0HfoR8RrgB8Dn8nMZ8/0cZk5k5lTmTk1MTHR7zDax/oESQPQ19U7EfEqqsDfm5m314ufjIhLMvN4RFwCnOh3kMWyPkFSw/q5eieAW4BDmfmNBZ/aD+ysP94J3NH78CRJTepnpv9u4G+BByLivnrZl4CvAj+KiOuAo8BH+xuiJKkpPYd+Zv4nEF0+vb3X55UkDY41DINihYKkEWQNwyBYoSBpRDnTHwQrFCSNKEN/EKxQkDSiDP1BsEJB0ogy9AfBCgVJI8rQHwQrFCSNKK/eGRQrFCSNIGf6klQQQ1+SCmLoS1JBDP1OrFCQ1FKeyF3MCgVJLeZMfzErFCS1mKG/mBUKklrM0F/MCgVJLWboL2aFgqQWM/QXs0JBUot59U4nVihIailn+pJUEENfkgpi6EtSQdoX+lYoSFJX7TqRa4WCJC2pXTN9KxQkaUntCn0rFCRpSe0KfSsUJGlJ7Qp9KxQkaUntCn0rFCRpSe26egesUJCkJbRrpi9JWpKhL0kFMfQlqSADCf2I+EBE/HdEPBIRNw7iNSRJK9d46EfE2cB3gauBNwEfj4g3LftAO3MkaeAGcfXO24FHMvNRgIj4N2AH8Nuujzh50s4cSVoFgzi8swl4bMH9Y/Wy7h5/3M4cSVoFgwj96LAsX7FSxK6ImI2IWZ5/vvMz2ZkjSY0aROgfAy5dcH8z8MTilTJzJjOnMnOKtWs7P5OdOZLUqMh8xSS8vyeMWAP8D7AdeBz4NfA3mflQt8ecG/HsX8F5seCXUMLpo3DkKTjZ6ADHx4XAU8MexAhwO7gN5rkdKgu3w9bMnFjJgxs/kZuZL0bE3wP/AZwNfH+pwAf4c+b58x9HxGxmTjU9rnHjdqi4HdwG89wOlX63w0C6dzLzZ8DPBvHckqTe+Y5cSSrIKIb+zLAHMCLcDhW3g9tgntuh0td2aPxEriRpdI3iTF+SNCCGviQVZKRCv8R2zoi4NCLujohDEfFQRFxfL98YEXdGxOH6dsOwx7oaIuLsiLg3In5a378sIg7W2+GHEdHlnXztERHrI2JfRDxc7xfvLG1/iIjP1j8PD0bEbRFxbgn7QkR8PyJORMSDC5Z1/N5H5Tt1Xt4fEVecyWuMTOj33M45/l4EPp+ZbwSuBD5Vf903AgcycxtwoL5fguuBQwvufw34Zr0d/ghcN5RRra5vAz/PzDcAb6XaHsXsDxGxCfg0MJWZb6Z6v8/HKGNf+BfgA4uWdfveXw1sq//tAr53Ji8wMqHPgnbOzHwemG/nbLXMPJ6Zv6k/fo7qB3wT1de+p15tD3DtcEa4eiJiM/BB4Ob6fgBXAfvqVVq/HSLifOC9wC0Amfl8Zj5NefvDGuDV9Tv81wHHKWBfyMxf8soWgm7f+x3ArVn5FbA+Ii5Z7jVGKfRX3s7ZMhExCVwOHAQuzszjUP1iAC4a3shWzbeALwCn6/sXAE9n5ov1/RL2idcBc8AP6sNcN0fEeRS0P2Tm48DXgaNUYf8McA/l7Qvzun3ve8rMUQr9M2rnbKuIeA3wY+AzmfnssMez2iLiQ8CJzLxn4eIOq7Z9n1gDXAF8LzMvB/5Eiw/ldFIfs94BXAa8FjiP6lDGYm3fF5bT08/HKIX+GbVztlFEvIoq8Pdm5u314ifn/6tW354Y1vhWybuBD0fE76kO7V1FNfNfX/8XH8rYJ44BxzLzYH1/H9UvgZL2h/cDv8vMucx8AbgdeBfl7Qvzun3ve8rMUQr9XwPb6jP0a6lO3Owf8pgGrj5ufQtwKDO/seBT+4Gd9cc7gTtWe2yrKTO/mJmbM3OS6nt/V2ZOA3cDH6lXK2E7/AF4LCJeXy/aTvVX50raH44CV0bEuvrnY34bFLUvLNDte78f+ER9Fc+VwDPzh4GWlJkj8w+4hqqW+X+B3cMezyp9ze+h+i/Z/cB99b9rqI5nHwAO17cbhz3WVdwm7wN+Wn/8OuC/gEeAfwfOGfb4VuHrfxswW+8TPwE2lLY/AF8BHgYeBP4VOKeEfQG4jeo8xgtUM/nrun3vqQ7vfLfOyweornZa9jWsYZCkgozS4R1J0oAZ+pJUEENfkgpi6EtSQQx9SSqIoS8BEfHputFy77DHIg2Sl2xKQEQ8DFydmb8b9likQXKmr+JFxD9TvfFnf0R8PiJ+UveT/yoi3lKvs7HTcmncGPoqXmb+HVVnyV8Dk8C9mfkW4EvArfVqX+myXBorhr70cu+hets/mXkXcEFE/OUSy6WxYuhLL9etrrbEmme1kKEvvdwvgWmAiHgf8FRWf9+g23JprHj1jgTUPf5TVH+16wdUf8DjFLArM++PiI2dlg9puFLPDH1JKoiHdySpIIa+JBXE0Jekghj6klQQQ1+SCmLoS1JBDH1JKsj/A1j2xjJu7B8FAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# linspace can be used to create an n-sized array of evenly spaced samples in some range:\n", "\n", "arr = np.linspace(0, 100, 20)\n", "print(arr)\n", "\n", "plt.axis([-1, 101, -1, 101])\n", "plt.xlabel('foo')\n", "lines, = plt.plot(arr, arr, 'ro')" ] }, { "cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[ 1. 1.12883789 1.27427499 1.43844989 1.62377674 1.83298071\n", " 2.06913808 2.33572147 2.6366509 2.97635144 3.35981829 3.79269019\n", " 4.2813324 4.83293024 5.45559478 6.15848211 6.95192796 7.8475997\n", " 8.8586679 10. ]\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXAAAAD4CAYAAAD1jb0+AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAO30lEQVR4nO3db4hdd53H8c8n/2im3ca2GaUmzdwKxa4061Yuu9WAiOODqq3pE6EyLUUW5smutrJQ6g5L8cHAPpAlebAsDLUacOgisaytiKtERfbBhr1JC2kbF6V2ptHY3Bqc6kZM3H73wdyJmZu5M5N7zj3n/M55vyDMzMmdOd9L0k9/+Z3f7/tzRAgAkJ4tZRcAABgOAQ4AiSLAASBRBDgAJIoAB4BEbSvyZrt3745Wq1XkLQEgeSdOnHgzIsb7rxca4K1WS51Op8hbAkDybC+sdZ0pFABIFAEOAIkiwAEgUQQ4ACSKAAeARBHgAJAoAhwAEkWAA0CiCHAASBQBDgCJIsABIFEEOAAkigAHgEQR4ACQKAIcABK1YYDbftr2OdsvXXHtZtvft/3T3sebRlsmAJRj/tS8Woda2vKlLWodamn+1HzZJV22mRH41yTd23ftCUnHIuIOScd6XwNArcyfmtf089NaWFpQKLSwtKDp56crE+IbBnhE/FjS+b7LByUd6X1+RNIDOdcFAKWbOTajC5curLp24dIFzRybKami1YadA39XRJyVpN7Hdw56oe1p2x3bnW63O+TtAKB4i0uL13S9aCN/iBkRcxHRjoj2+PhVZ3ICQGXt27Xvmq4XbdgAf8P2rZLU+3guv5IAoBpmJ2c1tn1s1bWx7WOanZwtqaLVhg3w5yQ90vv8EUnfyqccAKiOqf1Tmrt/ThO7JmRZE7smNHf/nKb2T5VdmiTJEbH+C+xnJH1E0m5Jb0h6UtK/S/qGpH2SFiV9OiL6H3Repd1uR6fTyVgyADSL7RMR0e6/vm2jb4yIzwz4rcnMVQEAhsZOTABIFAEOAIkiwAEgUQQ4ACSKAAeARBHgAJAoAhwAEkWAA0CiCHAASBQBDgCJIsABIFEEOAAkigAHgEQR4ACQKAIcQOXMn5pX61BLW760Ra1DrcqcAl81G/YDB4AizZ+a1/Tz05dPg19YWtD089OSVJmTcKqCETiASpk5NnM5vFdcuHRBM8dmSqqoughwAJWyuLR4TdebjAAHUCn7du27putNRoADqJTZyVmNbR9bdW1s+5hmJ2dLqqi6CHAAlTK1f0pz989pYteELGti14Tm7p/jAeYaHBGF3azdbken0ynsfgBQB7ZPRES7/zojcABIFAEOAIkiwAEgUQQ4ACSKAAeARBHgAJAoAhwAEkWAA0CiCHAASFSmALf9Bdsv237J9jO2r8urMADA+oYOcNt7JH1eUjsi7pK0VdKDeRUGAFhf1imUbZJ22t4maUzSL7OXBADYjKEDPCJ+IenLkhYlnZW0FBHf63+d7WnbHdudbrc7fKUACsfZlNWWZQrlJkkHJd0u6d2Srrf9UP/rImIuItoR0R4fHx++UgCFWjmbcmFpQaG4fDYlIV4dWaZQPibp5xHRjYhLkp6V9KF8ygJQNs6mrL4sAb4o6R7bY7YtaVLS6XzKAlA2zqasvixz4MclHZV0UtKp3s+ay6kuACXjbMrqy7QKJSKejIg7I+KuiHg4Iv6QV2EAysXZlNXHTkwAa+JsyurjTEwAqDjOxASAmiHAASBRBDgAJIoAB4BEEeAAkCgCHAASRYADQKIIcABIFAEOAIkiwIGEcMACrrSt7AIAbM7KAQsrPbpXDliQRH+ShmIEDiSCAxbQjwAHEsEBC+hHgAOJ4IAF9CPAgURwwAL6EeBAIjhgAf040AEAKo4DHQCgZghwAEgUAQ4AiSLAASBRBDgAJIoAB4BEEeAAkCgCHAASRYADOaJfN4pEP3AgJ/TrRtEYgQM5oV83ipYpwG2/w/ZR2z+xfdr2B/MqDEgN/bpRtKwj8MOSvhsRd0p6v6TT2UsC0kS/bhRt6AC3faOkD0v6iiRFxMWI+E1ehQGpoV83ipZlBP4eSV1JX7X9gu2nbF/f/yLb07Y7tjvdbjfD7YBqo183ijZ0P3DbbUn/JelARBy3fVjSWxHxj4O+h37gAHDtRtEP/IykMxFxvPf1UUkfyPDzAADXYOgAj4hfSXrd9nt7lyYlvZJLVQCADWXdyPM5SfO2d0h6VdJns5cEANiMTAEeES9KumpeBgAweuzEBIBEEeAAkCgCHI1F50Ckjm6EaCQ6B6IOGIGjkegciDogwNFIdA5EHRDgaCQ6B6IOCHA0Ep0DUQcEOBqJzoGog6G7EQ6DboQAcO1G0Y0QAFAiAhwAEkWAA0CiCHAASBQBjqTQvwT4E3qhIBn0LwFWYwSOZNC/BFiNAEcy6F8CrEaAIxn0LwFWI8CRDPqXAKsR4EgG/UuA1eiFAgAVRy8UAKgZAhwAEkWAA0CiCHAUhm3wQL7YSo9CsA0eyB8jcBSCbfBA/ghwFIJt8ED+CHAUgm3wQP4IcBSCbfBA/jIHuO2ttl+w/e08CkI9sQ0eyF8eq1AelXRa0o05/CzU2NT+KQIbyFGmEbjtvZI+KempfMoBAGxW1imUQ5Iel/T2oBfYnrbdsd3pdrsZbwcAWDF0gNu+T9K5iDix3usiYi4i2hHRHh8fH/Z2qAh2UwLVkWUO/ICkT9n+hKTrJN1o++sR8VA+paFq2E0JVMvQI/CI+GJE7I2IlqQHJf2A8K43dlMC1cI6cGwauymBasklwCPiRxFxXx4/C9XFbkqgWhiBY9PYTQlUCwGOTWM3JVAtHGoMABXHoca4Cmu6gbRxIk9DsaYbSB8j8IZiTTeQPgK8oVjTDaSPAG8o1nQD6SPAG4o13UD6CPCGYk03kD7WgQNAxbEOvAFY1w00C+vAa4J13UDzMAKvCdZ1A81DgNcE67qB5iHAa4J13UDzEOA1wbpuoHkI8JpgXTfQPKwDT8D8qXnNHJvR4tKi9u3ap9nJWYIZaJBB68BZRlhxLA8EMAhTKBXH8kAAgxDgFcfyQACDEOAVx/JAAIMQ4BXH8kAAgxDgJduoARXLAwEMwjLCEvWvMJGWR9cENIAr0U62glhhAiALArxErDABkAUBXiJWmADIggAfsfUeUrLCBEAWQwe47dts/9D2adsv2340z8LqYOUh5cLSgkJxeRv8SoizwgRAFkOvQrF9q6RbI+Kk7T+TdELSAxHxyqDvadoqlNahlhaWFq66PrFrQq899lrxBQFIUu6rUCLibESc7H3+W0mnJe0ZvsT64SElgFHKZQ7cdkvS3ZKOr/F707Y7tjvdbjeP2yWDh5QARilzgNu+QdI3JT0WEW/1/35EzEVEOyLa4+PjWW9XSYMeVPKQEsAoZeoHbnu7lsN7PiKezaektGymXzeHMQAYhSwPMS3piKTzEfHYZr6njg8xeVAJYNRGsZX+gKSHJX3U9ou9X5/I8POSxINKAGXJsgrlPyPCEfEXEfGXvV/fybO4Khk0z82DSgBlYSfmJqy3IYcHlQDKQoBvwnpdA9lNCaAsnEq/hvlT86tWjqz1kFL60zz31P4pAhtA4QjwPmstC7Ss0NWrdZjnBlAmplD6rDVdEgpZXnWNeW4AZSPA+wxa/hcK5rkBVErjp1D657tv3nmzfv37X1/1OjbmAKiaRgf4WvPd27ds146tO3Tx/y5efh3TJQCqqFEB3j/a/t3F3101333p7Uu6ZectumHHDfQvAVBpjQnwtUbbg5z//Xm9+fibRZUGAENpzEPMtVaXDMLyQAApqHWAX9m/ZL0R95WY7waQitoGeH//kkFu2XkLywMBJKm2c+CbmTIZ2z6mwx8/TGADSFJtR+Dr9eNmtA2gDmo7Ah/UhIoNOQDqorYjcPp0A6i72gY4fboB1N3QhxoPo46HGgPAqI3iUGMAQIkIcABIFAEOAIkiwAEgUQQ4ACSKAAeARBHgAJAoAhwAEkWAA0CiCHAASBQBDgCJIsABIFGZAtz2vbb/x/bPbD+RV1EAgI0NHeC2t0r6F0kfl/Q+SZ+x/b68CgMArC/LCPyvJP0sIl6NiIuS/k3SwXzKAgBsJMuRanskvX7F12ck/XX/i2xPS5ruffkH2y9luGeKdkt6s+wiCta099y09yvxnos2sdbFLAHuNa5ddTpERMxJmpMk2521mpLXGe+5/pr2fiXec1VkmUI5I+m2K77eK+mX2coBAGxWlgD/b0l32L7d9g5JD0p6Lp+yAAAbGXoKJSL+aPvvJP2HpK2Sno6Ilzf4trlh75cw3nP9Ne39SrznSij0UGMAQH7YiQkAiSLAASBRhQR407bc277N9g9tn7b9su1Hy66pKLa32n7B9rfLrqUItt9h+6jtn/T+vD9Ydk2jZvsLvb/XL9l+xvZ1ZdeUN9tP2z535b4V2zfb/r7tn/Y+3lRmjVIBAd7QLfd/lPT3EfHnku6R9LcNeM8rHpV0uuwiCnRY0ncj4k5J71fN37vtPZI+L6kdEXdpeQHDg+VWNRJfk3Rv37UnJB2LiDskHet9XaoiRuCN23IfEWcj4mTv899q+T/qPeVWNXq290r6pKSnyq6lCLZvlPRhSV+RpIi4GBG/KbeqQmyTtNP2NkljquH+j4j4saTzfZcPSjrS+/yIpAcKLWoNRQT4Wlvuax9mK2y3JN0t6Xi5lRTikKTHJb1ddiEFeY+krqSv9qaNnrJ9fdlFjVJE/ELSlyUtSjoraSkivlduVYV5V0SclZYHaZLeWXI9hQT4prbc15HtGyR9U9JjEfFW2fWMku37JJ2LiBNl11KgbZI+IOlfI+JuSf+rCvyzepR6874HJd0u6d2Srrf9ULlVNVcRAd7ILfe2t2s5vOcj4tmy6ynAAUmfsv2alqfJPmr76+WWNHJnJJ2JiJV/XR3VcqDX2cck/TwiuhFxSdKzkj5Uck1FecP2rZLU+3iu5HoKCfDGbbm3bS3Pi56OiH8uu54iRMQXI2JvRLS0/Gf8g4io9cgsIn4l6XXb7+1dmpT0SoklFWFR0j22x3p/zydV8we3V3hO0iO9zx+R9K0Sa5GUrRvhpgy55T51ByQ9LOmU7Rd71/4hIr5TYk0Yjc9Jmu8NTl6V9NmS6xmpiDhu+6ikk1pebfWCKrjFPCvbz0j6iKTdts9IelLSP0n6hu2/0fL/yD5dXoXL2EoPAIliJyYAJIoAB4BEEeAAkCgCHAASRYADQKIIcABIFAEOAIn6f0MjMvDzfWQhAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Likewise, logspace \n", "arr = np.logspace(0, 1, 20, base=10)\n", "print(arr)\n", "plt.axis([0, 11, 0, 11])\n", "lines, = plt.plot(arr, arr, 'go')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Creating a 2D array" ] }, { "cell_type": "code", "execution_count": 21, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[1, 2, 3],\n", " [4, 5, 6],\n", " [7, 8, 9]])" ] }, "execution_count": 21, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Create a 3x3 array\n", "\n", "import numpy as np\n", "\n", "example_array = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])\n", "example_array" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "As with everything else in Python, assigning an array to a variable does not copy its data. Instead, it creates an alias for the original data. For example, let's create an array of ones and assign it to a variable `first`, then assign the value of `first` to `second`: " ] }, { "cell_type": "code", "execution_count": 22, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[1. 1.]\n", " [1. 1.]]\n" ] } ], "source": [ "# Variable assignment does not *copy* data in array; it just creates a new pointer to that array\n", "first = numpy.ones((2, 2))\n", "second = first\n", "print(first)" ] }, { "cell_type": "code", "execution_count": 23, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[9. 1.]\n", " [1. 1.]]\n" ] } ], "source": [ "# Index assignment updates elements in an array\n", "# To index N-D arrays use a , to separate the indices into each dimension\n", "\n", "second[0, 0] = 9\n", "print(first)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "If we really want a copy of the array so that we can make changes without affecting the original data, we can use the copy method:" ] }, { "cell_type": "code", "execution_count": 24, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[1. 1.]\n", " [1. 1.]]\n", "\n", "first:\n", "[[1. 1.]\n", " [1. 1.]]\n", "\n", "second:\n", "[[9. 1.]\n", " [1. 1.]]\n" ] } ], "source": [ "# If we really want a copy we can use the .copy() method\n", "first = numpy.ones((2, 2))\n", "print(first)\n", "print()\n", "\n", "second = first.copy()\n", "second[0, 0] = 9\n", "print('first:')\n", "print(first)\n", "print()\n", "print('second:')\n", "print(second)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Arrays have properties as well as methods. We have already met dtype, which is the array's data type. Another is shape, which is a tuple of the array's size along each dimension: " ] }, { "cell_type": "code", "execution_count": 25, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[1. 1.]\n", " [1. 1.]]\n", "(2, 2)\n" ] } ], "source": [ "# Show the shape\n", "\n", "print(first)\n", "print(first.shape)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "It is possible to modify an array's shape, but be careful that the total number of elements in the array is still the same:" ] }, { "cell_type": "code", "execution_count": 26, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[1.],\n", " [1.],\n", " [1.],\n", " [1.]])" ] }, "execution_count": 26, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Modify the shape with a valid shape\n", "\n", "first.shape = (4, 1)\n", "first" ] }, { "cell_type": "code", "execution_count": 27, "metadata": {}, "outputs": [ { "ename": "ValueError", "evalue": "cannot reshape array of size 4 into shape (3,2)", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mValueError\u001b[0m Traceback (most recent call last)", "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[0;31m# Modify the shape with an invalid shape\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 3\u001b[0;31m \u001b[0mfirst\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mshape\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0;36m3\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m2\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 4\u001b[0m \u001b[0mfirst\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;31mValueError\u001b[0m: cannot reshape array of size 4 into shape (3,2)" ] } ], "source": [ "# Modify the shape with an invalid shape\n", "\n", "first.shape = (3, 2)\n", "first" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Notice that there are no parentheses after `shape`: it is a piece of data, not a method call. Also note that the tuple in `shape` is exactly what we pass into functions like `zeros` to create new arrays, which makes it easy to reproduce the shape of existing data:" ] }, { "cell_type": "code", "execution_count": 28, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[0.],\n", " [0.],\n", " [0.],\n", " [0.]])" ] }, "execution_count": 28, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# You can use .shape to define new arrays with that shape\n", "blank = np.zeros(first.shape)\n", "blank" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Other data members include `size`, which is the total number of elements in the array, and `nbytes` which is the total physical bytes of memory used by the array. As the code below shows, `size` is simply the product of the array's lengths along its dimensions, and `nbytes` is the product of `size` and the size of the data type:" ] }, { "cell_type": "code", "execution_count": 29, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "size: 84\n", "nbytes: 672\n" ] } ], "source": [ "# Arrays also have .size and .nbytes attributes\n", "block = numpy.zeros((4, 7, 3))\n", "print('size:', block.size)\n", "print('nbytes:', block.nbytes)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "There are also special methods for reshaping the array in common ways. We can rearrange the data in an array using the `transpose` method, which flips the array on all its axes. This does not actually move data around in memory. Instead, it creates an alias that appears to have the values stored differently. We also call this a *view* of the array:" ] }, { "cell_type": "code", "execution_count": 30, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[1 4]\n", " [2 5]\n", " [3 6]]\n", "\n", "The original array is unchanged:\n", "[[1 2 3]\n", " [4 5 6]]\n" ] } ], "source": [ "# .transpose() creates a *view* of the same data but transposed\n", "\n", "arr = numpy.array([[1, 2, 3],\n", " [4, 5, 6]])\n", "\n", "print(arr.transpose())\n", "print()\n", "print('The original array is unchanged:')\n", "print(arr)" ] }, { "cell_type": "code", "execution_count": 31, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[ 1 4]\n", " [ 2 42]\n", " [ 3 6]]\n", "\n", "[[ 1 2 3]\n", " [ 4 42 6]]\n" ] } ], "source": [ "# However, modifying the transposed array will modify the original array\n", "# (unless we explictly make a copy of it first)\n", "trans = arr.transpose()\n", "trans[1, 1] = 42\n", "print(trans)\n", "print()\n", "print(arr)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This is a common enough operation (especially on matrices) that `arr.T` can be used as shorthand for `arr.transpose()`. Note that this is not a function call with parentheses after it:" ] }, { "cell_type": "code", "execution_count": 32, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[ 1 4]\n", " [ 2 42]\n", " [ 3 6]]\n" ] } ], "source": [ "# We can also use arr.T\n", "print(arr.T)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The `ravel` method does something similar: it creates a one-dimensional view of the original data. As you'd expect, the result's shape has a single value, which is the number of elements we started with." ] }, { "cell_type": "code", "execution_count": 33, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[ 1 2 3 4 42 6]\n", "\n", "The original array is still unchanged:\n", "[[ 1 2 3]\n", " [ 4 42 6]]\n" ] } ], "source": [ "# .ravel() creates a 1-dimensional view of the original data\n", "print(arr.ravel())\n", "print()\n", "print('The original array is still unchanged:')\n", "print(arr)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "But even though the array returned by `ravel()` has a different shape it's worth emphasizing again that it is a *view* of the original array, and that modifying its contents will modify the original array too:" ] }, { "cell_type": "code", "execution_count": 34, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[ 1 2 3]\n", " [ 9 42 6]]\n" ] } ], "source": [ "# Remember, the result of .ravel() is just a view; updating the array returned by\n", "# .ravel() updates the original array data too\n", "arr.ravel()[3] = 9\n", "print(arr)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "What order do raveled values appear in? Let's start by thinking about a 2×4 array `A`. It looks two-dimensional, but the computer's memory is 1-dimensional: each location is identified by a single integer address. Any program that works with multi-dimensional data must therefore decide how to lay out those values. \n", "\n", "\n", "\n", "One possibility is *row-major order*, which concatenates the rows. This is what C uses, and since Python was originally written in C, it uses the same convention: \n", "\n", "\n", "\n", "In contrast, column-major order concatenates the columns. FORTRAN does this, and MATLAB follows along.\n", "\n", "\n", "\n", "There's no real difference in performance or usability, but the differences cause headaches when data has to be moved from one programming language to another. For example, if your Python code wants to call an eigenvalue function written in FORTRAN, you will probably have to rearrange the data, just as you have to be careful about 0-based versus 1-based indexing. Note that you cannot use the array's `transpose` method to do this, since, as explained earlier, it doesn't actually move data around.\n", "\n", "It's also possible that if your software is specifically tuned to operate on an array one row at a time or one column at a time it may be desireable to have the data arranged in memory so that it's accessed linearly. But this is the kind of performance tuning you won't need until and unless you *know* you need it. For the most part stick with row-major order (the default in Numpy) and don't worry about it." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "If we really want to change the physical size of the data, we have to use array.resize. This works in place, i.e., it modifies the array, rather than returning a new alias: " ] }, { "cell_type": "code", "execution_count": 35, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[0 1 2 3 4 5 6 7 8]\n", "\n", "[[0 1]\n", " [2 3]]\n" ] } ], "source": [ "# Resize a larger array to a smaller array\n", "\n", "block = numpy.arange(9)\n", "print(block)\n", "print()\n", "\n", "block.resize(2, 2)\n", "print(block)\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "As the example above shows, when we resize a 3×3 array to be 2×2, we get the values that were in the first two rows and columns. (And note that once again, the new dimensions are passed directly, rather than in a tuple.)\n", "\n", "If we enlarge the array by resizing, the new locations are assigned zero. Which locations are \"new\" is determined by the raveling order of the array; as the example below shows, the existing values are packed into the first part of memory, *not* into the upper left corner of the logical matrix: " ] }, { "cell_type": "code", "execution_count": 36, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[1. 1.]\n", " [1. 1.]]\n", "\n", "[[1. 1. 1.]\n", " [1. 0. 0.]\n", " [0. 0. 0.]]\n" ] } ], "source": [ "# Resize a smaller array to a larger array\n", "small = numpy.ones((2, 2))\n", "print(small)\n", "print()\n", "\n", "large = small.copy()\n", "large.resize(3, 3)\n", "print(large)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "It is, however, possible to fill the upper left corner by first allocating a new array of zeros and using slice assignment:" ] }, { "cell_type": "code", "execution_count": 37, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[1. 1. 0.]\n", " [1. 1. 0.]\n", " [0. 0. 0.]]\n" ] } ], "source": [ "# Fill a corner of an array from a smaller array\n", "large = numpy.zeros((3, 3))\n", "large[:2, :2] = small\n", "print(large)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Indexing Numpy Arrays\n", "\n", "Now we'll take a closer look at some of the ways we can index arrays. It may seem like a small topic at first, but clever indexing allows us to avoid writing loops, which reduces the size of our code and makes it more efficient.\n", "\n", "Arrays are subscripted by integers, just like lists and other sequences, so they can be sliced like other sequences as well. For example, if `block` is the array shown below, then `block[0:3, 0:2]` selects its first three rows and the first two columns: \n", "\n", "\n", "\n", "The comma syntax (`[X, Y]`) for indexing multi-dimensional arrays was lobbied for and eventually added to the Python language by the scientific community. Although it is possible to write `block[0:3][0:2]` this is comparatively inefficient, especially when doing many array indexing operations. (Exercise question: Why?)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "As with other sliceable types, it's possible to assign to slices (as we saw a bit earlier). For example, we can assign zero to columns 1 and 2 in row 1 of ``block in a single statement: " ] }, { "cell_type": "code", "execution_count": 38, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[ 10 0 0 40]\n", " [110 120 130 140]\n", " [210 220 230 240]]\n" ] } ], "source": [ "# Assign zero to a slice of the first row\n", "block = numpy.array([[10, 20, 30, 40],\n", " [110, 120, 130, 140],\n", " [210, 220, 230, 240]])\n", "block[0, 1:3] = 0\n", "print(block)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "And, as with most other operations (such as `transpose` and `ravel` like I saw earlier), slicing creates an alias rather than immediately copying data: " ] }, { "cell_type": "code", "execution_count": 39, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[0, 1],\n", " [2, 3],\n", " [4, 5]])" ] }, "execution_count": 39, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Create a 2D array\n", "original = numpy.arange(6).reshape((3, 2))\n", "original" ] }, { "cell_type": "code", "execution_count": 40, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[0, 1],\n", " [2, 3]])" ] }, "execution_count": 40, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Slice a corner of of the array\n", "slc = original[0:2, 0:2]\n", "slc" ] }, { "cell_type": "code", "execution_count": 41, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[0, 0],\n", " [0, 0],\n", " [4, 5]])" ] }, "execution_count": 41, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Set the values of the slice to zero--the original has changed\n", "slc[:, :] = 0\n", "original" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Notice in the example above how we used `slice[:, :]` to refer to all of the array's elements at once. All of Python's other slicing shortcuts work as well, so that expressions like `original[-2:, 1:]` behave sensibly." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Slicing on both sides of an assignment is a way to shift data along the axes. If `vector` is a one-dimensional array, then `vector[1:4]` selects locations 1, 2, and 3, while `vector[0:3]` selects locations 0, 1, and 2. Assigning the former to the latter therefore overwrites the lower three values with the upper three, leaving the uppermost value untouched:" ] }, { "cell_type": "code", "execution_count": 42, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([20, 30, 40, 40])" ] }, "execution_count": 42, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Assign a slice of a vector to another slice of the same vector\n", "vector = numpy.array([10, 20, 30, 40])\n", "vector[0:3] = vector[1:4]\n", "vector" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Compare this with loop-based code that shifts values up or down, and you'll see why most programmers prefer the vectorized programming model. " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We co do even more sophisticated things by using a list or an array as a subscript. For example, if `subscript` is a list containing 3, 1, and 2, then `vector[subscript]` creates a new array whose elements are selected from vector in the obvious way: \n", "\n", "" ] }, { "cell_type": "code", "execution_count": 43, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([ 0, 10, 20, 30])" ] }, "execution_count": 43, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Create the vector [0, 10, 20, 30] using arange() and *\n", "vector = numpy.arange(4) * 10\n", "vector" ] }, { "cell_type": "code", "execution_count": 44, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([30, 10, 20])" ] }, "execution_count": 44, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Index with a list of indices\n", "subscript = [3, 1, 2]\n", "vector[subscript]" ] }, { "cell_type": "code", "execution_count": 45, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([ 0, 10, 20, 30])" ] }, "execution_count": 45, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# It should be emphasized that such arbitrary subscripting returns a *copy* of the original array:\n", "sub = vector[subscript]\n", "sub[:] = 0\n", "vector" ] }, { "cell_type": "code", "execution_count": 46, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([42, 10, 20, 42])" ] }, "execution_count": 46, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# It is however possible to perform assignment into an array using an arbitrary subscript,\n", "# without copying\n", "subscript = [0, 3]\n", "vector[subscript] = 42\n", "vector" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Boolean indexing" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can use Boolean masking on the left side of assignment as well, though we have to be careful about its meaning. If we use a mask directly, elements are taken in order from the source on the right and assigned to elements corresponding to True values in the mask: " ] }, { "cell_type": "code", "execution_count": 47, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([0, 2])" ] }, "execution_count": 47, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Mask an array [0 1 2 3] with the mask [T F T F]\n", "a = numpy.array([0, 1, 2, 3])\n", "mask = [True, False, True, False]\n", "a[mask]" ] }, { "cell_type": "code", "execution_count": 48, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([101, 1, 103, 3])" ] }, "execution_count": 48, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Assign [101 102 103 104] using the mask\n", "a[mask] = numpy.array([101, 102, 103, 104])[mask]\n", "a" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Operators like `<` and `==` work the way we would expect with arrays, but there is one trick. Python does not allow objects to re-define the meaning of `and`, `or`, and `not`, since they are keywords. The expression `(vector <= 20) and (vector >= 20)` therefore produces an error message instead of selecting elements with exactly the value 20: " ] }, { "cell_type": "code", "execution_count": 49, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([ True, True, True, False])" ] }, "execution_count": 49, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Give the items <= 20 of [0 10 20 30]\n", "vector = numpy.array([0, 10, 20, 30])\n", "vector <= 20" ] }, { "cell_type": "code", "execution_count": 50, "metadata": {}, "outputs": [ { "ename": "ValueError", "evalue": "The truth value of an array with more than one element is ambiguous. Use a.any() or a.all()", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mValueError\u001b[0m Traceback (most recent call last)", "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[0;31m# Gives the items <= 20 and > 0\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 2\u001b[0;31m \u001b[0;34m(\u001b[0m\u001b[0mvector\u001b[0m \u001b[0;34m<=\u001b[0m \u001b[0;36m20\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;32mand\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0mvector\u001b[0m \u001b[0;34m>\u001b[0m \u001b[0;36m0\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m", "\u001b[0;31mValueError\u001b[0m: The truth value of an array with more than one element is ambiguous. Use a.any() or a.all()" ] } ], "source": [ "# Gives the items <= 20 and > 0\n", "(vector <= 20) and (vector > 0)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "One way around this is to use functions like `logical_and` and `logical_or`, which combine the elements' Boolean arrays like their namesakes: " ] }, { "cell_type": "code", "execution_count": 51, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([False, True, True, False])" ] }, "execution_count": 51, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Do the same expression with logical_and\n", "numpy.logical_and(vector <= 20, vector > 0)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Another is to use the bitwise and/or operators `&` and `|`:" ] }, { "cell_type": "code", "execution_count": 52, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([False, True, True, False])" ] }, "execution_count": 52, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Do the same expression with &\n", "(vector <= 20) & (vector > 0)" ] }, { "cell_type": "code", "execution_count": 53, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([1, 2, 3, 7])" ] }, "execution_count": 53, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Since comparison operators return an array of bools, a comparison expression\n", "# can also be used to index an array\n", "\n", "array1 = np.array([1, 1, 1, 2, 2, 2, 1])\n", "array2 = np.array([1, 2, 3, 4, 5, 6, 7])\n", "array2[array1 == 1]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "There are many other slick ways to index and get values from Numpy arrays, such as `array.where()`, `array.select()`, and `array.choose()`. Some of them can get very complex and difficult to understand, but have their uses--you can read more about them and other indexing tricks in the [Numpy documentation](https://docs.scipy.org/doc/numpy/reference/arrays.indexing.html), as well as its more advanced [indexing routines](https://docs.scipy.org/doc/numpy/reference/routines.indexing.html)." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Math\n", "### Arrays\n", "Math on arrays is vectorized and behaves more or less like \"most\" scientists would expect:" ] }, { "cell_type": "code", "execution_count": 54, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([3, 3, 3, 5, 5, 5, 3])" ] }, "execution_count": 54, "metadata": {}, "output_type": "execute_result" } ], "source": [ "array1 = np.array([1, 1, 1, 2, 2, 2, 1])\n", "array2 = np.array([1, 2, 3, 4, 5, 6, 7])\n", "\n", "array1 * 2 + 1" ] }, { "cell_type": "code", "execution_count": 55, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([ 1, 2, 3, 8, 10, 12, 7])" ] }, "execution_count": 55, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Multiply array1 and array2\n", "array1 * array2" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "But that doesn't mean they behave exactly like the matrices that mathematicians use. For example, let's create an array, and then multiply it by itself: " ] }, { "cell_type": "code", "execution_count": 56, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[ 1, 4],\n", " [ 9, 16]])" ] }, "execution_count": 56, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Create a 2x2 array and multiply it with itself\n", "arr = numpy.array([[1, 2], [3, 4]])\n", "arr * arr" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "NumPy does the operation elementwise, instead of doing \"real\" matrix multiplication. On the bright side, elementwise operation means that array addition works as you would expect: " ] }, { "cell_type": "code", "execution_count": 57, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[2, 4],\n", " [6, 8]])" ] }, "execution_count": 57, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Now add it with itself\n", "arr + arr" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "And since there's only one sensible way to interpret an expression like \"array plus one\", NumPy does the sensible thing there too. " ] }, { "cell_type": "code", "execution_count": 58, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[2, 3],\n", " [4, 5]])" ] }, "execution_count": 58, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Add 1 to all elements--note that this does not *modify* arr\n", "arr + 1" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Like other array-based libraries or languages, NumPy provides many useful tools for common arithmetic operations. For example, we can add up the values in our array with a single function call: " ] }, { "cell_type": "code", "execution_count": 59, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "10" ] }, "execution_count": 59, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Sum the array\n", "numpy.sum(arr)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can also calculate the partial sum along each axis by passing an extra argument into `sum`: " ] }, { "cell_type": "code", "execution_count": 60, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([4, 6])" ] }, "execution_count": 60, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Sum along the row-axis\n", "numpy.sum(arr, axis=0)" ] }, { "cell_type": "code", "execution_count": 61, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([3, 7])" ] }, "execution_count": 61, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Sum along the column-axis\n", "numpy.sum(arr, axis=1)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Matrices\n", "There also exists a special matrix class that behaves more naturally like a matrix. For example the `*` operator on two matrix objects performs matrix multiplication:" ] }, { "cell_type": "code", "execution_count": 62, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "matrix([[14],\n", " [32]])" ] }, "execution_count": 62, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Make a 2x3 and a 3x1 matrix/column vector and mupltiply the former by the latter\n", "matrix1 = np.matrix([[1, 2, 3], [4, 5, 6]])\n", "matrix2 = np.matrix([1, 2, 3]).transpose()\n", "matrix1 * matrix2" ] }, { "cell_type": "code", "execution_count": 63, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "matrix([[14, 32],\n", " [32, 77]])" ] }, "execution_count": 63, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Multiply a matrix by its transpose\n", "matrix1 * matrix1.transpose()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Various basic linear algebra operations are also possible on matrics; inverse, eigenvalues, determinants, etc.; many of these are available through the [`np.linalg`](https://docs.scipy.org/doc/numpy/reference/routines.linalg.html) module:" ] }, { "cell_type": "code", "execution_count": 64, "metadata": {}, "outputs": [], "source": [ "A = np.matrix([\n", " [2, 0, 0],\n", " [1, 1, 1],\n", " [0, 0, 2]\n", "])" ] }, { "cell_type": "code", "execution_count": 65, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "matrix([[ 0.5, 0. , 0. ],\n", " [-0.5, 1. , -0.5],\n", " [ 0. , 0. , 0.5]])" ] }, "execution_count": 65, "metadata": {}, "output_type": "execute_result" } ], "source": [ "np.linalg.inv(A)" ] }, { "cell_type": "code", "execution_count": 66, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "4.0" ] }, "execution_count": 66, "metadata": {}, "output_type": "execute_result" } ], "source": [ "np.linalg.det(A)" ] }, { "cell_type": "code", "execution_count": 67, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[1. 2. 2.]\n", "[[0. 0.70710678 0. ]\n", " [1. 0.70710678 0.70710678]\n", " [0. 0. 0.70710678]]\n" ] } ], "source": [ "eigenvals, eigenvects = np.linalg.eig(A)\n", "print(eigenvals)\n", "print(eigenvects)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Other functions are applied element-wise to the matrix, just as they are on plain arrays:" ] }, { "cell_type": "code", "execution_count": 68, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "matrix([[7.3890561 , 1. , 1. ],\n", " [2.71828183, 2.71828183, 2.71828183],\n", " [1. , 1. , 7.3890561 ]])" ] }, "execution_count": 68, "metadata": {}, "output_type": "execute_result" } ], "source": [ "np.exp(A)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Importing and Exporting Data\n", "The numpy function `genfromtxt()` is a powerful way to import text data.\n", "It can use different delimiters, skip header rows, control the type of imported data, give columns of data names, handle missing data, and a number of other useful goodies. See the [documentation](http://docs.scipy.org/doc/numpy/reference/generated/numpy.genfromtxt.html) for a full list of features of run `help(np.genfromtxt)` from the Python shell (after importing the module of course).\n", "\n", "The more recently added [`np.loadtxt`](https://docs.scipy.org/doc/numpy/reference/generated/numpy.loadtxt.html) is similar to `genfromtxt` but a bit simpler, and does not have options for handling missing data." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Basic Import and Export\n", "#### Import\n", "Basic imports using Numpy will treat all data as floats.\n", "If we're doing a basic import we'll typically want to skip the header row (since it's generally not composed of numbers)." ] }, { "cell_type": "code", "execution_count": 69, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "a, b, c, d\r\n", "1.0, 2.0, 4.0\r\n", "5.0, 6.0, 7.0\r\n", "1e2, 2e3, 4e5\r\n", "1e-2, 2e-3, 4.345e-5\r\n" ] } ], "source": [ "!cat data/examp_data.txt" ] }, { "cell_type": "code", "execution_count": 70, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[1.000e+00, 2.000e+00, 4.000e+00],\n", " [5.000e+00, 6.000e+00, 7.000e+00],\n", " [1.000e+02, 2.000e+03, 4.000e+05],\n", " [1.000e-02, 2.000e-03, 4.345e-05]])" ] }, "execution_count": 70, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Use genfromtxt() to create an array from the file\n", "data = np.genfromtxt('./data/examp_data.txt', delimiter=',', skip_header=1)\n", "data" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Export" ] }, { "cell_type": "code", "execution_count": 71, "metadata": {}, "outputs": [], "source": [ "# Use savetxt() to dump the output to a new file with a different delimiter\n", "np.savetxt('./data/examp_output.txt', data, delimiter=', ')" ] }, { "cell_type": "code", "execution_count": 72, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1.000000000000000000e+00, 2.000000000000000000e+00, 4.000000000000000000e+00\r\n", "5.000000000000000000e+00, 6.000000000000000000e+00, 7.000000000000000000e+00\r\n", "1.000000000000000000e+02, 2.000000000000000000e+03, 4.000000000000000000e+05\r\n", "1.000000000000000021e-02, 2.000000000000000042e-03, 4.344999999999999904e-05\r\n" ] } ], "source": [ "!cat data/examp_output.txt" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Importing Data Tables\n", "Lots of scientific data comes in the form of tables, with one row per observation, and one column per thing observed.\n", "Often the different columns to have different types (including text).\n", "The best way to work with this type of data is in a Structured Array." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Import\n", "To do this we let Numpy automatically detect the data types in each column using the optional argument ``dtype=None``.\n", "We can also use an existing header row as the names for the columns using the optional arugment ``Names=True``." ] }, { "cell_type": "code", "execution_count": 73, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "site,species,mass\r\n", "1,DS,125\r\n", "1,DM,70\r\n", "2,DM,55\r\n", "1,CB,40\r\n", "2,DS,110\r\n", "1,CB,45\r\n" ] } ], "source": [ "!cat data/examp_data_species_mass.txt" ] }, { "cell_type": "code", "execution_count": 74, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([(1, 'DS', 125), (1, 'DM', 70), (2, 'DM', 55), (1, 'CB', 40),\n", " (2, 'DS', 110), (1, 'CB', 45)],\n", " dtype=[('site', '" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "from scipy import special\n", "from matplotlib import pyplot as plt\n", "\n", "plt.xlim(0, 10)\n", "plt.ylim(-1, 1)\n", "x = np.linspace(0, 10, 1000)\n", "\n", "for v in range(6):\n", " y = special.jv(v, x)\n", " plt.plot(x, special.jv(v, x))\n", " \n", " # Add labels\n", " maximum = np.argmax(y)\n", " xpos = x[maximum] + 0.25 # trial and error\n", " ypos = y[maximum] + 0.05\n", " plt.text(xpos, ypos, f'$ J_{v}(x) $ ')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Numerical integration\n", "\n", "Functions can be integrated numerically over a given interval using, among other more specialized functions [`scipy.integrate.quad`](https://docs.scipy.org/doc/scipy/reference/generated/scipy.integrate.quad.html#scipy.integrate.quad). For example:\n", "\n", "$$\n", "\\int_0^{2\\pi} J_1(x)dx\n", "$$" ] }, { "cell_type": "code", "execution_count": 84, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "val: 0.7797230914600654\n", "error upper bound: 2.2568868939342564e-14\n" ] } ], "source": [ "from scipy.integrate import quad\n", "val, err = quad(lambda x: special.jv(1, x), 0, 2 * np.pi)\n", "print('val:', val)\n", "print('error upper bound:', err)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Solving ODEs\n", "\n", "Solve $ y(x) $ for\n", "\n", "$$\n", "x^2\\frac{d^2y}{dx^2} + x\\frac{dy}{dx} + (x^2 - 1)y = 0\n", "$$\n", "\n", "with initial values $ y(0) = 0, y'(0) = 0 $.\n", "\n", "using [`scipy.integrate.odeint`](https://docs.scipy.org/doc/scipy/reference/generated/scipy.integrate.odeint.html). `odeint()` solves first-order ODEs, so first we must convert to one by defining $ z = y' $ so that:\n", "\n", "$$\n", "\\frac{dz}{dx} = \\frac{y}{x^2} - \\frac{z}{x} - y\n", "$$\n", "\n", "We can then solve for $ y $ and $ z $ simulateously by packing them into a vector ``[y, z]``:" ] }, { "cell_type": "code", "execution_count": 85, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[ 2.22044605e-16 2.22044605e-16]\n", " [ 5.04406506e-02 4.98088155e-01]\n", " [ 1.00495655e-01 4.92369336e-01]\n", " ...\n", " [ 9.33283707e-02 -2.41613277e-01]\n", " [ 6.86195858e-02 -2.47203933e-01]\n", " [ 4.34727731e-02 -2.50283081e-01]]\n" ] } ], "source": [ "from scipy.integrate import odeint\n", "\n", "def func(yz, x):\n", " y, z = yz\n", " return [z, y/x**2 - z/x - y] # [dy/dx, dz/dx]\n", "\n", "# Solve over 100 points between 0 and 10\n", "# Use a small epsilon for the initial value to avoid divide by zero\n", "eps = np.finfo(float).eps\n", "x = np.linspace(eps, 10, 100)\n", "sol = odeint(func, [eps, eps], x)\n", "print(np.array2string(sol, threshold=10))" ] }, { "cell_type": "code", "execution_count": 86, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 86, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAtEAAAF3CAYAAABjZBdpAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nOzdeVxU1f/H8ddhUVDcd2XVcG1Bc19AEPc0tVXra5tfq2/pr7KyzEwrs83K1m+brWrWt3Kp1NxwX3BNzQVBUBAVcQNR1vP74w6GCMrAzNyB+TwfDx/KzJ173mN5+XDmc89RWmuEEEIIIYQQJedmdgAhhBBCCCHKGymihRBCCCGEsJIU0UIIIYQQQlhJimghhBBCCCGsJEW0EEIIIYQQVpIiWgghhBBCCCtJES2EEEIIIYSVpIgWQgghhBDCSlJEC6ellNqjlOppxfGNlVLpSqmLSqm1ZRx7glLqi1K8rq5Sap9Syqss45di3J+VUv0dOaYQouJTSmml1HVWHG/VddsRlFLxSqnIMrze6d6TcA5SRAunpbVuo7WOsuL4o1prH+ARG4z9mtZ6VCle+izwtdb6Yv4DSqlblVIblFKZSqmvy5qtGK8Dr9np3EIIUSLWXrfzKaUmK6W+t0Mka3NcUXCX9j2Jis/D7ABCVBRKKU/gfqBdoafOAm8BkUAVe4yttY5WSvkopTporaPtMYYQQggh/iEz0eKalFLBSqm1SqlzSqlkpdT4As/FK6XeUEodVUodVkr9q8BztZVS3yuljiulDiml/l3ovP9RSu23nHeNUqqh5fE3LG0ZeUV9BKeUelYpFaeUSlNK7VRKDbDhe73LMnZWwVkRpVSoUupgoWNnK6WeK/BQJ+Cc1jqx4HFa6yit9S/AqVLkqaOUOqKU6mr5up5SKkkp1bGIw6MAm/1dCCEqJqVUS8s196xSaqtSqnOB50KUUtstz71Q6HXFXtOvdt22fJ+YopQ6qJQ6pZR6pMBzPZRS6cAEIP/6m66Uql9M9vGWa2Ca5ftHL8vjrZRSUUqpM5b2i8FXef+Xtagopb5WSr2qlPoO8AcWWjI8WyB/ZEnGsRz7tFLqL8vf4VxHt/cJx5EiWpTEy8BOoCbQGlhX6PnuQEvgNuC/Sik/y+PfAdlAIMYs7GSlVHsApdRw4BngdqAW8BLgCaC1Hm9pyzhcTJ7TQH+gOvAC8KtSqm6Z36Ux9lzL2IVbI9YCVZRSIZb8lYGBwP8KHHMjsN8WOQrkSQUeAr5WSlUFPgU+1VpvLuLwvcBNthxfCFGxKKUUMBdYBNQB3sO4hnpbnpsNfAvUszxfULHX9BJctzthXCMfAt5USnlYXremwDV3rtbax/LrRBHZWwCPAx201tWAvkC85VPAhcCfQH1gDDDLcnyJaa3/Zck/yJLhzULjl3ScO4F+QJDlPd9vTQ5RfkgRLUoiD2gINNJan9ZaF75p7zOt9TlLG8EGoJ9SqhFGofuU1vqC1joWo+AcannNA8B0rfUurXWu1nqF1vpIScJorT/XWu/Xht+AMxjFvd1orfOAnzF+UADoDcRqrQvOTtcC0uww9p8Y3/DWAI2AqcUcmobxg44QQhQnCGiFcf3N0VrnF8adgKZAM+BDrXUW8G7+i0pwTb+WmVrrDOB3oBrG9xRr5QKVgdZKKU+tdbwlR2fAB3hda52ltV4B/AYML8UYV1PScd633KNzCqPoDrFxDuEkpIgWJfEscA7YYvk47o5Cz58o9OcGQP5s9CHLx15nMGYgGlge9wPiSxNGKfUvSxvHKct562GZxbazn/iniL7N8nVBpzG+OdjDp0BbjG9EucUcUx3jBwohhChOfeCM1jqzwGPHMIraepbnsi2PHy9wzLWu6ddyCsBSnANY3eJgmbR4ApgMnFBK/aCUagw0Bo5YJjvyJQBNrB3jGko6zrECf87AKLxFBSRFtLgmrXWS1voBjAvIm8BnhQ4peBGtj1FIHwEuAnW01jUtv3wKrHhxBOMjQasopQKAL4BH88+NUTiqAodlAe7WnrsE1gI1lVI3AoO5soj+C7Dq48OSUEq5A58D7wNT8nvHi9DKkkEIIYpzAuM6VrnAYw0xCub85/InJQpe2691TS+rvGsfAlrr2Vrr7kAAoIE3gKOAn1KqYE3jDyQVc5oMLr/Ju+A1VV9leGvHERWcFNHimpRSQ5VSjbXWGuMCU7hl4d9KqepKqQ5AF2CJ1joZWAW8rpSqqpTyVEp1VUrl9+x+DYxTSl2vlHJXSoUV6KW+mqqWDCcAd6XU01zZwrAfaKWUqlfEe2limU0fVsK3f4ll9uEXjGL2cKFWDoDNQA2llG+hMd0tN5a4WzJ75fcDFjpOK6UmFzH0BIzZof8D/gt8ZeldLKwn8IeVb0sI4VoOYVwjxymlPJRS92K0SGyyPBcDPGYppJ/Mf1EJrulldQxoWdS1MZ9SqoVSKsLyA8BF4AJGi8cm4DzwrCVXT2AQ8EMxp9oBjLBcm/sBYQWeO47R1lIUa8cRFZwU0aIk2gNbLXdQPw7cV+j5tcA+jALzMa11guXxezFmpg8CKRgzBvk3k8wG3rG85gxGn2+25aKWbhmr4F3SD1le9zcwHaNgTcZon7jsRhat9VbgG+Cg5bUFb47xxOj5q17Kv4ufMC64hWeh8z+m/Mryvgv6F8bF/jnLcxeAiQUPUErlz4qkFHq8I8bfef5sz2tAXctjhY9LK+aGQyGEAMAyGXI3xo3RqcDTwO1a6wzLc/dg3Ah3kkLXI4q5pl/rul1CP2K0DSYppRKLmgTBKPZft2Q7ZskywXLtHYzRs30S+BgYqbXeV8xY/4dR/J6xvN95BZ6bBky0tKw8XfBFpRhHVHDK+DcjROkopeKBUVrrZWZnsSWl1CtAE631g1a+ri7GDxVttdYXrHhdKMYNKn5a67NWhTVe/zPwpdZaZqKFEEIIB5DNVoQoxNKD3BVYYu1rtdYnMZb7s1YX4PPSFNCWcW+79lFCCCGEsBUpooUowDKzXg/YyJU3UNqN1voNR40lhBBCiLKTdg4hhBBCCCGsJDcWCiGEEEIIYSUpooUQQgghhLCS0/ZE161bVwcGBpodQwghSmXr1q0ntdZFLdNVYcl1WwhRXpXmmu20RXRgYCBbtmwxO4YQQpSKUirh2kdVLHLdFkKUV6W5Zks7hxBCCCGEEFaSIloIIYQQQggrSREthBBCCCGElWzSE62U6gfMANyBL7TWrxdxzJ3AZEADO7XWI2wxthBCCCFEeZednU1iYiIXL140O0qF5uXlha+vL56enmU+V5mLaMsWyR8BvYFEIFoptUBr/XeBY4KB54FuWuvTSqn6ZR1XCCFE6cjEhxDOJzExkWrVqhEYGIhSyuw4FZLWmtTUVBITEwkKCirz+WzRztEROKi1jtNaZwE/ALcWOubfwEda69MAWusTNhhXCCGElQpMfPQHWgPDlVKtCx1TcOKjDfCEw4MK4WIuXrxInTp1pIC2I6UUderUsdlsvy2K6CbAkQJfJ1oeK6g50FwptU4ptdEyCyKEEMLxZOJDCCclBbT92fLv2BZFdFFpdKGvPYBgoCcwHPhCKVXzihMpNVoptUUptSUlJcUG0YQQQhRi04kPuW4LIVyVLYroRMCvwNe+wNEijpmvtc7WWh8C9mMU1ZfRWn+mtW6vtW5fr55LbfQlhBCOYrOJD5DrthDCddmiiI4GgpVSQUqpSsDdwIJCx8wDwgGUUnUxZjnibDC2EEII69hs4kMI4TouXLhAWFgYubm5RT6flZVFaGgoOTk5Dk5mnjIX0VrrHOBxYAmwF/hRa71HKfWyUmqw5bAlQKpS6m9gJfCM1jq1rGMLIYSwmkx8CCGsNnPmTIYNG4a7u3uRz1eqVIlevXoxd+5cByczj002W9Fa/6G1bq61bqa1nmp5bJLWeoHlz1pr/ZTWurXW+gat9Q+2GLc8ys3LZXPSZjYnbeZY+jHydJ7ZkYQQLkQmPqyjtWZvyl7WJKzh0OlDZOZkmh1JCLvZtWsX3bp1u/T1tm3biIiIAGDWrFnceus/9yCHh4ezdOlSACZOnMjYsWMZMmQIs2bNcmxoE9lksxVxbVHxUczZNYd5++dx4vw/N7pXdq/MLc1v4bVer9G8TnMTEwohXIXW+g/gj0KPTSrwZw08ZfnlkjYnbea7nd/xW8xvxJ+Jv+y5GxvcyMQeE7mt9W24Kdn4V9jeE4ufYMexHTY9Z0jDEN7r995Vj2nTpg2xsbHk5ubi7u7OuHHjmD59OllZWcTFxREYGHjp2ClTpjBp0iROnDjB9u3bWbDA+EArOjraprmdmRTRdnYx5yJPLXmKT7Z8gk8lHwYED+DWFrfiU8mHw2cPcyD1AF/t+Ip5++bx8M0PM7nnZOpVlZtzhBDCDJk5mUxaOYm31r+Fl4cXkU0jea7bcwTVCiLpXBJHzh3hh90/cOf/7uTGBjfyeq/X6R/c3+zYQtiEm5sbbdq0Yc+ePcTExODv70+7du04evQoNWtefm9xaGgoWmveeecdoqKiLrV5VKpUibS0NKpVq2bGW3AoKaLtKPZULHf8dAfbj23nma7P8HL4y3h5eF1x3As9XuDlVS/z6dZP+ePgH6wYuYKgWmXfSUcIIUTJ7Tmxh3t+uYedx3cyut1o3u7zNtUqX1kIvNDjBebumcvkqMkMmD2At3u/zbiu40xILCqqa80Y21Pnzp1Zt24dH3/8MYsXLwbA29v7ig1Kdu3aRXJyMnXr1r2sYM7MzMTL68papyKSz6HsZNfxXdz82c3En4lnwd0LeLP3m0UW0AANfBrw0cCP2PDQBs5ePEvo16HEpMY4OLEQQriuncd20nVmV46mHWXB3Qv4dNCnRRbQAO5u7oy4YQS7Ht3FnW3u5OmlTzNuyTi5x0VUCJ07d2bixIkMHTqUJk2MJeRr1apFbm7upUI6OTmZe+65h/nz51O1alWWLFkCQGpqKvXq1cPT09O0/I4kRbQdHE8/zi1zbqFqpapsHb2VQS0Gleh1HZp0YOV9K7mYc5Gwr8PYm7LXzkmFEELEnY6j36x+VK9cnS2jt5T4ml3ZozJzbpvDmI5jeGfjO9w37z4ppEW517JlSypXrsz48eMve7xPnz6sXbuWjIwMhg0bxvTp02nVqhUvvvgikydPBmDlypUMGDDAhNTmkCLaxi7mXGTI3CGknE9hwd0LrG7LuKnhTay6fxUaTb9Z/Thz8YydkgohhDhx/gR9v+9LVm4WS+5dgn8Nf6te76bcmNFvBq+Ev8L3f33Pa2tes1NSIRxjxowZTJs2japVq172+OOPP84333xDlSpV2LBhA7179waM3ugNGzYAMHv2bEaPHu3wzGaRItqGtNY8tOAhNiZu5Pth33Nz45tLdZ7W9Voz7655HE07yuiFozFulBdCCGFLWblZ3DL7FpLOJfHb8N9oXa91qc6jlOKFHi9w7433MmnlJJYcXGLjpELYX2xsLC1btuTChQvcd999Vzzftm1bwsPDr7rZypAhQ2jRooW9ozoNKaJt6Nud3zJ712xeDX+VYa2GlelcnXw78Ur4K/z090/M3D7TRgmFEELke23Na0Qfjea7od/Rxa9Lmc6llOLTWz7lhgY3MOKXEVcsiyeEs2vWrBn79u3jyy+/LPaYBx988KqbrYwcOdJe8ZySFNE2kpqRytNLn6arX1ee7/G8Tc75bLdniQiKYOzisew7uc8m5xRCCAE7ju1g6pqp3HvjvdzW+jabnLOKZxV+vvNncvNyue3H28jKzbLJeYUQzkmKaBsZv2w8py+c5r8D/2uzxffdlBvfDf0Obw9vRv46Um5YEUIIG8jKzeL+efdTt0pdZvSbYdNzX1f7OmbeOpNtydt4b6N5y5QJIexPimgbWHt4LV9u/5KnujzFDQ1usOm5G1drzIx+M4g+Gs33f31v03MLIYQrem3Na+w8vpNPb/mU2t61bX7+Ya2GMbjFYKasmsLhs4dtfn4hhHOQIrqMsnOzeeS3R/Cv4c9LYS/ZZYzhNwynY5OOTFg+gfNZ5+0yhhBCuILYU7FMXTOVe264h8EtBtttnBn9ZqC15sklT9ptDCGEuaSILqOZ22eyJ2UP7/d7n6qVql77BaXgptx4p887JKUlMX3DdLuMIYQQrmBS1CQquVfird5v2XWcwJqBTAydyC97f2FRzCK7jiWEMIcU0WWQnZvNtLXT6Ozb2a4zGgDd/LtxR+s7eGPdGxxNO2rXsYQQoiLacWwHs3fN5v86/R+NqjWy+3jjuoyjeZ3mjFk0hos5F6/9AiFEuSJFdBl899d3JJxN4MXQF1FK2X281yNfJycvhxdWvGD3sYQQoqJ5YcUL1PKqxbPdnnXIeJU9KvN+v/eJPR3Ll9uKXzZMCFE+SRFdSjl5OUxdM5WbG91M/+v6O2TMprWa8niHx/l257fEnop1yJhCCFERrElYwx8xf/Bc9+eo6VXTYeP2adaH7v7deW3tazIbLUQp/frrr2ZHKJIU0aU0e9ds4k7HMSlskkNmofON6zoODzcP6Y0WQogS0lrz/PLnaeTTiMc7Pu7QsZVSTOk5haNpR/li2xcOHVsIW7lw4QJhYWHF7lYIxo6FoaGh5OTk2HTs2NhYli1bZnqOokgRXQq5eblMXTOVmxrcxKDmgxw6duNqjfnXjf/iqx1fceL8CYeOLYQQ5VFUfBTrjqzjxdAXqeJZxeHjhweG08O/B9PWTpPZaFEuzZw5k2HDhhW7WyEYOxb26tWLuXPn2nTsJUuWMGbMGNNzFEWK6FL4dd+vHEg9wMTQiQ6dhc73TNdnyMzJ5P1N7zt8bCGEKG/e3fgudavU5f6Q+00ZXynF5J6TZTZalAsRERGEhIQQEhKCl5cXP/30E7NmzeLWW2+9dEx4eDhLly4FYOLEiYwdOxaAIUOGMGvWrFKNu2vXLrp163bp623bthEREYGHhwctW7YEcEgOa3jYfYQK6KPojwisGcjQlkNNGb9F3RYMaTmEj6I/4rnuz+FTyceUHEII4exiUmP47cBvTAydiLent2k5wgPDCQ0IZdraaYxqNwovDy/Tsgjn98QTsGOHbc8ZEgLvlWATzRUrVgDwySefsHLlSgYPHsyYMWMIDAy8dMyUKVOYNGkSJ06cYPv27SxYsACA66+/nujo6MvO16NHD9LS0q4Y5+233yYyMvLS123atCE2Npbc3Fzc3d0ZN24c06dPp127doDRphEXF1fqHPYgRbSV9qbsJSo+imm9puHuVvzHCfY2vtt4ft33K59v/Zwnu8hi/kIIUZQZm2bg6e7Jfzr8x9QcSileCnuJXt/2YtZfs3io3UOm5hHiar799lsWLVrEzz//TEpKCjVrXn4zbmhoKFpr3nnnHaKioi61V7i7u1OpUiXS0tKoVq0aAGvWrCnRmG5ubrRp04Y9e/YQExODv7//pQIa4OTJk2XKYQ9SRFvp062f4unmyYNtHzQ1RyffToQFhPHOxncY02kMHm7yn1IIIQo6feE0X+34iuHXD6ehT0Oz4xAeGM6NDW7k/c3v82DbB01pBxTlQ0lmjO0lv31j/vz5eHp64u3tzcWLl/fy79q1i+TkZOrWrXtFkZqZmYmX1z+ftJR0Jhqgc+fOrFu3jo8//pjFixdf9lxZc9iD9ERbISM7g292fsPtrW+nftX6Zsfhyc5Pkngukd8P/G52FCGEcDqfb/ucjOwMnuzsHJ/WKaUY23Esfx3/i9UJq82OI8QVfvvtNz7++GN++eWXSwVorVq1yM3NvVTAJicnc8899zB//nyqVq3KkiVLLr0+NTWVevXq4enpeemxNWvWsGPHjit+FS6gwSiiJ06cyNChQ2nSpMllz5U1hz1IEW2FH3b/wJmLZ3ik/SNmRwFgYPOBNPJpxOfbPjc7ihBCOJXs3Gw+2PwB4YHh3NTwJrPjXDLihhHU9q7N+5vlxnDhfO677z4SExPp1q0bISEhfPmlsUlQnz59WLt2LRkZGQwbNozp06fTqlUrXnzxRSZPnnzp9StXrmTAgAGlHr9ly5ZUrlyZ8ePHF/m8o3KUlPQAWOGTLZ/Qpl4bevj3MDsKAB5uHjzU9iFeW/sah88exr+Gv9mRhBCChQuhUSNo3968DL/H/E7iuUQ+7P+heSGK4O3pzeh2o3lz/ZsknEkgoGaA2ZGEuCQ1NbXIxx9//HHeeecdIiMj2bBhw6XHQ0NDL/t69uzZTJs2rdTjz5gxg2nTplG1alVTc5SUzESX0JajW9hydAuPtH/EqfrYHmr3EFprZm6faXYUIYRAaxg2DDp0gOHDIS7OnBxf7fiKhj4NGdh8oDkBruLRDo+iUHwc/bHZUYQokbZt2xIeHn7NTU6GDBlCixYtrD5/bGwsLVu25MKFC9x3332m5bCWFNEl9PWOr/Hy8OJfN/7L7CiXCawZSN/r+vLl9i/JybP/7jxCCHE1OTnGrxtvhPnzoWVLeOsto7h2lGPpx/j9wO+MvHGkU9507V/Dn6Gthl7q2RaiPHjwwQevucnJyJEjS3XuZs2asW/fvkvtI2blsJYU0SWQnZvNj3t+ZFDzQdTwqmF2nCuMbjeaxHOJLD64+NoHCyGEHWVnG7+PGAExMTBoEDz7LIwbB3l5jsnw3c7vyNW5PND2AccMWApjO47l9MXTzN1t/13VhBD2IUV0CSyLW0ZKRgr33HCP2VGKdEvzW2jo05DPtn5mdhQhhIvLL6IrVYImTeCnn2DsWHj3Xbjvvn+etxetNV/t+Iouvl1oWbelfQcrg+7+3WlRpwVf7fjK7ChCiFKSIroEZu+eTU2vmvS7rp/ZUYrk6e7JAyEP8HvM7xxNO2p2HCGEC8vKMn7PX1nKzc1Y83bqVPj+e7jnHvu2dmxK2sTek3tNX8v/WpRS3B9yP2sOryEmNcbsOEKIUpAi+hoysjP4de+v3N7qdip7VDY7TrFG3jSSPJ3Hj3t+NDuKEMKF5c80F1yeVSmYMAFef92YmX7rLfuN/9X2r6jiWYU729xpv0FsZORNI3FTbny942uzowgnoR1584CLsuXfsRTR17Bg/wLOZ5/nnhuds5UjX8u6LWnXqB2zd802O4oQwoUVbOco7Nln4Y474PnnYfly24+dkZ3BnN1zuL317VSvXN32A9hY42qN6dusL9/s/IbcvOJXGxCuwcvLi9TUVCmk7UhrTWpqqs12MnS+25adzOxds2lSrQmhAaFmR7mmEdeP4OmlTxOTGkNwnWCz4wghXFDhdo6ClIKZM2HPHrjrLti6FQJsuEzygv0LSMtK4/6b7rfdSe3sgZAHuPN/d7Isbhl9r+trdhxhIl9fXxITE0lJSTE7SoXm5eWFr6+vTc4lRfRVpGaksujgIp7o9ARuyvkn7e+6/i6eWfoMc3bPYVLYJLPjCCFcUFHtHAX5+MCvvxrrSI8YAWvWGH3TtvDjnh9p5NOoXEx65BvcYjC1vWvz1Y6vpIh2cZ6engQFBZkdQ1jB+StDE/2892dy8nIYccMIs6OUiG91X8ICw5i9a7Z8HCSEMMXV2jnyNW8O778P69fDF1/YZty0zDT+iPmDO1rfgbtb8WvIOpvKHpW554Z7mLdvHqcvnDY7jhDCClJEX8Uve3+hWa1mhDQMMTtKiY24fgT7U/ez/dh2s6MIIVzQ1do5Cho5Enr2hPHj4fjxso+7YP8CMnMzy8UNhYU9EPIAmbmZzN0ja0YLUZ5IEV2MsxfPsuLQCoa2HOpU23xfy22tb8PTzVNuMBRCmOJa7Rz5lIJPPoHz542NWMrqx79/pEm1JnTx61L2kzlYSMMQWtVtxQ+7fzA7ihDCClJEF2PRwUVk52UzpOUQs6NYpbZ3bfoH92fO7jlyt7cQwuFK0s6Rr2VLeO45mDULli0r/ZhnL55l8cHF3NH6jnJx/0phSinuanMXqxNWy1r/QpQj5e9q4yDz9s2jftX6dPbtbHYUqw2/fjhH046y/sh6s6MIIVxMSds58k2YANddB2PGQE5O6cZcsH8BWblZ3HX9XaU7gRO46/q70Gh+2vOT2VGEECUkRXQRMnMy+SPmDwY3H1yublDJNzB4IJXcK/Hrvl/NjiKEcDElbefI5+UFb7wB+/bB7FJ2oc3dMxf/Gv50atKpdCdwAi3rtiSkYQg/7JGWDiHKCymiixAVH0VaVlq5a+XIV61yNXo37c0ve3+RVTqEEA5lTTtHvqFDoW1bmDLln9eX1OkLp/kz9k/uaH1Hubp/pSh3t7mbjYkbOXT6kNlRhBAlIEV0Eebtm0dVz6r0atrL7CilNqzVMBLOJrDj2A6zowghXIi17Rxg3GT46qsQFwdff23deAsPLCQ7L7tcrspRWP57+HHPjyYnEUKUhBTRheTpPObvn0//4P54edhmW0gzDGo+CDflJi0dQgiHsradI1///tClC7z8Mly8WPLXLdi/gMbVGtO+cXvrBnRCQbWC6NSkk7R0CFFOSBFdSHRSNMnpyQxpUT5bOfLVq1qPHv49+GXvL2ZHEUK4kNK0c8A/s9GJifD55yV7TWZOJktil1yaNKgI7r7+bnYc28H+k/vNjiKEuIaKcdWxofn75+Ph5sGA4AFmRymzYa2GsSdlDwdSD5gdRQjhIkrTzpEvIsLYgGXaNMjMvPbxUfFRpGelM7jFYOsHc1J3tL4DhZKWDiHKASmiC1l0cBHd/LpRy7uW2VHKLP/GyF/3SkuHEMIxStvOke/55yE5GeaWYPO+BfsXUMWzChFBEaUbzAk1qd6Ezr6dpRVPiHJAiugCktOS2XFsB/2u62d2FJvwr+FP+8bt5WIshHCY0rZz5OvdG9q0gXfegastLqS1ZsGBBfRt1rdc379SlKEth7L92HYSziSYHUUIcRVSRBewJHYJQIUposG4GG9K2kTSuSSzowghXEBZ2jnA6I1+8knYuROiooo/bsexHSSeS2RQ80GlG8iJ5X+KOG/fPJOTCCGuRoroAhYfXExDn4bc1OAms6PYzK0tbgXg95jfTU4ihHAFZW3nALjnHqhXz5iNLs6C/QtQKAY2H1j6gZxUcJ1g2tRrI58iCuHkpIi2yM3L5c/YP+l3Xb9yv2B/Qa3rtSagRoAU0UIIh7BFEe3lBf/5D/z2Gxwo5r7ohQcW0sWvC/Wr1j2BcbYAACAASURBVC/9QE5saMuhrDm8hpTzKWZHEUIUQ4poi+ij0Zy+eJp+zSpOKweAUoqBwQNZFreMizlWLL4qhBClkJUF7u7gVsbvLo8+CpUrw3vvXflc4rlEtiZvZXDzirMqR2FDWw0lT+ex8MBCs6MIIYohRbTF4oOLcVNuRDaNNDuKzQ1sPpCM7AxWJ6w2O4oQwkkopfoppfYrpQ4qpZ67ynG3K6W0UqpEu5lkZ5dtFjpfgwZGW8c338DZs5c/9/sB45O1QS0qXj90vrYN2+Jfw19aOoRwYlJEWyw6uIiOTTpSp0ods6PYXHhgOF4eXpe+8QghXJtSyh34COgPtAaGK6VaF3FcNWAssKmk587OLv3KHIU98ghkZMCcOZc/vjh2MQE1AmhVt5VtBnJCSimGtBjC0tilpGelmx1HCFEEKaKBkxkniU6Kpv91/c2OYhfent5EBEXwe8zv6KutGSWEcBUdgYNa6zitdRbwA3BrEce9ArwJlLgXLCvLNjPRAO3bw003wRdf/PNYdm42y+OW07dZ3wp1/0pRhrYaSmZuJosPLjY7ihCiCFJEA0tjl6LRFWppu8IGBg8k9nSs7F4ohABoAhwp8HWi5bFLlFJtAT+t9W9XO5FSarRSaotSaktKSorN2jmMc8OoUbB1K2zfbjy2MXEjaVlpFfp6na+7f3fqeNdhwf4FZkcRQhTBJkW0vXrrHGVJ7BJqe9fm5kY3mx3FbgYGG8tAySodQgigqCncSx9TKaXcgHeBcdc6kdb6M611e611+3r16tm0nQOMvmgvr39moxcfXIy7cq9QuxQWx8PNg37X9WPRwUXk5uWaHUcIUUiZi2h79tY5gtaa5YeW0yuoF+5u7mbHsZuAmgG0qddGimghBBgzz34FvvYFjhb4uhpwPRCllIoHOgMLSjIBYst2DoBateD222HWLKM/eknsErr4daGGVw3bDeLEBgYP5GTGSTYnbTY7ihCiEFvMRNutt84RDp46SOK5RJeY1RgYPJA1CWs4l3nO7ChCCHNFA8FKqSClVCXgbuBSz4DW+qzWuq7WOlBrHQhsBAZrrbdc68S2bOfIN2qUsULHzFnn2Jq8tcItRXo1fa/ri7tylwkQIZyQLYpom/XWmWHFoRUArlFENx9Idp5xU44QwnVprXOAx4ElwF7gR631HqXUy0qpMi2+bOt2DoDQUAgOhg8+MeZg+l7X17YDOLHa3rXp6tdVimghnJAtimib9dYVvkHFEVbEr6BJtSYE1w52yHhm6uLbBZ9KPvwZ+6fZUYQQJtNa/6G1bq61bqa1nmp5bJLW+oq72LTWPUsyCw22b+cA4wbDhx6CA9vrUyujA+0atbPtAE5uYPBAdhzbQdK5JLOjCCEKsEURbbPeusI3qNhbns5j5aGVRARFVPilkgA83T2JCIpgadxSs6MIISooe7RzAAwfkQcqD9+EZ3FTrrWw1C3NbwHkxnAhnI0trkR2662zt90ndpOSkUKvoF5mR3GY3k17E3s6lthTsWZHEUJUQPZo5wBI9dgJ/qtJ2dgbV1vuvnW91gTUCJAiWggnU+Yi2p69dfaW3w8dHhRuchLH6dOsD4DMRgsh7MIe7RyA0YZ24yyOJdRg2zbbn9+ZKaUYGDyQZXHLuJjjVPfmC+HSbPKZmL166+xtxaEVXFf7Ovxr+JsdxWGCawcTUCNA+qKFEHZhr3aO5YeW0zrsbypVMpa7czW3NL+FjOwMVsWvMjuKEMLCtRrLCsjJy2FVwioiAiv+qhwFKaXo06wPyw8tJycvx+w4QogKxh7tHJk5maw9vJbebTowYAD88APkutjeIz0De+Lt4c1vB5xukSshXJbLFtHbkrdxLvOcSyxtV1ifZn04l3lOFu8XQticPdo5NiZu5ELOBSKCIhgxApKTYeVK247h7Lw9vQkPCmdJ7BKzowghLFy2iM7vh+4Z2NPcICaICIpAoVgaK33RQgjbskc7x4pDK3BTboQFhHHLLVCtGsyebdsxyoO+zfoScyqGQ6cPmR1FCIELF9Er41fSpl4bGvg0MDuKw9X2rk2HJh34M076ooUQtmWPdo4V8Sto37g9Nbxq4O0Nt90GP/8MF13sHru+zYxNZmQ2Wgjn4JJFdE5eDusOr3PJWeh8fZr2YVPiJs5cPGN2FCFEBWLrdo70rHQ2Jm687P6V4cPh3Dn408XmAZrXaU5gzUApooVwEi5ZRG9P3s757PP08O9hdhTT9GnWh1ydy8pDLtZYKISwK1u3c6w9vJacvBx6Nf1nPf/wcKhVC/73P9uNUx4opejbrC/L45aTnZttdhwhXJ5LFtFrDq8BoEeA6xbRnXw7UdWz6qXecCGEsAVbt3Msj1tOJfdKdPXreukxT08YPBgWLDBmvl1J32Z9SctKY0PiBrOjCOHyXLKIXp2wmma1mtG4WmOzo5imknslegT0YEW8FNFCCNuxdTvHivgVdPHtQhXPKpc9fvvtcPYsrHCxS1hEUATuyp0lB6WlQwizuVwRnafzWHt4LaEBoWZHMV1EYAR/p/zNsfRjZkcRQlQQtmznOHXhFNuTt9MrqNcVz/XubazS4WotHTW8atDVryuLYxebHUUIl+dyRfTelL2kXkh16X7ofPlrZEtftBDCVmzZzhEVH4VGF7mef+XKMGgQzJsHOS62b1TfZn3ZlryNE+dPmB1FCJfmckV0fj+0zERDSMMQanrVlL5oIYRNaG38bquZ6JWHVlLVsyodmnQo8vnbboPUVFjlYjth973OWOpO1voXwlwuV0SvTlhNI59GNK3V1OwopnN3c6dnYE/pixZC2ISti+hVCavo6teVSu5FT2336wdVqhhrRruSdo3aUbdKXVnqTgiTuVQRrbVmdcJqQgNCUUqZHccpRARGEHc6jvgz8WZHEUKUc/lFtC3aOU5dOMWuE7sICwgr9pgqVWDAAPj1V8jNLfuY5YWbciOyaSTL4pah8//ShRAO51JFdPyZeJLSkqQfuoD8XkNp6RBClJUtZ6LXJJSs9e622+DYMdjgYiu+RQZFkpyezN8pf5sdRQiX5VJF9OqE1YD0QxfUul5r6letL0W0EKLMbFlEr05YTWX3ynRs0vGqx/XvDx4esHBh2ccsTyKbRgKwLG6ZyUmEcF0uVUSvObyGWl61aFO/jdlRnIZSioigCFYcWiEfCwohysSW7RyrD6+ms29nKntUvupxNWpAz56uV0QH1AwguHYwyw5JES2EWVyqiF57eC3d/LvhplzqbV9TRGAEyenJ7E/db3YUIUQ5ZquZ6HOZ59iWvK3EnxoOGgR790JMTNnGLW8im0YSFR8lW4ALYRKXqSZPZpxkf+p+uvl1MzuK0wkPCgekL1oIUTa2KqLXH1lPns676k2FBQ0aZPzuarPRkU0jSc9KZ3PSZrOjCOGSXKaI3pi4EYCufl1NTuJ8mtVqRpNqTViV4GKLrQohbCovz/i9rO0cq+JX4eHmQWffziU6PigIbrjB9Yro8MBwFEr6osup8+dh1iyYMAHuugs6dTKWbRw7Fj780PU+WSmPXKaIXn9kPR5uHrRv3N7sKE5HKUVYYBir4ldJX7QQotRsNRO9+vBqOjTuQNVKVUv8msGDYc0aOHWqbGOXJ7W8a9G+cXuWxsmmK+VJXByMGwe+vnDvvfDWW7BtG1SvDidPwldfwZgx0Ly5sYTj4sX//IAqnIvLFNEbEjcQ0jCEKp5VzI7ilHoG9OT4+eMcSD1gdhQhRDlliyI6IzuD6KRoq1dRGjTIWCt60aLSj10eRTaNZGPiRs5lnjM7iriGnByYMsUojt9/H/r2hdWr4cIFY9Z56VLYsgXOnYOEBOPY7duNFWi6d4fYWLPfgSjMJYro7NxsNidtpquvtHIUJyzQ6D2Mio8yN4gQotyyxeocGxM3kp2XbXUR3aEDNGjgei0dkU0jydW5l5ZwFc4pLg5CQ2HyZBg+3CiSf/gBevQwlmgsSCnw94dJk4zjvvwS/v4bbrrJ+LN8YOw8XKKI/uv4X2RkZ0g/9FUE1w6moU9D6YsWQpSaLWaiV8Wvwk25WX0TuJubMRu9aBFkZZV+/PKmq19XvDy8pC/aia1cCSEhRiE8ezZ89x00blyy11aqBA8+CLt2QceOMGqU0QKSLQuyOAWXKKLXH1kPyE2FV6OUIiwgjFUJ0hcthCgdWxTRa4+s5aYGN1HDq4bVrx00yPgofLULTcp6eXjRw7+HFNFOasUKGDgQ/Pxg505jFro0/Pxg2TJ49VWjEL/jDsjMtG1WYT3XKKIT1+Nb3Re/Gn5mR3FqYQFhHE07SuxpabwSQlivrO0c2bnZbEzcSHf/7qV6fa9extiu1hcdERTBnpQ9HE8/bnYUUcDy5UYB3bSpMRsdEFC287m5wQsvGCt3zJ8Pt94KGRm2ySpKxyWK6A1HNtDFt4vZMZxefl/0qnhp6RBCWK+sM9E7j+8kIzuj1Ov5V60KYWGuWUSD3NPiTNavh1tugeBgo4CuX992537sMaM3+s8/YcgQae0wU4UvopPOJZFwNkFaOUqgVd1W1KtST/qihRClUtYiet3hdQB08y/9plgDBhi7F8bHl/oU5U67Ru2oXrm6bJjlJI4dg9tvhyZNjNnoevVsP8aDD8Lnnxsrejz1lO3PL0rG49qHlG8bEjcA0g9dEpfWi5YiWghRCmVt51h7ZC0BNQLwre5b6gz9+8OTTxqz0Y8+WurTlCsebh6EBoSyMn6l2VFcXnY23HknnDljrO9sjwI630MPGTcrvvOOsdnQ6NH2G+tq8nQe25O3szRuKcsPLSc5LRmNcTFoXK0xEYER9Grai3aN2uHhVrHKzgo/E73+yHq8PLwIaRhidpRyISwgjMNnDxN/Jt7sKEKIcqYsM9Faa9YdXlemWWgw1uANCnLBlo7ACGJOxXDk7BGzo7i0Z54xNv354gu48Ub7j/fmm8Yuh489BqscPP+Vm5fL1zu+JmhGEO0/b8/zy5/nxPkTtKjbglZ1W9GqbiuOpx9nwooJdPqiE01nNOXzrZ+TnVtx+k8q1o8ERVh/ZD0dGnegknsZ96F1EWEB//RFB4YEmhtGCFGulKWIjj8TT3J6Mt39SndTYT6ljJaOr76CixfBy6tMpys38vuiV8avZORNI01O45p++QVmzID/+z8YMcIxY7q7w5w50LmzMQO+e7d9Z7/zLYpZxDNLn2FPyh46NO7A1IipRDaNpKFPwyuOPZ5+nBWHVjBj0wxG/zaaN9a9wWu9XuPONnfaP6idVeiZ6MycTLYf205n385mRyk32tRvQ23v2rJwvxDCavlbE5emnWPt4bVA2fqh8/Xvb6xasGZNmU9VbtzQ4AbqeNeRvmiTnDpltA+1a2ds4+1INWvC//5ntJCMHm3fzVhy8nJ4dumzDJg9gOy8bH664yc2jdrEvTfeW2QBDdDApwHDbxjOhoc2sHD4Qnwq+XDX/+7i3wv+zYXsC/YL6wAVuojecWwHWblZdGrSyewo5YabcqO7f3fWHHah7z5CCJsoy0z0uiPrqF65Om3qtSlzjvBwqFwZ/vijzKcqN9yUG+FB4aw4tELW+jfBU08ZhfTMmWVbJ720rr8epk6FefPg22/tM0ZqRir9Z/XnrfVv8Wj7R9n16C5ub307SqkSvV4pxS3Nb2Hr6K1M6D6BL7Z/QbeZ3Yg9VX6X1a3QRfSmpE0AdPKVItoaPfx7EHMqhmPpx8yOIoQoR8paRHf164q7m3uZc1SpAj17umZf9JFzR4g7HWd2FJfy55/wzTcwfryxNbdZnnzS2EZ8zBhju3BbOnL2CB2/6MjqhNV8MegLPh74canbZN3d3JnaayoLhy/k0JlDdPi8A9uSt9k2sINU+CK6cbXGZbrT2xX18O8B/PPxqhBClERpi+jTF06z+8TuUq8PXZT+/WH/fohzoXoyvy9aWjocJz3daKFo2RImTjQ3i7u7UcxrDQ88YLu2juPpx4n8LpKTGSdZdf8qHmr3kE3Oe0vzW9jy7y1Uq1yNyG8j2Z683SbndaSKXUQnbpJWjlJo16gdVTyrsCZBWjqEOS7mXGTJwSU88+czPLzwYUYvHM3ohaN5f9P7HEg9IB9XOymtjW/kblZ+Z8lfitTWRTQYy4y5iuZ1mtPIpxEr4qWIdpRXXoHDh43VOJzhJtagIJg+3djgZdassp/v9IXT9Pm+D4nnEvljxB82v8esWe1mrLxvJT6VfIj8LpIdx3bY9Pz2VmFX5ziZcZLY07H8u92/zY5S7ni6e9LZt7P0RQuH23lsJ5NXTWbJwSVcyLlAZffK1PKuBRjLKX2+7XMAgmoG8WDbB3mqy1NU8axiZmRRgNalbOU4vA4PNw86NulosyzBwRAYaGxG8Z//2Oy0Tk0pRURQBMvilqG1LnGvqiid+HhjNY6RI6Gb7X7+K7NRo4wdDZ9+GgYNgho1SneeC9kXGDB7APtO7uO34b/Z5KbfojSt1ZSo+6MI+zqMyG8j2fzvzTSt1dQuY9lahZ2J3py0GZB+6NLq4d+Dncd3ci7znNlRhAtIOZ/CI789QrvP2rEmYQ2j2o1i0T2LOD3+NMnjkkkel8yJZ04QOzaWjwd8THCdYF5c+SLNP2jOtzu/JU/nmf0WBEYRXZqVOdYdWUfbhm2pWqmqzbIoBX36wIoVrrUtcs/Anhw/f5z9qfvNjlLhTZhg/H/26qtmJ7mcmxt89BGcOAEvvVT684xZNIaNiRuZc9scejfrbbuARWhaqykrRq4gT+dx6w+3kp6VbtfxbKXCFtGbEjfhptxo37i92VHKpR7+PcjTeaw/st7sKKKCWxSziOYfNufL7V8ytuNYYsbE8H7/9+l3XT+8Pb0vO7ZpraY82uFRlty7hFX3r6JRtUbcN+8+Ir+NJDUj1aR3IPKVZiY6Jy+H6KPRdPHtYvM8ffrAuXOwebPNT+20wgPDAYiKjzI3SAW3ebOxPvO4ceDrhLddtW8PDz8MH34If/1l/eu/2fENX27/kgndJzCs1TDbByxCcJ1g5t4+l79T/mbkryPLxeRIxS2ikzbRpl4bfCr5mB2lXOrs2xkPNw/pixZ29d8t/2XQnEEE1gzkr0f+4t1+715q37iW0IBQNo3axOeDPmf9kfV0/KIje07ssXNicTWlKaL/Ov4XGdkZdPGzfREdEWHMyi1ZYvNTO62mtZriW91XtgC3I62NVon69Y0VOZzV1KnGGtKPPWbdTYa7ju/i0d8fJTwwnCnhU+wXsAi9m/Xm7d5v8+u+X3l1tZNN8RehQhbRWms2J22WmwrLoGqlqrRr1E76ooVd5Ok8nvnzGR79/VH6XteX1fevplW9Vlafx025MardKKLuj+J81nm6fNmFRTEutq6ZEylNO8eGI8ZNhfaYia5VCzp2NJYgcxVKKXoG9iQqPkpuwLWT+fONjXymTIFq1cxOU7zateG112DtWmP96JJIz0rn9p9up4ZXDWbfNhsPN8ffOvdE5ycYedNIXop6iWVxyxw+vjUqZBEdcyqG0xdPSz90GfXw78HmpM1k5mSaHUVUMM8ufZa3N7zNYx0eY/7d86lWuWzfiTr7dib639E0q92MIXOHsPKQzMKZoTQz0RsSN9DIpxH+NfztkqlPH4iONjbCcBU9A3py4vwJ9p3cZ3aUCicvz+gzbt7cuIHP2T34oLH83oQJkJNz7eNfWP4CMakxzLltTrE7ENqbUopPBn5CizoteGD+A5y5eMaUHCVRIYvoSzcVykx0mfTw70FmbibRR6PNjiIqkA82fcD0DdN5vMPjfND/A5vNdPjV8GP5yOUE1w5m8A+D2XJ0i03OK0qutEV0F78udltJok8fo/BZ4UKrvoUHSV+0vSxcaPQYv/ACeJSD9c08PIy2jn37jDWkr2bDkQ18sPkD/tPhP/QM7OmQfMWp4lmF74Z+R3JaMo//8bipWa6mQhbRmxI34VPJh9b1WpsdpVzr7t8dQPqihc3M3zef/1v8f9za4lbe6/eezQun2t61WXLvEup416H/rP7sPykrFDhSXp517Rwnzp8g7nScXVo58nXqBNWru1ZLR1DNIPyq+xGVEGV2lApFa2Nd6KZNYcQIs9OU3NChxr+Dl16CCxeKPiYzJ5NRC0fhW92Xab2mOTZgMTo06cCLoS8ya9csftrzk9lxilQxi+ikTbRv3N4m28e6sjpV6tC6XmvWHpGdC0XZbUvexvCfh9OhSQdm3zbbbv8+m1RvwtJ/LUWhGDB7gFN/FFjRWDsTbc9+6HweHtCrl1FEu0qLsPRF28fixbB1q9EaUR5mofMpBa+/DklJxmodRZm2dhp/p/zNJwM/KXN7nS1N6DGBDo078Mjvj3As/ZjZca5Q4YrozJxMdhzbIa0cNtLNrxvrj6wvF0vNCOeVkZ3BiJ9HUNu7NguHL7T7BinBdYKZd/c8Dp89zEMLHpJCwkGsLqITN+Dp5snNjW+2XyiMlo6EBIiJseswTqVnoNEXvffkXrOjVAj5s9D+/vCvf5mdxno9exq7eE6bBmfPXv7cvpP7eG3Nawy/fjgDmw80JV9xPN09+Xbot6RnpfP0n0+bHecKFa6I/uv4X2TnZdOhcQezo1QI3f27c+biGf5O+dvsKKIcG7dkHAdSD/Dt0G+pX7W+Q8bs6teVab2m8cveX/hwczHTL8KmrF2dY0PiBto2aouXh333S+7Tx/jdlVo68ntapS/aNlasgA0b4LnnSrehkDN49VU4fRo++ODyx5/+82m8Pb15r9975gS7hpZ1WzK+23hm7ZrldDeNV7giOv8muA5NpIi2hW5+xjafaw9LS4conQX7F/Dfrf9lXJdxRARFOHTscV3GMaj5IMb9OY7oJLlBNp9Sqp9Sar9S6qBS6rkinn9KKfW3UuovpdRypVRASc5rzUx0dm420Un22WSlsKZNISgIli+3+1BOI6hmEP41/KWItpE33oDGjY3VLsqrdu1g4EB4911ISzMeWxq7lN9jfmdij4kOm+Aojee7P09QzSD+88d/yMrNMjvOJRWyiK5ftT5+1f3MjlIhNK3VlIY+DVl3ZJ3ZUUQ5dCz9GA8teIiQhiG8GuH4hfOVUnw95GsaVWvE3T/fzfms8w7P4GyUUu7AR0B/oDUwXClV+C7s7UB7rfWNwP+AN0tybmuK6L+O/8WFnAsOKaLB6IteubJky3xVBNIXbTu7d8PSpTBmDFSubHaasnnxRWO5x08+gdy8XMb9OY7AmoGM6TTG7GhX5e3pzYcDPmTfyX1MXz/d7DiXVLwiOima9o3b2225JFejlKKbXzeZiRal8tSSp0jLTGPWsFlU9jDnu09t79p8N/Q74k7H8VLUS6ZkcDIdgYNa6zitdRbwA3BrwQO01iu11hmWLzcCJdrY2Jp2jg2JlpsK7bBTYVEiI41e0G3bHDKcUwgLCCMlI0X6osvovffA2xtGjzY7Sdl16gS9e8P06fDphu/YdWIXb0S+YfeWKlsYEDyAoS2H8srqV0g4k2B2HKCCFdHpWensPblX+qFtrLt/d+LPxJN0LsnsKKIciYqPYs7uOYzvNt705SZDA0J5+OaHeXfju7J+NDQBjhT4OtHyWHEeAkq0DaQ1M9EbEjfQuFpjh31qGGHpJFrm3Bug2VRYQBgAq+JXmZyk/EpJge+/h/vuM3YArAhefBFOnIDn3jSWl7yj9R1mRyqx9/q9h0YzYcUEs6MAFayI3pa8jTydJ0W0jeX3RUtLhyip7NxsHv/jcQJrBvJc9ytabk3xRuQbNKjagFELRpGdm212HDMV9TFdkZ/3K6XuBdoDbxV7MqVGK6W2KKW25OTklriI3pi4kc6+nR32qWG9enDTTa7VF920VlOaVGvCqgQpokvrv/+FzEwYO9bsJLbTowcEhiSQtuIRXgt9p1x9cu9fw58nOz/J7F2z2Xp0q9lxKlYRnT/D1L5xe5OTVCwhDUOo4lmFdYeliBYl81H0R+xJ2cN7fd/D29Pb7DgA1PCqwccDP2bn8Z1M3+A8PXUmSAQKTv/6AkcLH6SUigReAAZrrTOLO5nW+jOtdXutdXs3N/cStXOknE8h7nQcnZt0tjp8WfTqBevWFb/hREUjfdFlk5kJH38M/fpBq1Zmp7Gdc5nnOHnzE5DemIOrHPtv0BbGdxtP3Sp1eWbpM6b/f12hiujoo9H4VfejgU8Ds6NUKJ7unnRq0kk2XRElciz9GC9FvUS/6/oxuMVgs+NcZkjLIdzW6jamrJriND11JogGgpVSQUqpSsDdwIKCByil2gKfYhTQJ0p64pK2c2xK2gRAJ1/HrucfGWkURutcaD4gLCCM4+ePcyD1gNlRyp25c+HYMXjySbOT2NYHmz4g3XceLa4/z/Tpxk6j5UkNrxq8FPYSK+NX8kfMH6ZmqVhFdFK0LG1nJ938urHj2A7SMtPMjiKc3KSVk7iQfYH3+73vlB8TvtP3HQCeX/68yUnMobXOAR4HlgB7gR+11nuUUi8rpfJ/6nkL8AF+UkrtUEotKOZ0hc5dwiI6cRPuyp2bG9l3k5XCevQwdppzqb7oQEtftLR0WO3DD40Z6N69zU5iO+cyzzF9w3QGtRjEpOersm8f/GFuHVoqD9/8MMG1g3l22bPk5Jm35E6FKaJPXThF7OlY6Ye2k+7+3cnTeZdmkIQoSkxqDDO3z+SR9o8QXCfY7DhF8q/hz9NdnmbO7jlsTNxodhxTaK3/0Fo311o301pPtTw2SWu9wPLnSK11A611iOVXiT5SKOnqHJuSNnF9/eupWqlqWd6G1Xx8oEsX1+qLDq4dTEOfhlJEW2n7doiOhkceMbbNrig+2PQBpy+e5qWwl7jjDvDzg7ffNjuV9TzdPXk98nX+Tvmbb3d+a1oOmxTR9lq43xr5/dBSRNtHZ9/OKJQsdSeuavKqyVRyr8SEHs5x53RxxncfT0Ofhjy15CnTe+oqkpLMROfpPDYnbaazrzm9mL16wdatxlq5rkD6okvndTaBqwAAIABJREFU00/By6t8bvFdnEuz0M0HcXPjm/H0hCeegFWrjB8YypuhLYfSoXEHXln9imk3i5e5iLbnwv3WyN+N7ObGjv140FXU8KrBDQ1ukBU6RLF2n9jNnF1zGNtpLA19Gpod56p8KvkwNWIqGxI3MHfPXLPjVBglKaIPpB7gbOZZOjVxbD90vl69jJwrnWv3YLsKCwjjaNpRYk/Hmh2lXEhLg1mz4K67oFYts9PYzsfRH1+ahc43ahRUr14+Z6OVUkzpOYX4M/F8veNrUzLYYibabgv3WyP6aDTBtYOp6VXT1qcWFt38urExcSO5eblmRxFO6MWVL1KtcjWe7fas2VFK5L6b7iOkYQjjl43nYs5Fs+NUCCVp58hvoXH0TYX5OnWCqlVhxQpThjeFrBdtnTlzID0dHn7Y7CS2k5mTyYxNM+jTrM9lk43Vqxvv83//g0OHTAxYSv2u60enJp14dc2rpmwHbosi2mYL9xdcbzQlJcWqENFH5aZCe+vm1430rHR2ndhldhThZKKTopm3bx7juoyjtnf52JHA3c2dt3q/xeGzh/li2xdmx6kwrjUTvSlxE9UrV6dl3ZaOCVSIpyeEhrpWEd2ybkvqV61PVEKU2VHKhU8/hRtugM7lb/W3Ys3aNYtj6cd4puszVzw3dqzR9/3xxyYEKyOlFC+Hv8zhs4eZuX2mw8e3RRFts4X7C643Wq9evRIHSE5L5mjaUdo3kvWh7amrX1cAWS9aXOHl1S9Tx7sOT3R+wuwoVukV1IvQgFBeW/MaF7JdZPFgO7tmEZ20iQ6NO+CmzLuvPTwc9u2D5GTTIjiUUoqwgDBWxa+Svuhr2LLF2Br+4Ycrzg2FeTqPt9e/TUjDEHoF9brieV9fGDYMvvgCzp83IWAZ9W7am65+XZm6ZqrDP1W0xVXMpgv3l8bWZGPXGumHtq/AmoE08mnE+sT1ZkcRTmT3id38duA3xnYaS/XK1c2OYxWlFK+Ev0JyejKfbPnE7DgVwtXaOTKyM/jr+F+m9UPny98C3NX6oo+cO0L8mXizozi1zz6DKlXg3nvNTmI7i2IWsffkXp7u8nSxy46OGQNnzhi94OWNUoqXe75M4rlEvtz2pUPHtkURbbeF+0tq69GtKBRtG7a19alFAUopuvp1lZlocZm31r9FFc8qPNbhMbOjlEpoQCiRTSN5fe3rpGelmx2n3LvaTPS25G3k6lzT+qHzhYRAzZqu1dIRGhAKyHrRV5ORAT/8AHfeCTVqmJ3Gdt5a/xZ+1f24s82dxR7Tvbvx7+KDD4x7G8qbiKAIuvh24a31bzl0pY4yF9H2XLi/pLYmb6VF3RZUq1zNlqcVRejm142EswkknUsyO4pwAkfOHmH2rtmMajuKOlXqmB2n1F7u+TIpGSl8tPkjs6OUe1croi/dVGjyTLS7O/Ts6VpFdJv6bajtXZvVCavNjuK0fv3VWJnj/vvNTmI70UnRrEpYxROdn8DTvfh/nEoZs9G7d0NUlOPy2YpSigk9JpBwNoE5u+c4bFybNKXZa+H+ktqavNXhO1+5qm7+3QBYf0RaOgS8u/FdtNY81eUps6OUSRe/LvS/rj9vrn9TduUso6u1c2xK2kRgzUAa+DRwXKBiREQYqxHEx5udxDHclBuhAaEyE30VX38NQUHGzpYVxXub3qN65eqMajfqmscOHw516hiz0eXRwOCB3FD/Bl5f+zp52jF7mZf7HQuPpR/jaNpRKaIdJKRhCF4eXlJEC05dOMVnWz9j+A3DCahp8/2THG5yz8mX3pMovavNRG9O2kzHJh0dF+YqwsON312tLzrudByJ5xLNjvL/7N13fNXl3f/x15WdEPYKIySssPeSMAOyNMgGQaYoVtTW9m6t2ru9e7e21rvWUesoisowbAIRZA9BdtiEEPaeIcyQnev3xzfxRxGFJOec63zP+TwfjzwS4Jzv942G5JPrfK7P5XbOnLFOshw7FnxsXxlZLt2+xLykeYxvMf6h9qsEB8Ozz8LixXDqlAsCOphSitc6v0ZyajIJKQ5tePhRtv9U2XleNhW6UoBvAO1rtJdDVwQf7/iY9Jx0Xom2x1zoB2lfoz09avfgna3vkJXr0L3PXuXHiuiLty9y+sZp460chZo0gcqVvaul4/u+aJkX/QMzZli9wGPHmk7iOJ/u+pSc/Bwmt5v80M95/nnr/RSbriUMazKMOuXr8NeNf3XJJBr7F9EXZFOhq0XXjGb3xd3cybnz4AcLj5STl8OHOz6kT90+NKvazHQch3m106ucv3Wemftmmo5iWz/WzlF4qqy7rEQrZbV0rF1rz41UxdGiagvKBpaVvuh7aG21cnTtCnXqmE7jGDl5OXyS+Am96/amQaUGD/28WrUgNtYad5ft+rNLSszPx49Xol9hx/kdrDmxxun3s30RnXg+UTYVulinWp3Izc/9/pui8D7xh+K5cPsCL7V/yXQUh3q0zqO0rtaa/9v8f3IyZzH92Er09nPb8VW+brXg0aMHnD8PR46YTuIavj6+dK7VWfqi77F1q/U54EkbChenLObcrXO82O7FIj/3Zz+Dy5dh0SInBHOBcS3HERYaxt833/dIEoeyfREtmwpdr2PNjgDS0uHF/rX9X9QpX4e+9fqajuJQSile7fQqh68eZtEhm34HMexHi+jz22lapSmlAkq5NtBPKOyLXuP8BSu30S2iGylXU7h4+6LpKG7jyy+t2dBDh5pO4jj/2v4vIstF8lj9x4r83D59IDISPrbp6PwgvyBeav8SK4+tZP8l556wbOsiWjYVmlExpCINKzWUItpL7b24l42nNzK57WR8fXxNx3G4wY0GU69CPd787k053a0Y7tfOobV2q02FherVs05rs+NIr+LqFtkNQFo6CmRlwdy51ol9pT3kBe39l/bz7alveb7t88X6Gu3jY53YuH69dbKnHf2s7c8I8Q/hna3vOPU+ti6iZVOhOZ3CO7H17FaXjZER7uNf2/9FsF8wE1pNMB3FKXx9fHkl+hV2XtjJupNeNLrBQe63En007SjXM6+7XRGtlDUvev167+mLbhXWilL+pWRzYYHly62T+p56ynQSx/k48WMCfQOZ2Gpisa/x9NPWv+VPPnFgMBeqEFyBCS0n8NW+r7hw64LT7mPvIlo2FRoTHR5NWkYaKakppqMIF0rLSOOr/V8xuvloKgRXMB3Haca0GEPlkMq8t/U901Fs535F9PZz2wH32VR4t5gYq//z4EHTSVzD39efTrU6seG0rEQDxMVZU1p69jSdxDHSs9OZuW8mw5oMK9EBWFWqwJAhMG2adZKjHb38yMvk5ufyr+3/cto9bF9Ey6ZCM6LDowE5dMXbfLH7CzJyM2x7xPfDCvIL4vm2z7Pk8BKOXPWSXWcOcr92jh3ndxDiH0Ljyo1dH+gBCvuivaqlI6IbBy4fIPVOqukoRt26BQkJ1jHfPzXf3E7mJs3lVvYtJrWeVOJrPf+8tUo/d64DghlQr0I9BjYcyMeJH5Oene6Ue9i7iD4vmwpNiaoYRYXgClJEexGtNZ/s/IRO4Z1oEdbCdByne77d8/j7+vPPbf80HcVWfmwluk21Nvj5+Lk+0ANERlpjvbzt0BWAjac2Gk5i1qJFkJkJo0aZTuI4n+76lIaVGtK5VucSX6tLF2jYED791AHBDPl19K+5lnmNL/d86ZTr27aIvnT7EudunaN1tdamo3glH+VDx5od2XxWimhvseHUBo6mHWVSm5KvcNhBWGgYI5uO5Is9X3A987rpOLZxbxGdk5fDrgu73LKVA6y+6JgYayU630u2eLSt3pYgvyCvH3X31VfWD1EdO5pO4hhJl5PYcnYLz7R6BqVUia+nFDzzDGzeDElJDghoQHR4NB1qdOD9be87ZQ+XbYvoXRd2AchKtEGdwjtxKPUQV+9cNR1FuMBnuz+jTGAZhjb2oDlQD/DyIy+TnpPOZ7s+Mx3FNu5t59h/eT9ZeVluW0SDtbnw6lX7FgpFFegXSMeaHb26iL50CVavhpEjrWLRE3y661P8ffwZ28Jxxy6OHWv9YPyZjb8E/qLDLziSdoQVR1c4/Nq2L6JbhrU0nMR7FfZFbz271XAS4WzXM68z/+B8nmr2FCH+IabjuEzLsJZ0j+zOB9s/IDc/13QcW7h3JbpwU2G76u0MpHk43btb772tpWPvxb1cy7hmOooR8+ZBXp7ntHJk5mYyY98MBjUaROVSlR123cqVYdAgmD7dan2xoyGNh1AttBrvb3vf4de2bxF9cRf1KtSjbFBZ01G8Vrsa7fBVvjIv2gvE7Y8jMzeTZ1o/YzqKy73c4WVO3zgth688pPsV0ZVCKhFZLtJInocRGWm9edXmwshuaLTXfv2Oi4NmzaBpU9NJHGNh8kLSMtJ4tvWzDr/2s89CWhosXOjwS7tEgG8Ak9tNZsWxFRxKdezga/sW0Rd2ST+0YSH+IbSq1ko2F3qBz3Z9RquwVl75by42KpaIshF8tOMj01Fs4d52jsJDVhzRo+lMMTHw7bfe0xfdoUYHAnwDvHJe9OnTsGULPPmk6SSOM3X3VGqXq02P2j0cfu0ePaB2bXtvMJzUZhIBvgF8sO0Dh17XlkV0WkYaJ6+fpHWY931DdzfRNaPZfm47OXk5pqMIJ9l1YRe7L+72ylVosA5f+Vnbn7Hu5DoOXvGSYcIlcPdK9K2sWxy8ctCtWzkKxcRYq2379plO4hrB/sF0qNHBK/ui58+33g8fbjaHo5y6foq1J9YyoeUEfJTjyzofH5g40Xql5ohNJ35WKVWFUc1GMW3vNIduFLdlEb37wm4Ar1wVczfR4dFk5Gaw99Je01GEk3y26zOC/IIY1cxDmgeLYWKriQT4Bshq9EO4u4jedWEXGu3WmwoLFfZFe1VLR0Q3dl3Yxa2sW6ajuNTcudCqlXXsuyeYvnc6gEM3FN5rwgTw9YWpU512C6f7efufk56TztRdjvtL2LKILtxU2KqanFRoWqdanQA5dMVTZeZmMuvALIY0GkK5oHKm4xhTuVRlhjcZzvS9072u4Ciqu4voHed3AO69qbBQeDjUretdmwu7RnQlT+d5VV/0qVOwbZvnrEJrrfly75f0qN2DiHIRTrtP9erw2GPWCYa5Nt1j3apaK7rU6sJHiR+Rl5/nkGvas4i+uItaZWtRKaSS6Sher2aZmoSXCZci2kMtObyE65nXGddinOkoxr3Q7gVuZd9i5r6ZpqO4NZ+7vqvsOL+DiLIRDp0W4Ezdu8PGjd7TFx0dHo2fj59X9UUXtnIMG2Y2h6N8d/o7jl87zvgW451+r6efhosXYdkyp9/KaV5s/yLHrx1n+dHlDrmePYto2VToVqLDo71qJcObTN87neqlqztls4rddKjRgVZhrfhwx4dorU3HcUv37h3ccW4H7Wq4/yp0oe7d4do17+mLLhVQirbV23pVX/S8edC6tfWqgyf4cs+XhAaEMrjRYKff6/HHoWpV+Pxzp9/KaQY1HES10Gp8uONDh1zPdkX0raxbHL56WDYVupHo8GjO3jzLmRtnTEcRDnQl/QrLji7jqWZP4evjazqOcUopXmj3AklXkthwaoPpOG7p7iI69U4qJ66fsEUrR6Fu1mnYXtcXveP8DtKz001HcTpPa+VIz05n7sG5DG88nFIBpZx+P39/6/CVJUusw2rsyN/Xn+faPMeyo8s4mna0xNezXRFduIFNVqLdR6dw6Yv2RLMPzCY3P9epm1XsZmSzkZQNLMuUXVNMR3FLdxfRiecTAXv0Qxcq7Iv2tiI6Nz+XLWe3mI7idJ7WyrEgeQG3s28zvuV4l93z6aetnugZM1x2S4eb1GYSfj5+fLzj4xJfy3ZFdOGmQimi3Ufzqs0J8Q+RItrDTN83nVZhrWhaxUNOI3CAEP8QxjQfw/yD80m9k2o6jtu5u4jecW4HCkWb6m3MBSqG7t1hwwbv6YvuVKsTvsrXK/qi586FNm2gTh3TSRxj2t5p1Clfh861Orvsng0bQnS0NaXDrl1t1UpXY0ijIXy+5/MSvwJjyyI6LDSMaqWrmY4iCvj7+tO+Rnvpi/YgB68cJPF8oqxC38ekNpPIzsv+fqyU+P/u3VTYoFIDygSWMReoGLytL7pMYBlaV2vt8X3Rp0/D9u2eswp99uZZ1p1Yx5jmY1x+kNHEiXDoEGzd6tLbOtSL7V/keuZ14vbHleg6tiyiZRXa/XQK78Sei3u8oq/OG8zYOwNf5cvIpiNNR3E7zao2o2PNjkzZOUU2GN6j8Hu51pod53fYqpWjkDf2RXeP7M62c9vIyMkwHcVp4uOt94Odv//OJWbtn4VGM7r5aJffe/hwKFXK3jOjO4V3okXVFiXeKG6rIjojJ4ODVw7KpkI3FB0eTZ7O+34urLCvfJ3PV/u/onfd3lQNrWo6jlt6rs1zpFxN8fjVu6IqLKLP3TrHxdsXbVlEe2tfdHZeNlvP2nhp8QEWLoRmzaB+fdNJHGPm/pk8UvMR6lVw/YkxoaFWIT1nDqTbdN1MKcXkdpPZe2lviT7vbVVE77+8nzydJ4esuKFHaj4CyOZCT7Dp9CbO3DzDU82eMh3FbQ1rMszaYLhTNhjerbCI3nGu4JAVG423u5u39UV3rtUZH+XD+pPrTUdxikuXrPnfQ4aYTuIY+y7tY9+lfYxu5vpV6EITJsDt27BggbEIJTaq2ShKB5Tmo8Tin0RrqyJaNhW6rwrBFWhcubH0RXuAWQdmEewXzICGA0xHcVsh/iGMbTGWBckLZIPhXb4vos/vwM/Hj5ZhLc0GKiZv64suG1SWVmGtPPaVlcWLrU1wntLKMXPfTPx8/BjexNysvs6drWPT7TwzOjQglHEtxjE3aW6xv47bqojefWE35YPKE1HWeUdbiuKLrhnNljNbyNdesnzjgXLycph3cB79G/QnNCDUdBy39lyb58jOy2banmmmo7iNwiJ6+7ntNK/anCC/ILOBiskb+6K7RXRj69mtZOZmmo7icAsXWgVfUw8YNJSXn0fc/jj61utr9CRQpazV6G+/hePHjcUosefbPU92Xjaf7y7eTwP2KqIv7qZlWEuX70QVDyc6PJprmddISU0xHUUU05oTa0i9kyobCh9CkypNeKTmI0zdPVU2GBZQyuqpTzyfaMt+6ELh4dYYtHXrTCdxne6R3cnKy2Lb2W2mozjUtWuwZo3VyuEJpcO3p77l3K1zRls5Co0da03k+fJL00mKr3HlxnSL6MYniZ8U6/m2KaJz83PZf3k/rcKkH9pdRYdHA9IXbWezDsyibGBZ+tXrZzqKLTzT6hmSU5O94qCKh6EUHEs7xo2sG7YuosFq6di40Xv6ortEdEGhPK6lY8kS63AQT2rlKB1QmicaPGE6CjVrQu/eVhGdl2c6TfFNbjeZE9dPFOu5timiD6UeIjM3U/qh3VhUxSgqBleUvmibysjJID45nsGNBhPoF2g6ji2MaDqC0IBQpu6y8awnB1KK7yf0tK3e1nCakvG2vuhyQeVoGdbS4zYXLlhgFXtt7f3pCEBmbiYLkhcwpPEQgv2DTccBrJaOM2dg7VrTSYpvYMOBhIWGFeu5timid1/YDSCTOdyYUoro8Ggpom3qmyPfcCv7lrRyFEFoQCgjmoxgTtIcbmXdMh3HOKWsyRzBfsE0qdLEdJwSKeyL/tazFmZ/UreIbmw5u4Ws3CzTURwiPR1WrIBBg/7zICC7+ubIN9zMusmopqNMR/negAFQvry9NxgG+AbwSvQrxXqubT6tdl3YRbBfMA0qNjAdRfyE6PBoDl89LBMLbGjWgVlUKVWFmNoxpqPYyjOtnyE9J505SXNMRzFOKUi8kEiraq3w8/EzHadEatWy+qK9aXNh98juZOZmsv3cdtNRHGLFCsjM9JxWjrj9cW73NTowEEaNsg6zuXbNdJri+2XHXxbrebYpondf3E3zqs3x9fE1HUX8hE7hnQDpi7abW1m3WHpkKcMbD7d98eNqHWp0oEnlJny26zPTUYxTylrwsHs/dKHu3a2VaG/pi+4a0RWFYt1Jz9hRuWgRVKhgjWOzu5tZN1lyeAkjmoxwu6/REyZAVhbMnm06ievZoojWWrPn4h7ZVGgDbau3xd/Hn02npaXDTpYcXkJmbiYjmo4wHcV2lFJMbDWRbee2ceDyAdNxjMonjzs5d2zfD12osC96/37TSVyjfHB5j+mLzsmxNhX27w9+7lVzFkt8cjxZeVlu2W7XurV1GuQXX5hO4nq2KKJPXD/BjawbsqnQBoL9g2lTvY30RdvM3INzqV66+vcTVkTRjGkxBn8ff6/fYJincwA8ZiXaG+dFx0TGsPnMZtvPi9640foBaOBA00kcY9aBWUSWi/z+dGB3UjgzescOSEoynca1bFFEy6ZCe+kU3onE84kesznF093KusWyI8sY0mgIPsoWXxLcTqWQSgxoOICZ+2eSnZdtOo4xufm5lAksQ/2K9U1HcQhv7YvOysti69mtpqOUyKJFEBxsjWCzu0u3L7H6+GpGNR3ltudkjB5trfh722q0Lb5j7rqwC1/lS9MqHnDckBfoFN6JrLys749pF+5tyeElZOVlGT1C1hM83fJpUu+ksuTwEtNRjMnV2bSp1sajfhjr1s27+qK7RHTBR/nYuqVDa6uI7t0bQkJMpym5eQfnkafzGNnM/Vo5ClWuDLGxMGOG1UrjLWzxlW73xd00rtzYtkfIepvClgBp6bAHaeVwjN51e1OjdI1iHx/rCXLzczymlaOQt/VFlwsqR6uwVrbeXLh7tzW72JNaOZpVaeb2C4kTJsDly7BsmekkrmObIlpaOeyjamhV6pavK0W0DRS2cgxtNNSjVg9N8PXxZVyLcSw7uozzt86bjmOG0rSr4VlFtDfOi46JjGHr2a1k5GSYjlIsixZZc6FjY00nKblT10+x+cxmnmz6pOkoD9SvH1St6l0tHW7/XfPi7YtcvH2R1mGyqdBOOtXqxKbTm9Bam44ifsLXh78mKy+LYU2GmY7iESa0mkC+zmf63ummozyQUqqvUipFKXVUKfXqff48UCk1p+DPtymlIh981XyPmcxRKCICateGdfZdmC2ymNoxZOdl2/Y4+0WLoEsXqFTJdJKSm5s0F8AWRbS/P4wZY01FuXzZdBrXcPsiunBTYcuwloaTiKLoFN6JK3eucDTtqOko4ifMOzhPWjkcqF6FenSN6Or2LR1KKV/gQ6Af0BgYqZRqfM/DJgLXtNb1gHeBtx50XR8fHyLKRjg6rnHdu8OGDd7TF925Vmd8lS/rTtjvJ4fjx63WmwEDTCdxjDlJc2hXvR11ytcxHeWhTJgAubkwc6bpJK7h/kX0RSmi7UgOXXF/MpXDOZ5u+TRH0o6YjvEg7YGjWuvjWutsYDZwb9kxAJhW8PF8oKd6wGgAf18/t50eUBLdu0Namvf0RZcJLEOb6m1s2Re9eLH13hOK6CNXj7Dzwk5brEIXatwYOnSwjgH3hhei3f475+6Lu6lTvg5lg8qajiKKoFHlRpQLKid90W7smyPfWK0cjaWVw5GGNh5KaECo6RgPUgM4c9evzxb83n0fo7XOBW4AFX/qogF+/g6M6D66d7fee9Oou5jIGLaf2056drrpKEWSkABNm1qjCe1uTtIcANtNTnr6aWtedGKi6STO5/5F9IXdclKhDfkoHzrW7ChFtBtbkLyAKqWqSCuHg5UKKEXCkwmmYzzI/ZaL7103epjHoJSapJRKVEol+ivPnJFdOC/am/qiu0d2Jyc/x1Zfw9PSrENWPGEVGmD2gdl0qdWFmmVqmo5SJCNGWDO6vWGDoVsX0Tcyb3Ds2jEpom2qc63OHLxykLSMNNNRxD0ycjL45sg3DGo4CF8fX9NxPE5M7RjTER7kLBB+169rAveOFPn+MUopP6As8IN/zFrrKVrrtlrrtrXDKzsprnne2Bft5+Nnq77ob76BvDx44gnTSUruwOUDJF1JYkSTEaajFFnZsjBkCMTFQYY9B7w8NLcuovdd2gfISYV21blWZ0D6ot3RimMrSM9JZ0ijIaajCDN2APWVUrWVUgHAk8C9y+cJwLiCj4cCa7UXj9spnBe9b5/pJK4RGhDKIzUfYc2JNaajPLTFi6FaNWjrAQNi5hyYg4/yYWjjoaajFMuECXDjhjUpxZO5dRFduKlQVqLtqV31dvj7+PPd6e9MRxH3mH9wPuWDytM9srvpKMKAgh7nF4EVQDIwV2udpJT6k1KqcB1vKlBRKXUU+BXwgzF43qSwL9qbWjp6RPZg54WdXM+8bjrKA2VlwfLl0L+/NSPazrTWzE6aTY/aPagaWtV0nGLp3h0iI2HqVNNJnMutP9V2X9xNlVJVCAsNMx1FFEOwfzBtqreRItrNZOVm8fXhrxnQcAD+vp65EUw8mNb6G611lNa6rtb6LwW/9wetdULBx5la62Fa63pa6/Za6+NmE5sVHg5163rX5sKedXqSr/P59qT7nzSzfj3cvu0ZrRy7L+7maNpRW7ZyFPLxsVaj16yBEydMp3Ee9y6iCzYVeuLIJG/RObwzO87vIDM303QUUWDNiTXczLoprRxCFFFhX3RenukkrtGhRgeC/YJZe2Kt6SgPtHgxhIRAjx6mk5TcnANz8PPxY3CjwaajlMj48aCUZ28wdNsiWmtN0pUkaeWwuc61OpOdl03ieS+YdWMTCw4uoHRAaXrV6WU6ihC2EhMD16/D3r2mk7hGoF8gXSK6uH1ftNbWaLs+faypEHamtWbuwbn0qtOLCsEVTMcpkVq1rP8nX3zhuT94um0RnZGbQW5+rmwqtLnC8WnS0uEecvNzWZyymNioWAL9Ak3HEcJWunWz3ntTS0ePyB4kXUni0u1LpqP8qN274dw5qx/a7raf287J6ydt3cpxt4kT4exZWLXKdBLncNsi+k7OHUA2Fdpd5VKVaVipoRTRbmLjqY1czbgqrRxCFENa92BZAAAgAElEQVTNmlCvnncV0T3r9ARw65aOhASrbSA21nSSkpubNJcA3wAGNPSMYddPPAGVKsFnn5lO4hxuXUSHBoRSt0Jd01FECXUO78ymM5vI114yYNWNLUxeSJBfEH3r9TUdRQhbiomBb7/13Jen79UqrBXlgsq5fREdHQ2VbT6mPF/nM/fgXPrU7UO5oHKm4zhEQACMGWP9P7pyxXQax3NIEa2U6quUSlFKHVVK/WAMklIqUCk1p+DPtymlIh90zTs5d2hRtQU+ym3rfPGQOtfqzPXM6xy8ctB0FK+mtWZRyiL61O1DqYBSpuMIYUsxMXDzptVC4A18fXzpHtndbfuiz561/l94QivHljNbOHvzrMe0chSaOBFycmDGDNNJHK/EFapSyhf4EOgHNAZGKqUa3/OwicA1rXU94F3grQddNyM3Q1o5PEThoSvS0mFW4vlEzt48y6CGg0xHEcK2CudFr3XfhVmH6xHZgxPXT3DimvvNKvv6a+u9J4y2m5M0hyC/IJ5o4AF/mbs0aQKPPGK1dHjacU2OWOZtDxzVWh/XWmcDs4F7m3kGANMKPp4P9FQPmFuXn58vmwo9RJ3ydQgLDWPTmU2mo3i1+EPx+Cpf+jfwgCUbIQypVg0aNfKuQ1cK+6LdcTU6IcGa392woekkJZOXn8e8g/N4rP5jlA4sbTqOwz37LCQnw2YPO8DYEUV0DeDMXb8+W/B7931MwUlZN4CK915IKTVJKZWolEoE2VToKZRSdArvJCvRhi1MXkj3yO62H5skhGkxMbBxo/UStTdoVKkR1UtXZ/Xx1aaj/Ifbt61XBJ54wtpYaGcbT2/k4u2LHtfKUWjECChdGqZMMZ3EsRxRRN/vU/feBfuHeQxa6yla67Za67ZNqjShaZWmDogn3EHnWp05ef0kZ26cefCDhcMlX0km5WqKtHII4QA9ekB6OuzYYTqJayileLTOo6w5scatNoivXAnZ2R7SynFgDiH+ITxe/3HTUZyiVCkYPRrmzoVr10yncRxHFNFngfC7fl0TOP9jj1FK+QFlgbSfumiQX5AcSexBukZ0BayftoXrxR+KB2Bgw4GGkwhhf4Xzor2ppaNXnV6k3kllz8U9pqN8LyEBypWDTp1MJymZ3PxcFiQvoH9Uf4/e9P3ss5CZCTNnmk7iOI4ooncA9ZVStZVSAcCTQMI9j0kAxhV8PBRYq7WntZeLn9KiagtKB5Rm4ykpok1YmLyQDjU6UKPMvZ1WQoiiqlQJmjf3riK6Z22rL3rVMfc4NSMvD5YuhcceA3+br7etP7meK3eueGwrR6FWraBtW6ulw1MqwBIX0QU9zi8CK4BkYK7WOkkp9SelVOGLLFOBikqpo8CvgB+MwROezdfHl061OrHh9AbTUbzO6Run2XlhJ4MbDTYdRQiP0aMHbNoEWVmmk7hGtdLVaFqlKatPuEdf9NatkJrqGaPt5hyYQ2hAqFfM7580CQ4cgG3bTCdxDIcMYdZaf6O1jtJa19Va/6Xg9/6gtU4o+DhTaz1Ma11Pa91ea33cEfcV9tKlVhcOXjlI6p1U01G8yuJDiwFp5RDCkWJirJemt241ncR1etXpxcZTG8nIyTAdhYQE8PODfv1MJymZnLwcFh5ayIAGAwj2DzYdx+mefBJCQz1ng6GcZCJcprAvWqZ0uFb8oXgaV25MVMUo01GE8Bhdu4KPj3e1dPSq04usvCy3+Br+9ddWb3rZsqaTlMzq46tJy0jz+FaOQqVLw6hRMHu2Z2wwlCJauEy76u0I9A2UvmgXunrnKhtObZCpHEI4WLly0Lq1dxXRXSO64u/jz6rjZvuijx61Zg57RCtH0hzKBpald93epqO4zPPPQ0YGTJv24Me6OymihcsE+gXSoWYH6Yt2oSWHl5Cn86SVQwgn6NEDtmyBO3dMJ3GNUgGliA6PNl5EF55SaPciOis3i0WHFjGo0SAC/QJNx3GZli2hY0f4+GP7bzCUIlq4VNdaXdl9YTe3sm6ZjuIVFqUsomaZmrSp1sZ0FCE8To8e1oEr35nvbnCZXnV6sefiHq6kXzGWISEBmjaFOnWMRXCIFcdWcCPrhte0ctxt8mQ4fNg6LMfOpIgWLtUlogt5Oo8tZ7eYjuLx7uTcYcXRFQxsMBBl9+O8hHBDnTtb49XWuN9p2E7Tq24vAGOnF167Zp0WafdVaIDZB2ZTMbji9+MDvcnQoVCxInz0kekkJSNFtHCpjjU74qt8pS/aBVYeW0lGbgaDGkk/tBDOUKqU9bK0NxXRbaq1oUJwBVYcW2Hk/suWWTOi7X5K4Z2cOySkJDC08VCvPFguKAgmToTFi+HcOdNpik+KaOFSpQNL06paK+mLdoH4Q/GUDypPl1pdTEcRwmP17Am7dkHaT57B6zl8fXzpXbc3y48uN3IE+NdfQ5Uq0L69y2/tUEsPLyU9J90rWzkKPfcc5OfDp5+aTlJ8UkQLl+taqyvbzm4jK9dLTikwICcvh69TvqZ/g/5eucohhKv07Gltjlq/3nQS1+lbty+X0i+x79I+l943J8daiY6NtcYL2tnspNmEhYZ9P/rVG9WpY835njIFsrNNpykem38aCjvqGtGVrLwstp/bbjqKx9p4eiPXMq8xsIFM5RDCmdq3t9o6vKmlo3Ac2/Kjy116340b4cYN+7dy3My6ydLDSxnWeBi+Pr6m4xj1wgtw4QIsWGA6SfFIES1crktEFxSK9SfXm47iseKT4wn2C6ZPvT6mowjh0fz9rYNXvKmIrla6Gi3DWrq8iE5IsHppH33Upbd1uISUBLLysniy6ZOmoxjXty9ERcH775tOUjxSRAuXqxBcgRZhLVh/ar3pKB5Ja82ilEX0rtubEP8Q03GE8Hg9e0JKir03SBVV37p92XRmEzezbrrkflpbRfSjj1or/3Y2+8BsapWtxSM1HzEdxTgfH3jpJdi2zXqzGymihRHdI7qz5cwW6Yt2gl0XdnH25lk5pVAIF+lZMKHMm1aj+9brS25+LmtPuGbQb1ISnDhh/1aOtIw0VhxbwfDGw/FRUoIBjBsHZcrYczVa/g8KI7pFdiMjN4Md53eYjuJx4g/F46t8iY2KNR1FCK/QvDlUquRdRXTH8I6UDijNsiPLXHK/hATrfazNv6wtOLiA3PxcaeW4S+nS1ri7efPs92qOFNHCiK4RXaUv2kkWHVpE14iuVAypaDqKEF7BxwdiYqwi2u7HGD+sAN8AetbpyfJjy9Eu+EsvXmxt4qxWzem3cqq4A3FEVYyidbXWpqO4lRdftOZ/f/yx6SRFI0W0MKJCcAWaV20uRbSDHbl6hKQrSQxsKFM5hHClRx+1VtFSUkwncZ1+9fpx+sZpDqUecup9LlyA7dvt38px7uY5vj35LaOajpJTZO9Rp451CuW//w0ZGabTPDwpooUx3SO7s/nMZumLdqBFhxYBSBEthIv1sk7DZtUqszlcqU9da/rPsqPObelYssR6b/ciek7SHDSakc1Gmo7iln75S0hNhWnTTCd5eFJEC2O6R3aXvmgHiz8UT+tqralVtpbpKEJ4ldq1oW5dWLnSdBLXiSgXQZPKTVh6ZKlT75OQAJGR0LSpU2/jdLMOzKJNtTZEVYwyHcUtdesG7drB229brR12IEW0MEb6oh3rwq0LbD27VaZyCGFIr17WyYU5OaaTuE5sVCwbTm3gRuYNp1w/PR1Wr7ZWoe3cAXHk6hESzycysqmsQv8YpeCVV+DYMVi0yHSahyNFtDCmsC/621Pfmo7iERJSEtBoaeUQwpBeveD2bdi61XQS14mNiiU3P5eVx5yzBL96NWRmwoABTrm8y8w6MAuFYkTTEaajuLVBg6xXdN56yx6bdKWIFkZ1i+jGptObyM7LNh3F9hYeWki9CvVoUrmJ6ShCeKUePaxJHd7UF92xZkcqBFdgyZElTrn+4sVQtix06eKUy7uE1pq4/XF0jehKzTI1Tcdxa76+8Otfw44dsGGD6TQPJkW0MKqwL3rbWRseVeRGrmdeZ+2JtQxqOEh2fQthSLly1hg2byqifX18eaz+Y3xz5Bvy8h3byJqXB19/DY8/bh2vbld7Lu4h5WqKtHI8pHHjoHJl+L//M53kwaSIFkZ1i+yGQrns1CtPtfTwUnLzcxncaLDpKEJ4tV69rHFs166ZTuI6sfVjSb2TyrZzjl0M2bzZmtYw0OYdajP3zcTfx59hTYaZjmILwcHw85/DN9/A3r2m0/w0KaKFURWCK9C6WmvWnPCio76cYOGhhVQLrUb7Gu1NRxHCq/XqBfn5sG6d6SSu06deH3yVL0sOO7alY9EiCAiAvn0delmXysvPI+5AHI9HPU6F4Aqm49jGCy9YR4H/+c+mk/w0KaKFcT1r92Tr2a2kZ6ebjmJLGTkZLD+6nIENB+Kj5J+0ECY98giEhnpXS0e5oHJ0ieji0CJaa6sfukcP61hou1pzYg0Xb19kdLPRpqPYSvny1mr0ggWwf7/pND9OvuMK43rW6UlOfg7fnf7OdBRbWnlsJXdy7shoOyHcgL+/dQS4NxXRAP2j+rP/8n5OXT/lkOsdPGiNOvOEVo6ygWV5POpx01Fs55e/tH6AeuMN00l+nBTRwrhO4Z3w9/GXlo5iij8UT7mgcnSP7G46ihACq6Xj2DHrzVvERsUCOGw1unBOcP/+DrmcEenZ6SxMXsiwxsMI8gsyHcd2KlSAl16CefOsH6rckRTRwrhSAaXoGN5RiuhiyMnLISElgf5R/fH3tfH2dSE8SGEP74oVZnO4UlTFKKIqRrE4ZbFDrrd4MXToANWrO+RyRixOWUx6Tjqjm0srR3H96lcQEuK+vdFSRAu30LN2T3Zf2E1aRprpKLay4dQGrmVek1YOIdxIvXpQpw4sX246iWsNajiIdSfXlfjr+Llz1pxgux+wMnPfTGqVrUWXCBsPuTasYkV48UWYMweSkkyn+SEpooVb6Fm7JxotR4AX0YLkBYT4h9CnXh/TUYQQBZSyVqPXroWsLNNpXGdwo8Hk5ueWuKUjIcF6b+d+6Eu3L7Hy2EqeavaUbPguoV//2uqNfv1100l+SP7PCrfQrkY7SvmXYs1xael4WPk6n/hD8fSr148Q/xDTcYQQd+nbF9LTYdMm00lcp231ttQoXYP4Q/Eluk58PNSvDw0bOiiYAXH748jTedLK4QCVKsFvf2v9cPWdm80fkCJauIUA3wC6RnSVvugi2HxmMxdvX2RIoyGmowgh7hETY8049qaWDh/lw6CGg1h+dHmxR5ampVkztocMsVb07UhrzRd7vqBd9XY0rtzYdByP8PLLVn/8b35jjT90F1JEC7fRs3ZPUq6mcPbmWdNRbGFh8kICfANkdJIQbig0FLp0gWXLTCdxrcGNBpOZm8mKY8XbVblkCeTmwiAbb/PYc3EP+y/vZ3zL8aajeIyQEPjf/4WtW///5BZ3IEW0cBuP1nkUgNXHVxtO4v601ixMXkjvur0pE1jGdBwhxH307QsHDsBZL1oX6BLRhYrBFVmYvLBYz4+Ph5o1oW1bBwdzoS/2fEGAbwAjm440HcWjjB8PjRrBa69ZP2i5AymihdtoVrUZVUtVZeWxlaajuL2dF3Zy6sYpaeUQwo1546g7Px8/nmjwBEsOLyE7L7tIz01Pt9pfBg0CH5tWJ1m5WXy1/ysGNhxI+eDypuN4FD8/+NvfICUFPvrIdBqLTT9NhSfyUT70rtubVcdXka/zTcdxawsOLsBX+dI/ysYnEQjh4Zo0gRo1vKsvGqxRdzeybrDuxLoiPW/5csjMtHcrx9IjS0nLSGN8i/Gmo3ik/v2hd2/4/e/h4kXTaaSIFm6md93epN5JZdeFXaajuC2tNQuSFxBTO4aKIRVNxxFC/IjCUXerVrnPy8+u0KtuL0IDQlmQvKBIz4uPt+YCd7HxWOUv93xJtdBq9Krby3QUj6QUfPCB9cPWK6+YTiNFtHAzvev2BpCWjp9w4PIBjqQdkVYOIWygb1+4cQM2bzadxHWC/ILoH9WfBckLyMnLeajnZGfD119bB6z4+Tk5oJNcun2Jb458w9gWY/HzselfwgaioqwpHTNmwIYNZrNIES3cSpVSVWgV1qrYO7u9wfyD81EoBja08UkEwhilVAWl1Cql1JGC9z9o3FRKtVRKbVFKJSml9imlRpjI6gl697aKwqVLTSdxrSebPklaRtpDbxRfuxZu3rR3K8f0vdPJ03mMazHOdBSP9/rrEBEBkydDzsP9nOYUUkQLt9O7bm82n9nMzaybpqO4pXkH59E1oithoWGmowh7ehVYo7WuD6wp+PW97gBjtdZNgL7Ae0qpci7M6DHKlIFu3azRbd6kT90+lAsqx6wDsx7q8QsXWmMBH33UycGcRGvNp7s+pVN4JxpVbmQ6jscLCYH337eOAv+//zOXQ4po4Xb61O1Dbn5ukTeleIOky0kkpyYzvMlw01GEfQ0AphV8PA34wUsaWuvDWusjBR+fBy4DlV2W0MPExsLBg3D8uOkkrhPoF8jghoNZdGgRGTkZP/nY3FyrH7p/fwgKclFAB/v21LccSTvCpDaTTEfxGk88AcOHwx//CLt3m8kgRbRwO9Hh0YT4h0hf9H3MTZqLQjG40WDTUYR9VdVaXwAoeF/lpx6slGoPBADHfuTPJymlEpVSiVeuXHF4WE8QG2u998aWjlvZt1h29KdPnPn2W0hNhaFDXRTMCT7d9SnlgsoxrPEw01G8hlLWqLvKlWHsWGuzoatJES3cTqBfIDGRMdIXfQ+ttbRyiIeilFqtlDpwn7cBRbxONWAGMEHr+8+d1FpP0Vq31Vq3rVxZFqvvp149aNjQ+1o6YmrHUKVUlQe2dMybB6VKQb9+LgrmYFfvXGX+wfmMbjaaYP9g03G8SsWKMHWqdajRH/7g+vtLES3cUp+6fTh27RjH0u67+OWVkq5IK4d4OFrrR7XWTe/zthi4VFAcFxbJl+93DaVUGWAp8N9a662uS++ZYmNh/Xq4dct0Etfx8/FjWONhLDm8hFtZ9/+L5+Za/dCxsRBs0/pzxr4ZZOdlSyuHIf36wXPPwdtvwzoXd4FKES3cUp96fQAe+DKgN5mXNE9aOYQjJACF4wPGAYvvfYBSKgCIB6Zrree5MJvHio21xritfrhhFR5jZNORZOZmkpCScN8/37gRrlyBYTbtgtBaM2XnFDrU6ECzqs1Mx/Fab79tjb4bMQJOn3bdfaWIFm4pqmIU9SrUY+kRL2si/BFaa+YenEu3yG7SyiFK6m9AL6XUEaBXwa9RSrVVSn1W8JjhQFdgvFJqT8FbSzNxPUN0NJQr530tHR3DOxJeJpyv9n913z+fN8+atGDXVo7NZzaTnJosq9CGhYbC4sWQlQUDB8KdO665rxTRwm3F1o9l3Yl1pGenm45iXNKVJA6lHpJNK6LEtNZXtdY9tdb1C96nFfx+otb6mYKPZ2qt/bXWLe9622M2ub35+1sHryxdCvn37S73TD7KhzHNx7Di2ArO3zr/H3+Wl2e1cjz+uFVI29FHiR9RJrAMI5rIKHXTGjSAr76CPXvg2WdBa+ffU4po4bYej3qcrLws1pxYYzqKcTKVQwj7i42FS5dgxw7TSVxrfMvx5Ot8Zuyd8R+//9131n8Pu7ZyXLx9kXlJ85jQcgKlAkqZjiOw/o39+c8QFwd/+pPz7ydFtHBbXSO6EhoQypLDXvb65z201sw6MIuY2jHSyiGEjT32mHV64aJFppO4Vv2K9elcqzOf7/kcfdfy4Lx51mbCxx4zGK4E/p34b3Lyc3ih3Qumo4i7vP46jBtnzY/+61+dey8pooXbCvANoHfd3iw9svQ/vvB6m50XdnI07Sgjm440HUUIUQLly0P37lYLg7d9SZvQcgKHrx5my9ktgDWVY948a+WwlA0XcbPzsvlk5yf0q9eP+hXrm44j7qKUNfZu9Gj43e/grbecdy8pooVbi60fy/lb59lz0XvbMWftn4W/jz9DGg0xHUUIUUKDBsHhw5CcbDqJaw1rPIwQ/xC+2P0FYI0iu3wZRtp0bWDBwQVcvH2Rl9q/ZDqKuA9fX/jyS+vz69VXrRnSztiLIEW0cGv96ltbtr11Ske+zmdO0hz61utL+eDypuMIIUpoQMFxN/HxZnO4WunA0gxvMpzZSbNJz05n1iwoU8a+Uzk+2P4B9SvU/34cq3A/vr4wfTpMmGD1ScfGQlqaY+8hRbRwa2GhYbSr3s5r+6I3ntrIuVvnpJVDCA9RowZ06OB9RTRYLR23s28ze288CxfC4MEQFGQ6VdHtPL+TLWe38EK7F/BRUka5Mz8/q7Xj449hzRpo3Ro2b3bc9eX/vnB7sVGxbD+3ncvp9z1YzaPF7Y8jxD+EJxo8YTqKEMJBBg2CnTtdeyiEO+hSqwt1y9flnRnJ3Lhh31aOd7e+S2hAKONbjjcdRTwEpeBnP7MO9snPh06drIkwhw+X/NpSRAu393j9x9FovjnyjekoLpWdl8385PkMaDBAxicJ4UEGDbLee9uUDqUUk9pM4uDa5pSvmEuPHqYTFd3J6yeZfWA2k1pPomxQWdNxRBG0bw9JSdbUjmXLoHFj64TDuLjit3lIES3cXutqralZpibxh7zr9c9Vx1aRlpEmrRxCeJioKOsbuDe2dIyoPxEO96dah+/w8zOdpuje3fIuSilefuRl01FEMZQuDf/zP3DsGLzwgrXB9amnoHLl4l2vREW0UqqCUmqVUupIwfsf7HxSSrVUSm1RSiUppfYppeRYH1EkSikGNxzMiqMruJ1923Qcl4k7EEf5oPKycUUIDzRoEGzYAKmpppO41nerKkJOCMdqvMG1jGum4xTJ1TtX+Wz3ZzzV7CnCy4abjiNKoGpVeP99uHgRtm61ZksXR0lXol8F1mit6wNrCn59rzvAWK11E6Av8J5SqlwJ7yu8zOBGg8nKy2L50eWmo7jEzaybxCfHM6LJCAJ8A0zHEUI42ODBVn+mt61Gx8VBWPVssqqt5fPdn5uOUyT/2v4v7uTc4ZVOr5iOIhzEx8fa6PvnPxfz+SW8/wBgWsHH04CB9z5Aa31Ya32k4OPzwGWgmAvnwlt1rtWZSiGVWJi80HQUl1hwcAEZuRmMbTHWdBQhhBO0agX16sGcOaaTuM6lS7BiBUwYF0CXyM58uOND8vLzTMd6KOnZ6Xyw/QP6R/WnceXGpuMIN1HSIrqq1voCQMH7Kj/1YKVUeyAAOPYjfz5JKZWolEq8cuVKCaMJT+Lr48uABgNYcngJWblZpuM43fR906lfoT6P1HzEdBQhhBMoZW1qWrfOKi69QVwc5OXBmDHwUvuXOHH9hG02jH+++3OuZlyVVWjxHx5YRCulViulDtznbUBRbqSUqgbMACZore97bozWeorWuq3Wum3l4nZ5C481uNFgbmXfYs2JNaajONXJ6ydZf3I9Y1uMRSllOo4QwkmefNJq6Zg/33QS15g2Ddq1g0aNYGDDgdQoXYN/bv+n6VgPlJGTwd82/Y3OtTrTuVZn03GEG3lgEa21flRr3fQ+b4uBSwXFcWGRfN9BvkqpMsBS4L+11lsd+RcQ3qNn7Z6UDijt8S0dM/fNBGB089GGkwghnKlpU2jSxDtaOvbutd7GjbN+7e/rz0vtX2L18dXsOLfDbLgH+DjxY87fOs8bMW+YjiLcTEnbORKAgn8SjAMW3/sApVQAEA9M11rPK+H9hBcL9AskNiqWxSmLyc3PNR3HKbTWTN87nW4R3YgsF2k6jhDCyUaMsA6BOHvWdBLnmj4d/P2tv2+hye0mUz6oPG9sdN/i9Hb2bd787k0erfMo3SK7mY4j3ExJi+i/Ab2UUkeAXgW/RinVVin1WcFjhgNdgfFKqT0Fby1LeF/hpQY3GkzqnVS+O/2d6ShOse3cNo6kHWFci3EPfrAQwvYKi8q5c83mcKbcXPjqK3j8cahU6f//funA0rz8yMskpCSw9+JecwF/wj+3/ZPUO6myCi3uq0RFtNb6qta6p9a6fsH7tILfT9RaP1Pw8Uyttb/WuuVdb3scEV54n771+hLsF8zcJM/8jjNtzzSC/YIZ0niI6ShCCBeIioLWrT27pWPlSmvz5Lj7rA38vMPPKRNYhr9s/Ivrgz3A9czr/H3z3+kf1Z8ONTuYjiPckJxYKGwlNCCUJxo8wbyD88jJyzEdx6HSs9OJOxDHkMZDKBNYxnQcIYSLjBgB27fD8eOmkzjHtGlQsSI89tgP/6xcUDlebPci8w/OJ/lKsuvD/YS3N7/N9czr/CnmT6ajCDclRbSwnVHNRpF6J5XVx1ebjuJQ8w7O42bWTZ5t/azpKEIIFyps6fjqK7M5nCE1FRYtglGjIOBHzo36ZcdfEuIf4lar0Sevn+QfW/7Bk02fpGWYdKCK+5MiWthO33p9KR9UnrgDcaajONSnuz6lQcUGdKnVxXQUIYQLRURATIy1Yqu16TSONX06ZGfDpEk//phKIZV4od0LxO2PY/eF3a4L9xN+teJX+Cgf/t7r76ajCDcmRbSwnQDfAIY2Hkp8cjzp2emm4zhE0uUkNp/ZzDOtn5HZ0EJ4ofHj4dgx2LTJdBLH0RqmTIGOHa1xfj/ltS6vUTGkIi+veBlt+CeJVcdWEX8ont91+R01y9Q0mkW4NymihS2NajaK9Jx0vj78tekoDvHZrs/w9/GXqRxCeKkhQyA0FL780nQSx9m4EVJSfnoVulC5oHL8OebPbDi1gQXJC5wf7kdk52Xz8+U/p275uvxXx/8ylkPYgxTRwpa6RnSlRukaxO23f0tHZm4m0/dNZ2DDgVQuJSd1CuGNSpWCYcOsUXfpnvECG1OmQNmyMHz4wz3+mdbP0KxKM36z6jdk5GQ4N9yP+GDbBxxKPcR7fd8j0C/QSAZhH1JEC1vyUT6MbDqSZUeXcfXOVdNxSiQ+OZ60jDTZUCiElxs/Hm7dgvh400lK7upV6wHvfIoAABIESURBVDjz0aMhJOThnuPn48d7fd/j5PWTvLPlHecGvI+U1BR+v+73xEbFEhsV6/L7C/uRIlrY1qhmo8jNz2X+wfmmo5TIlF1TiCwXSc86PU1HEUIY1Lkz1KnjGS0dM2ZAVhY8W8S1gR61ezCo4SD++t1fOZp21Dnh7iMnL4fR8aMJ9g/m37H/dtl9hb1JES1sq2VYSxpXbsznez43HaXY9l/az/qT63m+7fP4KPnnKIQ38/GxDiRZuxZOnTKdpvi0hn//Gzp0gBYtiv78f/b7JwG+AYxeONpl5wG8seENEs8nMiV2CtVLV3fJPYX9yXdtYVtKKZ5p9Qzbz21n36V9puMUywfbPyDYL5hnWj9jOooQwg2MHWsVoXZejV61Cg4dguefL97za5apyb9j/822c9tcMjt669mt/GXjXxjbYqycFiuKRIpoYWtjWowhwDeAz3Z9ZjpKkV29c5WZ+2YyuvloKgRXMB1HCOEGIiOhb19rJTfHpoeyvvceVK0KTz5Z/GsMbzKcMc3H8OcNf2bLmS2OC3ePS7cvMXLBSGqUqcE/+/7TafcRnkmKaGFrlUIqMbjRYGbsm2FsN3dxTd09lYzcDF5q/5LpKEIIN/LCC3DhgnXSn92kpMCyZTB5MgSWcLjFB/0+ILxMOE8tfMopG8jv5NzhidlPcOn2JeYPm0/ZoLIOv4fwbFJEC9t7tvWzXM+8bnS2aFHl5ufy4Y4PiYmMoVnVZqbjCCHcSL9+1or0hx+aTlJ0//yndbz3c8+V/Fplg8oSNySO87fO81jcY9zOvl3yixbI1/mMiR/DjnM7iBsSR7sa7Rx2beE9pIgWttc9sjt1y9fl012fmo7y0BJSEjh94zQ/7/Bz01GEEG7G19dayf32WzhwwHSah3ftmtXLPWqU1c7hCNHh0cweOpvE84kMmTuE7LzsEl9Ta82vVvyKhckLeafPOwxsONABSYU3kiJa2J6P8uGZ1s+w4dQGUlJTTMd5KO9tfY+IshH0j+pvOooQwg09/TQEBcFHH5lO8vCmToU7d+AXv3DsdQc2HMiU2CmsPLaSsfFjSzSxIys3izHxY3h/2/v8osMv+EUHB4cVXkWKaOERxrccj5+PH1N2TjEd5YE2nd7ExtMb+UWHX+Dr42s6jhDCDVWsaG3Mmz4dbtwwnebBcnLggw+gWzdo2dLx15/YeiJvPfoWc5LmEDMthrM3zxb5GlfvXKXXjF58tf8r/trjr7zb512UUo4PK7yGFNHCI4SFhjGk0RA+2/0ZN7Numo7zk9787k0qBldkUptJpqMIIdzYCy9YR4BPm2Y6yYPNmAGnT8NvfuO8e7zS6RXiBsex5+IeWv27FSuPrXyo52mtWZi8kHaftmP7ue3MGjKL17q8JgW0KDEpooXH+HX0r7mZddOtx93tubiHpUeW8vIjL1MqoJTpOEIIN9a2rXWK4T/+4d7j7nJz4a9/hdat4bHHnHuvkc1GkjgpkbDQMPrM7EOfmX1YcngJefl5P3hsXn4eG05toMsXXRgydwjB/sGsHbeWJ5uWYPaeEHfxMx1ACEdpW70tXSO68v6293mp/Uv4+/qbjvQDb373JqUDSvNi+xdNRxFC2MDrr1uF6VdfwfjxptPc3+zZcOwYLFwIrljcbVipIdue2cY/Nv+DT3Z+Qv9Z/YkoG0Gzqs0IKxVGuaBy7Lu8jy1ntnAr+xZVS1VlSuwUJrSagJ+PlD3CcZTW2nSG+2rbtq1OTEw0HUPYzNcpX/PE7CeIGxzHyGYjTcf5DympKTT6sBG/7fRb3nz0TdNxhJMppXZqrduazuFK8nXb8bS2VngzMiApyZrc4U7y8qBpU/Dzg717raPLXSknL4dFhxYxc/9MTt84zaXbl0i9k0rDSg3pXKsznWt15okGTxAaEOraYMJ2ivM1W34kEx7l8ajHaVCxAW9veZsnmz7pVj1vb216i0C/QH7Z8ZemowghbEIpazV6+HCIj4ehQ00n+k8LFlhHfM+e7foCGsDf159hTYYxrMkw199ceD3piRYexUf58KuOv2LXhV18e+pb03G+l5KawvS905nUehJVSlUxHUcIYSODB0NUlNV37E4vHufnwxtvQIMG7lfcC+EKUkQLjzOm+Rgqh1Tmze/cp2XitTWvEewfzO+6/s50FCGEzfj6wquvwu7dsHy56TT/34wZsH8//OEP7tdmIoQrSBEtPE6wfzC/7fRbVh5byboT60zHYdPpTcQfiue3nX4rq9BCiGJ56imIiID//m9rBdi09HT43e+gXTtrnrUQ3kiKaOGRJrebTM0yNXltzWuY3DyrteY3q35DtdBq/PIR6YUWQhRPQAD85S+wa5c1qcO0f/wDzp2Dd94x0wsthDuQT33hkYL9g/ljtz+y7dw2Fh1aZCzHwuSFbDm7hT/F/EnmQgshSmTkSGt29GuvWcdrm3L+PLz1FgwZYs2xFsJbSREtPNa4luNoWKkhv1v7O3Lzc11+/8zcTF5d8yqNKzdmfMvxLr+/EMKz+PhYK7+FK8Cm/P731uEvb71lLoMQ7kCKaOGx/Hz8+EuPv5CcmsyMvTNcfv83NrzB0bSjvNP7HRnwL4RwiC5drGkdf/sbXLjg+vtv2QJffAEvvQR167r+/kK4EymihUcb1HAQ7Wu05/W1r3Mt45rL7rvv0j7e2vQWY5qPoU+9Pi67rxDC8731FmRnWxM7XCkjwzo1MTwc/vhH195bCHckRbTwaEopPnn8E66kX+FXK3/lknvm5ucyMWEi5YPK826fd11yTyGE96hXD155BaZPh6+/dt19f/97OHwYpk6F0qVdd18h3JUU0cLjtarWilc6vcKXe75kxdEVTr/f+1vfJ/F8Ih/0+4CKIRWdfj8hhPf5wx+gRQt49llITXX+/TZvtvqwn3sOHn3U+fcTwg6kiBZe4Q/d/kDDSg2ZtGQSt7JuOe0++y7t4/frfk//qP4MbzLcafcRQni3gABrJTotDX72M+eeZJieDhMmQK1a8Pe/O+8+QtiNFNHCKwT5BTH1iamcuXGG36z6jVPucfXOVQbOHkj54PJM6T8FpZRT7iNESSilKiilVimljhS8L/8Tjy2jlDqnlPqXKzOKh9O8OfzpT7BggfNmR+fnw9ixcPQofP65tHEIcTcpooXXiA6P5r86/hf/3vlvpuyc4tBr5+bnMnLBSM7dOsfC4QsJCw1z6PWFcKBXgTVa6/rAmoJf/5g/A9+6JJUoll//GqKjYdIk2L7d8df/wx9g4UJ4+23o0cPx1xfCzqSIFl7lzUffpF+9fkxeOpnVx1c77Lqvr3mdVcdX8fHjH9OhZgeHXVcIJxgATCv4eBow8H4PUkq1AaoCK12USxSDn59V5IaFQWystWLsKHFx1imJEyfCyy877rpCeAopooVX8fPxY/bQ2TSq3Iihc4eSfCW5RNfTWvPGhjf4++a/M7ntZJ5u9bSDkgrhNFW11hcACt5XufcBSikf4B/AA3uflFKTlFKJSqnEK1euODyseLCqVWH5cqv1ol8/cMT/hmXL4OmnoWtX+OgjkO40IX5IimjhdcoElmHJyCUE+QXRZ2Yf9l7cW6zr5Ot8/mvlf/H7db9nTPMxvNf3PQcnFaJ4lFKrlVIH7vM24CEvMRn4Rmt95kEP1FpP0Vq31Vq3rVy5csmCi2KLioKEBDh7Fnr3htOni3+tqVOhf39o0sTqtw4IcFxOITyJFNHCK0WUi2D56OXk63yiP49mwcEFRXp+Zm4mExMm8u7Wd/l5+5/z5cAv8ff1d1JaIYpGa/2o1rrpfd4WA5eUUtUACt5fvs8lOgIvKqVOAm8DY5VSf3PZX0AUS3Q0xMfD8ePQpg2sW1e052ttHaLyzDPWGLv166FSJWckFcIzSBEtvFbLsJbseHYHzas2Z+i8oby2+jVuZt184POWHl5K04+a8uWeL/ljtz/yXt/38FHyT0nYRgIwruDjccDiex+gtX5Ka11Lax0J/BqYrrV28fl4ojj69oUdO6ByZejVy+ppvn37wc/bssUqwv/3f61xdl9/LZM4hHgQ+c4vvFq10tVYP249T7d8mr9t+hvh74bz21W/5fSN0+i7Bq9eTr/MvKR5xMbFEjsrFn9ff1aNWcX/dP8fGWUn7OZvQC+l1BGgV8GvUUq1VUp9ZjSZcIioKNi2DQYNgv/+b+uY7tdfhxMnrL7pQlevwpIlMGKEVUCfPGm1ckydCv7ywpoQD6S0Mye0l0Dbtm11YmKi6RjCiySeT+Tvm//O/IPzydf5BPoGUq10NQJ8Azh89TAAZQPL8rsuv+MXj/yCAF9pFBQ/Tim1U2vd1nQOV5Kv2+5n61ZrPN3ChVa7hp8f1KhhFcmFkzxCQqxReb/5DYSGms0rhCnF+ZotRbQQ9zh+7TgJKQmcv3We87fOk56TTocaHYiJjKFN9Tb4+fiZjihsQIpo4U6OHYMVK+DMGestIwPatbNWoNu2tQppIbxZcb5mSzUgxD3qlK/Dy4/IUFQhhOeoWxcmTzadQgjPIj3RQgghhBBCFJEU0UIIIYQQQhSRFNFCCCGEEEIUkRTRQgghhBBCFJEU0UIIIYQQQhSRFNFCCCGEEEIUkRTRQgghhBBCFJEU0UIIIYQQQhSRFNFCCCGEEEIUkRTRQgghhBBCFJEU0UIIIYQQ4v+1dz+hUtVhGMe/D96itKJCilJBA7EkCEPCEiKyhVFkm6CgkGjZH4sgrE3bFhG1iCDMFBIjTEgi+oMF7STTIM0ksdBblkZU0sakp8WcC5eLmyP3nN9wfs9nMzNn9bzM3Oe+954zM9FSluiIiIiIiJayREdEREREtCTbpTOck6TTwOHSOXo2H/i9dIieZeY61DjzMtuXlg7Rp/R2FWqbFzJzLVp39kRXSWbBYdsrS4fok6S9mXn4MnMdJO0tnaGA9PbA1TYvZOZanE9n53KOiIiIiIiWskRHRERERLQ0zkv0m6UDFJCZ65CZ65CZ61DbzLXNC5m5Fq1nHts3FkZEREREjKtx/k90RERERMRYGsslWtJaSYclHZG0sXSerklaJOkLSYckHZS0oXSmvkiaI2m/pA9LZ+mDpMsl7ZD0ffN831o6U5ckPdO8pg9I2i7potKZuiBps6STkg5MO3alpM8k/dDcXlEyY5fS2ensoaqts6GO3p6tzh67JVrSHOB14G5gOfCQpOVlU3XuLPCs7RuAVcDjFcw8ZQNwqHSIHr0GfGz7euAmBjy7pAXAU8BK2zcCc4AHy6bqzBZg7YxjG4HdtpcCu5vHg5POTmcPXDWdDVX19hZmobPHbokGbgGO2D5q+wzwLrCucKZO2T5he19z/zSjH9IFZVN1T9JC4B5gU+ksfZB0GXA78BaA7TO2/yybqnMTwMWSJoC5wC+F83TC9pfAHzMOrwO2Nve3Avf3Gqo/6ex09iBV2tlQQW/PVmeP4xK9ADg+7fEkFZTTFEmLgRXAnrJJevEq8BzwX+kgPbkOOAW83ZwO3SRpXulQXbH9M/AycAw4Afxl+9OyqXp1te0TMFq6gKsK5+lKOjudPVRVdTZU39utO3scl2id41gVHyEi6RLgfeBp23+XztMlSfcCJ21/XTpLjyaAm4E3bK8A/mGgp/gBmuvJ1gFLgGuBeZIeLpsqOpDOTmcPVVWdDenttsZxiZ4EFk17vJABnkqYSdIFjMp4m+2dpfP0YDVwn6SfGJ3+vVPSO2UjdW4SmLQ99R+rHYwKeqjuAn60fcr2v8BO4LbCmfr0m6RrAJrbk4XzdCWdnc4eqto6G+ru7dadPY5L9FfAUklLJF3I6IL2XYUzdUqSGF1zdcj2K6Xz9MH287YX2l7M6Dn+3Pag/9q1/StwXNKy5tAa4LuCkbp2DFglaW7zGl/DwN+UM8MuYH1zfz3wQcEsXUpnVyCdDQy/s6Hu3m7d2ROdxjkPts9KegL4hNG7QjfbPlg4VtdWA48A30r6pjn2gu2PCmaKbjwJbGuWjaPAo4XzdMb2Hkk7gH2MPs1gPwP9FixJ24E7gPmSJoEXgZeA9yQ9xugX0wPlEnYnnZ3OHrhqOhvq6e3Z6ux8Y2FEREREREvjeDlHRERERMRYyxIdEREREdFSluiIiIiIiJayREdEREREtJQlOiIiIiKipSzREREREREtZYmOiIiIiGgpS3REREREREv/AytiJWEAiU3BAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Plot the solution against the known answer:\n", "_, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 6))\n", "\n", "ax1.plot(x, special.jv(1, x), c='g')\n", "ax1.set_xlim(0, 10)\n", "ax1.set_title(r'$ \\mathtt{special.jv(1, x)} $')\n", "\n", "ax2.plot(x, sol[:,0], c='g', label='$ y(x) $')\n", "ax2.plot(x, sol[:,1], c='b', label=\"$ z(x) = y'(x) $\")\n", "ax2.set_xlim(0, 10)\n", "ax2.set_title(r'$ \\mathtt{odeint} $ solution')\n", "ax2.legend()" ] } ], "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.8" } }, "nbformat": 4, "nbformat_minor": 1 }