{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Solving for eigenvalues and eigenvectors\n",
    "\n",
    "We saw in the previous notebook how to set up a matrix by integrating between pairs of basis functions.  In this notebook we will solve for the eigenvalues and eigenvectors of the Hamiltonian for a square well, exploring what happens when we use a finite basis set, and introduce changes to the potential in the square well."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Eigenbasis, adding potential\n",
    "\n",
    "We will start by using the eigenbasis (since it's simple) and will add small changes to the potential.  This will give us some insight into the effect of perturbations (and will feed directly into our studies of perturbation theory).  So we'll first set up the bits and pieces of code we'll need, some taken from the previous notebook."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "//anaconda/lib/python2.7/site-packages/matplotlib/font_manager.py:273: UserWarning: Matplotlib is building the font cache using fc-list. This may take a moment.\n",
      "  warnings.warn('Matplotlib is building the font cache using fc-list. This may take a moment.')\n"
     ]
    }
   ],
   "source": [
    "# Import libraries and set up in-line plotting.\n",
    "%matplotlib inline\n",
    "import matplotlib.pyplot as pl\n",
    "import numpy as np\n",
    "# This is a new library - linear algebra includes solving for eigenvalues & eigenvectors of matrices\n",
    "import numpy.linalg as la\n",
    "\n",
    "# Define the eigenbasis - normalisation needed elsewhere\n",
    "def square_well_eigenfunctions(n,width,norm,x):\n",
    "    \"\"\"The eigenbasis for a square well, running from 0 to a (width), sin(n pi x/a). \n",
    "    N.B. requires a normalisation factor, norm.\"\"\"\n",
    "    wavevector = np.pi*n/width\n",
    "    return norm*np.sin(wavevector*x)\n",
    "\n",
    "# We will also define the second derivative for kinetic energy (KE)\n",
    "def second_derivative_square_well_eigenfunctions(n,width,norm,x):\n",
    "    \"\"\"The second derivative of the eigenbasis for a square well, running from 0 to a, sin(n pi x/a)\"\"\"\n",
    "    wavevector = np.pi*n/width\n",
    "    return -wavevector*wavevector*norm*np.sin(wavevector*x)\n",
    "\n",
    "# Define the x-axis\n",
    "square_well_width = 1.0\n",
    "number_of_x_points = 101\n",
    "x = np.linspace(0.0,square_well_width,number_of_x_points)\n",
    "x_spacing = square_well_width/(number_of_x_points - 1)\n",
    "\n",
    "# Integrate two functions over the width of the well\n",
    "# NB this is a VERY simple integration routine: there are much better ways\n",
    "def integrate_functions(function1,function2,dx):\n",
    "    \"\"\"Integrate the product of two functions over defined x range with spacing dx\"\"\"\n",
    "    # We use the NumPy dot function here instead of iterating over array elements\n",
    "    integral = dx*np.dot(function1,function2)\n",
    "    return integral\n",
    "\n",
    "# These arrays will each hold an array of functions\n",
    "basis_functions = []\n",
    "second_derivative_basis_functions = []\n",
    "\n",
    "number_of_basis_functions = 10\n",
    "# Loop over first num_basis basis states, normalise and create an array\n",
    "# NB the basis_array will start from 0\n",
    "for n in range(1,number_of_basis_functions+1):\n",
    "    # Calculate A = <phi_n|phi_n>\n",
    "    integral = integrate_functions(square_well_eigenfunctions(n,square_well_width,1.0,x),square_well_eigenfunctions(n,square_well_width,1.0,x),x_spacing)\n",
    "    # Use 1/sqrt{A} as normalisation constant\n",
    "    normalisation = 1.0/np.sqrt(integral)\n",
    "    basis_functions.append(square_well_eigenfunctions(n,square_well_width,normalisation,x))\n",
    "    second_derivative_basis_functions.append(second_derivative_square_well_eigenfunctions(n,square_well_width,normalisation,x))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We have now set up all the mechanics that we need to create matrices for different Hamiltonians, except for one thing: the potential.  We have an implicit potential already in the infinite square well (we set $V=\\infty$ when $x=0$ or $x=a$).  If we create a potential function, then we can change this and create different Hamiltonian matrices.  So we'll do this next."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXkAAAEACAYAAABWLgY0AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAADh5JREFUeJzt3X+I5PV9x/HnS+8u4K8zURB6xmtrFEmpFdOaCxWcpgXv\nhGgJQn6USIWCpLVNKLSJJeXuj0LIf62V1hzYFAuikEC8NAYtrYtYoljjRRu9RFPxx8VeMXpH1Ety\n0Xf/mLnbdW9v57uzs7MzH58PGG5n57M7H77sPvd775nZTVUhSWrTSeu9AUnS2jHyktQwIy9JDTPy\nktQwIy9JDTPyktSwoZFPcm6S/0jyvSRPJPmzE6y7OcnTSfYmuWT8W5UkrdSGDmt+Afx5Ve1Nchrw\naJL7qmrf0QVJdgDnV9UFST4I3ApsW5stS5K6GnomX1X/W1V7B2+/BjwFbFm07Brg9sGah4HNSc4Z\n814lSSu0opl8kl8GLgEeXnTTFuCFBdf3c/wPAknShHWO/GBU81XgM4MzeknSlOsykyfJBvqB/5eq\nunuJJfuB9y64fu7gfYs/j78oR5JGUFUZ5eO6nsn/E/BkVf3dCW7fA1wHkGQbcLCqDiy1sKq8VLFz\n585138O0XDwWHguPxfKX1Rh6Jp/kt4E/AJ5I8hhQwF8BW/vNrt1VdU+Sq5I8A7wOXL+qXUmSxmJo\n5KvqP4GTO6y7cSw7kiSNja94XSe9Xm+9tzA1PBbzPBbzPBbjkdXOe1Z0Z0lN8v4kqQVJqDV+4FWS\nNIOMvCQ1zMhLUsOMvCQ1zMhLUsOMvCQ1zMhLUsOMvCQ1zMhLUsOMvCQ1zMhLUsOMvCQ1zMhLUsOM\nvCQ1zMhLUsOMvCQ1zMhLUsOMvCQ1zMhLUsOMvCQ1zMhLUsOMvCQ1zMhLUsOMvCQ1zMhLUsOMvCQ1\nzMhLUsOMvCQ1zMhLUsOMvCQ1zMhLUsOMvCQ1zMhLUsOMvCQ1zMhLUsOMvCQ1zMhLUsOMvCQ1zMhL\nUsOMvCQ1zMhLUsOMvCQ1zMhLUsOMvCQ1zMhLUsOGRj7JbUkOJHn8BLdfkeRgku8MLl8Y/zYlSaPY\n0GHNV4C/B25fZs0DVXX1eLYkSRqXoWfyVfUg8OqQZRnPdiRJ4zSumfyHkuxN8s0k7x/T55QkrVKX\ncc0wjwLnVdUbSXYAXwcuPNHiXbt2HXu71+vR6/XGsAVJasfc3Bxzc3Nj+VypquGLkq3AN6rq4g5r\nnwU+UFWvLHFbdbk/SdK8JFTVSGPxruOacIK5e5JzFrx9Gf0fHMcFXpI0eUPHNUnuAHrAWUmeB3YC\nm4Cqqt3AtUk+DRwBDgMfW7vtSpJWotO4Zmx35rhGklZsEuMaSdIMMvKS1DAjL0kNM/KS1DAjL0kN\nM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS\n1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAj\nL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kN\nM/KS1DAjL0kNM/KS1DAjL0kNGxr5JLclOZDk8WXW3Jzk6SR7k1wy3i1KkkbV5Uz+K8CVJ7oxyQ7g\n/Kq6ALgBuHVMe5MkrdLQyFfVg8Cryyy5Brh9sPZhYHOSc8azPUnSaoxjJr8FeGHB9f2D90mS1tmG\nSd/hrl27jr3d6/Xo9XqT3oIkTbW5uTnm5ubG8rlSVcMXJVuBb1TVxUvcditwf1XdNbi+D7iiqg4s\nsba63J8kaV4SqiqjfGzXcU0Gl6XsAa4bbGQbcHCpwEuSJm/ouCbJHUAPOCvJ88BOYBNQVbW7qu5J\nclWSZ4DXgevXcsOSpO46jWvGdmeOayRpxSYxrpEkzSAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS\n1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAj\nL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kN\nM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS\n1DAjL0kNM/KS1LBOkU+yPcm+JD9I8rklbr8iycEk3xlcvjD+rUqSVmrDsAVJTgJuAX4X+BHwSJK7\nq2rfoqUPVNXVa7BHSdKIupzJXwY8XVXPVdUR4E7gmiXWZaw7kySt2tAzeWAL8MKC6y/SD/9iH0qy\nF9gP/EVVPTmG/WmBL38ZXnpp6dsuvRSu9v9RkhbpEvkuHgXOq6o3kuwAvg5cuNTCXbt2HXu71+vR\n6/XGtIW2VcGNN8JNN8FJi/7/tX8/3HuvkZdaMTc3x9zc3Fg+V6pq+QXJNmBXVW0fXP88UFX1pWU+\n5lngA1X1yqL317D709LeeAPOPrv/72JPPgnXXtv/V1J7klBVI43Eu8zkHwHel2Rrkk3Ax4E9izZw\nzoK3L6P/w+MVNDYHD8LmzUvftnlz/3ZJWmzouKaq3kxyI3Af/R8Kt1XVU0lu6N9cu4Frk3waOAIc\nBj62lpt+Jzp0aPnIHzo02f1Img1DxzVjvTPHNSN76CH47Gf7/y5WBZs29Uc5GzdOfm+S1tZaj2s0\nBZYb1yRwxhmezUs6npGfEcuNa8CRjaSlGfkZcfAgnHnmiW8/80wffJV0PCM/IzyTlzQKIz8jjLyk\nURj5GeG4RtIojPyM8Exe0iiM/Iw4dGj4mbyRl7SYkZ8Ryz1PHvzVBpKWZuRnhOMaSaMw8jPCcY2k\nURj5GeG4RtIojPwMeOsteO01OP30E69xXCNpKUZ+BvzkJ3DaaXDyySde47hG0lKM/AwYNqoBxzWS\nlmbkZ8CwZ9bA/LjGX9cvaSEjPwOG/UoDgHe9qz/OOXx4MnuSNBuM/AzociYPPvgq6XhGfgYMe478\nUT74KmkxIz8DujzwCj74Kul4Rn4GOK6RNCojPwMc10galZGfAY5rJI3KyM8AxzWSRmXkZ4DjGkmj\nMvIzwHGNpFEZ+RnguEbSqIz8DOjyaw2gv8YzeUkLGfkZ4Jm8pFEZ+Sn385/DkSNwyinD1/rAq6TF\njPyUO3oWnwxf6wOvkhYz8lOu66gGHNdIOp6Rn3JdnyMPcMYZ/b8F+9Zba7snSbPDyE+5rs+Rh/4f\nDTnllP7fhJUkMPJTbyXjGnBkI+ntjPyUW8m4BnyGjaS3M/JTbiXjGvAZNpLezshPOcc1klbDyE+5\nrr/S4Ch/tYGkhYz8lPNMXtJqGPkp5wOvklbDyE85H3iVtBpGfso5rpG0GkZ+yjmukbQaRn7KOa6R\ntBpGfopVOa6RtDqdIp9ke5J9SX6Q5HMnWHNzkqeT7E1yyXi3+c50+DBs3AibNnX/GMc1khYaGvkk\nJwG3AFcCvwZ8IslFi9bsAM6vqguAG4Bb12CvTZmbmxu6ZqWjGpjNcU2XY/FO4bGY57EYjy5n8pcB\nT1fVc1V1BLgTuGbRmmuA2wGq6mFgc5JzxrrTxnT5Al7pqAZmc1zjN/M8j8U8j8V4dIn8FuCFBddf\nHLxvuTX7l1ijFVrprzQAOPVU+NnP+n8XVpI2TPoOP/KRSd/jdPr+9+HRR5df8/LLK498Au9+d/84\nb9w4+v4mqcuxeKfwWMzzWIxHqmr5Bck2YFdVbR9c/zxQVfWlBWtuBe6vqrsG1/cBV1TVgUWfa/k7\nkyQtqaoyysd1OZN/BHhfkq3AS8DHgU8sWrMH+BPgrsEPhYOLA7+aTUqSRjM08lX1ZpIbgfvoz/Bv\nq6qnktzQv7l2V9U9Sa5K8gzwOnD92m5bktTF0HGNJGl2rckrXn3x1LxhxyLJJ5N8d3B5MMmvr8c+\nJ6HL18Vg3W8lOZLko5Pc3yR1/B7pJXksyX8nuX/Se5yUDt8jZyTZM2jFE0n+cB22ueaS3JbkQJLH\nl1mz8m5W1Vgv9H9wPANsBTYCe4GLFq3ZAXxz8PYHgYfGvY9puHQ8FtuAzYO3t7+Tj8WCdf8O/Cvw\n0fXe9zp+XWwGvgdsGVw/e733vY7H4ibgi0ePA/BjYMN6730NjsXlwCXA4ye4faRursWZvC+emjf0\nWFTVQ1V19OVLD9Hu6wu6fF0A/CnwVeD/Jrm5CetyLD4JfK2q9gNU1csT3uOkdDkWBZw+ePt04MdV\n9YsJ7nEiqupB4NVllozUzbWIvC+emtflWCz0R8C31nRH62fosUjyS8DvV9U/Ai0/E6vL18WFwHuS\n3J/kkSSfmtjuJqvLsbgFeH+SHwHfBT4zob1Nm5G6OfEXQ2lpSX6H/rOSLl/vvayjvwUWzmRbDv0w\nG4BLgQ8DpwLfTvLtqnpmfbe1Lq4EHquqDyc5H/i3JBdX1WvrvbFZsBaR3w+ct+D6uYP3LV7z3iFr\nWtDlWJDkYmA3sL2qlvvv2izrcix+E7gzSejPXnckOVJVeya0x0npcixeBF6uqp8CP03yAPAb9OfX\nLelyLK4HvghQVT9M8ixwEfBfE9nh9Bipm2sxrjn24qkkm+i/eGrxN+ke4Do49oraJV881YChxyLJ\necDXgE9V1Q/XYY+TMvRYVNWvDi6/Qn8u/8cNBh66fY/cDVye5OQkp9B/oO2pCe9zEroci+eA3wMY\nzKAvBP5norucnHDi/8GO1M2xn8mXL546psuxAP4aeA/wD4Mz2CNVddn67XptdDwWb/uQiW9yQjp+\nj+xLci/wOPAmsLuqnlzHba+Jjl8XfwP884KnFv5lVb2yTlteM0nuAHrAWUmeB3YCm1hlN30xlCQ1\nzD//J0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1DAjL0kNM/KS1LD/B9yLB/TrDtUGAAAAAElFTkSu\nQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x108f10a10>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Define constants - here I use atomic units\n",
    "hbar = 1.0\n",
    "hbar_squared = hbar**2\n",
    "m_e = 1.0 # Mass of the electron\n",
    "two_m_e = 2.0*m_e\n",
    "# These are not needed but are easier\n",
    "pi_squared = np.pi**2\n",
    "square_well_width_squared = square_well_width**2\n",
    "\n",
    "# Define the potential in the square well\n",
    "def square_well_perturbing_potential(x,perturbing_potential,potential_magnitude):\n",
    "    \"\"\"Potential for a particle in a square well, expecting two arrays: x, V(x), and potential magnitude\"\"\"\n",
    "    # Zero the array\n",
    "    perturbing_potential[:] = 0.0\n",
    "    # Let's define a small bump in the well\n",
    "    bump_position = x.size/4\n",
    "    perturbing_potential[bump_position-1] = potential_magnitude\n",
    "    perturbing_potential[bump_position] = potential_magnitude\n",
    "    perturbing_potential[bump_position+1] = potential_magnitude\n",
    "    # Plot to ensure that we know what we're getting\n",
    "    pl.plot(x,perturbing_potential)\n",
    "    pl.ylim((0,2.0))\n",
    "    \n",
    "# Declare space for this potential (Vbump) and call the routine\n",
    "bump_potential = np.linspace(0.0,square_well_width,number_of_x_points)\n",
    "square_well_perturbing_potential(x,bump_potential,0.5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Solving for the eigenvalues\n",
    "\n",
    "Now that we have the potential, we need to build the matrix (remember that we have set $\\hbar = m = 1$):\n",
    "\n",
    "$H_{mn} = \\langle \\phi_m \\vert -\\frac{1}{2}\\nabla^2 + \\hat{V}\\vert\\phi_n\\rangle$\n",
    "\n",
    "By diagonalisation, we will find the eigenvalues.\n",
    "\n",
    "It's worth noting that we act with $\\hat{V}$ by multiplication in position representation:\n",
    "\n",
    "$$\\langle x \\vert \\hat{V} \\vert \\phi_m\\rangle = V(x) \\phi_m(x)$$\n",
    "\n",
    "Then we make matrix elements by integration in position representation; as we've defined a grid on which to represent the system, we'll just sum over the values on the grid (though there are much more accurate ways to do numerical integration).\n",
    "\n",
    "$$H_{ij} = \\langle \\phi_i \\vert \\hat{H} \\vert \\phi_j \\rangle = \\int dx \\phi_i^\\star(x) \\hat{H} \\phi_j(x)$$\n",
    "\n",
    "We will print out the matrix of _just_ the potential first: $\\langle \\phi_m \\vert \\hat{V} \\vert \\phi_n\\rangle$.  This will be useful when thinking about perturbation theory (the change in eigenvalues should be the same as the diagonal elements of the potential matrix, to first order in the potential)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Potential matrix elements:\n",
      "   0.015    0.021    0.015   -0.000   -0.015   -0.021   -0.015    0.000    0.015    0.021\n",
      "   0.021    0.030    0.021    0.000   -0.021   -0.030   -0.021   -0.000    0.021    0.029\n",
      "   0.015    0.021    0.015    0.000   -0.015   -0.021   -0.015   -0.000    0.014    0.020\n",
      "  -0.000    0.000    0.000    0.000    0.000   -0.000   -0.000   -0.001   -0.000    0.000\n",
      "  -0.015   -0.021   -0.015    0.000    0.015    0.021    0.014   -0.001   -0.015   -0.020\n",
      "  -0.021   -0.030   -0.021   -0.000    0.021    0.029    0.021    0.000   -0.020   -0.029\n",
      "  -0.015   -0.021   -0.015   -0.000    0.014    0.021    0.015    0.001   -0.014   -0.020\n",
      "   0.000   -0.000   -0.000   -0.001   -0.001    0.000    0.001    0.001    0.001   -0.000\n",
      "   0.015    0.021    0.014   -0.000   -0.015   -0.020   -0.014    0.001    0.015    0.020\n",
      "   0.021    0.029    0.020    0.000   -0.020   -0.029   -0.020   -0.000    0.020    0.028\n",
      "\n",
      "Full Hamiltonian\n",
      "   4.950    0.021    0.015   -0.000   -0.015   -0.021   -0.015    0.000    0.015    0.021\n",
      "   0.021   19.769    0.021   -0.000   -0.021   -0.030   -0.021   -0.000    0.021    0.029\n",
      "   0.015    0.021   44.428    0.000   -0.015   -0.021   -0.015   -0.000    0.014    0.020\n",
      "  -0.000   -0.000    0.000   78.957    0.000   -0.000   -0.000   -0.001   -0.000   -0.000\n",
      "  -0.015   -0.021   -0.015    0.000  123.385    0.021    0.014   -0.001   -0.015   -0.020\n",
      "  -0.021   -0.030   -0.021   -0.000    0.021  177.682    0.021   -0.000   -0.020   -0.029\n",
      "  -0.015   -0.021   -0.015   -0.000    0.014    0.021  241.820    0.001   -0.014   -0.020\n",
      "   0.000   -0.000   -0.000   -0.001   -0.001   -0.000    0.001  315.829    0.001    0.000\n",
      "   0.015    0.021    0.014   -0.000   -0.015   -0.020   -0.014    0.001  399.734    0.020\n",
      "   0.021    0.029    0.020   -0.000   -0.020   -0.029   -0.020    0.000    0.020  493.508\n"
     ]
    }
   ],
   "source": [
    "# Declare space for the matrix elements - simplest with the identity function\n",
    "H_matrix = np.eye(number_of_basis_functions)\n",
    "\n",
    "# Define a function to act on a basis function with the potential\n",
    "def add_potential_acting_on_state(H_on_phi,V,phi):\n",
    "    \"\"\"Add V(x)phi(x) onto an input array, H_on_phi\"\"\"\n",
    "    for i in range(V.size):\n",
    "        H_on_phi[i] = H_on_phi[i] + V[i] * phi[i]\n",
    "        \n",
    "print \"Potential matrix elements:\"\n",
    "# Loop over basis functions phi_n (the bra in the matrix element)\n",
    "# Calculate and output matrix elements of the potential\n",
    "\n",
    "for n in range(number_of_basis_functions):\n",
    "    # Loop over basis functions phi_m (the ket in the matrix element)\n",
    "    for m in range(number_of_basis_functions):\n",
    "        # Act with H on phi_m and store in H_phi_m\n",
    "        H_on_phi_m = np.zeros(number_of_x_points)\n",
    "        add_potential_acting_on_state(H_on_phi_m, bump_potential, basis_functions[m])\n",
    "        # Create matrix element by integrating\n",
    "        H_mn = integrate_functions(basis_functions[n],H_on_phi_m,x_spacing)\n",
    "        # The comma at the end prints without a new line; the %8.3f formats the number\n",
    "        print \"%8.3f\" % H_mn,\n",
    "    # This print puts in a new line when we have finished looping over m\n",
    "    print\n",
    "\n",
    "print\n",
    "print \"Full Hamiltonian\"\n",
    "# Loop over basis functions phi_n (the bra in the matrix element)\n",
    "# Calculate and store the matrix elements for the full Hamiltonian\n",
    "for n in range(number_of_basis_functions):\n",
    "    # Loop over basis functions phi_m (the ket in the matrix element)\n",
    "    for m in range(number_of_basis_functions):\n",
    "        # Act with H on phi_m and store in H_phi_m\n",
    "        # First the kinetic energy\n",
    "        H_on_phi_m = -(hbar_squared / two_m_e) * second_derivative_basis_functions[m] \n",
    "        # Now the potential\n",
    "        add_potential_acting_on_state(H_on_phi_m, bump_potential, basis_functions[m])\n",
    "        # Create matrix element by integrating\n",
    "        H_mn = integrate_functions(basis_functions[n], H_on_phi_m, x_spacing)\n",
    "        H_matrix[m,n] = H_mn\n",
    "        # The comma at the end prints without a new line; the %8.3f formats the number\n",
    "        print \"%8.3f\" % H_mn,\n",
    "    # This print puts in a new line when we have finished looping over m\n",
    "    print"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Notice that two things have changed compared to the perfect square well: first, the diagonal elements are _slightly_ larger; second, there are now off-diagonal elements.  Does it surprise you that these alternate (i.e. only in odd row and columns) ? Think about the symmetries of the system, particularly of the wavefunctions.\n",
    "\n",
    "What effect will this have on the eigenvalues and eigenvectors ? We'll diagonalise, print out the eigenvalues and plot the first few eigenvectors (as well as looking at their coefficients to get a rough idea of how much they've changed)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[   4.9498   19.7691   44.4282   78.9571  123.3851  177.6822  241.8203\n",
      "  315.8286  399.734   493.5083]\n",
      "[  1.0000e+00   1.4275e-03   3.7820e-04  -7.5022e-07   1.2595e-04\n",
      "  -1.2134e-04   6.2036e-05   3.5439e-07   3.7203e-05   4.2008e-05]\n",
      "[ -1.4279e-03   1.0000e+00   8.5641e-04   3.3707e-09   2.0269e-04\n",
      "  -1.8751e-04   9.3925e-05  -1.4036e-09   5.4284e-05   6.1205e-05]\n",
      "[ -3.7702e-04  -8.5703e-04   1.0000e+00   4.8328e-06   1.8593e-04\n",
      "  -1.5688e-04   7.5635e-05  -1.2205e-06   4.0245e-05   4.5583e-05]\n"
     ]
    }
   ],
   "source": [
    "# Solve using linalg module of numpy (which we've imported as la above)\n",
    "eigenvalues, eigenvectors = la.eigh(H_matrix)\n",
    "# This call above does the entire solution for the eigenvalues and eigenvectors !\n",
    "# Print results roughly, though apply precision of 4 to the printing\n",
    "np.set_printoptions(precision=4)\n",
    "print eigenvalues\n",
    "print eigenvectors[0]\n",
    "print eigenvectors[1]\n",
    "print eigenvectors[2]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can see that the eigenvalues look close to the perfect well values (we'll check them properly below).  The eigenvector coefficients show a single dominant value (corresponding to the unperturbed eigenvector), with very small contributions from other eigenvectors.  Now print the eigenvalues and calculate the change."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " Changed Original  Difference\n",
      "   4.950    4.935    0.015\n",
      "  19.769   19.739    0.030\n",
      "  44.428   44.413    0.015\n",
      "  78.957   78.957    0.000\n",
      " 123.385  123.370    0.015\n",
      " 177.682  177.653    0.029\n",
      " 241.820  241.805    0.015\n",
      " 315.829  315.827    0.001\n",
      " 399.734  399.719    0.015\n",
      " 493.508  493.480    0.028\n"
     ]
    }
   ],
   "source": [
    "# Now print out eigenvalues and the eigenvalues of the perfect square well, and the difference\n",
    "print \" Changed Original  Difference\"\n",
    "for i in range(number_of_basis_functions):\n",
    "    n = i+1\n",
    "    print \"%8.3f %8.3f %8.3f\" % (eigenvalues[i],n*n*np.pi*np.pi/2.0,eigenvalues[i] - n*n*np.pi*np.pi/2.0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compare the differences in eigenvalues due to the small potential we've added to the diagonal terms of the potential matrix.  How close is the agreement ? You might like to change the magnitude of the potential (higher up - it's passed as an argument to the potential function) and see how this agreement changes.\n",
    "\n",
    "Now we'll plot the eigenvectors (after building them) and look at the change with respect to the eigenvectors of the original system.  Remember that any function in this space can be written as:\n",
    "\n",
    "$$\\vert\\psi\\rangle = \\sum_i c_i \\vert \\phi_i \\rangle$$\n",
    "\n",
    "We'll use this to build the eigenfunctions of the perturbed system.  In this case, $c_i$ are the coefficients in the eigenvector of the matrix, and $\\vert\\phi_i\\rangle$ are the basis functions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAtAAAADSCAYAAABjAPe+AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzsnXd4FFUXh3+X3lsSEiAkEJHepSoI0qugFFGKggIWsEsR\nFERR6aA0CyIg0ouo4IfSJAQSIBBTCZBKQhJCKunJnu+PM4El7Ca72TLZ5L7PM08yM7ecmZ05c+69\n554riAgSiUQikUgkEonEMMqpLYBEIpFIJBKJRGJLSANaIpFIJBKJRCIxAmlASyQSiUQikUgkRiAN\naIlEIpFIJBKJxAikAS2RSCQSiUQikRiBNKAlEolEIpFIJBIjkAa0igghXhJC/KW2HLaAEOI5IUSE\nECJFCNHBgvX0EUJEWqp8c2ErckokpQkhxCYhxAKt/TeEEDGKXqorhHhKCBGs7D+rpqylBSFEfSHE\nv0KIZCHEChXleOi3l0iEjANtWYQQYQDqA8gFIAAQgJ+J6G015TIVIcTLAF4jot5Wqu8GgHeJ6A8L\n19MHwA4icrFkPaZiTTmFEK4AQgFUICKNpeuTSNRAS1fnAMgDEABgB4DvSceHUghRAUAKgG5E5Kcc\n+wfAYSJaby25zY0QYiuASCL6VG1ZAEAIsRBARyIaa2I5Jeq6rElZvnZLUkFtAcoABGA4EZ1SWxAz\nk98YKF5mIcoTUZ4RWVzBH7Ti1FVOGn4mkf9bC7UFkUgsyH1dLYSoCaAPgG8AdAcwTUd6JwCVAQRq\nHTNFTxmrE0s9QggBE+6pxDzIb6geiEhuFtzAPXf99Jx7GcBZrf1BAIIAJALYAOA0gGla56eBFcld\nAMcAuGid0wCYCSAYQAKA9crxSkp5rbXS2gNIB2Cv7I8AcEVJ5w6gnVZaZwAHAMQBuAP+oLQEkAHu\nqUkFkKCkrQVgu5I2FMCCAtfqDmA1gHgASwA8plxjkpJnl457VEmpIw/APQDXleOtAJxSZPYFMFIr\nz1YAGwH8qeR95P4DqAvgJwBRyv08qBzvAyASwPsAYpXzr2jlGwbAG0AygHAAi7TOuSq/wxTlXByA\nj7XOVwGwTfl9/AF8BO4VyD/fAMB+Jd9NALML5P1ZyesH4EMAEYU8d2sU+ZMB+ABoDaALgBgoI09K\nuucBXFX+7wbgopLnNoCVyvFw5f6ngnvcuhv4PL4Bfh6Tld/bDcA55ffeDe7RVv0dlZvciHTragBd\nlWe/tbK/VXmWH1f0UZ7yTvwD4Iayn64cqwjWiT8CiFb0yud4MPL7iE5Ujhut57XOT1fypih6oqNy\nXK9u0ZE/G0CmUsZvynG9+lZHGacAfAnAU3n3DwGoo3W+h6IHEsHfnT4F8n6h3Jc08PckG0CWIk8/\ncEN+nnK/7yi6RLv8Xlrlh4P1sc7r0iF7SwDHlXsfCGCc1rmt+b+Rsj9H+V1vAXhV+W3clHOVAKxU\n6r8N/h5VVs7p/caAdfBtPKyjnwPgo/xvtmsv7DeFjm8ogKHg71ZKvvxqv7Nqb6oLUNo3FG1A/6v8\nb68om1Fg3/S3FaUxTTk/Cqw0myvnPwZwTqssDYAjAGoCaAxWlIOUcz8C+Fwr7ZsAjir/d1Je4i7K\nyzlZkbmiUs9VRRFUUZTCkwVl1yp3O1hZVgMbk9cATNVKn6PUXU4p71cA85Xz98vWc680AJoq/1cA\ncB3AXOX/Z5SX+nHl/FZFKfTIL1tHeX8C2AX+wJUH0Fs53keRc5FyfChYkddWzj8NoI3yf1uwsntW\n2c83oL9Trqc9WGm1UM5/DVZYtQA0BBu2Eco5AeASgAVKvU3ASnKgVt4zAGoDaARWeDoNaHBD7CKA\nmsp+CwCOyv9+AAZrpT0Ido0BAA8AE5X/q4GHpvOvKw8PK3VDnsdDAKqDFXUmgL+VsmqCFfFktd9P\nucktf4MeXQ02RGYq/983ovS8F6EAntHaPwQ2RKqAdfwFANOVcwV1YmUD3yt9en4c2LDprOy7KWkK\n1S06rregoViovtWR/5QiRysAVcGG+w7lXCNwY2Gwst9f2bfTyhsGNmTLKfIWlOcdRVc1AH+nNgH4\nVes3SQEwXslbF0B7XdelQ+5qACLARqcA0AFspLbU8dsPARvPLZXfdofyLOQb0GsAHAbr6+oAfgOw\nVDlX1DfmOoD+WnLtBfCROa+9qN8Uj35DKyvXm//9rw2lcVaWN9UFKO0bWKGmgHsLEpW/ryrntA3o\nydBSlMqxCDwwoI9CMUaV/XLKS9dY2dcA6Kl1fg+AOcr//QHc0DrnjgeG0kYAnxWoNwhAb3BPQSyA\ncjqu6yEDWpEnC4qxqBybAeCkVvqwAmVsA7AZQCMD7qN2674XgOgC538F8Kny/1awn7m+spzAPum1\ndJzro9zXclrHYqEYkzrSrwGwSvk//4PaQOu8J4Dxyv83AQzQOvcqHhjQ3XXcn3kAtmjlHah1bjr0\nG9DPKL9hd2h93JVzcwD8ovxfT7nW+sr+abBStyuQJ/+6tO+JIc9jD63zl6B8BJT9lQBWq/1+yk1u\n+Rv0G9Dn8aChr8uALqerDLA/dSaUnkfl2IQidKIpev4v6OhZBvdqFqznvm7Rkb6gsVWovtWR/xSA\nL7X28xvQQtE/2wqk/wtKY1rJu7gIeQLwcCOlAbiHtZxyXQcMuS4d58cDOFPg2GYAn+j47bdAMYiV\n/cfw8DfqHpQOH2W/J4AQ5f9CvzHgUYp8vV9TKcvZnNde1G8KHd9QcMNmOpSOGbmRjMJhJUYRUT0i\nqqv83aIjTUNwq12bW1r/uwJYJ4RIEEIkgIeYCNyizydW6/90ADWU/08BqCqE6KpMCOsAbh3nl/tB\nfrlCiESw20ZDcO9FOBnm+2QPbslGaB0LLyBfwev7CPziewkhfIUQUw2oB9B9r4qqS5vGYLeTFD3n\n7xa45vv3UgjRXQhxUggRJ4RIAg+n2hfIr+93aIiHf1NtGV0ANCrwO8wHf4R15Q3Xd3HE/vbrwW5A\nsUKIzYpPJwD8AmCEEKIq+IPxLxHFKedeBfdWBwkhPIUQw/XVAcOexzit/zPw8H3JwIP7IpGUZBqB\nOz6MxRXcS3hb653ejIf1RUE9ZYqebwxuaOuSozDdUhSG6NuCaKcPB98He0WW8QVkeQrcqaErry5c\nARzSukcB4B5dR+i/B4bgCqBHAdleUsotSMF7cv9/IYQDuDf7spaMxwDYaaXX+40BG7LPCSEqgl3s\nLhNRvu4317UX5xs6BsBwAOFCiFNCiB4G1lVqkZMIrYMhk69uAygY9shZ6/9IAF8Q0S5jKycijRBi\nL1gZxAL4g4jStMpdSkRfPSI0vyAueiYQUIH9ePCL7Aru/YTyf5S+PIrhNkOp6ykA/wghzhBRSBGX\nFA1WFtq4gF1G9MmnTSSAekKIWoUY0frYCfYDH0xEOUKINXhYMRbGbfBvmn9/tCNoRIJ7KFroyZt/\nzfkTllwLq4g4CsB6IYQ9gH1gn+lFRBQthDgPVoaTwCMQ+Xlugp8RCCHGANgvhKgH3fcyAsV8HiUS\nW0EI0RVsbJwtRvZIcM+rHSldeDooeNyU9yoS3BOq63hhuqUomQzRtwXRTu8K/jbEK7JsJ6KZRtRf\nkPyR2fMFTyihPbsVs9xIAKeJaHAR6YAHujwfbV0eDzaI2xDRbQPKeggiChRChIPn27wINqjzMde1\nG/0NJaLLAEYLIcoDmA12LSnR0aosjeyBLjn8CaCtEOJZIUR5IcQsPNzy3QzgYyFEawAQQtQWQhgT\n1mcXgBfABpL2C/kDgNeFEN2UcqsLIYYJIaoD8AIriq+FENWEEJWFEE8q+WIBOCutZCgG9l4AS4UQ\nNZSe7vfAvmE6EUKMFULkt3iTwENghvR2ewJIF0LMEUJUEEL0BU+ENOijQ0Qx4B6BjUKIOkoZhobj\nqwEgUTGeu0ExOLUorLG0F8B8pc5GAN7SOucFIFW5pirKM9BGCNFFOb9PK68zgFn6KhFCdBFCdFPC\nbGWAP+La93UHeCi1LdgHOj/fRMXgBtgfn5R8d5S/2h/n72Da8yiRlFiEEDWFEPk6ZQcR6YsCofd9\nV/TMcQBrlPKEEMJNCPF0IVWb8l79COBDIURnJe9jQojGKFq3FCQW7D+djz59u7sQWSYJIVoKIaoB\n+AzAPqUR8QuAkUKIQUKIcoo8fYQQDQ28RoDv0ZdCCBflOh3Eg5jbOwH0V74t5YUQ9cSDdQMKXldB\n/gDQXAgxSbnOioou1dXw2AtgqtY1LoRicCrX+QOAtUpvNIQQjYQQg4y4xl/B/s69wbrf3Ndu1DdU\nuRcvKZ1O+RPKy3zEGGlAW4ffBQfWz98OFExARHfBk0BWgFuwLcF+o1nK+cPgiWS7BbsO/AeeyHC/\niIJFFijfC+x31QBsPOYfvwz2a1qvDAkFg33z8o3ikeAZ5xHgFvp4JetJ8ESwGCFE/lD92+CWdwiA\nf8G+tlsLuS9dAXgKIVLALiVvE1GYnrT3r4eIchS5hoHv1XqwD911PfdCF5PBftBBYOXyTiFptct7\nE8DnQohksNLcU0jagvtLwD3yoeAP6z48+H01YAXWUTkfB1bCtZS8n4F/g1Cwz+D2QuStpeRNUNLH\ng5+rfA6Be4UOElGm1vEhAPyV32MNgBeIKIuIMgAsBXBO8NBhN1OfR4mkhPK78m5HgN0cVkJ3CLt8\ninrOp4AnFAeA38d9eNhd4eHMJrxXRLQf/J7+qrzDhwDUM0C3FGQLgDbKu36wEH0brO86wI30beCe\nzkpQ9KviijAKPDnyDtht4EM8sEV06YmCx9aBJ+UdV34rDyg9r0QUqcj5Ifh+XwFP5n7kuh6phOge\neAL2BEXuaPBvUVlH2r/AI5GnwN/M/B7hLOXvXPBEzQvK73gcPDFUHwWvcTd4wvoJItJ2HzLLtRfz\nGzoZQKhyPTPwaOdRmcPkhVSU3rDt4N5SDYAfiOgbHem+wYPZpq8Q0VWTKi7lCCEE2Of1JSI6o7Y8\nEvMjhHgdbKQ+o0LdNwDMIKKT1q5bYn6EEEMArAUbIluIaJmONDp1cFF5hRAfgBtg9gU+5hLJIwgh\nToF77X9SWxZrIYRoCY6MVNnAOUOSUoA5eqBzwfEA24Bnmr6lPEz3EUIMBfAYET0OnnS12Qz1ljqU\nYa3aQojK4JBDAIc9kpQChBBOQognlaHcFgA+gJYLhRXlGANAI43n0oEQohy4B2kwgDYAXjRUBxeV\nV+kgGYhCJq1KJGURIcRoIUQlIURdAMsAHJHGc9nCZAOaiGLyezKUIZBAPDo7dxSUIWci8gRQWwih\na2ZrWacneBZtHHi26ygiyio8i8SGqAT2YctfeOEQOI6n1VB6hzaAXVEkpYNu4AWGwpWh2d1gnauN\nPh1cVN414Gg5EomhlBV3rZngb/V1PIjnLSlDmDUKhxCiCdjPyrPAqUZ4OCRKlHIsFpL7ENFnYF9X\nSSmEiCIAtFNZBqu7i0gsTkH9eguPzsbXlaZRYXmVyUmRROTLHmUSSdEQUT+1ZbAGRDRUbRkk6mI2\nA1oIUQO84tA7Sk90ccspK61XiURSCiEiW7A2C5VRcJzwj8HuG4XmkTpbIpHYMsXV2WaJwiE4XFb+\ncp2/6UgShYdjDjrj4fjAD2GOFWJsZVu0aJHqMshrltcsr9k8m0pE4eF4rLr0qz4drC/vY+Aln32E\nEKHK8ctCCJ2Lb6h93+WzLa9ZXq+85uJspmCuMHY/AQggonV6zh8Bh/PJX5wjiYik+4ZEIpGYzkUA\nzYQQrkKISuAwXEcKpNGng3XmJSI/InIiIjciagp27ehED1atlEgkkjKNyS4cgleQmwjAVwhxBTyB\n4GNwnFkiou+J6KjgxTlugEMoGbpks0QikUgKgYjyBC+8dBwPQtEFCiFmoggdrC+vrmpg2IqqEolE\nUiYw2YAmonMAyhuQTu/KaWWZvn37qi2C1ZHXXDYoi9esFsQLO7QocOy7Avs6dbCuvDrSFLaCW5mj\nLD7bZe2ay9r1AmXzmk3B5IVUzI0QgkqaTBKJRGIIQgiQbUwiNBtSZ0skElvFFJ0tl/KWSCQSiUQi\nkUiMQBrQEolEIpFIJBKJEUgDWiKRSCQSiUQiMQJpQEskEolEIpFIJEYgDWiJRCKRSCQSicQIpAEt\nkUgkEolEIpEYgTSgJRKJRCKRSCQSIzB5IRWJRCKRSCQSXcSnxyM6NRrx6fGIT49HRk4GcjW5yKM8\nVCxXEVUrVkW1itVQt0pd1K9eH/Wr10edKnUgRJkKpy6xQaQBLZFIJBKJxGTupt/FqbBTOBV6Cj6x\nPgiKD0Ie5cG5ljPsq9nDrqodqleqjvKiPMqL8silXKTnpCMtOw1JmUmIS4tDbFos8jR5aFq3KZrW\naYoWdi3QwakD2ju2Ryv7VqhYvqLalymRAJArEUokEonZkCsRSsoaSZlJ+NX3V2zz2YbAO4Ho5dIL\n/Zv2R+cGndHKoRUcqzsa3ZuckpWC0MRQ3Ey8iaD4IPjE+sAnxgdRqVHo3qg7nnZ9Gv2b9kcP5x4o\nX668ha5MUhYwRWdLA1oikUjMhDSgJWWFmwk3seTfJfgt6DcMbjYY0zpOwzNNn0Gl8pUsVmdiRiLO\nRZ7Dv+H/4vjN44hOjcbw5sMxttVYDG42GBXKyUF1iXGobkALIbYAGAEgloja6zjfB8BvAEKUQweJ\n6As9ZUllLJFIbBK1DGghxBAAa8ETw7cQ0TIdab4BMBRAGoBXiOhqYXmFEEsAjAKgARCr5InRUa7U\n2WWImHsx+PzM59jjvwfvdH8Hb3Z9E3bV7FSRJSwpDEeuHcEuv12ISI7Ayx1exmudX4NbXTdV5JHY\nHiXBgO4F4B6A7YUY0B8Q0bMGlCWVsUQisUnUMKCFEOUABAPoDyAawEUAE4goSCvNUACziGi4EKI7\ngHVE1KOwvEKIGkR0T8k/G0BrInpDR/1SZ5cBiAjbfLbho78/wuT2k/Fx749hX81ebbHu4x/njy1X\ntmC7z3YMfGwg5j01Dx2cOqgtlqSEY4rONksYOyJyB5BYRLIyNawpkUgkVqIbgOtEFE5EOQB2g3uO\ntRkFYDsAEJEngNpCCMfC8uYbzwrVwT3RkjLI3fS7GLdvHFafX40TU05g9eDVJcp4BoA29dtg9eDV\nCHknBJ2dOmPozqEYtXsUrsVfU1s0SSnFmnGgewohrgoh/hRCtLZivRKJRFKaaQQgUmv/lnLMkDSF\n5hVCfCGEiADwEoBPzSizxEa4GnMVnb7rBNfarvCa7oX2jo8MMpcoalWuhY+e+ggh74TgaZen0Wtr\nL7x97G3cTb+rtmiSUoa1PO4vA3AhonRlKPEwgOb6Ei9evPj+/3379kXfvn0tLZ9EIpEYzenTp3H6\n9Gm1xSgOBo0IEtFCAAuFEHMBzAawWFc6qbNLJ3/f/BsTD07ExuEbMbb1WLXFMYoqFarggyc/wJQO\nU7Do9CK02dgG64asw/g242WM6TKMOXW22aJwCCFcAfyuywdaR9pQAE8QUYKOc9KfTiKR2CQq+UD3\nALCYiIYo+/MAkPZEQiHEZgCniGiPsh8EoA+ApkXlVY43BnCUiNrpqF/q7FLIdp/t+Ojvj7B/3H70\ndu2ttjgm43nLE1N/m4qW9i2xcfhGONVwUlskSQlAdR/ofDmgp1dD8bXL/78b2HB/xHiWSCQSidFc\nBNBMCOEqhKgEYAKAIwXSHAEwBbhvcCcRUWxheYUQzbTyjwYQaNnLkJQUtvtsx4KTC3D65dOlwngG\ngO7O3eE90xst7Vui03ed8L8b/1NbJImNY64oHL8C6AvADhzuaBGASuCejO+FEG8BeANADoAMAO8p\nE1l0lVUiejPS8vLwT2IiAtPS7h9zqlQJQ+3s4FipQJzLnBzA3x+4dAkIDgYyM4GsLKBcOcDNDXj8\ncaBtW6BZMxgLESE2LRZhSWGISI7ArZRbSM1KRXpOOjJyMx5aCtWphhMa12qMxrUb47G6j5XqFZsy\nMoA7d4C7d4HERCAhAUhNBe7d4y0jg3+GzEwgNxfIy+NNCP5ZypUDKlYEKlcGqlQBqlUDatTgrVYt\noG5doE4dwN4eqF+f05iTHI0G7snJ8EpNRf7zXqtCBQypVw9uVauat7KShEYD3LwJXL/OW1wcUKkS\n3+B69YDOnfldqVxZbUmLhcph7NbhQSi6r4UQM6HoYCXNegBDwGHsphKRt768yvH9YFc7DYBwAK8T\n0W0ddZcInS0xD0euHcHMP2bi5JSTaOXQSm1xLMKZsDOYeHAipnacisV9F8vFWMowqoexMydqK+Nj\nd+9iY3Q0ziQloVvNmuhcs+b9bvqbmZn4OyEBratXx2R7e0y/fBkVtm0D/v4baNIE6NIFaNUKqFqV\nDYDc3AfGwpUrQM2awHPPAePGsaGgg4ycDLhHuONU2Clcvn0ZV25fgYY0cKvrBpfaLnCu5YzalWuj\nWsVqqFKhCnI0OcjIyUBaThpi7sUgMiUSEckRiEqJQkv7lujo1PH+ylCudVytdh9NITMTCA0FwsIe\nbLduAVFRvMXGcvukfn3Azo6N3bp12fCtWZON4KpV2SarXBmoUAEoX543gG24vDxu92RlcX3p6UBa\nGhvhyclAUhIb5vHxbONVqQI4OQENGwKNGgGNGwOurry5uQFNmxpm811JTcXKyEgcS0hAs6pV0bt2\nbVRU/PHicnJw9O5d2FesiHH16+N9Z2fUrFAKFgbIyQFOnAAOHwZ++41vVIsW3KB0cgKys/mHiIsD\nvL2BGzeATp2AiROBCRPYsLYR5EIqElvmdNhpjN83HkcnHkWXhl3UFseixN6LxaRDk5CnycO+cftU\ni2UtURdpQJuB8MxMvHvjBvzS0rDI1RUj7OxQp+KjPbhZKSk4vW0blgO4a2+PTRoNeo4YAdSuXXgF\nRNxDfegQ8OuvwGOPAQsXAn37Ii79Dg4EHMCBwAO4cOsCOjp1RP+m/dG1UVd0cuqEhjUbGj3pIT0n\nHX5xfvC+7Y0z4WdwMvQkalWuhZHNR2J8m/Ho3qi76hMpEhIAX18gIIC3wEBua8TGAi4ubJg2acJG\nqrMzG66NGrHNVasW9yZbAyI2qmNigOhoNuIjI4HwcDb0Q0OBiAigQQMebGjdmttRbdoA7dpxT3Zy\nbi4+CQ3Fnrg4fOzqirEODmikw+LWEMErJQUbo6NxKikJqx97DGMdHFT/rYpFZiawdSuwbBn/aGPG\nAKNH800qjPR04MwZYNs24K+/gKFDgUWLgJYtrSO3CUgDWmKrBN4JRJ+f+2DP2D14pukzaotjFfI0\neZh/Yj4OBR3CHy/+gRb2LdQWSWJlpAFtIj9ER2N+SAjednbGnMaNUaW8juGcnBzghx+Azz8H+vcH\nffop9tSujQ9u3sRz9vZY26wZKpQz0KU8Jwe5v2xHxpJPEF4xAzOH5sLl6REY33o8+rv1R63Ktcx7\ngQA0pIFvrC8OBR3CHv89SM9Jx+T2kzHjiRlwqe1i9vq0IWKj89Il4PJl7mT87z82Stu2ZUOzVSve\nmjdng9nWOl5zctiIvnaNGwKBgezV4+cHVH0qESmzAtEmxQ7v1XRDvycqomHDosv8NykJbwYHw6VK\nFfzaqpXOBl2JhAj46Sfg00+5J3nBAqBnz+KVlZQEfPcdsHIlj94sWsStqBKKNKAltkhKVgq6/dAN\nc5+ai6mdpqotjtXZemUr5p2Yh53P78QAtwFqiyOxItKALiZEhEVhYfg1NhZH27dH82rVdCcMCgIm\nTeKuxBUr2ChQSM7NxYsBASgHYE+bNqiuy/jWIi4tDhu8NuAH7x/Qom4zLAtvga6bDkO8/gYbGlbw\n/SQi+Mb54kfvH7HTdyd6ufTCu93fRd8mfc3S05mZycbyuXOApydw4QK7TDzxBG+dOwMdO7KhbGib\nw1b5JSYW7wbfwOyk1sj2rIvLl/neVKkCdO8O9OgB9OrF90TXT5+j0eDDmzdxMikJf7Vvr7PXukRx\n8yYwYwa3jjZvZrcmc5CYCHz9NbBlC/doT5tmvSEII5AGtMTW0JAGY/aOgVN1J2wasUltcVTjTNgZ\njN8/HmsGr8FL7V5SWxyJlZAGdDHI1WjwxvXruJKaiqPt26N+wYmBAPekbd7MPWmffw7MnKnzo52j\n0WD6tWsISk/HH+3awV5HWVEpUVjhsQLbfbZjfJvxmN1tNtrUb8Mno6OBN99k38/9+606VJ2WnYad\nvjux0mMl7KvZY0HvBRj2+DCjDOmMDMDDAzh1ikfer1zh3uQnn+SOxx492FgugfaORVkdGYm1t27h\naLt2aFujxv3jREBICDcuzp/nhkZwMDcu+vblrWfPBxMXiQgrIiOxMSoKx9q3R6vq1VW5nkIh4p7i\nhQuBefOAd9+1zDCCnx/7Rru5Ad9/Dzg4mL8OE5AGtMTW+OrsVzgSfASnXz6NyhVKeAPdwvjH+WPI\nziGY8+QczO4+W21xJFZAGtBGQkSYFBiI+JwcHGjTBjV0fegzMoCpU9mo3bmTJz0VUebHoaH4PT4e\n5zp3Rm2lzISMBCz9dym2Xt2KVzq+gg+f/BANa+oYvydif9F583j4e8QIc1yqweRp8rA/YD+Wnl2K\nqhWrYuXAlXrDF2k07Ibx99/A8ePAxYtA+/bAM8+w8dejB0/mK8ssi4jAtpgY/NW+PVwMCOGRksLG\n9OnTvPn5cQ/1wIG8dewI/BIXg7khITjTsaP+0RI1yMoCZs3iCzh4kP1wLF3fJ5/wXILDh83Xy20G\npAEtsSXcI9wxbt84XJp+CY1qlVzXKGsSlhSGQTsG4cW2L2Jx38W2Of9EYjDSgDaSRaGh+F9CAk51\n7Iiqulwu4uKAUaN4BtvWrUbFMJsVHIzgjAwcat0C31/ahK/cv8JzLZ/DZ898Zljg9gsXgLFjgddf\nZ5cOK7+8GtJgl+8uLDi5AB2dOmLFwBV43O5xJCbyfK5jx/ivnR0waBAbd336SINZm/1xcXjv5k1c\n6Ny52C4XKSlsSP/zD/C//7FHxODBQOUx0TjpEAmvLp1RryT4RMfEAM8/z5MEt22z7oNw+DAwfTrP\nTRg92nr1FoI0oCW2QmpWKjp+1xFrBq/Bsy2eVVucEkVcWhwG7RiEgW4DsXzgcmlEl2KkAW0EO2Nj\nsSAkBJ6wwx4PAAAgAElEQVRPPPFoPGeAZ4ENGwa89BLw2WdGO+nmajTodfEsgqP+Rdf0c1g9aNUD\nVw1DiY4GRo7krtxvv1XFUTgzNxOLj32D9VeWwz50Fu4enoe+vapg2DAOitCkidVFsgm8UlIw3NcX\nx9u3RyczGpMhIdx4OXoU+Puxm6jeOQULUzpg7KhycFUrOmFYGNC/PzBlCvcIq+HQfukSN3bffx/4\n4APr118AaUBLbIWZv89EriYXW0ZtUVuUEklCRgIG7RiEXi69sGbwGmlEl1KkAW0g55KT8ZyfH052\n6PCQT+p9goLYIFiyBHj1VaPLT8xIxJy/5+DPkNOo8MQmvNekBd5r3Lh4wiYnsyHfogX3sBUxOdFc\nBAWxG/bBgxx7uf9ztxDe+h3E4T98/+x36Ne0n1XksEUiMzPRw9sbm5o3x7P29harJ+UeYeB5PySE\nV0Ti/BZwdRF4/nkOL25p74n7BAfz8MOcOcBbb1mpUj1ERHD3/EsvsSGvItKAltgCfwb/iVnHZsHn\ndR+LRH0qLSRlJmHIL0PwRIMnsH7YemlEl0JKylLeJZqknBy8FBCArS1b6jaer10DBgwAvvyyWMbz\n8ZvH0X5ze1QqXwmBr1/C2a698HVEBC6lpBRP4Nq12cE4PJwjgOTmFq8cA7hxA1i6lP2Y+/dnD5Y1\na4Dbt4Fd3znD450DWDt0DV4+/DLePvY20nPSLSaLrZJHhImBgZjVqJFFjWcAqFVD4GS/1qjcPgUr\nrsZi1Sr2pOjTB+jQgX/LmzctKIC/Pzu8L1qkvvEMcNDwU6fYJ/qLL9SWRiIp0SRkJGDGHzPw86if\npfFcBHWq1MHxycdxJeYK3j72NmRDUaJNmemBnhQQgDoVKmC9ri6669eBfv2453mqcTEw07LTMOfv\nOfg9+Hf8NOqnh2JI7omLw6LQUHh36YJqxe1BzszkIerGjbkn2kwt4NhYYPdutjnCwrj38oUXgKee\n0j8Sn5iRiNnHZsMrygs7ntuB7s7dzSJLaWB5RAT+vHsXJzt2RHkr9VL43LuHAT4+uNi5M5pUrYq8\nPI7osXs3jyK4uXHAihde4FUbzUJICPD008Dy5dzjW5K4fZsN+ylTgI8/VkUE2QMtKenM+H0GKpWv\nhPXD1qstis2QnJmMATsGoLdLb6watEr2RJciTNLZRFSiNhbJvOyOjaUWFy5QWm7uoydjYoiaNiX6\n7jujy/WP86fWG1rTSwdeosSMRJ1pJvr701vXrhld9kOkphJ17Ur08ccmFZORQbR7N9GwYUS1axNN\nnkz0119EOTnGlbPffz/VX1GflrsvpzxNnkkylQaupKSQvbs7hWVkWL3uZeHh9LS3N+VqNA8dz84m\nOnqUaNIk/q1HjCDau5efgWITG0vUrBnRxo2mCW1JoqNZxmK8z+ZA0V+q61FrbpbQ2RLL4B7uTg1X\nNaSkjCS1RbE5EtITqNPmTjT377mkKaBvJbaLKTpbdeX7iEBmVsa3MjPJwd2dLiYnP3oyLY2oWzei\nTz81utytV7aS/XJ72uK9pdCXKTE7mxp7eNCx+Hij63iIuDii5s2J1q0zKptGQ3TpEtGbbxLVq0c0\nYADRjh1E9+6ZJk5YYhj1+LEHDd85nO6k3TGtMBsmIzeX2nh60rbbt1WpP1ejod7e3rQiPFxvmtRU\nom3biPr1I7KzI5o1i+jyZSMrSkkheuKJYr0rVuf6dSJHR25BWBlpQEtKKtm52dR2Y1va47dHbVFs\nlvi0eGq3sR0tPrVYbVEkZkJ1AxrAFgCxAP4rJM03AK4DuAqgYyHpzHZjNBoNDffxoc9CQx89mZtL\nNGoU0ZQpbGUaSGZOJk0/Mp1arW9FfrF+BuU5kZBAjc6doxRju3oLEhZG1KgR0eHDRSZNTCRav56o\nQweiJk2IliwhKsTGKhbZudn04f8+JNc1rnQ52liLrHQw/+ZNet7XV9UeidD0dLJ3d6egtLSi04YS\nLVpE5OpK1Lkz0aZNRLralg+Rm0s0ZAjRzJlGvSuq4uFB5OBQjJaCaahlQAMYAiAIQDCAuXrS6NTB\n+vICWA4gUEl/AEAtPeVa5mZKzMrXZ7+mIb8Mkb2nJhKTGkPNv21OK86tUFsUiRkoCQZ0LwAd9RnQ\nAIYC+FP5vzuAC4WUZbYb89udO9TS05Oy8nS4GXzwAXfJZWUZXF50SjT1/LEnPbf7OUrJTDFKlpcD\nAujDGzeMyqMTT08ie3siP93Gu5cX0Suv8LD9+PFEf/9NpOvyzclev71kv9yedv6307IVlTCC0tLI\n7uxZisrMVFsUWhURQYOuXjX445iXx+47Y8YQ1alD9NprhdiaH31E1L+/8b4+arN/Pzc4b92yWpVq\nGNDgyeA3ALgCqKgYvC0LpNGpgwvLC2AAgHLK/18D+EpP/Ra7nxLzEJYYRnbL7Ohmwk21RSkVRCZH\nUtO1TWmjVwl2Z5MYhOoGNMsA10IM6M0AXtDaDwTgqCetWW5Kem4uNT1/nv6+e/fRk7t3s9+zrnN6\nuBh1kZxXO9Nnpz8rlt9vTFYW2bu7k7+pvhNEPB7frNl9+dPTibZs4RH2Jk2IvvqK3VWtiU+MD7mt\nc6OPjn9UJvyiNRoNDbp6lVZGRKgtChERZeflURtPTzoQF2d03tu3iZYuJWrcmKh7d6Lt24nutwl+\n+YXIzY3IVBcktfjiC6KePY1qKJuCSgZ0DwDHtPbnFeyF1qeDDcmrHB8NYIee+i1wJyXm5IV9L9Ci\nU4vUFqNUcTPhJjmvdqZtV7epLYrEBEzR2dYKY9cIQKTWfpRyzGIsi4jAEzVrYkC9eg+f8PfnZYcP\nHAAKntPD4aDDGLpzKL4d+i0+7fMpygnjb5tjpUr4xNUVs69fz//oFJ8pU4ARI5AxegI+npMLFxe+\nnCVLOCTdvHlmjLpgIO0d28PrNS94Rnli3L5xpT7U3aH4eNzKysLbjUrG8rcVy5XD+scfx3s3biAt\nL8+ovE5OHLQiNJT/7tgBuLoCm167jLy33+UV/+zsLCS5hZk/H3Bw4IVWSi8F9estPKpf9aUxJC8A\nTANwzGRJJVbHPcIdHpEemPPUHLVFKVW41XXD8UnHMfefuTgQcEBtcSQqUEFtAXSxePHi+//37dsX\nffv2NSp/SEYG1kdF4UqXLg+fSE7mZYdXrQI6dSqyHCLCOs91WOGxAscmHkOXhl2KzFMYbzZsiC23\nb2PfnTsYX0wLlwjw8AC+iViBWecHoU/255h2/jM0a2aSaGbBrpodjk86jtd+fw3PbHsGRyYcgWMN\nR7XFMjvpeXl478YNbGvZEhXVWH1PD33r1sWTtWvjy/BwLHVzMzp/+fLAs8/yFuyVhHoDxmKaZjPy\nlrXDu+8CBV8nm6BcOWD7dqBrV24ZTJ5s1uJPnz6N06dPm7VMK2Fw2CYhxAIAOUT0q740pupsiWXQ\nkAbv/vUuvh7wNapVrKa2OKWOVg6tcPSloxiycwiqVayGoY8PVVskSRGYVWcXt+u64AbjXDiCYEEX\njud8fenLsLCHD2o0RGPHEr3+ukFl5Gny6N1j71KbDW0oLDGs6AwG8m9iIjl7eFC6rpB6hZCTw54n\nXbqw98a33xKlBkcTOTkRnTxpNvnMgUajoUWnFlHTtU3p+t3raotjdhaFhNAEf3+1xdDJrcxMsjt7\nlm6mpxe/EI2GnaNnz6bERKKVK4lcXIh69eL5q5b2qbcIvr48d8DX16LVQD0Xjr+09g1x4QjCAxcO\nvXkBvALgHIDKhdRviVspMQM/ef9EPX/sKScOWhiPCA9yWO5Ap0JPqS2KxEhM0dnmVOJNAPjqOTcM\nDyaw9IAFJxGeT0qixh4elFHQQP3xR6L27Q0KhJudm02TDk6ip7Y8RQnpCSbJo4vnfH0LDTumTWoq\n0dq1HDWhd28dBsxff/FEKWs7PRvA95e+pwYrG5B3tLfaopiN2Kwsqnf2LIWaYqBamEUhITQpIKD4\nBWzcSNSpk5Yj9MMNuObNOcyyCmGvTeOnn4jatuVJAxZCJQO6PB5MBKwEngjYqkAanTq4sLzg6Bz+\nAOyKqN9St1NiAimZKdRgZQPyuuWltihlgpMhJ8lhuQN5RHioLYrECFQ3oAH8CiAaQBaACABTAcwE\nMEMrzXpFUfsA6FxIWcW+ERqNhvp4e9OW6OiHTwQFFRq5Qpu07DQatnMYDd85nNKyiw4LVhz8790j\nB3d3SszO1psmNpZo4UIWe+xYogsXCilw7lyioUNLZNfggYAD5LDcgU6GlKxe8uIyOziY3g4OVluM\nQknJySFHd3fySU01PvPVq/zQ6blGjYbo9GlejMfJiSesJtnKmgwaDYemmT3bYlWoYUDTA2P3GjhM\n3TzlmEE6WFde5fh1AOEAvJVto566LXU7JSbw6clPadLBSWqLUaY4dv0YOSx3oEtRl9QWRWIgqhvQ\n5txMUcbH4uOppacn5WgbkllZHJ5iw4Yi8ydnJlPvn3rT5IOTKTtXv3FrDqYGBtLHNx8NKRQWxgtd\n1K3L3ibXDfGAyM4m6tHD6EVWrEV+y/zP4D/VFsUkQtLTqd7ZsxRrpYgOprA2MpJG/PefcZnS04la\nteIwHAbg40M0cSIvzvLxx7zWT4knIYH9Uf74wyLFq2VAq7lJA7rkEZ0STfWW1TOr+6HEMA4FHiLH\nFY7kE+OjtigSAzBFZ5ecGVAmoiHC/NBQLG3aFBW0J3Z9+inQsCHwxhuF5k/MSMTAHQPRxqENfh79\nMyqWr2hReRc3aYJN0dG4nZUFAAgOBqZOBTp3BqpXBwICgE2bYNjkwIoVeaLUkiXAtWsWlbs4PNP0\nGRx58Qim/jYVhwIPqS1OsVkUFoZZjRqhfqVKaotSJK83bAjfe/fgnpRkeKaFC4G2bYFJkwxK3r49\n8MsvgJcXkJAAtGgBvPsuEB1dTKGtQd26LPRrrwGxsWpLI5FYhMWnF2Nax2lwreOqtihljtEtR+Ob\nod9g8C+D4Rfnp7Y4Egsi2AAvOQghqDgy7Y6Nxepbt+DZuTOEUCaYX7gAjB4N/PdfoXHd4tPjMXDH\nQPR17YvVg1c/yG9h3r9xA9HxGoh1zXHiBDB7NkfYq1u3mAVu2MCG9LlzQIWSF2DF+7Y3hu0chjWD\n1+DFdi+qLY5R+N67hwE+PrjevTtqlcB7q4ttMTH4IToaZzt1KvqZ/vdf4MUXAR8fwN6+WPXdvg2s\nXAls3QpMmMDhFF1cilWU5Zk/n1ut+/cDZnzfhRAgIusokBJCcXW2xDIE3gnE0z8/jeBZwahbtbgf\nE4mp7PLdhQ+Of4ATU06glUOrItNnZHBHREICBwxLTeUtMxPIyQGys1lVVajAW9WqQO3avNnbA40a\nAdVkoBWjMUVnlwoDOo8Irb28sOHxxx/Efc7M5FB1n30GjB+vN298ejz6b++PYc2G4cv+X1rNeP7v\nP2DB8mz8OcELc290wcevVkHNmiYWqtEAgwYB/fpxQN8SiF+cHwbtGISVg1bipXYvqS2OwYzx88OT\ntWvjg8aN1RbFYPKI0O7iRaxp1gyDC4t5npoKdOgArFsHjBxpcr1xccDq1cAPP/CrN39+CTSks7J4\nuOeTT9jaNxPSgJaozejdo9HbpTc+ePIDwzJkZ/P7kJPD+3XrmrVRWZb55b9fMOfvudjy9D+olNIK\nkZFAVBSP0kVFATExPBAWFwfk5nK4/bp1gTp1gJo1gRo12FCuWJE3IThdTg6Qns6GdnIycOcOl1m1\nKtC0KdCyJW/t2gHdurFxLdFNmTeg98TF4Ztbt+Cu3dM2dy4QEgLs26c33930u+i/vT+GNhtqNePZ\n15dt+nPngI8+Am4Nv4lskYf1zZubp4LISDYMTpzgMfYSiF+cHwbuGIg1g9dgQlvzGS+Wwj8tDf2u\nXkVIjx6oXr682uIYxa+xsdgcHY1/C4t7/vrr/BH96Sez1h0fzyHXv/+eDekFCwBnZ7NWYRoXLwIj\nRnBr1tE88cqlAV1CyM3lRbN8fdnCSElhI9HRkV36Gjdm66KiZV31rI17hDsmHZyEoFlBqFKhyqMJ\nUlL423DqFBAYCAQFsQVXqRLfCyK+Ty4ugJsb0KsXd8h06VLq7pU5uXePb2dwMG83brD5ERoK3HXe\nAeo/F538j6NVvbZwduZHsEEDXsTK0ZG3GjVMa7cQAXfvcr1BQbz5+LCaq1AB6NmTf8r+/dndTraR\nmDJtQGuI0OHSJSx3c8PQ/NXSDHDdSMhIwIDtAzDAbQCWDVhmceM5IIAN5zNn2HB+4w0ebonNzkYr\nLy/4d+2KBpUrm6eyH38EvvuO70MJNfh8Y30xcMdAfDP0G4xvo3+EoCQwMSAA7apXxzxX2/MnzNVo\n0OriRWxp0QJP16nzaIIzZ9jn2c+PxwItQHw8sGIFP5aTJnGPtJOTRaoynvnzed7AgQNm+aJIA1pF\n4uKAXbuAgwcBb29urXXowCvO1qrFRmJcHHfVhYYCYWFAjx5Anz7cwjNXJ4ZKEBF6be2FmU/MxJQO\nUx6cuHcP2LMH2LmTramePYEBA7gB0aoVG8va84bS0oCICLYE//0XOHmS79fYscC0aZy/DFtfUVHA\npUvA5cvAlSusOmNj2Sht0QJ4/HHe3Ny4N7hBA2CP/y68f/x9HJt4DB2dOlpVXiIgPBxwd+ef8sQJ\nPvbss8CoUfz428C0Hothks4u7uxDS20wckb3obg46nzx4oNA8VlZRK1bE+3apTdPcmYydf2+K733\n13sWDzB/4wbRpElEDg5EX3/NcZ0L8nZwML1vULgNA9FoiPr2JVq92nxlWgCfGB9yXOFIhwMPqy2K\nXoLT0sje3Z2Sc3LUFqXYbImOpgFXrz56IiOD6PHHObi4FYiJIXr3XaJ69Tjy4t27Vqm2cDIzWV/s\n3WuW4iCjcFiff/8lGjmSqHZtosmTiY4eJUpMLDpfQgLR778Tvf02kaMjUbduROvXE6WkWF5mC3Ak\n6Ai13diWcvOUNRD8/YmmTyeqU4fo2WeJDh0iuneveIVHRxMtW0bUogVRy5a8roINRCMyBxERRD/8\nwI9WkyYcdWjIEKIFC4gOHOCIn4asi7bPfx/VX1GfPG95Wl7oQtBoiAICOARpjx6sj199lejECcOu\no7Rhis5WXfk+IpARylij0dATFy/SQe34WV98wUFq9RjGadlp1Pun3vT6769b1HiOjCSaMYNfts8+\nI0pOLiRtRgbVPXuW4sypkIKDufKQEPOVaQEuRl0kh+UOdPzGcbVF0cm0wEBaHBqqthgmkZWXRy4e\nHnS+YMDmjz/mIONWRvvdWLJEd6PSqri7EzVsaJjRVQTSgLYiXl5EgwezVfPDD6Y9SDk5RMeOEY0b\nxw/mnDlEt26ZT1YLk5uXS202tKEjQUeILl4keu45ovr1+eMTFWW+ijQaolOniAYM4NW9Nm/mMKql\nCI2Gb+GcObzukp0d0Ysv8uJRAQF6TQuD+P3a7yVuXYSICKIVK3jtrIYNiebPNzB8bimhzBrQx+Lj\nqY2nJ+XlP9HXr/PTrsfgyczJpME7BtPkg5MpT2OZRUfi44k++IBbdXPm8L4hzAwK0hkX2iS++opo\n0CDT3ngrcDb8LDksd6Cz4WfVFuUhwjIyqN7Zs3S3FHwgNty69XBc6KtXeVjk9m3VZLp+neill3hB\nlm+/VblDa+ZMDrxuItKAtgIJCURTp/LXfuNG8z84ISFE77zDwfjffNMmDOmfr/xMkz/rSJqRI4mc\nnXlNgDTLLAR2Hw8PNqRbtuTGh41z4wb3Kj/2GFGzZvz/+fPm75U9HXqaHJY70G9Bv5m3YDPg60v0\n3nv8aejbl1efLe0DDWXWgO7l7U07Y2J4R6Mh6t+faOVKnWlz83Lp+T3P05g9Yygnz/zD8ffucee3\nnR3RG28Y3+jPX6QjyZyuAtnZRB06EO3YYb4yLcTxG8fJYblDiVr2+61r12jujRtqi2EWMnJzqcG5\nc3Q1NZW/CN26ca9dCeDKFR40atqUaOdOlRbUTEggatCA6Nw5k4qRBrSF+e03okaN2LC1tKvFnTtE\nH33EvSHvvVdiVwrKDA6kA09UpyyHekRr1rBrlrXQaIiOHGGLc/hwInN3AlmYnBx2wxg4kBdgff99\nosuXLd/ndDHqIjmtdKKtV7ZatqJikpXFXm3PPMMdHAsW2EQ7sliUSQP6fFISNTl//sGqgzt28BiE\nDgNUo9HQa7+9RgO2D6DMnEyDyjeU7GyiTZu4M2TCBNOGPl7096cV4eHmE46Im9ANGphleNrSHAg4\nQA1XNaTrd9UfP4rPzqa6Z8/S7UzzPi9q8lVYGE0OCOBh16eeKnFLv58+zXZ9p05E//ufCgLs3k3U\npo1JXS7SgLYQGRnsz/vYY/ygWJPbt3l52Hx/PNV9jhTi44nee4/Sa1WjX8Y0V9d3OzOTRzzt7Hju\nTQl3ps3IYDXo5kbUsyfRL79Yt91BRBR0J4iarm1KS04vsfhcLFMICHiwOvILL/DAQ2miTBrQY/38\naG1kJO8kJnIzyctLZ9r5/8ynrt93pdQs8yk+jYbo4EGi5s254/vSJdPLvJSSQs4eHpRlbsNm+nSi\n2bPNW6aF+P7S99R0bVOKTolWVY7PQ0NpWmCgqjKYm4TsbKp75gxFtmjB63CXQDQaon37uENr0CD2\nNLFq5UOG8GSpYiINaAsQEkLUuTPR+PHqGok3bz7wOVq3jo1GNbh3j+jLL4ns7Slz5mvU5lN7+i/m\nv6LzWYPgYKI+fbglHBCgtjSPkJlJ9M033Kc0fDjRWZW9Bm+n3qbO33Wm1357zSIj4+YkOZlo7Vpu\nw3bvzv0NNjy3/j5lzoC+mZ5OdmfPUkr+rzd7Nvsw6mC1x2pqub4l3Um7U2S5huLhQfTkk0Tt2xP9\n9Zd5h3ueuXKFdpjbLzU+nieUXL5s3nItxFdnv6J2G9tRYoY6veYZubnk6O5OfsWdsV6CeWfDBvqo\nhLhuFEZ2NvtF169P9MorVhw+zJ9Hkd84NxJpQJuZ48f5IVi7tuTM5bhyhWjECCIXF45GYa05EvnW\nn5MTdwVeu0YLTyyklw+9bJ36DSUvj4dl7e3ZpaQEjHTl5XEvc9OmREOH8k9YUkjJTKEhvwyhQTsG\nqfbNM4bcXA7c9PTTRI0b8wTEgvPTbQnVDWgAQwAEAQgGMFfH+T4AkgB4K9vCQsoq8oJnBwfTvHxf\nqytXWMHqmK23y3cXOa92pvAk87hF3LjBQQucnYl+/tkyo1R/xsdTR+2wfObixx+52VgClFlRaDQa\neufYO9Rnax/KyLHyuBoR/RAVRcNKaA+tSZw9S6EdOlC9f/+1mbB8yck8K7xePaJPPrHS6PnChdzb\nWQzUMqCL0sFKmm8AXAdwFUDHovICGAvAD0AegM6F1F2se1UkW7ZweLkzZyxTvqmcO0fUrx9Ho/j2\nW8tN2ktJYSulYUOeLODN80Rup96mesvqUVhimGXqNZUbN9hVrG9fInO7JhrBxYvcId6tm/W9fwwl\nJy+HZh+dTS3Xt6Qbd21n3s2lSzwoU7cuhyi1xYBVqhrQAMoBuAHAFUBFRTm3LJCmD4AjBpZX6MXe\nVXxTozIzuUfiqac4vkwBToScIIflDmYZ2kpI4MkFdnZES5dadnJznkZDrT096Z+EBDMXnMfOXjbQ\n+0hElKfJo/H7xtO4veMsFjFFd70aaunpSSfNff/VJieHqF07ot276QU/P1oVEaG2REYRHk40cSIP\nvf74o4VdLNPSODTaP/8YnVUNA9pAHTwUwJ/K/90BXCgqL4AWAB4HcNKqBrRGw7OW3NyIgoLMW7Yl\nOH+eaPRo7siZN4/o2jXzlOvvzyGd7Oy4x9n74QnWb/7xJr3313vmqctS5Oayb7SDA9GePVat+u5d\nDpfp6Ej000820XdEG702kuMKRzoVekptUYwiIuLBfNtx44guXFBbIsNR24DuAeCY1v68gj0gigH9\nu4HlFXqxS8PC6JV839Rt24i6dHnka+oT42OWWIvZ2ezq5uDAHiL5AT8szZboaBpiiR7Qy5dZm9jA\nhEIiooycDOqztQ+9c+wdq02y+P3OnYcX5iktrF/PU6o1GrqYnEyNPTwo2xa+KAXw8uI2c8eORCct\nGUr10CGiVq2MnlCokgFtiA7eDOAFrf1AAI4G5j1lNQM6N5do2jTuLoyNNV+51iAwkOjDD9mQ7t2b\n3U4CAw13PdFoiPz82O2hRw9uLc6bxz25Bbh+9zrZLbMzq2uiRfHy4kWbXnnFKn7s+/fz7XvrLZv5\n3N3nn5v/kOMKR1pzfo3NfYdSUvjxdXVlPX3wYImfT6q6AT0GwPda+5MAfFMgTR8A8Urvxp8AWhdS\nnt4LzcrLowbnztF/qak8ttugAZHnw6v6RCZHkvNqZ9rlq38lwqLQaIj++IMnCA4axLERrUlmXh45\nnTtH/pbwwZ0+ncdabITEjERqs6ENrfawzqqKz1y5Qr9Yq6VkLeLjuRWoFQf6aW9v2mWj15k/0bBJ\nE+74s0ikwfwJhatWGZVNJQPaEB38O4Antfb/BtDZwLzWMaCzsth1pn//khPpojhkZ7OT6GuvsZNo\n48a8sMmcOUTff88P7/79HD9t82buunvuOXbRaNKEdfQffxQ6Q2vs3rG09N+lVrwoM5CaykveNWvG\nfhUWIDaW3SxbtOD1kWyV0MRQ6ri5I008MJHSsi0cz9sC5OTwgEPXrvxzr19f/EUwLY0pOrsCrMNl\nAC5ElC6EGArgMIDm+hIvXrz4/v99+/ZF3759AQD779xBq2rV0K5GDWDePGDQIKBbt/tpU7JSMPzX\n4ZjdbTYmtJ1QLEH9/YH33+e149esAYYOBUTxVkkvNpXLlcPMBg3wbVQUNjXXe5uKx9KlQOvWwPTp\n/LeEU6dKHRydeBRPbnkSLrVdMKb1GIvV5XvvHoLS0zHOwcFidajCJ58A48cD7drdP/SOszNWRkZi\ngqOjioIVDyGAsWOBESOAtWuB7t2BadOAhQuBWrXMWMmaNUDv3sCkSUD9+jqTnT59GqdPnzZTpVbF\nrPMlU98AACAASURBVFpNn842mIwMYNw4oFw54I8/gCpVzCmedalYERg1ijciIDgY+O8/4Pp1wMMD\nSEnh40SAnR3w2GPAhAnAsmVAs2ZFfnDORZyD5y1PbB+93UoXZCZq1AB+/BHYuxcYNgz48EPeypUz\nS/FHjgAzZwJTpgA7dtj2I9SkThOcm3YOM/+YiR4/9sD+8fvR3M7MtoAFqVCBPznjxgHu7sDq1cDi\nxcCMGcCsWUCDBurJZladXVzLO38DDwH+pbX/yBCgjjyhAOrpOae3pdD90iU6fOcOdznVq/fQaiXZ\nudk0aMegYi/RHR/Pwz0ODjzypvbic9GZmVTn7FlKsIQga9Zw5HgbGh7yjvYmh+UO5BFhuSCU04OC\naIktzoIojKtXeUj57t2HDucoy3t7FbbGvI0QHc0L0zk58bwzs3qmvPceO1IaCNRz4ShUB+NRF44g\nPHDhKCqvZXug09NZH02YoL7iLeFoNBrq/kN32n51u9qimEZYGFGvXjzaYGKIndRU7uxv2lT9sHTm\nRqPR0OaLm8l+uT3t8bOuD7m5uX6dbay6dVlfW3tkXx+m6GxzKO/yeDAJpRLYTaNVgTSOWv93AxBW\nSHk6L/JCcjI1PX+ecjUaHrf98sv75/IXShm2c5jRsRRzcngCtYMD/7iGLr1tDSb6+9NKS0z2ys5m\n/87Dh81ftgX5M/hPclrpZJFZyvHZ2VTn7FmKKU3rlmo0HJN140adp5eHh/PCKqUELy+eJ/vEEyYv\nKPiAxESeN2Bg3CuVDGhDdPAwPJhE2AMPJhEakvcUgCcKqb84d5bJyCAaPJjoxRdLR1BZC7Pbdzd1\n/q6zVSdWW4ycHF6YxtGR5xwUgytX2NXy5ZfZq7O0cjn6Mrmtc6M3/niD0rPT1RbHJOLjedVmJyf2\nkjtxQt2+PFUNaK4fQwBcA4dImqccmwlghvL/W+BwSFcAeADoXkhZOi9yor8/Rw745x9uamotG7Tc\nfTl12NSBUjKNm5xw4gRR27bcCC4prSFtPJOTqUl+o8Hc/PUXOyfZmMG4wWsDtVzfkhLSzRslY1l4\nOE0pRcYkEbGPZdu2eg2Tu0qjoTSttqjRcLzXRo2IJk16aJCq+GzezEFPDXgP1TCgyQAdrOyvV4xl\nH+0eZV15leOjAUQCyABwG1qTDQvUXbz7mpnJQXlfeEEazwaQmZNJTdY2MXlyfInDw4O/6dOnG+z7\nrtGwX62DA9HOnRaWr4SQlJFEL+x7gdpubEt+sX5qi2MyGRkcUallS14nac8edSYcqm5Am3PTpYzz\n3RkSMzM5FNf+/ffPHQg4QI1WNaLIZMMXPggNJRozhudrHDxYsj0Zul+6RIfi4ixT+NChvOyqjfHO\nsXeo37Z+lJ1rnuHefHeGS2qucmZuMjM5DNjx44UmmxEURItLm9sK8Xd4/nyOAPb11yYuGpebS9Sh\nA9HevUUmVcuAVnMrlgGdnU307LM840sazwaxzH0Zjfx1pNpiWIbkZO5GbtasyBhoSUlEzz/PRtf1\n69YRr6Sg0Whoi/cWsl9uT5subrK5KB26yMsjOnKEF6d77DHur7Dmsuql3oD+NCSE3rh2jWcwa/UE\ned3yIvvl9nQpyrB1tNPTiRYv5o/qkiW8X9L5NSaGnrHUskkBAbxaVEnyWzGA3LxcGvHrCHr1t1fN\nokD2x8XRUzaySqPBrFjBa9UWgd+9e+R07pz5l48vIVy/zovGPf440dGjJhR04gT3khVhiUsD2gBy\nc9nfecQI6fNsIJHJkWS3zI6u3y3lFuO+fTxnY9Einc+GtzcbWW++qd5K6iWBwDuB1Pm7zjTy15EU\nd89CHWwqcPYsrxXUsCH37VkjckepNqCzlJBuAbGxHLZOCX8TmRxJjVY1okOBRftOaTTc0+zqykG+\nVVwUyWjyQ/f5Wiqs01tvEc2aZZmyLUhqVip12NSBVpxbYXJZfby9abetxZwtjLg4biXmx0svgv5X\nrtBOGw1pZyh//smdW6NGEYWEFLOQkSO5YVII0oAuAo2GZ3z162fdbiYbZ8L+CbTgxAK1xbAOt26x\nc+wTT/BiMsSPzfffc3/PruJHqC1VZOVm0dy/51KDlQ3oz+A/1RbHrHh78+BU/fo8gmjJqJal2oD+\nNSaG+l25witTTZ5MRET3su5Rp82daJn7siJvTmAgx3Ju3drCCy9YkEX5PfCW4M4d1ko26P8bnhRO\nDVc1pN+Cfit2Gb6pqdTg3DmbXFREL2++STR7tsHJD8bF0ZOlrQdeB5mZvJKonR13cBk9AhUUxO/K\nHf2LV0gDuhA0Go5q0qOHbcd5tjInQ06SyxoXm4wHXGw0Gh7Lt7en7K9W0tQpudS6tcF9AmWKU6Gn\nyHWNK804MoNSs0rXe+Xvz1MkHB2Jli+3TI90qTagn7p8mQ4EBHDYuogIytPk0ejdo+mVw68UOnyf\nmspx6+3teSjAlkcKbyk+4MmW8hVctYqHU22QC5EXyH65PfnEFG/lxjevXaNPi90lWQIphltOTl4e\nOXt40JXS5ANeCOHhPAeiaVP2vTOKWbMKHbGRBnQhLF3Kc1gSzDsBuDSTnZtNbTa0of3++4tOXAoJ\nO3mTLlbvQ8F2PejeJWk96yMpI4leOfwKua1zI/dwG15BRg9+ftwj3bAh0YYN5o19UGoN6KupqeTs\n4UE5U6YQLVxIRETz/5lPvX/qTZk5uh2gNBqi3buJnJ2Jpkwhun3b+BtaEhnr50ffRho+UdIoMjPZ\nmjhxwjLlW5hdvrvIZY0LxaQa54aQnJNDdc+epVulyZlu5EiilSuNzvZ5aChNDwqygEAll+PHOQTW\n8OFGrGaYP2KjpytMGtB62LSJJ7VGRxedVnKfFedW0MDtA0vFZDFj+f13jrLx7bo80ny7noeOli2T\nk04L4VDgIWqwsgF98L8PbD7cnS4uX+aol25u7Mpjjtei1BrQM4KCaImHB/s+p6TQDp8d1HRtU71O\n8/7+7FrXoUPpC6h+KiGBWnl6Wk6R7tlD1KmTmVehsB6fnvyUev7YU2/DShfrb92isX62Hw7oPqdO\nGTTRTRcxWVkc6caWh2qKQWYm0VdfGenWsWIFR5DQgTSgdbB3L3cdWWTN9dLLtfhrZLfMziJx70sy\nubncX+bsXCCe+82b/IHv0oXIp3gjjmWBO2l3aPy+8dRyfUu6EFl4RBNb5eRJdpHv3t30JdtLpQGd\nmB+jduRIos2b6XzkebJfbk++sY8GbE5JIfroI+4Y+vbb0tlA1Wg01MbTk05aavhTo+GncbttrnCV\np8mjMXvG0JRDUwxqZGg0Gmrt6UmnSstwcl4ea5Tdu4tdxIv+/rTWUqMcJZx8tw43N6I//igicUYG\nx8A8c+aRU9KALsA//3A34tWrhd9TyUPk5uXSU1ueonUX1qktilW5c4fnLPXpQ6RzXrNGQ/TDD/yx\nX7iwbIfiKII9fnvIcYUjffi/D0tlb3ReHtGOHUSNG3NwiLCw4pVjis42zyL0FmBbbCyGZGbC6fp1\nRI4dhDF7x+CnZ39C2/pt76chAvbuBVq3BuLiAD8/Xme9QgUVBbcQQgi82agRNkRFWaoCYNUqYMEC\nICPDMnVYkHKiHLaN3gbfWF+sOr+qyPRnkpIAAH3q1LG0aNZh1y5+8MePL3YRbzVqhI1RUdCwUVSm\ncHEB9u8HNm4E3n0XGD0aCA/Xk7hKFWDpUuDDDwGNxqpy2hTe3sCLLwL79gEdOqgtjU2x3ms9hBCY\n1W2W2qJYjYsXgS5d+FH55x/A0VFHIvH/9s47PMpia+C/SQgl1BACoYQWepGiVFGCSBMBFUG4oiKf\noKBy1esVFAsqUhURURTL1asoCFylqAgoERBC7y2BEEISQoAUEkjf8/0xAQOkb8sm83ue99l3d6ec\neXf37Hlnzpyj4IknYN8+OHgQ2reHLVscLqsrMLz1cA6OP0j4pXDaf9qeLeEl6zq5ucGoUXD8OLRp\nA7feClOnwpUrDhSiqJa3vQ5ALBaLNAsKkk2DBknKih+k46cdb4q4ceyYyN136z0pmzYV7c7D1bjk\nCJ/doUOvS5PuaoTHh0ud9+rImuN5TyMOO3RIFkREOEgqO5OcLFK/vtV+SxaLRdrt2CHrLl60kWCu\nSXKyjhNfvbr+KeS4YeXqjP8NMbUwM9CakBDteve//+V6nQ05E3IxRLxneUvwhWBni+IQLBaRTz/V\nCxUrVhSy4vLl2j3oySd1hhVDjqw4skJqv1tbJqyZIAkpJTPn+enTeia6YcPCbQ63Rmc7XfneJBDI\nhthYabN2rWTeeacM/2GYPPK/R64tyycl/Z1h7P33S6a7Rl6MP35c3rBn1IjgYH1x7ZX90AFsDd8q\nPrN95HDM4Rzfj0pJES97RjVxNLNni9x/v02a+iQyUu4vjnntnUBoqN5g2KJFLvtrN27U2jrbDa0x\noEWvvTdurK0iQ6FISU+R2xbdJvO2zXO2KA7h8mWdgLB1az0pViTi4kTGjROpW7eQFnjpIvZKrPzf\nyv8Tv7l+supYYcMPuQ7r1+vN4YMHF8ytwxqdrXT94oNSSobu28ddc+dSvnsFFrnvI3B0IOXcy7Ny\npV5e7d4d3n0X6tRxtrSO52BSEv0PHCCsa1c83OzkgfPss3p95IMP7NO+A/h639e8velttj+xHW9P\n7+veezssjMjUVD5p3txJ0tmQixehRQu9jGmD8SRlZFA/KIiDnTpRt1w5Gwjo2ojAqlXwz39qvfPe\ne1C7drYCgwZBQAD861+AdrUSEeUUYZ2EUkqu/Y8kJurrMXgwvPGG1W2fv3yegzEHOX7hOMcvHudU\n/Cmik6KJToomNjmWDEsGmZZMlFJUr1Ad7wre+FbypbVPa9rWakt73/Z08O2Au5u71bI4gqd/fpqz\nSWdZMXwFSpXsr9HJkzB0KLRuDYsWQcWKVja4aROMGwctW8KCBVC3rk3kLGn8ceoPnlrzFG1qtmH+\ngPnUq1LP2SLZnNRUmDMH5s2DV16BiRNzd+21RmcXSwO62rp1bP1oJn17BrPjiR0kn6/NxIn6B/fR\nR3DXXc6W0rn02LOHF/z8eMDHxz4dnD+vlVBQEDRpYp8+HMCL615kb/Re1j68Fg93DwAyLBYabd/O\nmrZtaVepkpMltAEvvAApKdp510Y8HRyMj4cHUxs1slmbrs7lyzBtGnz+Obz+Oowfn6WQjxzRBuPx\n4+DlVboN6LQ0fUPRsCF88on2Vy0Escmx7IjcwfaI7ew6u4u9Z/eSlJZE21ptaeHdguY1muPv5Y9v\nJV98K/lSvUJ1PNw9KONWBotYiE2O5cKVC0QlRnEo5hCHYg6xK2oXZ5POcnfjuxnQZAD3t7ifquWr\n2udCWMl3B7/jjcA32DV2V7GV0Vb89JO2dd94AyZMKPRXJXdSU2H6dK0P33wTnnpKTwYZriMlI4WZ\nW2ayYMcCXrnjFZ7t/Oy1/8iSREiI/grEx8Nnn0HHjjeXsUpnF3XqOvsB9AeOAcHApFzKzAdCgH1A\n+zzakrGTJ8ttk7xky6kd8vbb2qNgxgzbBs92Zb6Njpa77b2r/Z13RIYPt28fdiYjM0P6f9tfnv3l\n76x8P50/L91KSta90FDtqGvjNNwHExOlTknLzmgjDh8WCQjQER+DrkaIGjtWhwES57lwWKODc6sL\neAHrgOPAb0DVXNrVPuGjRul10wK4RqVlpMnuqN3y8Y6P5dEfH5VmHzaTytMrS8BXATJp/SRZfni5\nnIw9aZOwnREJEfKfvf+R+5fcL1VnVJURy0fI2pC1kmkpPt/vwzGHpcbsGrLvbMmOVpKern8q9etn\n+/3Yg0OHRLp318fhnF35DDpUYt9v+krLBS1l/cn1zhbHLlgsIl99pdOC//vf2m0oO9bobFsobjfg\nBNAA8MhSzi1uKDMA+DnrvAsQlEd7MnXYLTJ58WJp2lTkvvuKHp6kpJKSmSk+W7ZI8I3fBFty+bL2\nKbOrlrM/cclx0uzDZvLZ7s9ERKTfvn3y35KSXWfECJE337RL0z327JEVLuwHb08sFh0+yddX712K\nPRSpb2TCwpxiQFujg/OqC8wCXso6nwTMzKV/nfa1W7eb/51Eh5gMvhAsiw8slufXPi+3f3G7VHyn\norT6qJU8/tPjsmjXIjkQfUAyMjNs9yHlwoXLF+SjHR9Jh086SPMPm8snOz9xeors8PhwaTivoXy9\n72unymFvoqJE7rxTh6k7f94BHWZminz8sQ5598YbJuRdLlgsFvnp6E/SaF4juW/JfXL8wnFni2QX\nzp3TacGbNNHbV67ibAO6K/BrtueTb5wBAT4BHsr2/ChQK5f2pOv/PSuNGulMRIaceenECXkhJMS+\nnXzxhcgdd9gm3Y8TOXb+mPjM9pHvT2ySGlu2SHKG/f+o7c7OnTrKQWKiXZpf7IhVDhcnLk5kwgSR\nWrVE9g1+TSyPPOIsA7rIOjivulmz0rWyzn2BY7n0r3daXrggiamJEnQmSD7b/Zk88/Mz0uPLHlJ5\nemVp8H4DGbp0qMzYPEM2nNwg8cnOjZhgsVgk8FSgDP5+sNScU1Omb5oul1Icn8r+XNI5af5hc3lv\n63sO79uR/PmnDpbxxhs6UYpDOXNGZMgQkZYtb8jMYshOcnqyzNw8U2rMriHj14wvdGZfV2HVKp2k\n58knRRISrDOgrfaBVkoNBfqJyLis56OAziIyMVuZ1cAMEdma9XxD1szGnhzakymvZTDlZXcqVLBK\ntBJNaHIynXfv5ky3blRwt9MGmcxMHZRzxgzt2+jCrD2xlgd3beDhWx7m01YdnC2OdYjojQAjR2pH\nQjuQarFQf9s2NnfoQDNPT7v0UVLYuRP+NS6RyfRj4L5tiIN9oIuog9ejZ5Ub5VZXKRUnIl7Z2ogV\nkeo59C9dv32JEOK4lJaET8Va+Faug2+lOvhWrk2tirUp71EBASwiCHD1Xyen/x+lFAquHW5K4XbD\no3v2c9DPczjPXjZ7W9n7iLgUwdLDS9h7di+Dmw/m3maDqODhSdY/K5YseTNFyMx6npF1nilCRi5H\netZx7bnFct37VzLTWB/6O7Ur16O5T6vr2rdke8ztH/rGa5N9vGWyHR7ZHj2UwsPNjbJZ52WzzrM/\n3vi6h1K6fNZ59nbK3NDX1et/9VAC8z9QvD8XvvwCBvTT1/3q53yVq9+Dq98NyRp39mtgueHaZL9e\nmbm8dt3xxx9kvv8+mb16kTluHJkVKuRfL6u/G8/z6t9yw2s3yn31UeT6z9iSbdzZr8O1a5TH9yD7\n7+bqa1e/Hzee5/R7yv6YkpHM5tOb2Ht2N51q38qdDe6garnKf/92bvgdueXSz7UjB7lulDf7OLKP\n87rzXK6F3PiY7TeTm66JT0pn+YoI0sIOs3fRk0XW2cUy5UgZt7eZNUufBwQEEBAQ4FR5iiONK1Sg\nU5Uq/HD+PI/5+tqnE3d3mD1bJ4wYMMClM9QENOoDEe5sCZrElaY/4enhwkbhr79CdDSMGWO3Lsq5\nufG4ry+fRkXxngtvJLU3gYGBBAYG0r5PPMNOemkHCNegKH8Yuc62RK2LpI5nBfzLeFO/6200bNb1\nmkGngHSRa3++7lmGkwK44U/0agfZjYn0rNmeG42SvIyf7Oe5GWJ/G2puVGzwDzrUHsIfsSH8sH0Z\njb0a07haIzzc3K/90btnM1LLXDUS4W/jMZth6Q54uLlRPsvIvGrEXj2PT47lwx0L6OLdjDGtA669\n55bVj8r2eKOhcZXshtdVwy0z27ivGvHXjHmL5ZpRn5Z1npiZee08zWIhNatcWrYy2c+vtpP9BiF7\nXzcanilpgrQHt8XCYMDyZ/5fspsMPq43/q5+BlevlXu289xuqNyVwt3PD/c5c3APC8N99WrcmzTB\n3dv72k1H9huQG28ErrWfQ5/X3cxl1VVZj2WV0udZ8lz3ud5wfp3heYMReu3a3LDb8sYbj2vnORjl\nuRnt134jWZ9d+TIV6O3fj451u7EtYhvv7lxE65pt6VC7A5XLVc3xd5T9piCnm4Hsct0oL9ysWHIy\nrG86z3YtVA6PCriSkUxCSjyXUhJIjY2hzPmLpO46TOrRU5TLgFQP6yYfbWERRQL1sz2vl/XajWX8\n8ilzjalTp9pArJLP+Dp1mH76tP0MaNCG87vvwn/+A2PH2q8fO7Ps/Hlu96pJLS9fHl/5OEuGLnHN\nMFGZmfDSSzBzpt1vaJ6sU4fOe/YwrVEj+61yuDgBAQHc1v02Oi4ehsc9z8P/fnGGGNbo4LJ51I1W\nStUSkXNKKV8gJjcB0td/y5erdCY516YbR88fZeqfU/lz55+80O0Fxt82nsrlKtushz/D/mTETyOY\ndPsk/tnln66ph/LhyBF44AEdoOaDDyCniJgi4pyxd+4M69frjIZ33AHvvw/e3vnXK420bEdUYhTz\ngubx5c+T6VG/B892fpZejXrhpopPdJPY5FgOnDvAgXMHOBRziMPRB3E/eIg7w4S+keVpf/IylK9A\nUqdb8OjRE6//m0iZ9h3Bw8O672BRfT+uHoA7f29CKYueg2l5Q5l7+HsDS1fy2URoKBgZFov4bd0q\ney7Z2Xdv507twJaUZN9+7Ei33bvlx5gYSU5Pls6fdZZpf05ztkhF44svRHr0cJhf+oD9++WrkrLp\n0g5kWjJl6NKh0nztJ/JGaKizfKCLrIPzqoveRHjVHzrPTYQ//qgzyf32m90utcM5eO6gjFw+UmrM\nriFTN06VmCTrNtWmpKfIrC2zpOacmrLuxDobSVn8WLFC79v78ktnS5IPiYkizz2ndwMvWeLye33s\nTVJqknyy8xO5ZeEt4jfXTyavnywHzx20SaScgnJ1Q/Kyw8tkyu9TZODigVJvbj2pPL2yPDCroyx+\n+k45cVd7SfOqIulN/XWCncWLRcLDc23TGp1tkzjQSqn+wAfo1ZYvRGSmUurJLMEWZZVZgA6XdBl4\nXHLwf84qJ7aQqbQwLSyM8NRUFtk7KcjIkdCqFbz2mn37sQP7k5K49+BBTnXpQhk3N6ISo+jyeRcW\nDFjAkBZDnC1ewblyRSdLWbYMunZ1SJerL1zgndOnCbr1Vof052q8Gfgmv4RuJKTZ2xzq1Jm65cs7\n3AcarNPBOdXNer068AN65vo0MFxE4nPoW0SELVt0Yow5c+DRR+09YscRfDGY2X/NZsXRFQxsOpDx\nt42nu1/3As9ciQhLDy/lld9foXXN1szrNw//6v52ltrxWCw6Rvo338CKFS60GrFtm56NbtJEx482\nCVjy5cC5A3x74FuWHFqCu5s7/f37069JP7rW64pvJetXxEWEmMsxHLtwjEMxhzhw7gAHYw5yMOYg\n3hW8aefbjg7ebegdWY62O05T9Y+/ULGx0K8f9Omj9wjVK1iCmBKXSKW4yVSciU5NpeXOnYR17UpV\ney7ph4bqpa/Dh6FWLfv1YweeOn6cOuXK8XrDhtde2xG5g4HfDeSPR/+gba22zhOuMEyfDnv3agPa\nQWSK0DgoiJ/atKFDZdstY5cEVhxZwfO/Pc/TQ35j55VMlrdpU7oTqaCX7gcM0IlmJk2yYYKMYkBs\ncixf7/uaRXsWkZiayMCmAxnYbCAdfDtQt0rd65a0My2ZbIvYxurjq1l5fCWVylbi3b7vEtAwwHkD\nsCOXLsE//qETUS5bBjVrOluiQpI9Acu0adpd0SRgyRcR4cj5I6w9sZZ1oevYGbmTimUr0rF2R5p4\nNaFhtYbUr1ofrwpeVCpbiYoeFRGEtMw00jLTiE+J58KVC5y/fJ6ISxGEJYQRFh9G8MVg3JQbLWq0\noI1PG9rWakvbmm25pXITvDbt0Hdov/wC/v46wMHAgdChQ5E+M2NAl3KGHz7MnVWr8kwB77iKjB2y\n3tmbSxkZNAgK4nCnTtS5wRFv8YHFvLbxNXaM3UENzxpOkrCAxMToFQAnZId85/RpTqek2H+Vw4XY\nH72fu7+5m7UPr+XxSAvvN2lC79KeiTCLyEi45x7o0QPmz9d7kUsawReDWRO8hl9CfuHohaPEJsfi\nV8UPQYhLjiMhNYE2Ndtwb9N7GdR8ELfVua1Y+YzaklOntA1zxx368/Zw5YR2Bw9q47lcOZ1f3Oi8\nQiEinIo/xd6zewmNCyUsPozwS+EkpCSQmJZIUloS7sodD3cPPNw88KrgRQ3PGnhX8Mavih8NqzWk\nQbUGNK3eFJ+KWZmWr1yBn3+GH36AdeugUye91DV4sE1WC4wBXcrZGBfHMyEhHOrUyb6bMi5ehBYt\nYPNm/egCfBwZyR9xcSxv0ybH9ydvmExQRBDrH1lfvFOZPvOMvrueP9/hXTtslcNFiLkcQ+fPOjPr\n7lnU8+vPmGPHONa5s945bwxoABIS9H9cxYrw/fdQ0iMhXkm/wun407i7ueNV3otq5asVb31iI7Zs\ngWHDYMoUePrpErLikJkJH30Eb70Fzz4LkyfnvAvSYD/S0uC337Ty+OUXvfo9fDjcdx/UsO1klzU6\nu2TeEpcyAqpVwwJsSkiwb0fe3joCxOTJ9u3HRogIC6OiGJ/HXeo7d71D5XKVefbXZ3OMSVssCA6G\nJUu0g6ET8C1Xjn5eXvw3Otop/RcnUjNSGfrDUEbdMoqH2jzEwshInqpTp0RGUrCGqlX1/161atCr\nl15AKcl4enjS0qclzbyb4VPRp1QYz0uX6kgbX32l7+9LzE/A3R0mTtTucvv26VwIGzc6W6qSj8UC\ngYE6t0Ht2jqEbo8eEBKiZ56feMLmxrO1GAO6BKCUYkKdOnwcmWtkQNvx7LNasWzebP++rGRzQgLp\nFgt3VauWaxl3N3cWP7CYLeFb+GjnRw6UrhC8/LKOxe1E5TGhbl0WRkUV35sMByAijP95PD6ePrzV\n6y1i0tL4OTaW0fYMI+nClC2rjau+ffWe16NHnS2RwRaIwNy5WiVt2KD3bZVI/Pzgxx+1ITd6tHby\njopytlQlCxHYtUt/merXh+ef1y6KV22MCRPAx8fZUuaKMaBLCI/6+rIuLo6zqan27ah8eb3Z1FhJ\nTQAAIABJREFU4sUX9R1jMebjyEgm1K2b7+xglXJVWDVyFdM2TWP9yfUOkq6A/PWXTnX3z386VYw7\nqlbFXSk2xt8UhKHUMC9oHnvO7uG/9/8XN+XG52fPMrRGDbxc2unTvigFb7+tF0969oQ//nC2RAZr\nsFj0Vpgvv4StW+GWW5wtkQMYPFjvjm3cWA94zhy96dBQNERg/3549VVo2lRH+KpQQc8y792rV7nr\n18+/nWKAMaBLCFXLlOEhHx8WnT1r/85GjtSadOlS+/dVRM6mpvJbXFyBk8w09mrM0geXMurHUQRf\nDLazdAVERP9bvfMOzs5rf3WV4yNHrHIUQ34N+ZU5W+ewcoSOqJBhsfBJVBQTTMirAjF6tFYXI0dq\n48vgemRk6OSnO3dq32c/v/zrlBgqVtTRObZtg02boGVLvamtFK/IFQqLRX9xpkzRGzOHDNE3IUuW\naBfFt9/Wm+RdDGNAlyCerluXRVFRpNt7ZtjNTWcnfPllHZWjGPLZ2bM85ONTqE1vPRv2ZFqvaQz6\nfhBxyXF2lK6ALF2qN7Q8/LCzJQFgVK1abIyPJ6KYfub24nDMYR776TGWDVtGg2oNAPg5Npa65crR\n0YT2KzC9esGff8KMGfDvf+uvtsE1SE2FESPg7Fm9tysPr7iSTdOmsHo1fPEFzJoF3bppPxZjSN/M\n5cv6Wk2YoGeUH31U/+gXL9ahW+bM0cHCXdh53hjQJYi2lSrhX6ECKy9csH9nPXtC+/ZOiQqRH+kW\nC4uKODs49taxDGw6kGHLhpGemW4H6QpISoq+QXn33WITj7RymTI8XKuWY1Y5ignnL59n0PeDeK/v\ne9xe//Zrr38UGcnTdeo4UTLXpEULHYlx1y64/34dN9hQvElJ0cEPLBZYtUpPxpZ6evXSM6oTJ+rw\nIz176o2GpdmQzsiA7dth5kydzMTXV6dJb9RI32QcParf69TJpY3m7JgwdiWMpTExfBIVxcb27e3f\nWXAwdO+ufxjFyNF/xfnzzIuIYHOHDkWqn2nJZPCSwdSvUp+PB37snAgLc+boddKVKx3fdx4cvXyZ\nu/bv53TXrpQtJoa9vUjNSOXub+7mjvp3ML339GuvB1+5wh179xLerRvlbrgGJoxdwUhL05Ebtm3T\nRlmjRnYSzmAVqanaeK5aFb79FkwUyxzIyNDh1t5+W0/Nv/CCjuFY0vdGXLig74S3bdPH9u3QoIHO\nAtirl350gRU6EwfacI00i4UGQUFsaNeO1o6YKpg4UU9NLFhg/74KyF379jGudm1GWJEx8VLqJW7/\n8nbGdhzLxC4TbShdAbhwQfvYbdlSLAP52+L6FndEhDGrxpCQksDy4cuvS4LxXEgIFdzdmdG48U31\njAFdcES02njnHe0KGRBge9kMRSctTduB5ctr+9AYz/mQmQlr1ugQJadOaYfx0aMhWwZclyQ1FU6c\n0FmIDx3SyWb27IH4eJ39r3t3fXTtWuzCzBUEY0AbrmPqqVOcS09nYbNm9u/s4kVt7G3cCK1b27+/\nfDiUlETfAwcIs8EMaVh8GN2/6M5ngz5jYLOBNpKwAEyYoP+tiqF7DMCP58/z7pkz/NWxo7NFsRsz\nNs9g+dHlbBq9iYpl/74RTczIoGFQEHtuu40G5cvfVM8Y0IXn99+1m/9rr+mvfglZ3XVp0tN13gql\n9FaMkj6ZanP27tW7Zb//XkfuGDZMp2u0d7bgoiACsbEQHg5hYdr4Dw2Fkyf1KnNkpJ5Zbt0a2rTR\nR8eOOipJCViFNAa04TquZo4L7dLFMSG2PvhAZ01Yu9bp/35PHj9OvXLleM1Gd/1BEUEM+n4QGx7Z\nQDvfdjZpM08OHdJLX8eOQfXq9u+vCGRYLDTZvp3lrVtzW5UqzhbH5vxw+AdeXPciQU8EUafy9X7O\nH0VGsjGPzJbGgC4aJ09qV4EuXXQSOJP4zXlYLPD443D+PPz0k47nbSgiqak6DfWPP+r/yIYNoXdv\nnfe8e3ednMyeiOhJrqgobQhHROjHM2euP8qU0Rv9GjbUhvHVo3lz/VoJvoNymgGtlPIClgINgDBg\nuIjclA5PKRUGJAAWIF1EOufRpjGgbcAjR4/SrmJFXnREPMX0dGjbVi9d3XOP/fvLhdj0dPy3b+dY\n587UsqHWz8ugsikiOuvE4ME6YU0xZk54OAcvX+a/LVs6WxSbcvWGaf0j62nve/0+AosIrXbsYFHz\n5tyZSxgCRxvQhdDB/YF56I3jX4jIrLzqK6WqA8uBTsB/RCRXPyZb6eykJG24nTkDK1aAiRDoHP79\nbx1+fv16s2HQpmRk6Av755/aPS8oSE+StGqlj4YNoU4dfXh56YtfsaI2XkX0kZ6uo1tcuQKXLkFc\nnD4uXtTpPmNi4Nw5HS7l7FmIjtYhUOvW1e3Wq/f34ef392OpDaviXAN6FnBRRGYrpSYBXiJyU55n\npVQocKuI5BsbzBjQtmHnpUsMO3yYk1274u6IWeGff4Z//Uv7RznpbnV2eDhHLl/mKzsYddM3T2fF\n0RX8OfpPKpWtZPP2AR3y56WX4MCBYn/Hf/Vm5WinTviWkOnC0LhQenzZg0WDFnFvs3tven/txYu8\nfOoUe269NdeNpU4woPPVwUopNyAY6A1EATuBESJyLLf6SilPoD3QBmjjCAMatI0wYwZ8/LF2Hbj9\n9vzrGGzHnDnw9dc61HExXQArOWRmapeJI0e0f3F4uDZ6o6K0f3FSkj6uxntUSi8HeHrqo3JlbWh7\neemZ7Fq1oGZNfdSurQ9fX13WkCvONKCPAT1F5JxSyhcIFJEWOZQ7BdwmIhcL0KYxoG1Etz17mOTn\nx32OiJAhAgMG6MMJWfMyLBb8t2/nf23acKsddv6KCGNXjyU6KZqfRvxEGTcb76hJS9O+ZR98oK+h\nC/Dk8ePUKVeON1x9kwwQmxxL9y+682znZ3m689M5lrnnwAGG+fjweO3aubbjBAM6Xx2slOoKvCEi\nA7KeTwZERGblV18p9Rh68sMhBvRVfv1V779680148kmne4aVCr7/HiZP1hkGzey/obRgjc621gO8\npoicAxCRaKBmLuUEWK+U2qmUGmtln4YCMrFuXeY7KnOcUtqFY9o0vYzkYFZevIhfuXJ2MZ5B/8gW\nDlxIWmYaE3+diM1v8ubPhyZNXMZ4Bni2bl0+iYoirZindM+P1IxU7l96P/c2uzdX4zn4yhV2JSYy\nsmZuKs5pFEQH1wXOZHsekfUaQK0C6nCHMmCAXu1esADGjjWZk+3NX3/peY+ffzbGs8FQUPKdRlNK\nrQeyx6tSaIP41RyK52ZV3C4iZ5VSPmhD+qiIbMmtz6lTp147DwgIIMDENyoSQ318ePHkSQ4kJXFL\nJTu5HWSnVSt45BF45RX4/HP795eN+RERTLTzDmcPdw+WD19Ojy978O7Wd/n37f+2TcNRUTrA/LZt\ntmnPQbSpVInWnp78EBPDqAKmTC9uWMTCmFVj8PH0YXaf2bmW+zAykrG1a1Pe3f261wMDAwkMDLSr\njDbSwQWlSPXtobObNNE/idGjdYi7FSu0G6fBtpw8CQ8+CN98oxfBDIaSjC11trUuHEeBgGzLfxtF\nJE8HVKXUG0CiiMzN5X3jwmFDpp8+TUhyMv9pcZNnjX1ISNBh7Vau1BmHHMDOS5d48PBhTnTpgocD\nwupEXIqg+xfdmXn3TP7R9h/WN/jII3ozx4wZ1rflYNZcuMDrYWHszsMvuDgzaf0ktpzZwoZHNlDB\no0KOZWLT02myfTsHO3Wibj7+3k5w4chXB2e5cEwVkf5Zz7O7cORZ31kuHNkRgenTtV/0smU6eIHB\nNsTH6/C9zz0HTz3lbGkMBsfjTBeOVcDorPPHgJvSpimlPJVSlbLOKwJ9gUNW9msoIE/VqcPKCxeI\nctQaaNWq+t/umWd0PCQH8N6ZMzxXr55DjGeAelXq8cvDv/Dc2uf4PfR36xrbsgUCA2HKFJvI5mju\n8fYm2WJhY3y8s0UpNPO3z2dV8CpWjViVq/EM8GlUFIO9vfM1np1EvjoYvWmwiVKqgVKqLDAiq15B\n6zv1zkgp/fNYtEiHuvvyS2dKU3LIzISRI6FfP2M8GwxFwVqLYxbQRyl1HL3DeyaAUqq2UmpNVpla\nwBal1F4gCFgtIuus7NdQQKp7eDCqVi0+dJQvNMCjj+oA6199ZfeuwpKT2RAXxxN5bOyyB21qtmHZ\nsGWMXDGSfdH7itZIZqYOVzdnDjjCxcYOuCnFv+rV490zZ/IvXIxYfmQ5s/+azdqH1+LtmXss1lSL\nhQ8jI/mXn58DpSsU+epgEckEngHWAYeBJSJyNK/6WW2cAt4DHlNKhSulHLSMlTMDB+roELNmaX/d\njAxnSuP6TJmifcvffdfZkhgMrolJpFIKCE1OpvPu3Zzq2pXKjsrHunu3/sc7csSu8ZCeP3ECD6WY\n7e9vtz7yYtnhZTz323Nsfnwzjb1uTu2cJwsWwPLlOoujC7o/XCUlM5OGWenj27jAjcAfp/5gxPIR\nrHtk3U2xnm/kP2fPsiQmht/aFSyJjkmkYn/i42HECH3/uWxZqQ5hW2S+/15vVdm50yWzLxsMNsOZ\nLhwGF6BxhQr08vLiy+hox3V66606fenkm8KC24y49HS+jo5mohO3jQ9rPYxX73iVvt/0JTqpENc3\nKkrH6Fq40KWNZ4Dy7u48U7cucyMinC1KvuyO2s2I5SNYNmxZvsaziPDemTO8WHxnn0sl1arpaBFt\n2kC3bnDihLMlci3274eJE3WWQWM8GwxFxxjQpYQX/fyYFxFBhiNDjk2bpv/ptm61S/OLzp7lXm9v\n6pUvb5f2C8r4TuN5rN1j9P+2P/EpBfQFfuEFGDdOb7gsAYyvW5efLlzgbDGONxZ8MZh7v7+XRYMW\n0bNhz3zL/xYbi7tS3O3l5QDpDIXB3R3ef1+7cvTooV07DPkTHw9Dh+pw8wVcVDEYDLlgDOhSQpcq\nVahXrhzLz593XKdVq+rY0E89pVOQ2pCUzEzmR0QUG9/UV+98lZ4NejLo+0FcTrucd+HffoMdO1x2\n42BOeHt48I+aNR0Xd7yQnEk4Q99v+jL9runc1+K+AtWZfeYM//Lzc8noIqWFp57S4dcefFC7JRhy\nx2KBxx7TMbb/YYPgQQZDaccY0KWIl+vX553wcCyO9DEfPlynFJ03z6bNfhUdTftKlWhXTHxulVK8\n3/99/L38uX/p/aRm5DITm5wMTz8NH31U4lKsvujnx6KoKOJsfLNkLdFJ0fT+b2/+2eWfPN7h8QLV\n2ZqQwKmUlOKYOMVwA336wO+/w6RJeoOh2UKTM3PmwLlz8N57zpbEYCgZGAO6FDGgenXKKsXqi/lm\nVLcdSukArrNm6Yj9NiDdYmFmeDivNmhgk/ZshZty4/PBn1OlXBVGrBhBemYOhuRbb0HHji6VcbCg\nNKxQgSE1ajg24ks+xCbH0vebvoy6ZRTPd3u+wPWmnT7N5Pr1HRYa0WAdbdtqT7HvvtP3p5mZzpao\neBEYqOcwli2DsmWdLY3BUDIw/w6lCKUUrzZowLTTp22fijov/P3h5Zd1Tl4b9Lv43DmaVKhAt6pV\nbSCcbSnjVobvhn5HSkYKj/30GJmWbP/ku3frILYffug8Ae3M5Pr1+TAyksRiEGMsISWB/t/2p59/\nP16787UC19udmMiBpCRGu2h2xdJKvXqweTMcO6ajdBRjd3yHcu4cPPwwfP01FBOPN4OhRGAM6FLG\nkBo1SM7MZF1cnGM7fu45uHwZPvvMqmYyRZheDGefs1PWvSz/G/4/Yi7H8PjKx7URnZYGY8bo9dNa\ntfJvxEVp5ulJHy8vFkZFOVWOS6mX6L+4P53rdmZ2n9mF8mN+5/Rp/l2/PuXM7LPLUaUK/PKLvk8f\nOBASE50tkXPJzNT+zmPGQN++zpbGYChZmH+IUoabUkxp0IC3w8IcOwvt7q5nX6dMASuSbiyLiaGm\nhwc9i3nw1woeFVg1chURlyJ4YvUTWGbM0NM/Dz/sbNHsziv16zP3zBmuOGkdPTE1kQGLB9DBtwMf\nDviwUMbzoaQktiYkMNbBiXkMtqN8eVi6VC983XUXONJjrbjx9tt68+DUqc6WxGAoeRgDuhQyvGZN\nYtLTCXR0+uXWrXUA0iefLJIrh0WEd7Jmn10hMoKnhyerR66GAwe5PHcGmR9/5PIxnwtCm0qV6F61\nKoucMAudmJrIPd/dQxufNiy4Z0GhvyfvhIfzvJ8fnu7udpLQ4Ajc3eGTT6B3bwgIAEeGwC8u/P67\nTn/+3Xf6ehgMBttiDOhSiLtSvN6gAVNOnXLsLDToxCoxMVqzF5IlMTFUcnennx0zG9qailKGz5el\nsvDBRjy2awoZFuf7BjuCNxo2ZGZ4uEN9oeNT4unzTR9a+7Rm4b0LcVOFU2/7EhPZGBfHhDp17CSh\nwZEoBTNmwEMPwZ13WrXw5XKcPQuPPKJD/JnFFIPBPhgDupQyslYtkjIzHRuRA8DDQ2v1V1+F4OAC\nV0uzWHjt1ClmNm7sErPP13j1Vdz9m/Dswt1cuHKBkStG5hydo4TRrlIlent5Mc9B2QkvXLnAXV/f\nRbd63Vg4sPDGM8CUU6d4pUEDx6W7N9gdpbSqmTBBG9GnTjlbIvtz1e953Dg9A28wGOyDMaBLKe5K\nMb1RI14JDSXT0bPQLVvCG2/oKZICxgz+7OxZmnl6Fnvf5+sIDITFi2HRIiqU9WTliJWkZabxwA8P\nkJye7Gzp7M5bjRrxQUQEF9LS7NrP2cSz9Pq6F/2b9Gduv7lFusHaFB/P4cuXedLMPpdInnsOXnpJ\nu3PYKJpmseWtt8DNDV4reOAZg8FQBKwyoJVSDyqlDimlMpVSHfMo118pdUwpFayUmmRNnwbbMdDb\nm2plyrD43DnHd/700+DlpdN950NSRgbTTp9meqNGDhDMRsTH67Rfn38OPj4AlCtTjuXDllO1XFX6\nfduv4Gm/XRT/ChV4qGZNZoSH262PkIsh3P7l7fyjzT+Y3nt6kYxnEeHl0FDeatTI5SJvKKW8lFLr\nlFLHlVK/KaVyjO2Ymw7Orb5S6m6l1C6l1H6l1E6lVC9HjclejB+v9zD36gUhIc6Wxj6sX69VzuLF\nxu/ZYLA31v5bHATuB/7MrYBSyg1YAPQDWgMjlVItrOzXYAOUUsxs3JjXT50i1WJxdOc6KseiRbBx\nY55FP4iMJKBaNTpUruwg4axEBJ54AgYNgnvuue4tD3cP/nv/f+ng24GArwKITirZu5tea9CAr6Kj\nCU9JsXnbe87uoedXPXm5x8u8fMfLRW5nzcWLXMrM5GHXDC84GdggIs2BP4CbLkQ+Oji3+ueBe0Wk\nHTAa+Maeg3AU48bpiBS9esHx486WxrZERMCjj8K334IJYW4w2B+rDGgROS4iIUBe0z6dgRAROS0i\n6cASYIg1/RpsR49q1WhbqRILnJE9rk4dHd1/1Khct8nHpKXx/pkzvNWwoWNls4b58yEsLNecuW7K\njXn95/Fgqwfp/kV3jl045lj5HIhvuXKMr1OHV23sfLr+5Hr6f9ufBfcsYOytY4vcTrrFwsuhobzT\nqBHuruRb/zdDgK+zzr8G7suhTF46OMf6IrJfRKKzzg8D5ZVSHvYZgmMZM0a7OfTuXahtGMWa9HQY\nPlwHOerl8msFBoNr4Ij1yrpA9v3PEVmvGYoJ7/r7MzM8nLPOSN3Vt6+erR05Msf8u5NCQxnt60tT\nT0/Hy1YUgoLgnXfghx+gXLlciymlePXOV3m95+v0/Konf4bluojj8kyqX58/4uLYmpBgk/Y+3/M5\no34cxfLhy3mg5QNWtfVRZCS1y5VjkLe3TWRzAjVF5BxAlsFbM4cyeengWvnVV0o9COzJMr5LBNmN\n6JLgzjFpElSvrh8NBoNjyHe7uVJqPZB9bVMBAkwRkdX2EGpqtqjvAQEBBAQE2KMbQxbNPT35P19f\nXgoN5ZuWLR0vwOuvQ79+emNhNp/obQkJrIuN5Vjnzo6XqShcvKhjZn32GTRuXKAqo9uPxq+KH8OW\nDWNuv7mMumWUnYV0PJXLlGGOvz9Ph4Sws2NHyhTRz9giFl774zWWHl7K5sc308y7mVVynU1N5Z3w\ncDa3b1/kyC6BgYEEBgZaJUd+5KGDX82huLU7gq+rr5RqDcwA+uRVyRV19pgx+vGuu+CPP6BpU+fK\nU1SWLYMff4Tdu/XmQYPBkDs21dkiYvUBbAQ65vJeV2BttueTgUl5tCUGx5OYni71tm6VP+PinCPA\nuXMifn4iP/wgIiIZFot02LlTFkdHO0eewpKaKtKrl8hLLxWp+qFzh6TRvEby0rqXJCMzw8bCOR+L\nxSI99+yRD8+cKVL9SymX5L4l98ntX9wuMUkxNpFp1JEj8tKJEzZp6ypZ+ssmerUgB3AUPYsM4Asc\nzaFMrjo4r/pAPeA40DUfGWx6DR3NZ5+J1KsnEhLibEkKz/79IjVqiOze7WxJDAbXxBqdbcv71dym\ncHYCTZRSDZRSZYERwCob9muwAZXKlGGuvz/PhISQ4egNhQA1a8LKlTpg665dfBoVRRV3d0bWzGlF\nupghoqOKVKwI06cXqYnWNVuzY+wOdp3dxb3f30tccpyNhXQuSik+ataMN0+fJqaQYe1OxJ6g6xdd\nqelZkz8e+wOfij5Wy7M5Pp7A+Hhea9DA6raczCr0Jj+Ax4CVOZTJSwfnWF8pVQ1Ygza0g+wieTHh\niSf0Ithdd8GJE86WpuDExsL998O8edAx1xhYBoPBXlgbxu4+pdQZ9AzHGqXUr1mv11ZKrQEQkUzg\nGWAdcBhYIiJHrRPbYA8e9PGhpocHHzhjQyFAhw7w2WecHT2aqaGhLGja1DWSprz/PuzYYXXO3Bqe\nNfht1G80925O5887sy96nw2FdD6tK1bk0Vq1eLEQgXhXH1/N7V/ezsTOE/l00KeUdS9rtRypFgsT\nQkJ4z9+fSq6fNGUW0EcpdRzoDcyEQungHOsDTwP+wOtKqb1KqT1KqRqOGpSjGTtWJ1y56y7X8InO\nzNTbRu67Dx5+2NnSGAylEyWOTqKRD0opKW4ylSZOJifTdc8eAtu3p3XFig7vX0S4Z+VKOm/fzpsv\nvwxVqjhchkLx00969nnbNqhf32bNfn/weyaunci0XtMYd+s417iRKABJGRl02L2bWY0b84BP7jPJ\naZlpvLzhZZYfXc73Q7+nu193m8kw6eRJgpOT+V/r1ja/rkopRKRkfFgFpCTp7M8+gzffhN9/h+bN\nnS1N7jz3HBw6BGvXguvfAxoMzsManW22HBiuw79CBWY2bsyoo0dJc4IrxydRUVz08+PVxES49164\ncsXhMhSYDRt0YNmVK21qPAOMbDuSv8b8xce7PmbkipElxqWjUpkyfNOiBROCg3ON+hIaF8qd/7mT\nkNgQ9ozbY1PjeVN8PN+cO8eiZs1KzE2JwXaMHav3Md91Fxw54mxpcuaDD3TClOXLjfFsMDgTY0Ab\nbmKMry/1y5XjjbAwh/Z7/MoVXg8L45uWLfGYPx8aNYIHHgBnhNfLjy1b9BrqihVw22126aKZdzOC\n/i8IH08f2n3SjvUn19ulH0fTtWpVxtWpw/8dP072mUsRYdHuRXT5vAsPtX6IlSNW4u1pu/ByCRkZ\nPHbsGJ81b45PWetdQQwlk9GjYdYsHeJu715nS3M9P/0Es2fDL79AtWrOlsZgKN0YFw5DjsSkpdF+\n1y6WtGrFnQ7Q1GkWCz327mW0ry8T6maFqM3I0GHhRGDJEiguRs+uXTrD4Lff6jjWDmD9yfWMWTWG\nIc2HMKP3DCqXc5GsjLmQbrHQbc8extSuzYS6dYm8FMm4NeOITormm/u/oZVPK5v2JyKMPnaMCm5u\nfGLHtXnjwlFyWLFCp/9evhzuvNPZ0mgvsSFDtPFsp3t2g6HUYVw4DDanZtmyfNm8OSOOHCEsOdmu\nfYkI44ODqVO2LOPr1Pn7jTJl9MY8i0WnxU5KsqscBeL337Xx/PnnDjOeAfr49+HAUwe4nH6ZVh+3\n4qdjPzmsb3vg4ebGd61aMTUsjGe3fk67T9rRqU4ngv4vyObGM8D8yEh2Jybyrr+/zds2lEyGDtXq\nZ+hQWLPGubLs3KmN56+/NsazwVBcMDPQhjyZHxHBoqgo/urYkap2cribFR7O0pgYNrVvn3NUhIwM\neOopOHAAfv4Z8th8Zld++AGeeUZnLujZ0zkyAIFhgTy15ima12jOe33fo0n1Jk6TxRq2ndnGo39+\nwOnao1jVoiH967axSz9rLlxgXHAw2zp2pEH58nbp4ypmBrrksWMHDB6sMxeOG+f4/vfsgQED9D37\noEGO799gKMmYGWiD3Xi2bl16VqvG8MOH7RIfesX58yyIjGR127a5hxQrU0Zvj+/bF3r0gKMOjoIo\nAu++C88/r3fvONF4BghoGMD+p/bTtW5Xun7elefXPk9scqxTZSoMJ2NPMmzZMIYvH84bt9zLhy06\nMDHiMrHpts8UvT8piTHHj/NjmzZ2N54NJZPOnWHzZnjvPa0CMjMd1/eePXrB65NPjPFsMBQ3jAFt\nyBOlFB80aYKbUjx+/DjpNjSi18fG8lRwMCvbtKFuuXL5CaK3x0+erB0Sv/nGZnLkSVyczlbwww+w\ndSu0a+eYfvOhXJlyvHzHyxx5+gipmak0X9CcNwPfLNbROsLiw5jw8wS6fN6FDr4dOP7McUbdMoon\n69ZlsLc39x48SJwNjeijly9z78GDLGjalC7FPRyioVjTtCkEBelFsMGDIT7e/n2uWQP9+sHChVoF\nGQyG4oUxoA35UsbNjeWtWxObns6QQ4e4bIMpmMXnzjHq6FH+17o1HSsXYkPc449rP+Rp03QKsUuX\nrJYlV7Zu1Sm+GjXSUTeKYda6mhVr8vHAj9k6ZiunE07T5MMmvPL7K0ReclIynBw4cv4IY1aO4dZF\nt1KtfDWOPn2UV+54BU8Pz2tlZvv706VKFe7Yu5eIlBSr+9yWkECvffuY1qgRw10hm6XbRcM/AAAO\nJElEQVSh2OPlpeMu+/vrnE9//WW/vubP1+4iq1cb49lgKK4YH2hDgUm3WBgXHMyRy5dZ07ZtkUKB\niQhzIyL4ICKCX9q2pU2lSkUTJjFRZxNYu1bHnHr4YT1LbQvOndMz3evWwUcf6XRfLkJYfBhz/prD\n94e+p3fj3ky4bQI9G/bETTn2XjktM40fj/7Iwl0LCb4YzPjbxvNM52fwquCVax0R4d0zZ1gQGcmv\nt9xCqyIm8ll94QJjjh/nvy1aMMDbdmHwCoLxgS4drFqlDdwJE+CVV2wXjzk+XruJbN+ut3s0amSb\ndg0GQ85Yo7ONAW0oFCLC62FhfHn2LO/6+zOiZs0CJ6QITU5mYkgI4amp/Ny2LX628EkNCtKZAD09\n4bXX4O67wa2IxuLFi7BoEcydq4PBvvZa8c+EmAuJqYl8e+BbFu5aSFxKHA+1foiHWj/ErXVutZsx\nnZ6Zzh+n/mDp4aWsPL6SdrXaMaHTBIY0H4KHu0eB2/kmOprnT5zgpfr1ea5ePcoW8POMTU9nyqlT\nrLxwgZ/atKGzEz47Y0CXHiIjYcwY/fj++9Cnj3Xt/fqrNsrvvVfHei7MwpzBYCgaxoA2OJxtCQlM\nCAnBq0wZ5vj707FSpVwN6dj0dD6KjOSDiAhe9PPjBT+/AhtFBSIzU8dknjtXJ1155hm97nk1nnRe\nZGToGFFffKEDvw4ZApMmQcuWtpPPyRyKOcTSQ0v54cgPxKfE07tRb/o07kPXel1p5t0Mdzf3IrWb\nmpHK4fOH2Xx6MxtObWDT6U20rNGSh1o/xLDWw6hXpV6RZT5x5QrPnjjB6ZQU5vr706d6ddxz+X5d\nyczk+5gYpoSG8qCPD9MaNaKaR8ENdltiDOjShYiejf7Xv6BVK5gyRW86LMxi2M6dehFt1y6thnr3\ntp+8BoPheowBbXAKGRYLH0dF8X5EBACDvb3pVqXKNUMnIjWVVRcvsjsxkYHe3sxs3Ni+kRBE9Hb5\nhQu1+4WPj87J26QJ1KqlnycnQ0wMREfr+FSbN4Ofn07YMm4clHB/2bD4MDaEbmBD6AZ2Re0iOima\ntrXa0rR6U+pXrY9fFT+8Knjh6eFJRY+KZEomV9KvkJyezLnL5ziTcIYzl85w9MJRQi6G4F/dn271\nunF347vp1bAXPhVtF2JQRPjxwgXeCgsjMi2NgdWr09vLi/JZN1+xGRn8cvEiG+Pj6VKlCjMbN+ZW\nJ0/bGQO6dJKaCh9/rD2+KlbUKcHvuUe7YORkTJ85A5s26dB0oaHabeOJJ6CoHm0Gg6FoOM2AVko9\nCEwFWgKdRGRPLuXCgATAAqSLSOc82iz1ytjVEBEOX77MqosX2Zst2Yl3mTIM9Pamt5cXnu5Fm+Us\nMhYL7N8PgYFw+rQ2mmNitKtHzZramO7QAQICSrzRnBcJKQnsP7ef0LhQwhPCCU8IJyE1gSvpV7ic\ndhl3N3c8PTypUKYCPp4++FX1w6+KH81rNKdNzTaUL+OY0HCnU1JYfeECWxISuLqF1dPNjb7VqzOg\nenWqO2nG+UYcbUArpbyApUADIAwYLiIJOZTrD8xDbxz/QkRm5VVfKdUJWJStiTdFJMfsPUZn/43F\nolXOF1/oe/PERK1mqleH9HR9HD2qc0LdcYdO0jJ8OBSTr6/BUOpwpgHdHG0Ufwq8mIcBHQrcKiL5\nxtgyythgMLgqTjCgZwEXRWS2UmoS4CUik28o4wYEA72BKGAnMEJEjuVWXylVHkgTEYtSyhfYD9QW\nkZviWBqdnTsxMbB3rzakPTz0ZsPGjaFFC9vteTYYDEXHGp1t1d5hETmeJUB+nStMyDyDwWCwNUOA\nq5l9vgYCgck3lOkMhIjIaQCl1JKsesdyqy8i2WMJVkBPlBgKSc2aOpazwWAoeTjKqBVgvVJqp1Jq\nrIP6NBgMhpJOTRE5ByAi0UBO/kh1gTPZnkdkvQZQK7f6SqnOSqlD6Nnnp3KafTYYDIbSSr4z0Eqp\n9UCt7C+hDeIpIrK6gP3cLiJnlVI+aEP6qIhsya3w1KlTr50HBAQQEBBQwG4MBoPBcQQGBhIYGGjX\nPvLQwa/mUNxaX4pr9UVkB9Amy1Xvv0qpX0UkLadKRmcbDAZXwJY62yZROJRSG4F/5eYDfUPZN4BE\nEZmby/vGn85gMLgkTvCBPgoEiMi5LF/ljSLS8oYyXYGpItI/6/lkQERkVkHqZ9X5Hfh3Tjre6GyD\nweCqWKOzbenCkaMASilPpVSlrPOKQF/gkA37NRgMhtLKKmB01vljwMocyuwEmiilGiilygIjsurl\nWl8p1VAp5Z513gBojo7SYTAYDAasNKCVUvcppc4AXYE1Sqlfs16vrZRak1WsFrBFKbUXCAJWi8g6\na/o1GAwGAwCzgD5KqePoKBsz4XodLCKZwDPAOuAwsEREjuZVH+gB7FdK7QFWAONFJNZBYzIYDIZi\nj0mkYjAYDDbCJFIxGAwG16G4uHAYDAaDwWAwGAwlHmNAGwwGg8FgMBgMhcAY0AaDwWAwGAwGQyEw\nBrTBYDAYDAaDwVAIjAFtMBgMBoPBYDAUAmNAGwwGg8FgMBgMhcAY0AaDwWAwGAwGQyEwBrTBYDAY\nDAaDwVAIjAFtMBgMBoPBYDAUAmNAGwwGg8FgMBgMhcAY0AaDwWAwGAwGQyEwBrTBYDAYDAaDwVAI\nrDKglVKzlVJHlVL7lFIrlFJVcinXXyl1TCkVrJSaZE2fJY3AwEBni+BwzJhLB6VxzI5GKeWllFqn\nlDqulPpNKVU1l3I56uD86iul6iulEpVSL9h7LK5Eafxul7Yxl7bxQukcszVYOwO9DmgtIu2BEODl\nGwsopdyABUA/oDUwUinVwsp+Swyl8Qtrxlw6KI1jdgKTgQ0i0hz4g8Lr4Pzqvwf8YifZXZbS+N0u\nbWMubeOF0jlma7DKgBaRDSJiyXoaBNTLoVhnIERETotIOrAEGGJNvwaDwWAAtC79Ouv8a+C+HMrk\npYNzra+UGgKEAoftILfBYDC4NLb0gR4D/JrD63WBM9meR2S9ZjAYDAbrqCki5wBEJBqomUOZvHRw\nrRvq1wJQSlUCXgLeBJR9RDcYDAbXRYlI3gWUWk+WUr36EiDAFBFZnVVmCtBRRIbmUH8o0E9ExmU9\nHwV0FpGJufSXt0AGg8FQjBERmxqceejgV4GvRKR6trIXRcT7hvq56mClVJyIeN1YXyk1B9guIsuV\nUm8ASSLyXi7yGZ1tMBhclqLq7DIFaLhPXu8rpUYD9wB35VIkEqif7Xm9rNdy68/MdhgMBkMWeelg\npdQ5pVQtETmnlPIFYnIolpcOjs6lfhdgqFJqNuAFZCqlkkXk4xzkMzrbYDCUOqyNwtEf+DcwWERS\ncym2E2iilGqglCoLjABWWdOvwWAwGACtS0dnnT8GrMyhTF46OMf6InKniDQWkcbAPGB6TsazwWAw\nlFas9YH+EKgErFdK7VFKfQyglKqtlFoDICKZwDPoiB2HgSUictTKfg0Gg8EAs4A+SqnjQG9gJhRK\nB+dY32AwGAx5k68PtMFgMBgMBoPBYPgbp2QiLEhiFaXUfKVUSFaSlvaOltHW5DdmpdQ/lFL7s44t\nSqm2zpDTlhQ0gY5SqpNSKl0p9YAj5bMHBfxuByil9iqlDimlNjpaRltTgO92FaXUqqzf8sGsfRMu\ni1Lqiyzf4wN5lClV+iurTKkas9HZRme7KkZn51im8PpLRBx6oI32E0ADwAPYB7S4ocwA4Oes8y5A\nkKPldMKYuwJVs877l4YxZyv3O7AGeMDZcjvgc66KXkavm/W8hrPldsCYXwZmXB0vcBEo42zZrRhz\nD6A9cCCX90uj/iqNYzY62+hslzuMzs7x/SLpL2fMQBckscoQ4L8AIrIdqKqUqoXrku+YRSRIRBKy\nngbh+rGyC5pA51lgOTlHD3A1CjLmfwArRCQSQEQuOFhGW1OQMQtQOeu8MnBRRDIcKKNNEZEtQFwe\nRUqd/qIUjtnobKOzXRSjs2+mSPrLGQZ0QRKr3FgmMocyrkRhk8k8Qc5JaVyJfMeslKoD3CciCykZ\nyRoK8jk3A6orpTYqpXYqpR5xmHT2oSBjXgC0UkpFAfuBfzpINmdRGvVXaRxzdozOdk2MzjY6G4qo\nv/KNA21wLEqpXsDj6CWHks48ILv/VUlQyPlRBuiIjpteEdimlNomIiecK5Zd6QfsFZG7lFL+6Kg9\nt4hIkrMFMxisxejsEo/R2UZn54gzDOiCJFaJBPzyKeNKFCiZjFLqFmAR0F9E8lpucAUKMubbgCVK\nKYX2sxqglEoXEVeNE16QMUcAF0QkBUhRSm0C2qF90lyRgoz5cWAGgIicVEqdAloAuxwioeMpjfqr\nNI7Z6Gyjs10Ro7Nvpkj6yxkuHAVJrLIKeBRAKdUViBeRc44V06bkO2alVH1gBfCIiJx0goy2Jt8x\nS1aiBhFphPapm+DCihgK9t1eCfRQSrkrpTzRGxZcOS56QcZ8GrgbIMuvrBkQ6lApbY8i99m3Uqe/\nKIVjNjrb6GwXxejsmymS/nL4DLSIZCqlrgb1dwO+EJGjSqkn9duySER+UUrdo5Q6AVxG3w25LAUZ\nM/AaUB34OOvuPl1EOjtPauso4Jivq+JwIW1MAb/bx5RSvwEHgExgkYgccaLYVlHAz3ka8FW2EEIv\niUisk0S2GqXUd0AA4K2UCgfeAMpSivVXaRwzRmcbne2CGJ1tO51tEqkYDAaDwWAwGAyFwCmJVAwG\ng8FgMBgMBlfFGNAGg8FgMBgMBkMhMAa0wWAwGAwGg8FQCIwBbTAYDAaDwWAwFAJjQBsMBoPBYDAY\nDIXAGNAGg8FgMBgMBkMhMAa0wWAwGAwGg8FQCP4fLMYuhCsoUG4AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10b179310>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Define a figure to take two plots\n",
    "fig = pl.figure(figsize=[12,3])\n",
    "# Add subplots: number in y, x, index number\n",
    "ax = fig.add_subplot(121,autoscale_on=False,xlim=(0,1),ylim=(-2,2))\n",
    "ax.set_title(\"Eigenvectors for changed system\")\n",
    "ax2 = fig.add_subplot(122,autoscale_on=False,xlim=(0,1),ylim=(-0.004,0.004))\n",
    "ax2.set_title(\"Difference to perfect eigenvectors\")\n",
    "for m in range(4): # Plot the first four states\n",
    "    psi = np.zeros(number_of_x_points)\n",
    "    for i in range(number_of_basis_functions):\n",
    "        psi = psi+eigenvectors[i,m]*basis_functions[i]\n",
    "    if eigenvectors[m,m] < 0:  # This is just to ensure that psi and the basis function have the same phase\n",
    "        psi = -psi\n",
    "    ax.plot(x,psi)\n",
    "    # Now subtract the unperturbed eigenvector to see the change from the perturbation\n",
    "    psi = psi - basis_functions[m]\n",
    "    ax2.plot(x,psi)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Notice how the eigenvectors for this changed system are *very* close to the original system.  This is a perfect example of a good system to study with perturbation theory.  The changes in the eigenvectors can be related to matrix elements of the potential and the eigenvalues of the unperturbed system, but we won't do this just yet.\n",
    "\n",
    "In the next notebook, we will explore a more complex perturbation (a truncated quantum harmonic oscillator, or QHO) and the effects of varying the magnitude of the perturbation.  As well as improving your understanding of matrix mechanics, this should give a useful insight into the limitations of perturbation theory and finite basis sets."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}