{
 "metadata": {
  "name": "",
  "signature": "sha256:a5146b953e6ce4107eb6016b9e4a4a34b2153b7ea3de0a799042c5e879c6ce01"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "<img src=\"https://raw.githubusercontent.com/mar-one/ACM-Python-Tutorials-KAUST-2015/master/images/ACM_chapter_logo.png\" width=\"250px\" style=\"float: left;\">\n",
      "<img src=\"https://raw.githubusercontent.com/mar-one/ACM-Python-Tutorials-KAUST-2015/master/images/KAUST_logo.png\" width=\"250px\" style=\"float: right;\">\n",
      "\n",
      "<h1 align=\"center\">Introduction to NumPy & SciPy - the core numerical scientific packages</h1>\n",
      "<p> <br style=\"clear: both;\"/></p>\n",
      "<p> <br style=\"clear: both;\"/></p>\n",
      "<p> <br style=\"clear: both;\"/></p>"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "<p style=\"color: #999;\">This tutorial is prepared by ACM Student Chapter of King Abdullah University of Science and Technology (KAUST). <br/>\n",
      "Material has been adapted from the following:</p>\n",
      "\n",
      "<p style=\"color: #999;\"> - David Ketcheson's [Introduction to NumPy and Matplotlib](http://nbviewer.ipython.org/github/ketch/AMCS252/blob/master/2_Introduction_to_Numpy%20and_Matplotlib.ipynb)</p>\n",
      "<p style=\"color: #999;\">- Software Carpentry's [Numerical analysis with NumPy](http://nbviewer.ipython.org/github/geocarpentry/2014-01-30-mit/blob/gh-pages/lessons/vignettes/Numerical%20analysis%20with%20NumPy.ipynb)</p>\n",
      "<p style=\"color: #999;\">- J.R. Johansson's [Numpy - multidimensional data arrays](http://nbviewer.ipython.org/github/jrjohansson/scientific-python-lectures/blob/master/Lecture-2-Numpy.ipynb) and [SciPy - Library of scientific algorithms for Python](http://nbviewer.ipython.org/github/jrjohansson/scientific-python-lectures/blob/master/Lecture-3-Scipy.ipynb#SciPy---Library-of-scientific-algorithms-for-Python)</p>\n",
      "<p style=\"color: #999;\">- Official [Tentative NumPy tutorial](http://www.scipy.org/Tentative_NumPy_Tutorial)</p>\n",
      "\n",
      "**Prerequisites:** Introduction to Python\n",
      "<p><br/></p>"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## NumPy"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Python includes a package for numerical computation called numpy, which is an essential tool for scientific computing. \n",
      "It is based on multidimensional arrays (vectors and matrices) and provides a large library of functions and operators that work efficiently on these objects.\n",
      "In this way, any numerical algorithm is expressed as operations on arrays and in many cases runs almost as quickly as the equivalent C code.\n",
      "\n",
      "Like MATLAB, `NumPy` relies on BLAS and LAPACK for efficient linear algebra computations.\n",
      "\n",
      "`NumPy` is used in:</br>\n",
      "Numerical Analysis, Linear algebra, Solution of differential equations, Statistical analysis and many more...\n"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Why to use NumPy?\n",
      "\n",
      "Consider the following list and suppose we want to perform some simple algebraic operations:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "price = [5.99, 10.25, 2.0, 40.99, 5.60, 63.49]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 2
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "price*2"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 3,
       "text": [
        "[5.99, 10.25, 2.0, 40.99, 5.6, 63.49, 5.99, 10.25, 2.0, 40.99, 5.6, 63.49]"
       ]
      }
     ],
     "prompt_number": 3
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "price+0.5"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "ename": "TypeError",
       "evalue": "can only concatenate list (not \"float\") to list",
       "output_type": "pyerr",
       "traceback": [
        "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m\n\u001b[0;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
        "\u001b[0;32m<ipython-input-4-ba6ebb983a2c>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mprice\u001b[0m\u001b[0;34m+\u001b[0m\u001b[0;36m0.5\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
        "\u001b[0;31mTypeError\u001b[0m: can only concatenate list (not \"float\") to list"
       ]
      }
     ],
     "prompt_number": 4
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Clearly, if we wanted to double the value of each element and add half, this cannot be done easily. Instead, we would have to use a loop:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "for i in range(len(price)):\n",
      "    price[i] = price[i]*2\n",
      "    price[i] = price[i]+0.5\n",
      "print price"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[12.48, 21.0, 4.5, 82.48, 11.7, 127.48]\n"
       ]
      }
     ],
     "prompt_number": 5
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "But, there should be a more convinient way..."
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Getting started\n",
      "\n",
      "We start, by importing the NumPy module. Import statements like this are the typical way of getting access to functions in Python. The most general way to import a module is the following:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from numpy import *"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 6
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "The above imports all functions included in the package and can call them simply by their name."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "sin(3*pi/2)*exp(2)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 7,
       "text": [
        "-7.3890560989306504"
       ]
      }
     ],
     "prompt_number": 7
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Instead we can import the numpy module and tell Python that we want to refer to numpy by the short abbreviation *np*. This is recommended since we keep track of which module the functions belong to."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import numpy as np"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 8
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Arrays\n",
      "\n",
      "NumPy's main object is the homogeneous multidimensional array. \n",
      "It is a table of elements - all of the same type - indexed by a tuple of positive integers.\n",
      "\n",
      "NumPy's array class is called *ndarray* and the package includes basic utility functions for creating and manipulate arrays in a similar fashion to MATLAB.\n",
      "For example we can create a one-dimensional array by:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "price = np.array([5.99, 10.25, 2.0, 40.99, 5.60, 63.49])\n",
      "print type(price)\n",
      "price = price*2+0.5\n",
      "print price"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "<type 'numpy.ndarray'>\n",
        "[  12.48   21.      4.5    82.48   11.7   127.48]\n"
       ]
      }
     ],
     "prompt_number": 9
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "There are other ways to create one-dimensional arrays:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "x = np.linspace(-1,1,5)\n",
      "print x"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[-1.  -0.5  0.   0.5  1. ]\n"
       ]
      }
     ],
     "prompt_number": 10
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "y = np.arange(0,1,0.2)\n",
      "print y"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[ 0.   0.2  0.4  0.6  0.8]\n"
       ]
      }
     ],
     "prompt_number": 11
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "---\n",
      "### Exercise 7\n",
      "\n",
      "Change the inputs of the *linspace* and *arange* function until you understand exactly what they do. Using both ways, create a vector of starting from -1.5 to 4.5 inclusive with a step of 0.5. "
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# type your solution here"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 12
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "##### Solution"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "#%load solutions/1D_array.py"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 13
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We can always find out what a function does by typing its name, followed by a question mark:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "np.arange?"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 14
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Resize the help window (to get it out of your way) by dragging the divider. Double click on the divider to close it."
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Multidimensional arrays\n",
      "\n",
      "We can create multidimensional arrays, like the following:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = np.array([[1,2.4,-13],[4.1,5,0],[7.2,8,9]])\n",
      "print \"A = \", A\n",
      "print type(A)\n",
      "print A.shape\n",
      "print A.dtype"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "A =  [[  1.    2.4 -13. ]\n",
        " [  4.1   5.    0. ]\n",
        " [  7.2   8.    9. ]]\n",
        "<type 'numpy.ndarray'>\n",
        "(3, 3)\n",
        "float64\n"
       ]
      }
     ],
     "prompt_number": 15
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "To see all different attributes of the array object type its name followed by a fullstop and press *tap*. Note that some attributes are basically functions and require parentheses (with or without arguments), e.g. *conjugate* and others don't, e.g. *size*:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = np.array([[1,2+4j,3j],[4j,5,6-5j],[7j,8,9+3j]]) # dtype=complex\n",
      "print A.conjugate()\n",
      "print A.size"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[[ 1.-0.j  2.-4.j  0.-3.j]\n",
        " [ 0.-4.j  5.-0.j  6.+5.j]\n",
        " [ 0.-7.j  8.-0.j  9.-3.j]]\n",
        "9\n"
       ]
      }
     ],
     "prompt_number": 16
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We can also define the shape of the array:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "np.random.uniform(0,1,size=(5,5)) # 5 x 5 matrix with elements from a uniform distribution"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 17,
       "text": [
        "array([[ 0.18894904,  0.4016313 ,  0.30167486,  0.04452608,  0.29880108],\n",
        "       [ 0.85862765,  0.18898635,  0.46421297,  0.6303044 ,  0.56185042],\n",
        "       [ 0.78056999,  0.82705839,  0.63247392,  0.07986216,  0.64999473],\n",
        "       [ 0.07852346,  0.27204531,  0.85521149,  0.94891731,  0.7810846 ],\n",
        "       [ 0.47172365,  0.57284027,  0.57618836,  0.30150177,  0.54945194]])"
       ]
      }
     ],
     "prompt_number": 17
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "#### Special arrays\n",
      "\n",
      "The function *zeros* creates an array full of zeros, the function *ones* creates an array full of ones and *eye* an identity array. Function *diag* is used for a diagonal array with specified diagonal elements. A matrix with random entries is generated by the *random* function."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print np.zeros((3,4))\n",
      "print np.ones((2,3))\n",
      "print np.eye(3)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[[ 0.  0.  0.  0.]\n",
        " [ 0.  0.  0.  0.]\n",
        " [ 0.  0.  0.  0.]]\n",
        "[[ 1.  1.  1.]\n",
        " [ 1.  1.  1.]]\n",
        "[[ 1.  0.  0.]\n",
        " [ 0.  1.  0.]\n",
        " [ 0.  0.  1.]]\n"
       ]
      }
     ],
     "prompt_number": 18
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print np.diag(x,0)\n",
      "print np.diag([2,3,4])\n",
      "print np.random.random((2,2))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[[-1.   0.   0.   0.   0. ]\n",
        " [ 0.  -0.5  0.   0.   0. ]\n",
        " [ 0.   0.   0.   0.   0. ]\n",
        " [ 0.   0.   0.   0.5  0. ]\n",
        " [ 0.   0.   0.   0.   1. ]]\n",
        "[[2 0 0]\n",
        " [0 3 0]\n",
        " [0 0 4]]\n",
        "[[ 0.89073408  0.92386659]\n",
        " [ 0.43499535  0.1673178 ]]\n"
       ]
      }
     ],
     "prompt_number": 19
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Notice that there are two pairs of parentheses used when calling the *zeros*, *ones* and *random* command, because we are passing as argument (e.g. the Python tuple *(3,4)*), which is the shape of the resulting array."
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Reshape\n",
      "\n",
      "Manipulation of arrays can be done in many ways:  "
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = np.arange(1,10)\n",
      "print \"A = \", A\n",
      "A = A.reshape((3,3))\n",
      "print A"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "A =  [1 2 3 4 5 6 7 8 9]\n",
        "[[1 2 3]\n",
        " [4 5 6]\n",
        " [7 8 9]]\n"
       ]
      }
     ],
     "prompt_number": 20
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "B = arange(15).reshape((3,5))\n",
      "print \"B = \", B\n",
      "print B.ravel() # flatten the array"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "B =  [[ 0  1  2  3  4]\n",
        " [ 5  6  7  8  9]\n",
        " [10 11 12 13 14]]\n",
        "[ 0  1  2  3  4  5  6  7  8  9 10 11 12 13 14]\n"
       ]
      }
     ],
     "prompt_number": 21
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Indexing\n",
      "\n",
      "We can slice numpy arrays just as other programming languages, but the indexing is a little different. See the examples below:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print A\n",
      "print A[0,0]\n",
      "print A[1:3,0:2]\n",
      "print A[1:3,[0,1]]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[[1 2 3]\n",
        " [4 5 6]\n",
        " [7 8 9]]\n",
        "1\n",
        "[[4 5]\n",
        " [7 8]]\n",
        "[[4 5]\n",
        " [7 8]]\n"
       ]
      }
     ],
     "prompt_number": 22
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Note that the fundamental difference from MATLAB is that the indexing starts from zero and is exclusive of the second index. Hence, in the above example 1:3 refers to the second and third row and 0:2 to the first and second column."
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We can index in some clever ways too:  "
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print A"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[[1 2 3]\n",
        " [4 5 6]\n",
        " [7 8 9]]\n"
       ]
      }
     ],
     "prompt_number": 23
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print A[0,:] # row selection\n",
      "print A[:,0] # column selection"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[1 2 3]\n",
        "[1 4 7]\n"
       ]
      }
     ],
     "prompt_number": 24
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print A[:,1:] # all rows, from 2nd column until the end\n",
      "print A[-1,:] # last row, all columns"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[[2 3]\n",
        " [5 6]\n",
        " [8 9]]\n",
        "[7 8 9]\n"
       ]
      }
     ],
     "prompt_number": 25
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We can change array's elements:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A[:,2] = [.3,.4,9.12]\n",
      "print A\n",
      "print A.dtype"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[[1 2 0]\n",
        " [4 5 0]\n",
        " [7 8 9]]\n",
        "int64\n"
       ]
      }
     ],
     "prompt_number": 26
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Notice that the elements of A remain integers. That is because A has type *int*. The following does the desired change:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = A.astype(float)\n",
      "A[:,2] = [.3,.4,9.12]\n",
      "print A"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[[ 1.    2.    0.3 ]\n",
        " [ 4.    5.    0.4 ]\n",
        " [ 7.    8.    9.12]]\n"
       ]
      }
     ],
     "prompt_number": 27
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "There is an easy way to reverse a vector:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "a = arange(10)\n",
      "print a\n",
      "print a[: :-1] # reverse a"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[0 1 2 3 4 5 6 7 8 9]\n",
        "[9 8 7 6 5 4 3 2 1 0]\n"
       ]
      }
     ],
     "prompt_number": 28
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We can also index an array using boolean arrays:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "b = a>5\n",
      "print b\n",
      "a[b] = 0  # All elements of 'a' higher than 5 become 0\n",
      "print a"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[False False False False False False  True  True  True  True]\n",
        "[0 1 2 3 4 5 0 0 0 0]\n"
       ]
      }
     ],
     "prompt_number": 29
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Array and matrix multiplication"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print a**2\n",
      "print 10*sin(np.pi/a[1:5]) # avoiding division with zero!\n",
      "print A*A"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[ 0  1  4  9 16 25  0  0  0  0]\n",
        "[  1.22464680e-15   1.00000000e+01   8.66025404e+00   7.07106781e+00]\n",
        "[[  1.       4.       0.09  ]\n",
        " [ 16.      25.       0.16  ]\n",
        " [ 49.      64.      83.1744]]\n"
       ]
      }
     ],
     "prompt_number": 30
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "By default, NumPy uses component-wise multiplication. For a matrix-matrix (or matrix-vector) multiplication, we use the *dot* function:  "
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "np.dot(A,A)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 31,
       "text": [
        "array([[  11.1   ,   14.4   ,    3.836 ],\n",
        "       [  26.8   ,   36.2   ,    6.848 ],\n",
        "       [ 102.84  ,  126.96  ,   88.4744]])"
       ]
      }
     ],
     "prompt_number": 31
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Other basic operations"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = array(A,dtype='int')\n",
      "print A"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[[1 2 0]\n",
        " [4 5 0]\n",
        " [7 8 9]]\n"
       ]
      }
     ],
     "prompt_number": 32
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print A.sum() # sum of all elements\n",
      "print A.min()\n",
      "print A.max()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "36\n",
        "0\n",
        "9\n"
       ]
      }
     ],
     "prompt_number": 33
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print A.sum(axis=0)     # sum of each column\n",
      "print A.min(axis=1)     # min of each row\n",
      "print A.cumsum(axis=1)  # cumulative sum along each row"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[12 15  9]\n",
        "[0 0 7]\n",
        "[[ 1  3  3]\n",
        " [ 4  9  9]\n",
        " [ 7 15 24]]\n"
       ]
      }
     ],
     "prompt_number": 34
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "In Python *axis* refers to the dimensions of arrays; *axis=0* refers to the colums and *axis=1* to the rows."
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Stacking and splitting\n",
      "\n",
      "The *vstack* and *hstack* commands can be used to stack arrays and *hsplit* to split."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = floor(10*random.random((2,2)))\n",
      "print A\n",
      "B = floor(10*random.random((2,2)))\n",
      "print B"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[[ 4.  0.]\n",
        " [ 6.  7.]]\n",
        "[[ 1.  2.]\n",
        " [ 5.  3.]]\n"
       ]
      }
     ],
     "prompt_number": 35
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print np.vstack((A,B))\n",
      "print np.hstack((B,A))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[[ 4.  0.]\n",
        " [ 6.  7.]\n",
        " [ 1.  2.]\n",
        " [ 5.  3.]]\n",
        "[[ 1.  2.  4.  0.]\n",
        " [ 5.  3.  6.  7.]]\n"
       ]
      }
     ],
     "prompt_number": 36
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = floor(10*random.random((2,12)))\n",
      "print A"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[[ 6.  7.  5.  2.  3.  6.  2.  4.  7.  5.  2.  9.]\n",
        " [ 8.  4.  4.  5.  7.  0.  5.  3.  3.  3.  6.  3.]]\n"
       ]
      }
     ],
     "prompt_number": 37
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print np.hsplit(A,3)       # Split A into 3 arrays\n",
      "print np.hsplit(A,(3,4))   # Split A after the third and the fourth column"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[array([[ 6.,  7.,  5.,  2.],\n",
        "       [ 8.,  4.,  4.,  5.]]), array([[ 3.,  6.,  2.,  4.],\n",
        "       [ 7.,  0.,  5.,  3.]]), array([[ 7.,  5.,  2.,  9.],\n",
        "       [ 3.,  3.,  6.,  3.]])]\n",
        "[array([[ 6.,  7.,  5.],\n",
        "       [ 8.,  4.,  4.]]), array([[ 2.],\n",
        "       [ 5.]]), array([[ 3.,  6.,  2.,  4.,  7.,  5.,  2.,  9.],\n",
        "       [ 7.,  0.,  5.,  3.,  3.,  3.,  6.,  3.]])]\n"
       ]
      }
     ],
     "prompt_number": 38
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "---\n",
      "### Exercise 8\n",
      "\n",
      "1. Create a vector with values ranging from 10 to 99\n",
      "2. Find indices of non-zero elements from [1,2,0,0,4,0]. Hint: use <em>np.nonzero</em>\n",
      "3. Declare a 8x8 matrix and fill it with a checkerboard pattern of zeros and ones, without a loop.\n",
      "4. Create a 5x3 matrix whose $(i,j)$ entry is equal to $i\\times j^2$ without using a loop. Hint: use <em>np.fromfunction</em>.\n",
      "5. Create a 4x4 array with numbers from a normal distribution with mean = 3 and standard deviation = 1.5. Hint: use <em>np.random.normal</em>. Sum the values of each row and store the result in a new array."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# type your solution here"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 39
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "##### Solution"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "#%load solutions/array_operations.py"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 40
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Now let's work with some real data. Ensure that the <em>stockholm_temp.txt</em> in the same directory with this notebook. \n",
      "This file contains the [daily average temperatures](http://bolin.su.se/data/stockholm/homogenized_daily_mean_temperatures.php)\n",
      "according to observations of Stockholm for the years 1756 - 2012.\n",
      "The following code checks if the file exists, opens it, prints the first lines and stores the data in an array."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import os\n",
      "filename = 'data/stockholm_temp.txt'\n",
      "if os.path.exists(filename):\n",
      "    # look at the first 3 lines\n",
      "    with open(filename,'r') as f:\n",
      "        print '\\n'.join(f.readlines()[:3])\n",
      "    # alternative look at the first head lines\n",
      "    !head 'data/stockholm_temp.txt'\n",
      "    \n",
      "    data = np.loadtxt(filename) # alternative we can use getfromtxt command\n",
      "    \n",
      "else:\n",
      "    print 'Weather data does not exist, please download \"stockholm_temp.txt\" from Github.'"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "1.756000000000000000e+03 1.000000000000000000e+00 1.000000000000000000e+00 -8.699999999999999289e+00\n",
        "\n",
        "1.756000000000000000e+03 1.000000000000000000e+00 2.000000000000000000e+00 -9.199999999999999289e+00\n",
        "\n",
        "1.756000000000000000e+03 1.000000000000000000e+00 3.000000000000000000e+00 -8.599999999999999645e+00\n",
        "\n",
        "1.756000000000000000e+03 1.000000000000000000e+00 1.000000000000000000e+00 -8.699999999999999289e+00\r\n",
        "1.756000000000000000e+03 1.000000000000000000e+00 2.000000000000000000e+00 -9.199999999999999289e+00\r\n",
        "1.756000000000000000e+03 1.000000000000000000e+00 3.000000000000000000e+00 -8.599999999999999645e+00\r\n",
        "1.756000000000000000e+03 1.000000000000000000e+00 4.000000000000000000e+00 -7.700000000000000178e+00\r\n",
        "1.756000000000000000e+03 1.000000000000000000e+00 5.000000000000000000e+00 -7.200000000000000178e+00\r\n",
        "1.756000000000000000e+03 1.000000000000000000e+00 6.000000000000000000e+00 -1.600000000000000089e+00\r\n",
        "1.756000000000000000e+03 1.000000000000000000e+00 7.000000000000000000e+00 6.999999999999999556e-01\r\n",
        "1.756000000000000000e+03 1.000000000000000000e+00 8.000000000000000000e+00 1.300000000000000044e+00\r\n",
        "1.756000000000000000e+03 1.000000000000000000e+00 9.000000000000000000e+00 2.399999999999999911e+00\r\n",
        "1.756000000000000000e+03 1.000000000000000000e+00 1.000000000000000000e+01 8.000000000000000444e-01\r\n"
       ]
      }
     ],
     "prompt_number": 41
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Since the data are in the form\n",
      "*Year\n",
      "- Month\n",
      "- Day\n",
      "- Temperature*\n",
      "we would like to express them in a better way. We define a `dtype`, that is the type for each column: int, int, int, float:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "dt = np.dtype([('Year', 'int16'), ('Month', 'int8'), ('Day', 'int8'), ('Temp', 'float64')])\n",
      "data = np.loadtxt(filename,dtype=dt)\n",
      "data[:10] # first 10 entries of our data"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 42,
       "text": [
        "array([(1756, 1, 1, -8.7), (1756, 1, 2, -9.2), (1756, 1, 3, -8.6),\n",
        "       (1756, 1, 4, -7.7), (1756, 1, 5, -7.2), (1756, 1, 6, -1.6),\n",
        "       (1756, 1, 7, 0.7), (1756, 1, 8, 1.3), (1756, 1, 9, 2.4),\n",
        "       (1756, 1, 10, 0.8)], \n",
        "      dtype=[('Year', '<i2'), ('Month', 'i1'), ('Day', 'i1'), ('Temp', '<f8')])"
       ]
      }
     ],
     "prompt_number": 42
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We can now look for the example only on the temperatures:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "data['Temp']"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 43,
       "text": [
        "array([-8.7, -9.2, -8.6, ...,  0.2,  2.6,  2.2])"
       ]
      }
     ],
     "prompt_number": 43
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "or find all leap years:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "years = np.unique(data['Year']) # making an array with each year appearing only once\n",
      "leap_years = np.zeros((1,),dtype='int16') # initializing array for leap years\n",
      "for i in range(len(years)):\n",
      "    if np.fmod(years[i],4)==0 and ((np.fmod(years[i],100)==0 and np.fmod(years[i],400)==0) or np.fmod(years[i],100)!=0):\n",
      "        leap_years = np.append(leap_years,years[i]) # adding a leap year\n",
      "leap_years = leap_years[1:] # removing dummy first value\n",
      "print leap_years"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[1756 1760 1764 1768 1772 1776 1780 1784 1788 1792 1796 1804 1808 1812 1816\n",
        " 1820 1824 1828 1832 1836 1840 1844 1848 1852 1856 1860 1864 1868 1872 1876\n",
        " 1880 1884 1888 1892 1896 1904 1908 1912 1916 1920 1924 1928 1932 1936 1940\n",
        " 1944 1948 1952 1956 1960 1964 1968 1972 1976 1980 1984 1988 1992 1996 2000\n",
        " 2004 2008 2012]\n"
       ]
      }
     ],
     "prompt_number": 44
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "---\n",
      "### Exercise 9\n",
      "\n",
      "Make an array with the temperatures of the 5th day of every month."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# type your solution here"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 45
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "##### Solution"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "#%load solutions/temperatures.py"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 46
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Linear algebra\n",
      "\n",
      "Handling *ndarrays* is really powerful. Of the many things one can do with NumPy's arrays is linear algebra:  "
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "x = np.random.random(10)\n",
      "print x\n",
      "np.linalg.norm(x,np.inf) # maximum norm of a vector"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[ 0.05834554  0.05282837  0.39608145  0.37972821  0.59102374  0.38928025\n",
        "  0.00254778  0.4471012   0.22718756  0.32243128]\n"
       ]
      },
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 47,
       "text": [
        "0.59102373985998502"
       ]
      }
     ],
     "prompt_number": 47
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We can use this module to solve a linear system of equations $Ax=b$:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = np.array([[0,2],[8,0]])\n",
      "b = np.array([1,2])\n",
      "print A\n",
      "print b"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[[0 2]\n",
        " [8 0]]\n",
        "[1 2]\n"
       ]
      }
     ],
     "prompt_number": 48
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "x = np.linalg.solve(A,b)\n",
      "print x"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[ 0.25  0.5 ]\n"
       ]
      }
     ],
     "prompt_number": 49
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We can also solve eigenvalue problems:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "lamda,V = np.linalg.eig(A)\n",
      "print lamda\n",
      "print V"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[ 4. -4.]\n",
        "[[ 0.4472136  -0.4472136 ]\n",
        " [ 0.89442719  0.89442719]]\n"
       ]
      }
     ],
     "prompt_number": 50
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Notice we have put two variables on the left of the equals sign, to assign the outputs of *eig()* to two different variables. Here *lamda* are the eigenvalues and the columns of *V* are the eigenvectors."
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "---\n",
      "### Exercise 10\n",
      "\n",
      "Create a diagonal matrix using the *lamda* variable and then multiply from the left with V and from the right with the inverse of V (use *np.linalg.inv()*). What do you find?"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# type your solution here"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 51
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "##### Solution"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "#%load solutions/linear_algebra.py"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 52
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "<p></br></p>\n",
      "<p></br></p>\n",
      "<p></br></p>"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## SciPy"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "`NumPy` is complemented by `SciPy`; a library that adds more MATLAB-like functionality.\n",
      "\n",
      "`SciPy` is a collection of mathematical algorithms and convenience functions that provides high-level commands and classes for Nowdays, everything from parallel programming to web and data-base subroutines and classes have been made available to Python All of this power is also accessible through the the mathematical libraries in SciPy `NumPy` and `SciPy`.\n",
      "\n",
      "`SciPy` is used in:</br>\n",
      "Numerical Analysis, Optimization, Linear algebra, Image and signal processing, Statistical analysis, Data processing, Data transformation and query and many more..."
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Getting started"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "As with NumPy, we can access all packages in SciPy:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from scipy import *"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 53
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Special functions\n",
      "\n",
      "SciPy provides implementations of a large set of special functions, usually used in mathematics, physics and engineering. \n",
      "Available functions include airy, elliptic, bessel, gamma, beta, hypergeometric, parabolic cylinder, mathieu, spheroidal wave, struve, and kelvin.\n",
      "\n",
      "For a list of functions see SciPy's documention at http://docs.scipy.org/doc/scipy/reference/special.html#module-scipy.special.\n",
      "or try"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from scipy import special\n",
      "#help(special)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 54
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "For example we can look at the Bessel functions.\n",
      "\n",
      "Bessel functions are a family of solutions to Bessel\u2019s differential equation with real or complex order alpha:\n",
      "\\begin{align}\n",
      "    x^2 \\frac{d^2 y}{d x^2} + x \\frac{d y}{d x} + (x^2 \u2212\\alpha^2)y = 0\n",
      "\\end{align}"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from scipy.special import jn, jn_zeros\n",
      "\n",
      "x = 0.0\n",
      "\n",
      "n = 0    # order\n",
      "# Bessel function of first kind, zero order \n",
      "print \"J_%d(%f) = %f\" % (n, x, jn(n, x))\n",
      "\n",
      "n = 1    # order\n",
      "# Bessel function of first kind, first order \n",
      "print \"J_%d(%f) = %f\" % (n, x, jn(n, x))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "J_0(0.000000) = 1.000000\n",
        "J_1(0.000000) = 0.000000\n"
       ]
      }
     ],
     "prompt_number": 55
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "The next two lines allow plotting inside the ipython notebook. More during the next sessions..."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# The next two\n",
      "%pylab inline --no-import-all\n",
      "import matplotlib.pyplot as plt"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "Populating the interactive namespace from numpy and matplotlib\n"
       ]
      }
     ],
     "prompt_number": 56
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We can now plot the Bessel functions of first kind for different values of order:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "x = linspace(0, 12, 100)\n",
      "\n",
      "fig, ax = plt.subplots()\n",
      "ax.plot(x,zeros(len(x)),'--k')\n",
      "for n in range(4):\n",
      "    ax.plot(x, jn(n, x), label=r\"$J_%d(x)$\" % n)\n",
      "ax.legend();"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAD7CAYAAABpJS8eAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzsnXlcVNUbxp/JJZdUkEVc0yT3JTVNW9wyyzAtW0zNyjbN\nXDLLTEtbLXHL3FL7aVnumOa+awqIILIIioigILKq7AzL3Of3xxETZWCGuXcG4Xw/n/nkcM89553b\nzHPPfc973ldHEhKJRCIp/9xnawMkEolEYh2k4EskEkkFQQq+RCKRVBCk4EskEkkFQQq+RCKRVBCk\n4EskEkkFobKtDShAp9PJ+FCJRCIxE5I6U9uWqRk+yXL5mjlzps1tkJ9Pfj75+crfy1zKlOBLJBKJ\nRDuk4EskEkkFQQq+Fejdu7etTdAU+fnubeTnqzjoSuMH0gKdTseyYotEIpHcC+h0OtCMRdsyE6Uj\nkUgqLjqdyZpVYVFjQiwFXyKRlAnkE75x1LohSh++RCKRVBCk4EskEkkFQQq+RCKRVBCk4EskEkkF\nQQq+RCKRaERUVFSJbeLi4pCVlWUFa6TgSyQSiSZERkbCx8enxHZOTk5wd3e3gkVS8CUSiUQTli9f\njmHDhpXYrnLlynBzc8OaNWs0t0kKvkQikZjB+vXr4ebmhm7dusHLy6vINkFBQWjUqJHJfXbt2hUH\nDx5Uy0SjSMGXSCQSMxg2bBhq166NSZMm4Yknniiyzc6dO9G3b1+z+nVyckJERIQaJhpFCr5EIpGY\nAUkcOXKkWEH38/NDmzZtzOq3Y8eO8Pf3t9S8YpGpFSQSicQMgoKC4OjoiHr16hltk5WVdVc6hO3b\nt6NSpUo4fvw42rdvj71792L69Olo1aoVAMDe3h7h4eGa2i4FXyKRlHnUyq2mRrqeQ4cOoV+/fsW2\nMRgMhd5HR0ejTZs2cHV1xYwZMzB16lTUqVMHTZo0udWmevXqyM3NtdzAYrBI8HU63SoAbgASSbY3\n0uYXAAMAZAF4m2SAJWNKJJKKR1nKq3bw4EGMHTsWN27cwIoVK+Ds7IwOHTqgS5cut9pUrlxYWguE\nPSEhAbVq1YKdnR0GDhxYqE1qairq1q2rqe2W+vBXA3jO2EGdTvc8AFeSDwP4AMAyC8eTSCQSq5Kd\nnY0GDRogMjIS8fHxCA0NRd++fbF69Wr06dMHI0eOxPz58wud4+LigoyMjFvvw8LCEBQUhN27d6Nn\nz54AxMLu7cTFxcHV1VXTz2KR4JM8DuBGMU0GAfjjZtuTAOx0Op1xx5dEIpGUMapWrYqxY8ciODgY\n8+bNw969e1GzZk1ERUWhfv36qFy5Mq5fv17onF69esHX1/fW+/3792Pnzp0gCb1ej61bt8LZ2bnQ\nOYGBgUajftRCax9+QwAxt72/AqARgASNx5VIJBJVqFSpEr788ksAwIsvvnjr74qioFKlSgDuzlc/\nZMgQzJ0791Ykz4QJE4odQ6/Xo3bt2qhWrZqapt+FNRZt71xuMeqNc3D4Gvb2wIMPAh991BtDhvTW\n1jKJRCIpJS1btkRCQgLq1q2L2rVrFzpmZ2cHR0dHJCcnw9HRscS+NmzYgNGjR5fY7ujRozh69Ghp\nTba8pq1Op2sKYEdRi7Y6ne5XAEdJbrj5PgxAL5J3zfB1Oh3PniUiIoDt2wEPD+Dpp4Hx44FevSwy\nUSKRlHFu1ma1tRlmce3aNaxatQp16tRB+/bt0aNHj0LHSeK3337D+++/X2w/MTExOH36NAYPHmy0\njbHrY25NW60F/3kA40g+r9PpugP4mWR3I/0UKmKelgasXQu4uwvBX7AAsLe3yFSJRFJGuRcF35qo\nJfgWLdrqdLr1ALwBtNTpdDE6ne4dnU43WqfTjQYAkrsBROp0uggAywGMNbXv2rWBDz8EzpwBHngA\naN8e2LXLEmslEomkYmPxDF8t7pzh38nRo8CoUcAbbwDffqveRgyJRGJ75Ay/eMqMS0ctShJ8AEhK\nAp57DujWDVi8GLi5QC6RSO5xpOAXT5lw6VgbJyfgyBEgLAwYPhzQeBeyRCKRlCvuKcEHhG9/zx4h\n9q++CtyRskIikUgkRrjnBB8AqlUDNm0CMjKAyZNtbY1EIpHcG9yTgg8AVaoAW7YA+/cLf75EIpFI\niueeTo9sZydCNR9/HGjWDHBzs7VFEolEUna5p6J0jOHjA7zwAuDtDTz8sMqGSSQSzZFROsVTIaN0\njNG9OzBzJjBiBJCXZ2trJBKJRBAVFVVim7i4OGRlZVnBmnIi+ADw0UcibHPmTFtbIpFIJEBkZCR8\nfHxKbOfk5AR3d3crWFSOBF+nA1avBn7/XezKlUgkEluyfPlyDBs2rMR2lStXhpubG9asWaO5TeVG\n8AHA2Rn43/+AN98E7qhHIJFIJKqwfv16uLm5oVu3bvDy8iqyTVBQEBo1amRyn127dsXBgwfVMtEo\n5UrwAWDAAGDQIOCzz2xtiUQiKY8MGzYMtWvXxqRJk4xWqNq5c+et4iem4uTkhIiICDVMNEq5E3wA\nmDVLxOcfP25rSyQSSXmDJI4cOVKsoPv5+aFNmzZm9duxY0f4+/tbal6x3NNx+MaoXRv4+WdgzBgg\nIACoWtXWFkkkkvJCUFAQHB0dUa+e8fLcWVlZd5U93L59OypVqoTjx4+jffv22Lt3L6ZPn45WrVoB\nAOzt7REeHq6p7eVS8AFgyBBg1Spg/nxg6lRbWyORSCxB9406+dA50/JY/0OHDqFfv37FtjHckeQr\nOjoabdq0gaurK2bMmIGpU6eiTp06aNKkya021atXR67GGSHLreDrdCLlQteuwNChYieuRCK5N1FD\nqNXi4MGDGDt2LFJTU3Ho0CGcP38eX3zxRaE2lSsXltYCYU9ISECtWrVgZ2eHgQMHFmqTmpqKunXr\namp7ufThF9CsGfDpp0AJBeMlEonEKNnZ2WjQoAEiIyMRHx+P0NBQPP3006hTpw66dOlS5KzcxcUF\nGRkZt96HhYUhKCgIu3fvRs+ePQGIhd3biYuLg6urq6afpVwLPgB88glw9ixw+LCtLZFIJPciVatW\nxdixYxEcHIx58+Zh7969qFGjRrHn9OrVC76+vrfe79+/Hzt37gRJ6PV6bN26Fc7OzoXOCQwMNBr1\noxbl1qVTQNWqwI8/ijBNPz/gvnJ/i5NIJGpSqVIlfPnllwCAF1980aRzhgwZgrlz596K5JlQgptB\nr9ejdu3aqFatmmXGlkCFkL9XXwUqVwbWr7e1JRKJpDxhLOGbnZ0dHB0dkZycbFI/GzZswOjRo9U0\nrUgqhODrdMDcucD06YBeb2trJBJJeSAjIwNbtmyBv78/QkJC7jo+ceJEbN26tcR+YmJiYG9vj5Yt\nW2phZiHKRXpkU3nxReCJJ+QuXImkrCHTIxePWumRK5Tgnz8PPPmk+K/G0U8SicQMpOAXj8yHXwpa\nthSz/PnzbW2JRCKRWJ8KNcMHgKgo4NFHgfBwwMFB8+EkEokJyBl+8cgZfilp1kykXViwwNaWSCQS\niXWpcDN8ALh0CejSRc7yJZKygpzhF4+c4VtA06bAyy/LWb5EIqlYVMgZPiBn+RJJWULO8ItHzvAt\npGCW//PPtrZEIpFIrEOFneEDQEQE0KOHiNx54AGrDi2RSG5DzvCLR87wVcDVFejTB1i50taWSCSS\n8khUVFSJbeLi4pCVlWUFayq44APA55+Lxdu8PFtbIpFIyhORkZHw8fEpsZ2TkxPc3d2tYJEUfHTp\nArRoITNpSiQSdVm+fDmGDRtWYrvKlSvDzc0Na9as0dymCi/4ADBlCuDuDiiKrS2RSCRlnfXr18PN\nzQ3dunWDl5dXkW2CgoLQqFEjk/vs2rUrDh48qJaJRpGCD+CZZ0ShlD17bG2JRCIp6wwbNgy1a9fG\npEmTjFao2rlz563iJ6bi5OSEiIgINUw0ihR8iHz5U6YAs2fb2hKJRFLWIYkjR44UK+h+fn5o06aN\nWf127NgR/v7+lppXLFLwb/LKK8Dly8CpU7a2RCKRlGWCgoLg6OiIevXqGW2TlZUFna5wtOT27dux\na9cuTJ06FWvXrsXIkSMRFhZ267i9vT2uXLmimd2AFPxbVK4MjB8PLFxoa0skEsld6HTqvFTg0KFD\n6NevX7FtDAZDoffR0dFo06YN3NzccODAAbi5uWHo0KFo0qTJrTbVq1dHbm6uKjYao9wXMTeHd98F\nHnoIiIsD6te3tTUSieQWZWhT1sGDBzF27FhERETgzJkzCA4OxgsvvIDOnTvfalO5cmFpLRD2hIQE\n1KpVC3Z2dhg4cGChNqmpqaircWUmOcO/DXt7YPhwYNkyW1sikUjKCtnZ2WjQoAEiIyMRHx+P0NBQ\n9O3bFzt27EDDhg3xySefYO7cuYXOcXFxQUZGxq33YWFhCAoKwu7du9GzZ08AYmH3duLi4uDq6qrp\nZ5GCfwcTJgDLl8ti5xKJRFC1alWMHTsWwcHBmDdvHvbu3YuaNWti0qRJ6NatG2JiYtCsWbNC5/Tq\n1Qu+vr633u/fvx87d+4ESej1emzduhXOzs6FzgkMDDQa9aMWFTqXjjHc3ERitXfesbUlEknF4F7O\npfPDDz9g0qRJqFGjxq2/paSkYO7cufj+++9N6kOv12PatGmYb6T+qsyloyETJ4osmvfo908ikViJ\n7du3Y8KECYiNjS30dzs7Ozg6OiI5OdmkfjZs2IDRo0drYWIhpOAXwTPPAAYDcPSorS2RSCRlla1b\nt+K7777DkCFDsGnTpruOT5w4EVu3bi2xn5iYGNjb26Nly5ZamFkI6dIxwtKlwOHDgIeHrS2RSMo/\n97JLxxqo5dKRgm+E9HTgwQeB4GDAlJQYChXEpMYgLiMOde6vA/vq9qhbvS6qVqqqvbESyT2OFPzi\nkYJvBcaNE+UPv/nm7mMkcTL2JP4I/APeV7xx4doF2Fe3R8NaDZGWk4Yb+hu4kX0D7eu1xwDXARjg\nOgCPNXoMle+TWx8kkjuRgl88UvCtQGio8OdfuiSSqwFAVl4Wlvguwf8C/geCeLvj2+jfvD9aOLRA\nrftrFTo/z5CHE1dOYM+FPdh1YRfSc9Mx/anpeLPjm3LmL5HchhT84pGCbyX69AE+/BB49VViW9g2\nTNo3Cd0adsOk7pPQvVH3u/JlFIdXtBe++fcbhF8Lx1c9v8KoTqNwn06um0skUvCLRwq+ldi8GZj7\n2yXYjRyNK2lXsHjAYvRp1seiPr1jvDF5/2TUrFITv7/4OxrVNj1vtkRSHpGCXzwyDt9K2D1yBKce\n6Y421foicHSgxWIPAI83fhzHRx1Hrwd7ocuKLtgculkFSyUSiaR45AzfCCSxyHcRZh2fhf4Za1Ez\n4WlNcuz4xvrijb/fwPMPP4/5z86XLh5JhUTO8IunzLh0dDrdcwB+BlAJwG8kZ99xvDeAfwBE3vzT\nFpJ37TcuS4KvUMGHOz+ET6wPtg3dhqpZzdC+PRAdDTzwgPrjpehTMHjDYLg84II1L67B/ZXvV38Q\niaQMIwW/eMqES0en01UCsBjAcwDaABim0+laF9H0X5Kdbr5MSy5hIxQqGL1jNM4ln4PXO15oZt8M\nDRsCPXtqV+jcrpod9r2xDwbFgOfWPodUfao2A0kkkgqNpf6DbgAiSF4imQdgA4DBRbRTp/KAxpDE\n2F1jcS75HHYN34UHqv43nR89WmTR1Ipqlath4ysb0c6pHfqu6StFXyIpB0RFRZXYJi4uDllZWVaw\nxnLBbwgg5rb3V27+7XYI4HGdThek0+l263Q68wo9WgmSGL9nPIISgrB7xO67Yur79weSkgAtS05W\nuq8SfhnwC3o06oFBGwYhOy9bu8EkEommREZGwsfHp8R2Tk5OcHd3t4JFlle8MsXpdhpAY5JZOp1u\nAIBtAFoU1fDrr7++9e/evXujd+/eFppnOu5e7vCO8caRt46g9v217zpeqRLw/vvAihXazvR1Oh1+\nGfALRm4didc8XsPfr/2NKpWqaDegRCLRhOXLl2P27NkltqtcuTLc3NywZs0avPnmm8W2PXr0KI5a\nkNXRokVbnU7XHcDXJJ+7+f4LAMqdC7d3nBMFoAvJ63f83WaLtjvO78CYXWNw8r2TxcbEX70KtG0r\nFm9r1TLaTBXyDHl4aeNLsKtmhzUvrSk5eic7G4iJAWJjhaEZGWKFuWZNoE4doF07wMlJW6MlklJy\nLy3arl+/Hn/99ReSkpKwYMGCIouWBAUF4dixYxg/frzJ/b755ptYs2ZNkcfUWrS1dIZ/CsDDOp2u\nKYCrAIYCGHaHQfUAJJKkTqfrBnGTuX5nR7YiNDEU72x/BzuG7ShxA1SDBmLn7bp1wqevJVUqVcHm\nVzej35/98N2/32Fm75mFG+TmAkeOiBzO//4LBAUJAwteDzwAZGUJ4b9xQ2SBq10b6NwZ6NcPeO01\n4I6KOxKJpGSGDRuG7du344033jBaoWrnzp148cUXzerXyckJERER2pY5JGnRC8AAAOcBRAD44ubf\nRgMYffPfHwEIARAIwBtAdyP90NokZyaz+cLm/CPwD5PP2buX7NSJVBQNDbuNuPQ4NprfiNvObRN/\nOHuWnDyZdHYme/QgZ8wgDx0iMzOL70hRyIgIcsMGcsQIsk4dsn9/8q+/yLw87T+IRFIMtvj9lxZF\nUVivXj3Gx8cbbTN48GAqZorEH3/8wQ0bNhR5zNj1ufl3k/W6wm68IokX1r+Alg4tMe/ZeSafpyhA\n8+Yi5cKjj2po4G34xvpixo/94XGmDR4IjwLeekvUX2xR5FKIaWRlATt2iMT/sbHAV18BI0YAlWU2\nT4n1uZdcOoGBgXjjjTcQEhJitE3//v2xf//+Qn/bvn07KlWqhOPHj6N9+/bYu3cvpk+fjlatWgEA\nduzYgfDwcEyePPmu/sqKS+eeZYnfEiRkJuDvoX+bdd599wHvvgv89puVBD8wEN2+/BYe/pXwU++L\n+HRnCOxqqeCLr1EDGDpUvI4eFTmgv/sO+OUX4PnnLe9fIlERnUrl56hCIMihQ4fQr1+/YtsYDIZC\n76Ojo9GmTRu4urpixowZmDp1KurUqYMmTZrcalO9enXk5uZabF+xmPM4oOULVnykO5Nwhg6zHRie\nHF6q869cIe3tyYwMlQ27nfR0cuJEsl49cuFCMjub43aN40sbXjL7UdFk9u0jH3qIHD6cTEzUZgyJ\npAis+fu3lOeee47bt2/npUuXuGnTJs6aNYunTp0q1KZ///5FnhsfH89evXoVeczDw4O//vprkceM\nXR+Y6dKpcIlb9Pl6DN8yHLP7zcbDDg+Xqo+GDYEnngCKKGOpDrt3i6ialBSRlH/CBKBaNcztPxdR\nKVFYeXqlNuP27y8Wd11cgPbtgS1btBlHIrmHyM7ORoMGDRAZGYn4+HiEhobi6aefhpeXFxwcHPDw\nww8jPDy80DkuLi7IyMi49T4sLAxBQUHYvXs3evbsCUAs7N5OXFyctgu2qIAunakHp6KFQwu80+kd\ni/p5/31g9mxg1CiVDANEaOWECcChQ8DKlaL6ym3cX/l+rH95PZ5a/RSeavIUWjsVlcXCQmrWBObN\n+8/dc+IE8NNPRfr2b+TlwSs1FQEZGYjS6xGZnY3onBxkKwrySeQpCmpWqoT6Vaui/v33o8n996Nr\nrVp4rHZttKxRA/eZUUtAIrEVVatWxdixYxEcHAwvLy/s3bsXNWrUwPDhwxEVFYX9+/fj22+/LXRO\nr1694Ovri759+wIA9u/fj/T0dNSvXx96vR5bt25Fw4aF96gGBgbivffe0/SzVKhFW69oL7y6+VWE\njA1B3ep1LeorPx9o0gQ4cEDE5ltMWJgIlWzXTuzsKibQf6X/SizxW4KT753UNtHatWvA8OEiBHTj\nRihOTjiRlobNiYk4nJKCKL0ej9Wqha61a6N5tWpoVr06Hrz/ftSoVAlVdDpU1umQYTAgLjcXcbm5\niMzOhm96Ok6mpeFGfj4G1K2L15yc8FzduqhWqZJ2n0NS5rmXFm3vxMfHB9u3b8esWbNu/S0lJQVz\n587F99+bljpMr9dj2rRpmD9/fpHH1Vq0tbnvvuAFjX14+jw9Wy9uzU0hm1Trc9o0ctIkFTpat450\ndCRXrDAp3lNRFA7ZOIQf7/lYhcFLID+fYT/8wE+mTGGjf/9l25Mn+U1UFH1TU5lrMJS626t6PZdd\nucI+AQG0O36co86d45n0dBUNl9xLaP3714IpU6YwNDSUp0+f5pAhQ+46vmDBAiYlJZnU1+rVqxkW\nFmb0uLHrAzN9+DYX+luGaPw//OsjX/OFdS+ouuAZESF0Wq8vZQcGAzl9OtmsGRkUZNap17KuscG8\nBjx26VgpBy/BNEXh7uRkPhcURGdPT37xzz8M6dKF9PFRfaz4nBz+cOkSXby8OCAoiIevX1d9DEnZ\n5l4UfG9vb/7zzz+cNWsWQ0JC7jquKApXrFhRYj/R0dHctm1bsW3UEvwK4dI5l3QOT61+CgGjA9C4\nTmNV+376aeHPf/11M0/MyhLx9HFxwNatpUp78E/YP/j0wKcIGhOEGlVqmH1+UZDE9mvXMCMqCvfp\ndJjYsCFed3YWLpedO8Wixdq1YoFXZfQGA9YmJuKn6Gi4Vq+Oec2bo03NmqqPIyl73MsuHWtQZgqg\nqIVWgq9QQc/VPfF6u9cxrts41ftfvx5YtUr48k0mMREYOBBo2VIE9N9fej/8sC3D0LBWQ8ztP7fU\nfRRw4Pp1TIuKQp6i4PtmzeDm4HB3kXZPT+Dll8U6g5lbx00lV1GwNDYWP0RHY5izM75p2hT2VWQC\nufKMFPzikT58E/nf6f/xsZWPMd+Qr0n/2dmkgwMZGWniCZcukS1akF99pUp+hqTMJLrMdeGJmBOl\n7iMiK4svBAfT1ceHmxISaCjJrlOnSCcncs+eUo9pCkk5ORxz/jwbeXtz77Vrmo4lsS1a/f7LC8au\nD6QP/z9S9al0metCv1g/1fu+nfHjhX6XyNmzZOPG5M8/qzr+xpCNbL24NfV55i0mZOXn86vISDoc\nP85Zly5Rb84irJeXWMA4csQ8Y0vBwevX2cTbm6PDwpgu8/6US6TgF49agl+uXTpTDkxBclYyVg1e\npWq/dxIUJDw0ly6JvPlF4u8PuLkB7u5ACTmvzYUkXtr4EjrX74wZvWaYdI5nSgreOX8eHR94APOb\nN0fjatXMH/jwYRGrv2MH0L27+eebQWp+Pj6OiIBXaiq2tmuHttK3X66QLp3ikS6dEjiffJ4Osx0Y\nlx6nar/GePTRYjwcvr4iu+XWrZqNfznlMh1mO/DCtQvFtsvIz+eE8HDW9/Li32qkT9i1S6R/OH/e\n8r5M4I+4ODp6enJTQoJVxpNYB7V//+UNY9cHMrWCYPL+yfj8ic/h8oCLVcYrSKh2F76+Ymb/22+a\nLXICQJM6TfD5E59j3O5xRmdKfmlpeOTUKVzPz8eZrl3xkhoFUZ5/XiRdc3MDkpMt768E3nRxwb4O\nHTAlMhKfX7wIg5wVSiSmY87dQcsXVLzD77mwh66/uJrt07aElBSRYr7QpNnHRyxu7thhFRty83PZ\nbmm7uzaX5SsKZ126RGctZ8ZTp5KPPy5Wsa1AUk4O+wQEcMiZM8zK12ZBXmI9IMqlylcxL2PXjRXZ\nh29QDOj4a0fMenoWBrUcpIJlpvPWW0DHjsAnnwA4fRp47jkRszlwoNVs8Iz2xOser+PsR2dR+/7a\nuJqTgxHnzoEk/mzdunS+elNQFJGGARAlwe7T/uExR1HwdlgYruTkYHu7djJ0U1LhMNeHX+5cOn8G\n/wn76vZ4ocULVh+7wK3DkFDh6li+3KpiDwBPNnkS/Zv3x9dHv8bB69fxqL8/+trZ4dAjj2gn9oAQ\n+N9/By5fFsnWrMD9992Hta1bo1utWngyIAAxer1VxpVI7lXK1Qxfn69Hi0UtsOGVDXi88eMqWWY6\nJNC/2QXszOqD+392/2/Ga2XiMxLR/J+ZqPngMGxo2x597e2tN/jVq0DXruLON2CA1YadFxODxbGx\nONKxI5pWr261cSUSW1KhZ/hLfJegc/3ONhF7ANDFRMMj9RlsaPm1zcT+Rl4eRkXGo34jN7SJ/hl9\n7Oysa0CDBsDGjcDbbwMXL1pt2MmNG2Nyo0boHRiIqOxsq40rkdxLlBvBT9GnYLbXbMx6elbJjbUg\nKQno3x+6jyfi45D3cFvtA6sRnJGBrv7+aF2jBs480R8JN85hZ/jOkk9UmyefBGbMAF56CcjMtNqw\n4xo1wpQmTdA7MBAX7yHRJ4ErV4AzZ4Djx4Fdu4Bjx4CICKtePkkFoNy4dKYdmoaEjAT8b/D/VLTK\nRNLSgL59xSLt999j0CARgfmOZTVWzGJjYiLGXbiAha6uGF6vHgBgX8Q+jNszDiEfhmibN78oSJFo\nLT8f+PNPwIrFTpZfvYpZly/jWKdOeFDLdYtSQgIBAcDevYC3N+DjI+rLODoCdnZA7driK3X1qsit\n5+AA9OkjXk8/DTz4oK0/gaSsUCE3XsWlx7Hu7LqMTokudR+lJjub7NOHHD36Vm6cbdtEhKI1MCgK\nv7h4kU1PnGBAWtpdxweuG8jZnrOtY8ydZGaSbdqQq1ZZfeifY2L4sI8PE3JyrD62Ma5dI3/5hezY\nkWzaVNRS8PAgY2ONn6Mo5IULolTC66+LbBZ9+pCbNpFl6KNJbAQqYi6dj/d8zIl7Jpb6/FKTn0++\n9BL52mvi3zfJzSVdXETqHC1JzcvjwOBg9jp9mklGfv3hyeFW3XF8FyEhQqVCQ60+9IzISD7i58cU\nG+ffSUsTuZbs7YVoHzwoSiGUhpwccsMGsndv8R1zd7fa1gdJGcRcwb/nXTpX06+i3dJ2CB0bivq1\n6mtgmRFIYMwYIDJS5Im/I8Xx1KmAwQDMmaPN8BFZWRgcEoJednZY6OqKKsXEvX+y7xNk5mZi+QvL\ntTGmJFatAubPF7uOa6iTt98USGJiRAQCMjKwr0MH1LByGcX8fGDFCrERuX9/8d8mTdTrPzQU+Oor\nseXj++9FnIAW2x+upF3Bvoh9CIwPRGRKJC5ev4i4jDjooEOl+yqhyn1V0LhOY7RybIWWDi3RpX4X\n9G7aG9UulvmnAAAgAElEQVSryGgpralw+fAn7JmAyvdVxvxni64FqRkzZ4rVtSNHiqw/Gx4O9OwJ\nxMQAau8HOnLjBoadPYuZTZviwzsKIRfF9ezraLW4FY68dQRtndUowGsmJDByJFC9uijObkUUEm+e\nO4csRcHmtm1RyUprCRcuiI9co4aoCd+pk3ZjeXoCn34q9r79/jvQpo3lfV64dgG/nf4NuyN2Iy49\nDs80fwbdGnRD87rN0dy+ORrUagAAMNCAPEMeLqdexvnk8ziXfA7eMd4IjA9Er6a9MKjFILze7nXU\nut94jWZJ6alQgn8l7Qo6LOuAcx+dQ70H6mlkWREsXQosWAB4eQHOzkab9ewJfPwxMGSIekMvv3oV\nM6KisL5NG7Pi63/2+Rn7L+7H7hG71TPGHNLTheq5u6t7QUwgV1HwbHAwOj3wAOa7umo6Fim2IEyb\nJuYEH31knfVqUtxLp00TT5eTJhWTudUIChXsjdiLRb6L4H/VH6MeGYUhrYfg0QaPotJ95nV2I/sG\n9l3cB4+zHjgcdRgjO4zEuG7j8LDDw+YZZQokEB0tniBjYsRKd1wckJcnJhnVqgH29kCrVuJu2KpV\nkZO0e5EKJfgf7foINarUwJz+GvlNisLDA5g4UcTPPfRQsU3XrBEh6bt2WT5svqLgk4sXceDGDexo\n1w6uZrpGcg25aLu0LZY+vxTPNH/GcoNKg48PMHiwCFFp0MCqQ9/Iy8MTAQH4sEEDjG/USJMxMjPF\n9oOICFEFUo2ZtrlERgobAJHhwtSPeuzyMXyy7xMYaMDExybi9Xavo1pldSKcYlJjsOzUMvx2+jf0\nadYHPz79Ix6yL/63UyIZGcDff4sfl5eX8J917w40bQrUrw+4uABVqwJ6PZCdDVy/DoSFAWfPisfv\n1q3Fbng3N+DRR62SCkQLKkyUTnRKNO1/smdChhXT5B4+LJKhBQSY1Dwzk6xbl4yJsWzY67m5fCYw\nkM8GBvJGbm6p+9lydgvbL22vWfUvk5g5k3zmmdKvWlpAVFYW63t58Z+kJNX7jo0lO3cm33rLgqL2\nKmEwkLNmkfXrk8dKqHF/8fpFDtk4hE0WNOG64HVUVKjCZoyMnAx+e/Rb1p1dl5P2TuK1LDOrmCmK\nWPEePlxkKhw4kFy9moyIMK96XG6uKNzz2WciiqxJE/L778k49QIbFEVhRFYWV129yikREXwlJISd\n/fzYyNub9b286OzpyXqenuzg68sBQUF8PyyMC2Ni6JeaylwzfhuoKFE6H+36iJ/t/8yscywiIECI\n/eHDZp324Yfkd9+VftjwzEy28PHhxPBw5lkokoqi8KlVT/E3/98s6sci8vLIxx4jFy60yfC+qal0\n8vRkUHq6an0GBIhCZj/8oErVStXYu1eUYViy5G67FEXhMr9ldJjtwO///Z5ZuVlWsys+PZ6jd4xm\nvTn17srsWiSKQm7fTnbrRrZuTS5adEdaWgs5dYp87z3Szk6EUZUyoizPYODO5GS+dfYsG98U9mGh\noZx16RI3JCTQLzWVl7OzeUWvZ5xez1i9nqfT0rg9KYlLr1zhB2FhbOfry5r//sunAwL4a2wsE0uI\nva0Qgn817Srtf7JnfHq8yedYxMWLZIMG5ObNZp/q7y9irkuj1QeuXaOzpyeXFxeobSY+MT5sOK8h\nM3MzVevTbC5cEKGaISE2GX5dfDybnjhR4o/JFI4cER9l40bL7dKCCxfIdu3IMWP+ixyOS4/j82uf\nZ5flXXgu6ZzNbDsRc4ItF7Xka5tfY2KGEQE/coTs1ElsXti8Wdsnwxs3yNmzxcTuzTdNLlR9LiOD\nn1y4wHqenuzh789FMTEMy8ws9dNSSl4etyQmcmhICGsfO8b+gYHcmpjI/CL6qxCC/+m+Tzl+93iT\n21tEXBzZvDm5dGmpu+jUiTxwwPT2iqJwUUwM63l68sj166Ue1xivbHqFs47NUr1fs1i+nOzSRTxe\n24AvLl5kz9OnmWOBgBw6JMT+0CEVDdOAtDTy6afJIUPIvWGH6TLXhdMPTWdOvu13bmXlZnHyvsl0\nmevCfRH7/jsQF0e+8YZ4dNq82bqPTikp5IwZwh/7+efCN1sEgenpfCUkhE6envzi4kWGGWlnCZn5\n+VwXH89up06x+YkTXBQTw4zb9vyUe8FPzkym/U/21tlVm5JCPvII+c03FnWzZAk5dKhpbXMMBn4Q\nFsa2J0/yYpY2j9kFm7GSMtX3ZZuMopDPPmvxtS0tBkXhoOBgfhAWVqqZ2MGDYiJ49KgGxmmAXk92\nHr2EVaY5858zZe8OdSTqCOvPrc/vjn5Lw4rl4k46ZQqpouvNbOLihIvnoYfIff/djM5nZnJQcDDr\ne3lxXnQ0062wsU9RFHqmpPClM2dY38uLv8bGMs9gKP+C/9Xhr/jeP++ZfcHMJjub7NWLHDfO4tnF\njRtijamktcKEnBw+dfo0XwgOZqrGX6KPdn1km93JtxMTI1Tz9GmbDJ+Wl8e2J0+a7TI7dEiYXdKC\naFkhNz+XY3aMYevFrfnG+Ag+8giZnGxrq+7m6uVQHu7iwItNajHF97itzfmP3bvJpk2ZMWoUp549\nS4fjx+l++bLNKq35paayT0AAW508Wb4FPyU7hQ6zHRhxLaJUF8pkcnPJwYPJYcNU8xmOHEnOn2/8\n+Km0NDbx9uZXkZE0WOHxNSEjgXVn1+XF6xc1H6tY1qwRTmYbhbaEZ2bSydOTJ1JSTGp/+vS9NbPP\nzM3ks38+y+fXPs9UfSoVRXgpOnYsY6J/7BjZpAnzx4/j5H8+YotFLbT/nZvBjitX2HjnTr4xaxZj\nT5ywtTlUFIW7k5PLt+DPOjaLI7aMKNUFMhmDgRwxgnz+eVWzU/37rwgwKErL/4qPp6OnJz3UjDww\ngW+OfsOhm030NWmFooib6xdf2MyEf5KS2Mjbm/El/P+OjBRr91u2WMkwC0nJTuGTq57km1vfZJ7h\nvyfGMif6ixeT9eoVqv281HcpXea60Dva24aGkel5eXw/LIxNT5wQ62lbtojQp59+sklo8Z2UW8HP\nys1ivTn1GJKgYWSHoohwhl69SJX954pCtmxJenr+97ccg4Hjw8PZ/MQJBtvAV5mek876c+vzVOwp\nq49diLg48SM6ZTs7voqMZK/Tp43GQCcmkg8/LLTpXiAxI5Gdfu3EcbvG0aDc/ZnKhOjn5oq45TZt\nRCTcHewK30Undyd6hHrYwDjyREoKm584wVHnzhV2sV6+LNLhvviiWBG3IeVW8Jf4LuGg9YNKfWFK\nRFHEIlHXrmRqqiZDzJsnor1IMlav5xP+/nwhONiizVSWstR3Kfut6Wez8W/x559k+/Y2y/mbrygc\nEBTETy5cuOtYdjbZvTs5bZoNDCsFCRkJbL24Nacfml7sgrSikJ9+KsLbrT7fuH6d7NtXPEkX83s7\nffU0G8xrwJX+K61mWkGUnLOnJ7cYe+rOySE/+IBs21bEvtqIcin4eYY8Nvu5Gb2ivSy6OMUyY4bw\nJWs43UlKEou32y5fZwMvL35/6ZJV/PXFkZufS9dfXLk/Yr9N7aCiiJ2TM2bYzIRrubl80Nubf9/2\nI1cUsXv21VfL1qYqY1zLusYOyzrwq8NfmdReUch33xVhm1ZbRomNFb+1iRMLpRU3RnhyOB9c8CDn\nexezCKYSmfn5HHn2LDv4+jLClKf8ZcvE0+nBg5rbVhTlUvDXn1nPp1Y9ZdGFKZZvvhGPlQnapmnI\nVxS2mxPF2vu9eOCamdvKNWRjyEZ2Xt65yEd/q3LlilgRDQy0mQknb+7ELfixz58v3B4ZGTYzyWRS\n9ansuqIrP9n7iVmhpvn55Msvi5fmgScXLpDNmoncD2bYGJ0SzRaLWvDrI19rlv4hJjubnfz8ODw0\nlJnmXIijR4Xor1mjiV3FUe4EX1EUdlzWkTvP77T44hTJ99+TrVqR8dru2o3T69k3IICPHA5gqyf0\nZWq2aFAMfHTFo1x/Zr2tTRHVsTp1stmGLJJcfOUKH/Hz4479+XRxIS9dspkpJpOZm8mnVj3FMTvG\nlEoQ9XqyXz+RYUCz72ZQkFj1Xr68VKfHp8ezw7IOnHpgquqiH5CWxkbe3nS/fLl0fZ89K7bUWzm/\nRrkT/L0X9rLd0nba3NW//55s0YK8elX9vm/jn6Qkunh5cUZkJPMMCl1dyTIQ2VWIQ5GH2Hxhc9vv\nvizYkPXjjzY0QeHAk6Gs9mXYPRF+mW/I56D1gzhiywiLntLS00UCuG+/VdG4AgICRCSOhTkokjKT\n2H5p+xLXJ8xh77VrdPL05CZLn/CvXhWTldGjrfCoJCh3gt/7995cE6jyo5KiiDDANm2KLyhqIRn5\n+Rx9M6TL87Y4b3d3ctQozYYtNf3/7M/FJ8tAGMqlS2KnZViYTYbX68lHHs+j8x4f/qliBkUtUBSF\nH+78kE//8bQqN+u4ODFR/f13FYwrIDBQiL2HOtE2iRmJbLe0HWcctny95/e4ONbz9Cz0+7SItDSx\nGD10qFUCEMqV4J+8cpJNFjRhbr6Kj/cGAzl+vJjKaJAmtwDPlBQ+7OPDkWfP3lVTNSFBJOa7cUOz\n4UvF6aun6TLXhek5NtzOXsCiReQTT9gk1vmjj0TemcC0dDp6evJcGXbg/3T8J3ZY1oEp2SoJFoV3\nwtnZvPxPRgkKEsV3N5mQFdMMEjIS2GZJG35ztPSpORbGxLCxt7f6/3+zs8lBg0QEkkbpUQooV4L/\n6qZXueDEAlUuDEmRmnfUKBFDq9Yd/Q4y8/P58YULrO/lZTykiyJFh40yBBfL8C3DLfoRqYbBIATf\nyoHvGzeK1CkFN+MVsbFs5+tr3iKelVgXvI6N5zfmldQrqvd97JhYPz9zxoJOzp4VYq9RKtH49Hi2\nWNTCbI1QFIXfRUXR1ceHl7SqAJ+bKzZw9uypWZg3WY4E/+L1i3SY7cA0vUobG9LSyOeeIwcM0Czk\n4sC1a3T18eHw0FAml7DoWNzOW1tScN2Npqu1JufOCdeOlVZNw8OFyPn7//c3RVE4PDSU79nIvWQM\nnxgfOro7Mjg+WLMx/vpLuHdKFc9w+bLIdKmqb6iIYVIu88EFD5pc40FRFE6JiGA7X1/GaR2HajAI\nf36PHppNMMuN4H+06yNOPTBVnasSGyuyXr7/vpjlq0ysXs/XQ0P5oLc3t5voJlIUsYRQFhcFx+8e\nzwm7J9jaDMGsWWIRV+M7Y06OyNZc1ANFWl4eW/j48C+NI7lMJSY1hg3nNeQ/Yf9oPtaMGWLTmVme\niYQEEQyxQMWn82IITw5ng3kNuDGk+CcJRVE4+cIFdvbzK3FCphqKInyE3bpp4sMtF4KfnJlMu5/s\neDVNheiZwEDywQc1CZfSGwycGx1Nh+PHOe3iRbMf+3/5xfS0ydakILFa5HXTCkBoSm6uuFlrHOM8\ndSr5wgvGvyKB6cKff16DnOfmkJmbyS7Lu1itnoGiiO/o66+b+PNJTRXrY19+qblttxMUH0TnOc7c\nc2FPkccVReGnERHs5OfHa9YO+VUUscns0UfFDmMVKReC/+3RbzlqmwphLH/+KVwC69WNL1cUhRsT\nEtjsxAm6BQWVetEnJUUs3paRiWMhvj7yNYdvGW5rMwT+/mIVUaMLdeSIqP9aUlTesitX2NHXl9k2\n8ucrisLXPV7n8C3DNa09eydZWWKWP3NmCQ1zcsSW3TFjbOKr9Ir2opO7E31ifAr9vcCN84gtxP4/\nI8hPPhE3QxVF/54X/IIkaaGJpasrSVJ88caPF5WqgtXzcSqKwn3XrrG7vz8f8fPjQRX+x733nnj4\nKGsUJFbzv+pfcmNrMGUK+dprqnd7/bpwNe8pemJYCEVR+EpICMeeP6+6HaYwx2sOuyzvYtX6swXE\nxYla30aDbRRF5AAfPNhqMehFsSt8F+vNqceziWdv/e2ryEh28PW1nhvHGIpCfvyxyNelkk//nhf8\n5aeW022tW+mvQHi4KJLt5qaaz0xRFO69KfStTp7k2vh41XLgnDolfkhlMAik7CRWI8U009WV3LZN\ntS4VRdxDJpixXHEjN5fNTpzgZo3TcNzJ4cjDrDenHi/duGTVcW/n9GnxwOxf1Bxg+nTxu7Oxy4sk\n1wSuYeP5jRmdEk33y5fZ6uRJJtgoKd9dKIooqtS9uyrRO/e04BsUA1ssasHDkYfN/+SKIlbcHByE\nc1yF+G29wcDf4+L4iJ8fW588yfXx8UUWEraUrl0LpQIvM+Tm57LFohaFa43akqNHyYYNVbuRr10r\nFs7NDZX2vZlvR6sSlHcSnRJNl7kuPHBRjcB4y/DwEE9EhTanL18ubsZWrudQHPO859Fl7Qd80NuL\nMVqFXpaWgjTsTzxhcZrSe1rwt4dtZ+flnc33T54/Tz7zjJhhqPC4HZmVxRmRkXTx8mL/wEDuTk7W\nNKvlH3+Q/ftr1r1FeIR6sOOyjrZPrFbAmDEivaOFxMaKEMzSpuBfEB3NrqdOWVQE3RT0eXp2W9mN\nPx3/SdNxzOHbb0XQSXY2yb17xS5aG6YILoo1cXGsdXg3u/wx0CYusBIxGMT3uE+fUm/OWr/+Hhf8\nnqt7cl3wOtM/cUoKOXmyeM6cO9eikMuUvDz+ERfH3gEBdPT05LjwcIZYaYelXi/WJMtYqDdJ4c7q\n/lt39dNblJbUVDHFtGAbqKKI7RglLkIW24cogj5JY6H7cOeHHLJxiFUXaUuiIHLn84EhVMpgcd/t\nSUms5+nJM+lpHLFlBAetH1So4leZIT9fbM569lmzc1Nv2CD2tFld8AE8ByAMwAUAnxtp88vN40EA\nOhlpw8bzG5uWRiEtTeStrVdPrHqWMnrjql7PlbGxHBAUxFrHjnFgcDA3JyRQb4Pt/NOnC9deWeTY\npWNssqAJs/PKyKPxzaLSpX0cXrlSnYSc13Jz2cTbm/9olKJjbfBauv7iqmraBLXIjEpgTNVm3DH0\nT1ubUohjN27Q0dOTJ2/6x3Pyc/jMmmf4/vb3y9RN8xZ5eeQrr4iYYBPXGTZtEmIfHGxlwQdQCUAE\ngKYAqgAIBND6jjbPA9h989+PAfAx0hfnec8r/pPGxoqAaQcHMcUICDDpAhUQp9dza2IiJ4SHs83J\nk7Q7fpyvhoRwQ0JC4RJmNiAmhrS313QXtkW8uOHFMuVW4MiRIrbZTKKixAOhRSkDbsMrJYXOnp68\nrLKf+GziWTq6OzIgzrzvuFXIziZ79GDKhC/ZoEHZWX8KTE+ns6cn999RayJNn8Yuy7vw6yNf28iy\nEsjNFYL/yisleim2bBHz3IKSEdYW/B4A9t72fiqAqXe0+RXA0NvehwGoV0RfTNUXoXZxcWIxtlcv\nUS5q3Lgi61/ejt5gYEhGBjclJPDrqCgOOXOGjb29aX/8OPsHBnLWpUv0TU3VZAHWEl57Taw3l0XO\nJ58vOykXSPLaNRE8f3uR4BJQFBEmrnbmZffLl9nd3181f35GTgbbLmlr1bJ+JqMowg3x6qukwUBv\nb3EDDdGw1LQpXMzKYgMvL240Ej0Vnx7P5gubc/mp0uXi1xy9Xrh2RowwGrK3fbtw/Z4+/d/fzBV8\nnTindOh0ulcAPEvy/Zvv3wDwGMnxt7XZAeBHkt433x+86frxv6MvXjl2DExNBaOioJw5A0NoKAzX\nryP36aeR8+yzyOneHZmVKiHdYECGwYCU/Hwk5eUhOS8Pibm5iMnJQXRODq7n5aFZtWpoXbMm2tSo\ngbY1a6JbrVpoXr06dDpdqT+v1nh6Au++C5w7B9x3n62tuZsJeyZAoYLFzy+2tSmCrVuBzz8HgoKA\n6tVLbL5iBfDbb4C3N1C5snpmKCQGh4SgRfXqmOfqanF/b297GwoV/PHiH2Xv+/rTT4CHB3DsGFCj\nBgBgzRrgm2+AkycBR0frm5SYm4snAgIwqVEjjG3Y0Gi7iOsR6Lm6J5a5LcPgVoOtaKGJZGUBbm5A\n8+biy3qbCOzbB4wcCezaBXTt+t8pOp0OJE3/kphzd7jzBeBlACtve/8GgEV3tNkB4Inb3h8E0LmI\nvlj9lVdY7ZVXeP8rr7DaV1+xyrp1rLtrF9ucPMlH/PzY7dQp9g0I4KDgYLbfsIGYOJF46y3ixReJ\n3r2JNm046Ycfipy5z5w5kwDues00snJnq/YuLleL3ARUFuxPykyio7sjw5LCyoQ9JMWe/8mTS2x/\n+XJhV47a9hTUw912059f6v4fAfERiCpl7/s5CGBqrVqiFOUdPP74cQJHCFSx6vfhi2+/5aOnTvHL\nyEiT2r8/8306uTvx+OXjNr+eRbZPTycff5wnu3a9rV0fAgkEevCtt97izJkzb71gZZdOdxR26XyB\nOxZuIVw6r9/23qhLRyIq/A0YYGsrjOPu6c5B6wfZ2oz/SEoSK1hexgvcFxTR+u47bU05kZJCJ09P\nRpYyzC40MZSO7o4MSbCxf6QogoLEHfPkySIP5+cLN/T771svq0KOwcD+gYF8LyzMrAXZfRH76DzH\nWdNMoxaRkiI253z8MY8fU+jkZDzJorUFvzKAixCLtlVR8qJtdxSzaCsR62EuLmSoBZkltCQ7L5tN\nf25aus1xWuHhQbZsaTSe2ZplcudHR7OLn5/Z+XYyczPZdklb/u/0/zSyzAISE0VU1LriQ6bT0sh2\n7ciff9beJIOi8I2zZ/lCcDDzSrF2si54HRvNb2TTncvFcv06M1p24sLqn3P/PuM3M6sKvhgPAwCc\nh4jW+eLm30YDGH1bm8U3jwehCHcOpeAX4ttvVdlbpBmbQzezw7IOzDeUoXwQQ4eK5FR3ULDBqiCq\nQWsK8u2MNnNTxTvb3uEbf79R9kIHc3LIp54ip00zqXlUlFhL37lTW7M+i4hgD39/iwrTLPRZyBaL\nWpSdQITbCAggWzomM7Vpe/Krr4y2s7rgq/WSgv8fSUkiRLMsZtEkhaj1Wt2Lv/r9amtT/iMpSSjN\nbc++iiJyeRXze9GE1Jv58/8wsR7uX0F/seWilmWjtOTtKIrw0QwebFaqkhMnxE02KEgbsxZER7PV\nyZOqJEObdnAauyzvUnSEoI0ICRFP+R4eFE9Xbdsa3SUoBb+cMGaM9YXKHALiAlhvTj1ez1I3v7dF\n7NhBNmsmfAsUuxHbtDF7E6MqnLmZPz+4hM1h55PP09HdkUHxGqmjJSxaJHw0aeZXnduwQSQFvKpC\nSYtC/SYksJG3t2qlCRVF4ZgdY9hrda8ykYIhLIxs0EBUG7tFQoIQ/Rkz7logkYJfTjh/XsySykDy\nQaOM3jGaH+/52NZmFObdd8n33ru1luvjU/IpWvFnXBxdfXx4w8hMNDsvmx2XdeQyv2VWtswEDhwQ\nO3zuiH4xh+++EzU/1MpQcvD6dTp7ejLIwoRjd2JQDBzmMYxua91M2+mvERERZKNGYs3pLhISxM33\nyy8Lib4U/HLEoEHksjKoBQUkZiTS0d2xUO5xm5OaSjZtyjm9dhTl0rc6E8LD+XxQUJHJ98buHMtX\nN71a9vz24eFih4+F9TcVhXz7bXLgQMsri/qnpdHJ05P/alAmkBSZYd3WunGYxzCbrE1FRYnCfL8W\n5yVNTBSiP3XqLdGXgl+O+Pdf8uGHVcn0rBkLTixg/z/7lynR8v7xKOMr1WdmlHVz1hdFrsHAXqdP\n3xUnvjl0Mx9a+FDZy5Nz44aIeFqxQpXucnNFJtjRo0sfrhmRlcX6Xl7conH65azcLPb5vQ9HbRtl\n1eywly4JT6RJu+yTkkTI2YQJpKJIwS9PKIpIQ+vhYWtLjJObn8u2S9pyc+hmW5tCUoQwN2pEXho2\nVRTBKQM3ooScHDbx9ubfNwUr4loEHd0d6XvF18aW3UFentiwYE5FGBNISxNliUtT2S0+J4fNT5zg\nsiI2e2lBRk4Gn1r1FD/Y/oFVJjGXLwuxNyuU9cYNskcP8r33pOCXN7ZtE2Uwy4BuGeXYpWNsNL8R\n0/TmL+6pzQcfiBdzcoQDedEiW5tEkvS7WTTFP+U6Oy/vzIU+C21t0t1MmCDqSmiQSPDqVeGyWL3a\n9HNu5Oayo68vv46KUt2e4kjTp7HHbz04btc4TUX/8mXyoYfIBQtKcXJ6OtmnjxT88obBINx2ptRc\ntSVvbn2Tn+771KY2HD4sZve3yoWGh6ubGtNC/oyLY+1DO/n8RusWITeJJUvIVq1UqyZWFGFhYiF9\n69aS22bm5/PJ06c5PjzcJtcqJTuFXVd05fjd4zVx71y6JMR+/nwLOsnKkoJfHlm7lnzySVtbUTzx\n6fF0cnfimQTbiGtGhvgB3ZWqd/Vqcce0UjnC4vAI9WCdDZPZ/ZSvTWouGGXfPhGRExGh+VD+/iL6\n7OBB421yDQa6BQVxRGioppXmSiIlO4U9fuvB97e/r6roX7ggnnYWqvCQJwW/HJKXRzZvLhZxyzJL\nfJfwqVVP2WRGNmkSOXx4EQcKyjONHm11m27nwrULdHJ3ok/MSb505gzfPneubMzyQ0OFAluxatW/\n/4ohi0rLk68oHB4aSregIOaWgZtiek46e//emyP/HqlK1axz58RT6HKVsjRLwS+nrFxZduveFpBv\nyGfXFV2tnsfdy0u4CowWnkpNJVu0uGM3i/XIys1ix2UdufjkYpJkRn4+H/Hz46xLl2xizy3i48WK\noTmOdZXYuVNEft5eU1hRFL4XFsbeAQHMsiBlgtpk5may/5/9+fLGly2q+hYUJDaD//GHerZJwS+n\n5OSImYFvGQvsuJOg+CA6ujsyNi3WKuNlZQkt31xSkFBBtkcbZKV79593+brH64Vm9Ff0ejbx9uYa\nE9MvqE5GhsjIOGOGbcanCEgoEH1FUTghPJzd/f2ZZuPqc0Whz9Pztc2vsdfqXryRbf46R8FTzcaN\n6tolBb8c88svItKwrPPV4a84eP1gq7gsPvtMFF8yif/9T+RasFJxepJcdXoVWy9uXWSenNCMjCJL\n8mlOfr7Y1ffmmzYP/9q2jXRyVjjK6yI7+fkZ3ZVcFjAoBn6852O2W9qOMakxJp+3bZsQ+wMH1LdJ\nCjrE/+gAABTmSURBVH45Rq8nGzcmvb1tbUnx6PP0bL24NTeFbNJ0HB8fsdZopKrd3SgK+dZbwtlv\nBaELjAuko7sjQxONP1Ucu3GDTp6eDChFvppSoSiiTGjfviYXzdbWHIWv7rnISmtOcqen7e0pCUVR\nOMdrDhvPb8zTV0+X2H7lSuFu9PPTxh4p+OWclSvJPn1sbUXJeEd702WuC5MzkzXpPztbTNbXrzfz\nxMxMsbHB3V0TuwpIzkxms5+bcf2Zkg3cnJDABl5eDLdG4qQffxRRSxqGX5qKoiicdvEi2/v6cv2e\nHDo5kbt22doq09gcupmO7o5cF1x0jQCDgfz8cxE5dv68dnZIwS/n5OaSrq7Fh7WVFSbumcihm4dq\n4tr57DPy5ZdLOVGPjharZxqpS54hj/3W9ONn+z8z+ZyVsbFsomIWyKIHWSkKmVhp12pxKIrCqTfF\nPvHmk8aJE+KJbc0aGxtnIoFxgWz2czNO2T+lUP6djAzypZdEGQGjgQQqIQW/ArB2LfnYYzZ3v5ZI\nZm4mWy9uzb+C1I2OOX5cPCZblFrF01M4Vs+dU82uAqbsn8J+a/qZHca3MCaGzU+c4BUt8jlv2SJu\ncuHh6vdtJgZF4ccXLrCDry+T7nArnT0r0irPmlX2v9+keJLr+0df9v69N2NSY3j5sniAfOst66Tl\nloJfASjYfbt9u60tKRn/q/50cndSrZRcerrYk2DKbs0S+e03kZ0uWT2304YzG9j056aldmX9dPky\nW508yXg1/euHDokIJX9/9fosJXkGA98+d449/P153cgCbUyMyIrx+utlOz14AfmGfH7373e0+96Z\ndXp4cM4c692szBX8+yC557jvPuC774Dp0wGDwdbWFE/n+p0xucdkvLXtLRgUy439/HPgiSeAF19U\nwbh33wUGDwYGDQKysy3uzi/WD+P2jMPWoVvhUMOhVH183qQJhjs7o2dAAGL0eottwvHjwNChwObN\nQOfOlvdnAXqDAa+ePYu4nBwc6NgR9lWqFNmuUSPg2DGgcmXgySeB6GgrG2ourIScA1+iypZ/UHPw\nFJx9+B3c0F+3tVVFY87dQcsX5AzfLBRFpFtQKYutpuQb8tlzdU/+dPwni/rZt0/sRVB1vdFgIEeM\nEGX8LIj/jk6JZsN5Dbnt3DZVzJoXHc2mJ07wgiVTXC8v7eIBzeRabi57BwTwtZAQ5pi4g1ZRyHnz\nhF9flSc6DYiIIHv2FIEU8fFkqj6VY3eOZb059bg6YLXmocmQLp2Kw6lTwpedUsZSqhfFpRuX6DzH\nmZ6XPUt1fmKiKP126JDKhpEiPLFfv1InbU/PSecjvz7COV5zVDVreWwsG3p58UxpKjydPCnEvgxk\n3QvPzGQLHx9OvnCB+aW4vl5eItrlvfeES68sYDCQixeTDg7ipnTnxmC/WD8+uuJRPrnqSfrEaFd2\nTQp+BWPUKBGxci+wK3wXG85ryLh083aXKorYJ6Tp50xNFYUlpk83S/TzDfkctH4Q3/3nXU1mc+vi\n4+lk7uas48eF2N+VSc76HLtxg/U8PflrrGU7r9PSxHe9eXPbR6gFB5O9epHdu4sMoMbIN+RzxakV\nbDy/MZ/76zmeiDmhui1S8CsYV6+SdeuKDHz3AjMOz2Cv1b3MimD59VcR+aD5PqGCEnJ31A01hqIo\nHL1jNPv+0Zc5+doZ9+9N0VxqSjjl7t1C7Pfv18weU1AUhYtiYujk6cl9Ku4k3rZNRJYOGWJRud1S\nkZhIjhkjLu8vv9w9qzeGPk/PX/1+ZZMFTdhrdS+uCVzDzFx1VqOl4FdAZs0iX3zR1laYRr4hn8/+\n+azJufPPnhUBJhpETxaNGaI/4/AMdl7e2SqFXy5kZrKljw8nhIcbzyK5YYNITmPjrdhpeXkcGhLC\njr6+lq1BGCErSxRIr1uXnDJFTHq0JCGB/Oor8T38+GPy+vXS9ZOTn8PNoZs54K8BtP/Jnu/+8y49\nQj1KFdGVZ8ij/1V/swVfRyG2Nken07Gs2HKvodcDrVsDK1YAzzxja2tK5lrWNTy68lF83+d7jOgw\nwmi77Gyge3fgo4+ADz6wooFJSUDfvsALLwA//ADodHc1Weq3FAt8FsDrHS8413S2ilk38vIw4tw5\n3MjPx7rWrdGsenVxgAQWLgTmzAH27AE6dLCKPUURlJGBoaGheMrODr+4uqJ6pUqajXXlCjB7NrB2\nLfDSS8Cnn4rfgVqcPSsu66ZNItBp8mTg4YfV6Ts2LRYbQzfiYORBeEZ7onnd5ujk0gkP2T+E5vbN\n0aBWA1S+rzIq3SeuX3xGPKJToxGTGoOghCD4XPFB4zqNcfajsyB59xfUCFLwywl79ghhPHMGqFnT\n1taUTEhiCPr+0RcbX9mIPs36FNnmvfeArCzxgy5Cc7UlKQlwcwPatBF30qpVbx1af2Y9PjvwGY6P\nOo5m9s2sapZCYsGVK5gdHY0lDz+MV+3sgHHjAB8fYMcO4MEHrWpPAbmKgh+jo7EkNhbzmjfHSBcX\nq4197RqwdCmwZIkI6Xz5ZfFq0cK8fkggLAzw8BBRrMnJ4js4bhzgrOE9Pc+QB7+rfghNDEXkjUhE\npkTiavpVGBQDFCogiHo166Fx7cZoUqcJWju1xhONn4BDDQfodDop+BWVESOA+vWBuXNtbYlpHIk6\ngqEeQ3H4rcNo59yu0LE//gB+/BHw8wNq1bKRgZmZ4qKmpwNbtgB2dlh3Zh0+3f8p9o/cf5fN1uRU\nWhpeDwnBY6dOYZ63N1xWrrTZhQpIT8eosDA0vP9+LG/RAo2qVbOJHfn5YtuBhwfw999AlSpAp05i\n+0GbNoCdHVCnjrhMGRniRpGcDFy4IL5nfn5AtWriaeG114DHHxd7XsoyUvArMElJQPv2wM6dwKOP\n2toa01gbvBbTDk/DiXdPoEGtBgCAkBCgTx/gyBGgne00VWAwAJ98Ahw8iO0/vYMx4fNwYOQBtHVu\na1u7PD2R+c47+Pbzz7GqZUvMbNoUHzZsiEpWfBSKy8nBjEuX8E9yMuY0b44369WDzuqPYkWjKEBU\nFHD6NBAQAJw/D6SkAKmp4v79wAOAoyPg4AA0awZ07SpeDRva2nLzMFfwbb5YW/CCXLRVhT//JDt0\nEEnW7hV+PP4jWy9uzbj0OKamki1bkr//bmurCnNi+ttMrnkfryz+0baG5OWRX38tdiPdzK0RkpHB\nnqdPs4OvLz0SEzWvA3sjN5dfR0Wx7vHjnBIRUaZz2Jd3IBdtKzYkMGCA2JL+5Ze2tsZ0vvv3O6w9\nsxYNDx5Ci/oNsWyZrS0SkMQc7zlY5LsIRx9ZiOZjp4up4JIl1nehnDsHvP++8DusWQM0aFDIzu3X\nruGHy5eRaTBg2oMP4hUnJ9yvok8iPCsLv1y5gnWJiRjo4IBvmzZF04KFY4lNkC4dCaKjhSZt2wb0\n6GFra0yn57TZOGX4DSFTDuMhh8a2Ngf5Sj4m7JkAz2hP7P5/e/ceHFV9BXD8e5pAKiSI4RGCpCJa\nHmIRgRJG6oAIrWMVGZEqvmiZqWhtDXWmCGVU8IEMHdROMk4By0MNSovW1yAQKVEcnCAaHoIZTC0G\ngklMwSQbmgSS0z9+twzoJpDsJjfZez4zmbl3c/fuucnO2d/+7u93fndsoF+3fq5fPyMDNm92N0um\nTWv9O8pVVa540qpV8Mgj7u58I4lcVck5dozFRUXsDoWY2quXq83TvTvfa0GcRTU1vFFezj/Ky/m0\nupp7UlP5zYUX0jchIdKrMlFgCd8ALtnPnu36Ly+4wO9ozm71anjiCZjxl6dZtS+Lt6a/5Ws/eVVt\nFbe/djs1J2tYP20953///DMPeO89eOABSE52Y/daYyhkbS28+CIsXOiGiS5ZAikp5/z0opoaXikr\nI7u0lKLaWtKTkhjTrRujkpJITUigd6dO9OrUiQagqr6eqpMnOVJXx+5QiN2hEDurqjhcW8sNPXow\npWdPfpac3KrDLE3zWcI3p2RkwKFDboBJO7mXFtYHH8DNN7scOmQIvLTnJR7c9CDLb1zOlMHRKIvZ\nPPlf5XPr+lu5pv81ZF2fRae48FUdOXnSDdlcsADS090n7IQJkf+xKyth2TJ49ln3QfLww27ISARK\n6+rIq6zkw8pKdoVClNbVUVZXR9mJE8SJkBQXR1JcHCmdOzOsa1euSExkeGIiP05KIr69D1UJMLtp\na06pqXElCbKy/I6kcbt2ucmhmzad+fiOwzs07ek0fXTro1rfcG7VFSPV0NCgmXmZ2mtJr3NamvCU\n6mpXtnToUDdL96mnVPPzm1eIraJCde1at4xXt26q06e7P44xTcBu2prTFRa6xuGrr8LVV/sdzZkK\nC2HcOHjmGTfu+dtKQiVM+/s04iSOZTcsY1DPQa0Xy9FCMjZmUBIqYd0t67g0+dLmn0TVjSV9/XU3\nEy4Uchc4YIAb+5eW5lr/tbXu59Ah2LPHzZb7/HP3D5o61dXo79GyevomWKxLx3xHTg7cdRfk5sLg\nwX5H4xQXu/w2d27TZRPqG+rJ2pHF4+8/TkZ6Bg/95CE6x3Vu/AnNVF1XzaJti1j28TLmjJ3D7DGz\no3f+wkI3A/bgQTcovKjIJfyEBPeTmuq6bIYNcxMOOsIUadOuWMI3Ya1eDY89Btu3QxvOeg+ruNjV\n/Ln7bpfwz0VRRRH3b7ifgvICt4LWFTM4r1PLhwQe/e9Rnv/keTJ3ZDLuonEsmbTk1MQvYzoKS/im\nUQsXunIrublupqEfCgtdsr/vPpgzp3nPVVXe//J9ln64lLziPO4deS+3XHYLl/e+/JxmeNbV15F3\nOI/svdms27eOGwfeSEZ6BiP7jmzh1RjjL0v4plGqMGuWqwL49tuutkhb2rPHTQpbsMDNH4pEQXkB\nWTuyeKfwHY6fOM7EARMZ0WcEKYkp9O7am24J3Sg/Xk5pqJQjVUfYfng7277cxsAeA5k8aDL3jLyH\nPok+f9UxJkKW8E2TGhpcaZitW2HjRteN3BY2b3b3ETIzw9+gjcQXx74g51857P96P2XHyygNlVJZ\nW0nPLj1JSUyhT9c+jOo7igkXT2jx4uLGtEeW8M1ZqcKiRbBypUvEl1zSeq/V0OAmiS5fDmvXukEr\nxpjoaG7Cj2/NYEz7JALz57tqgWPHurIwU6dG/3XKy+HOO91CJjt3tt23CWNMeDaFLsBmzYI33oB5\n89yImYqK6Jy3vt5NFB06FIYPhy1bLNkb0x5Ywg+49HRXb6drVzccfPVqqKtr+fm2bXOF27KzYdMm\nWLwY4u17pDHtgvXhm1Nyc10BswMH3PqdM2eeWwXgigrXP79ihVtk4skn4bbb2nf9HmNigd20NRH7\n6CO3OPTGja5bZtw4t5h4UpJb2jU+3k0ezc+HXbtc//ykSW7G7LXXtv9l4YyJFZbwTdTU1EBenqti\n+fHH7uZrba3r8klLc+uFDh8Oo0db6Rdj/GAJ3xhjAqK5Cd++fBtjTEBYwjfGmICwhG+MMQFhCd8Y\nYwLCEr4xxgSEJXxjjAmIFk96F5FkYB1wEXAQ+IWqfhPmuINAJVAPnFDV0S19TWOMMS0XSQt/LpCj\nqgOBLd5+OAqMV9UrLdkbY4x/Ikn4k4E13vYaYEoTx1pVFWOM8VkkCT9FVUu97VIgpZHjFHhXRHaK\nSIQL2xljjGmpJvvwRSQHCLfw5/zTd1RVRaSxughjVfUrEekF5IhIgapuC3fgggULTm2PHz+e8ePH\nNxWeMcYESm5uLrm5uS1+fotr6YhIAa5vvkREUoGtqjr4LM95FAip6tIwv7NaOsYY0wxtWUvnTWCG\ntz0DeD1MMF1EJMnb7gr8FNgbwWsaY4xpoUha+MnA34AfcNqwTBHpC6xQ1Z+LyADgNe8p8UC2qj7V\nyPmshW+MMc1g5ZGNMSYgrDyyMcaYsCzhG2NMQFjCN8aYgLCEb4wxAWEJvw1EMlGiI7Dr69js+oLD\nEn4biPU3nF1fx2bXFxyW8I0xJiAs4RtjTEC0q4lXfsdgjDEdTYecaWuMMaZ1WZeOMcYEhCV8Y4wJ\nCN8TvohcJyIFIvK5iDzkdzzRJCJpIrJVRPaJyKci8oDfMbUGEYkTkXwRecvvWKJJRLqLyHoR+UxE\n9ovIGL9jiiYRmee9N/eKyFoRSfA7pkiIyEoRKRWRvac9liwiOSJyQEQ2i0h3P2OMRCPX9yfv/blb\nRF4TkfObOoevCV9E4oAs4DrgMmC6iAzxM6YoOwH8XlWHAmOA+2Ps+v4vA9iPW84ylvwZ2KCqQ4Bh\nwGc+xxM1ItIf+DUwQlV/BMQBt/kZUxSswuWS080FclR1ILDF2++owl3fZmCoql4BHADmNXUCv1v4\no4FCVT2oqieAV4CbfI4palS1RFV3edshXMLo629U0SUi/YDrgeeJocXqvZbS1aq6EkBVT6pqhc9h\nRVMlrkHSRUTigS5Asb8hRcZbOvXYtx6eDKzxttcAU9o0qCgKd32qmqOqDd5uHtCvqXP4nfAvBA6d\ntn/YeyzmeC2qK3H/lFjyDPAHoOFsB3YwFwNfi8gqEflERFaISBe/g4oWVT0KLAWKgCPAN6r6rr9R\ntYoUVS31tkuBFD+DaWUzgQ1NHeB3wo+1LoCwRCQRWA9keC39mCAiNwBlqppPDLXuPfHACOA5VR0B\nVNOxuwPOICKXALOB/rhvnYkicoevQbUyb4WlmMw5IjIfqFPVtU0d53fCLwbSTttPw7XyY4aIdAJe\nBV5S1e+s+9vBXQVMFpF/Ay8DE0TkBZ9jipbDwGFV/cjbX4/7AIgVo4DtqvofVT2JW4r0Kp9jag2l\nItIHQERSgTKf44k6Efklrlv1rB/Yfif8ncAPRaS/iHQGbsUtjh4TRESAvwL7VfVZv+OJNlX9o6qm\nqerFuBt+/1TVu/2OKxpUtQQ4JCIDvYcmAvt8DCnaCoAxInKe9z6diLvxHmveBGZ42zOAmGp0ich1\nuC7Vm1S15mzH+5rwvZbFb4FNuDfbOlWNmZEQwFjgTuAab9hivvcPilWx9nX5d0C2iOzGjdJZ5HM8\nUaOqu4EXcI2uPd7Dy/2LKHIi8jKwHRgkIodE5FfAYmCSiBwAJnj7HVKY65sJZAKJQI6XX55r8hxW\nWsEYY4LB7y4dY4wxbcQSvjHGBIQlfGOMCQhL+MYYExCW8I0xJiAs4RtjTEBYwjfGmICwhG+MMQHx\nP9YqKxhgYoJdAAAAAElFTkSuQmCC\n",
       "text": [
        "<matplotlib.figure.Figure at 0x106e73350>"
       ]
      }
     ],
     "prompt_number": 57
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "and check for their zeros..."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# zeros of Bessel functions\n",
      "n = 0 # order\n",
      "m = 4 # number of roots to compute\n",
      "jn_zeros(n, m)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 58,
       "text": [
        "array([  2.40482556,   5.52007811,   8.65372791,  11.79153444])"
       ]
      }
     ],
     "prompt_number": 58
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Integration\n",
      "\n",
      "The numerical evaluation of the definite integral of a function is called numerical quadrature, or simply quadature. \n",
      "SciPy provides a series of functions for different kind of quadrature, for example the quad, dblquad and tplquad for single, double and triple integrals, respectively."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from scipy.integrate import quad"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 59
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "xl = 0  # the lower limit of x\n",
      "xu = 10 # the upper limit of x\n",
      "\n",
      "val, abserr = quad(lambda x: jn(0,x), xl, xu)\n",
      "\n",
      "print \"integral value =\", val, \", absolute error =\", abserr "
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "integral value = 1.06701130396 , absolute error = 7.43478946065e-14\n"
       ]
      }
     ],
     "prompt_number": 60
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We can pass extra arguments to integrand function, by using the `args` keyword argument.\n",
      "For example we can use `args=(3,)` for the third order Bessel function."
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "---\n",
      "### Exercise 11\n",
      "\n",
      "Create a function that returns the Bessel function of first kinf for different orders and then find the integral."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# type your solution here"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 61
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "##### Solution"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "#%load solutions/bessel_integral.py"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 62
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Interpolation"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "The `interp1d` function accepts given arrays describing points X and Y data.\n",
      "It returns a function that can be evaluated everywhere in the domain for an arbitrary value of $x \\in X$ and the image it is corresponding interpolated value. Options for interpoland kinds are linear, quadratic and cubic spline interpolation."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from scipy.interpolate import interp1d"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 63
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "x = np.linspace(0, 20, 20)\n",
      "y = lambda x: jn(5,x) # sampled on given values x\n",
      "f1 = interp1d(x, y(x)) # linear interpolation\n",
      "f2 = interp1d(x, y(x), kind='quadratic')\n",
      "f3 = interp1d(x, y(x), kind='cubic')"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 64
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# new points where we want to evaluate the interpoland\n",
      "xnew = np.linspace(0, 20, 100)\n",
      "\n",
      "fig, ax = plt.subplots(figsize=(10,4))\n",
      "plt.plot(x,y(x),'o',xnew,f1(xnew),'-', xnew, f2(xnew),'--',xnew, f3(xnew),'--',xnew,y(xnew),'k')\n",
      "plt.legend(['data', 'linear', 'quadratic','cubic','true'], loc='best')\n",
      "plt.show()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAlwAAAEACAYAAABmuKihAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzs3XdcVuX/x/HXYQgoiiAORIZ7pSXuCU5cuc2RW3Ml0rfS\nVByYmqVmirtylJblIBeJGoo4UDAtzS0uFLeIILKv3x+aP3cg49zA5/l43I+4ua9znfftCe4P51zn\nujSlFEIIIYQQIvMY6R1ACCGEECKnk4JLCCGEECKTScElhBBCCJHJpOASQgghhMhkUnAJIYQQQmQy\nKbiEEEIIITJZugsuTdNaapp2StO0s5qmffaadjU1TUvSNK1TevcphBBCCJGdpKvg0jTNGJgPtAQq\nAT00Tav4inZfAf6Alp59CiGEEEJkN+k9w1ULOKeUuqiUSgR+Adq/pJ0HsA64lc79CSGEEEJkO+kt\nuOyB8KeeX3n8vSc0TbPnURG26PG3ZGp7IYQQQuQq6S24UlM8zQHGqEdrCGnIJUUhhBBC5DIm6dz+\nKuDw1HMHHp3lelp14BdN0wBsgVaapiUqpTY93UjTNDnzJYQQQohsQymV6pNIWnoWr9Y0zQQ4DTQF\nIoAQoIdS6uQr2i8HNiulfF/ympKFtEVqeHt74+3trXcMkU3I/y8iteT/FZEWmqalqeBK1xkupVSS\npmkjgG2AMbBUKXVS07Qhj19fkp7+hRBCCCFygvReUkQptRXY+tz3XlpoKaX6p3d/QgghhBDZjcw0\nL7IdNzc3vSOIbET+fxGpJf+viMyUrjFcGUnGcInsyM8vCB+f7cTHm2BmlsTIkS1o06aR3rGEEEJk\nsiwdwyVEbubnF4Sn5zbCwqY9+V5YmBdAjii6Ht9ZLAyE/EEqRPYmZ7iEeEMtWnixY8dQ7AoFY2Vz\nHBPTK2B0C2USQdmitgzu148WXbtibGysd9Q38vivN71jCORYCGGI0nqGSwouIdLo/v37/PTTT4wa\nNZkHD1IoYG2BKlIQU7MCGJvYYJk/iTyWGg9PneL+nTu4urrSpEkTmjVrRqVKlfSOn2ryIW845FgI\nYXjSWnDJoHkhUilg714GDx6Mk5MTf/zxBxaNLWHCHe53sifapi93L63i1uGNlDN34dSWLVw6d46T\nJ0/SrVs3jh07RvPmzenerRtXf/1V77cihBAii0nBJcR/8N++Havy5WndoQNOTk4cPHKQBx0f4NSk\nEM7rRsLuiVD4BAypjtmI4th1ucXNBzcBKFasGD169OC7777j7NmzFC1UCKcPPmBokyYkX7ig8zsT\nQgiRVaTgEuIVkhMT6d6vH627dsW9Zk2iz5+n94jedNnahZIFS3LAcx/zZ3bEvfQ+XKNK0PzoYCbU\n+xSjwolUW1KNnRd2PtNf3rx5mbtwIYu3bOGHW7co3LAhWz/9FBITdXqHOUe/fv2YMGGC3jGEEOKV\nZAyXEC9xZPduWgwbxsOEBNbNmUPLtm35M+JP2v3Sjk/qfsL/6vzvtXfx7QjbQZ8NfRhWYxheDb0w\nNnp24HxsUhI9pk1j85w5uJcqhV9AAEYFC2b220qT7DRuqH///jg4OPD555+/tp2bmxu9e/dm4MCB\nWZQsY2SnYyFEbiFjuIRIp9WrV9OkXTsqV6nCzWPHaNm2LRtPbaTlTy2Z32o+H9f9+D+nTGheujl/\nDv6TgAsBtPyp5ZNLjP/Ka2LCxkmT2Pf33+yNi6Pfxx+TkpKSmW8rw/j5BeHuPh43N2/c3cfj5xek\nSx/PS01BIlNdCCF0o5QyiMejKELoJykpSQ0bNkxVqFBBHTp0SCmlVEpKivom+BtV/OviKuRKSJr7\nTExOVF4BXsr+a3sVeCHwpW2ioqJUw4YNVf/+/VVSUlK63kNGetnP5JYtu1Xp0uMUqCeP0qXHqS1b\ndqe634zo4/Dhw6patWoqf/78qlu3bqp79+5q/PjxKjIyUrVp00YVLlxYWVtbq7Zt26orV64opZQa\nN26cMjY2Vubm5srS0lJ5eHgopZQaOXKkcnBwUAUKFFDVq1dXe/bsSXWOrCK/H4UwPI9/LlNf56Sl\ncWY+5BeK0FNCQoLq0aOHcnNzU/fv31dKPSqWRviNUJUXVFYXIy+mq3//s/6q6MyiauruqSo5JfmF\n12NiYpSbm5vq3bu3wRRdL/uZbNHC65lC6d+Hu/v4VPeb3j7i4+OVo6OjmjNnjkpKSlLr1q1Tpqam\nasKECerOnTvK19dXPXz4UEVHR6uuXbuqDh06PNnWzc1NLV269Jn+Vq1ape7evauSk5PV119/rYoV\nK6bi4+NT/X6ygvx+FMLwpLXgkkuKIteLjYykU6dOREdH8/vvv5M/f36i46Np/0t7Tt85zb4B+3Aq\n6JSufbiXcefQ4EP4h/nT+qfW3Hpw65nX8+XLh5+fHxEREfTu3ZukpKR07S+zxMe/fHGKuLjUT+6a\n3j4OHDhAUlISnp6eGBsb07lzZ2rWrAmAjY0NHTt2xNzcHEtLS8aNG8fu3buf2V49d+nx/fffx9ra\nGiMjIz7++GPi4+M5ffp0qt+PEEKkhhRcIle7sXs3jrVqEQP4+vpiYWHB1ftXabSiEcUti+PX0w8r\nc6sM2VeJAiXY1XcX1YpVw+VbF/Zc2vPM63nz5mXz5s3cuXmT5pUqEX/1aobsNyOZmb28EDQ3T86y\nPiIiIrC3t3/me05Ojwrihw8fMmTIEJydnbGyssLV1ZWoqKhniqznx3HNmjWLSpUqUbBgQaytrYmK\niuL27dupfj9CCJEaUnCJXOviunWU7dcP67Jl2erri6mpKX9d/4s6S+vQvXJ3vn33W0yNTTN0nyZG\nJkxvNp0lbZfQZW0Xvtz7JSnq/wfLW1hYsH7TJo7ky0fDzp3BwAbSjxzZgtKlvZ75XunS4/DwaJ5l\nfdjZ2XH1uWL00qVLKKWYNWsWZ86cISQkhKioKHbv3v30sIUXiq09e/Ywc+ZM1q5dy71794iMjMTK\nykruCBRCZDhZvFrkShfWraPyRx9RrkED/vz5Z4yNjPA740e/jf1Y2HohXSt3zdT9ty7bmtAPQum+\nrjt7Lu/hhw4/YJvXFgDLvHk56OdHZRcXJg4dyufffpupWdLi30W5582bQFycMebmyXh4tEzTYt3p\n7aNevXqYmJjg4+PDsGHD2Lx5M6GhoTRp0oSYmBgsLCywsrLi7t27TJ48+ZltixYtSlhY2JPn0dHR\nmJiYYGtrS0JCAl9++SX3799P9XsRQohUS8uAr8x8IINCRRaJOnBAWZUsqd7u1k0lJz8awD7/4HxV\nbFYxtf/y/izNkpCUoD7d9qlymO2g9l7a+8xrC9auVUYFCqj9a9ZkaaZ/GfLP5KFDh164S3HChAkq\nIiJCubm5KUtLS1W+fHm1ZMkSZWRk9OQ4BwcHq3Llyilra2vl6empkpOT1YABA1SBAgWUnZ2dmjFj\nhipZsqQKCAjQ+R0+y5CPhRC5FWkcNC8Tn4pcJSUlhbZNm3IqJYXTO3dipMGn2z/FP8wfv55+lLIu\npUuuzac3M2jzIEbVG8UndT95cumr1ZAh7Nu2jeuhoeQtXDhLM8lkm4ZDjoUQhietE59KwSVylc8+\n+4x9+/bxxx9/kGyUTE/fnkTHR7P+vfVYW1jrmu3SvUt0W9eNIvmKsKLDCmwsbEhKSaFMw4Y0c3Hh\n+3nzsjSPfMgbDjkWQhgemWleiFdYvHgxv/32Gxs3biQyMRLXFa7YWNjg38tf92ILwKmgE0H9gyhr\nUxaXJS4cvHIQEyMjDm/ezI5Nm9i4caPeEYUQQrwhOcMlcgU/Pz8GDRrE3r17ibWMpe3qtnzg8gFe\nDb0McrmXjac2MnjLYMbUH8NHdT4iODiYDh06EBoa+mQKhMwmZ1UMhxwLIQyPXFIU4jm7N22iy8CB\nbN68mSjbKHr/1pu5LefSo0oPvaO91oXIC3Rb143i+YuzvP1yvpv3HRs2bGDPnj0YG6d+otE3JR/y\nhkOOhRCGRy4pCvGUMxs20HTAAMbNnMnfJn/Td0NffLv5GnyxBVDSuiR7B+zFycoJl29dcO3hiqmp\nKd8uWqR3NCGEEGmU7jNcmqa1BOYAxsD3Sqmvnnu9PfA5kPL4MUoptfMl/cgZLpGhEm/epHijRji/\n8w6NBzqy8fRG/Hr6UcamjN7R0mz9ifUM8xvGBza9+cZjEY3r9udBcmHMzJIYObJFmubBSi05q2I4\n5FgIYXiy9JKipmnGwGmgGXAVCAV6KKVOPtUmn1LqweOvqwC/KaVe+MSTgktkKKVo2qwZodev03RC\nJe48vMFv3X6jUN5Ceid7Y2F3w6i/pAHGvo4Ql4eIfx4tDVS6tBdz57pneNElH/KGQ46FEIYnqy8p\n1gLOKaUuKqUSgV+A9k83+LfYeswSkEXKRKabMWYMgUeOUGFgEcxMjNnRe0e2LrYAStuUxv4vN+65\nRHHtyj+UclgGQFjYNObN26FzOiGEEK+T3qV97IHwp55fAWo/30jTtA7AdMAOaJHOfQrxAj+/IHx8\nthMfb4JR8m32HltFuZ4VaVC5GrNazMJIyxnDFfPfKE9sMpSqEsW18GlAb8CUuLjMH0RvKJydnVm6\ndCkhISGcP3+e7777Tu9IQgjxn9JbcKXqHLdSagOwQdO0hsBKoPzL2nl7ez/52s3NDTc3t3TGE7mB\nn18Qnp7bCAubBiSCaW2Ma8dTv7oLs91n6x0vQ5mZJcG2bzg/5C3Mw60o85Yn5/5ZiLl5st7Rsoym\naWiaxtixY/WOIoTIRQIDAwkMDHzj7dNbcF0FHJ567sCjs1wvpZTao2maiaZphZRSd55//emCS4jU\n8vHZ/rjYAgr0gkLHSY5cxpW1p2Cgvtky2siRLQgLm0dY4DQKuazi/NY1ODub4+HRQe9ouUZKSgpG\nRjnjjKkQIvWePxE0efLkNG2f3t8ah4CymqY5a5qWB+gGbHq6gaZppbXHM0tqmuYC8LJiS4g3FR//\n+O+Gop8D6yD5Fzj+fo68zNamTSPmznWnhW040U4XKWCfl6pVwzPlLkVDppTC29ub3r17A3Dx4kWM\njIz48ccfcXJyonDhwnzxxRfPtP/yyy8pU6YMtra2dOvWjcjIyCevd+3aFTs7OwoWLIirqysnTpx4\n8lq/fv0YNmwYrVu3xtLSMl1/4Qohcq90FVxKqSRgBLANOAH8qpQ6qWnaEE3Thjxu1hk4pmnaEWAu\n0D09+xTieWZmSWB2CqIng+XXcLkjQI69zNamTSO2+U9l35itGL8fy779uzh27JjesbLcy1YI2Ldv\nH2fOnCEgIIDPP/+c06dPA+Dj48OmTZsICgri2rVrWFtb8+GHHz7Zrk2bNpw7d45bt27h4uLC+++/\n/0y/q1evZsKECcTExFC/fv3MfWNCiBxJZpoX2d7mLYH0GtmL+0ZFIexPAEqXHsfcuS1z/JmfcQHj\n2L56O4XCC7Ft27YM7fu/piLQJmfMkkhqUtp+7kuWLMn333/P3r17OXfuHCtXruTixYuUKlWKK1eu\nULx4cQBq167NJ598wnvvvUfFihVZsGABTZo0AeDatWs4OTkRFxf3wuXBe/fuYWNjQ1RUFPnz56df\nv34ArFixIt3v9U3JtBBCGJ60TguR3jFcQujOZ+N3JKYY07JMfR6W8MbcPBkPj5xfbAGMbzSeX/7+\nhXD/cPbv30+9evWybN9pLZSyQrFixZ58nTdvXmJiYgC4dOkSHTt2fKa4MjEx4caNGxQpUgQvLy/W\nrVvHrVu3nrS5ffs2+fPnR9M0SpQokbVvRAiR40jBJbK1XwK2ErBmC/MneTL848/1jpPl8prmZU6b\nOYw6MJzRgwax96mxR+L/OTo6snz5curWrfvCaytXrmTTpk0EBATg5OT05AyXnFESQmQkudVGZFsP\n4uPpN/hDmtSpkiuLrX+1LdcWoxrmHLh9m++//VbvOFkmLQXR0KFDGTduHJcvXwbg1q1bbNr06P6e\nmJgYzMzMsLGx4cGDB4wbN+6N9yOEEK8iBZfItuq83xVLU1M2/LxO7yi6MtKMGNJoBDUavsVYHx+9\n42SJf+fienrg/MsG0f/L09OTdu3a0aJFCwoUKEDdunUJCQkBoE+fPjg5OWFvb89bb71F3bp1X+j3\ndX0LIURqyKB5kS35BfjRqUM3fpk1no5DxugdR3f34u5RdkZpIucbs2TKFAYOGfLfG/0HGahtOORY\nCGF4snTx6owkBZdIrejoaOzK2tF8SDN+m7xB7zgGY+iWoZxcFc7JYxe5efx4uvuTD3nDIcdCCMOT\n1YtXC5HlugzugomTCasn/KJ3FIMyotYIzlc9yt1bt1i2erXecYQQQjxFCi6Rraz/fT1/bP2D9cvW\nY25irnccg/JWkbcoY1+Gth90ZNXy5XrHEUII8RQpuES2ER0dTZ/+fegyqgtNKzfVO45B8qjlwU2n\nY5w/c4b9+/frHUcIIcRjMoZLZBuVWjXmzq1/uHzgCmYmZnrHMUhJKUmUmluKvsl9CdkRkq7Z52Xc\nkOGQYyGE4ZExXCJHGrd8CWdCjvLrRxOl2HoNEyMThtUYRrhzOKdPn+bAgQN6RxJCCIGc4RLZQHhk\nJM4VK/BB7Yos3hiodxyDd+vBLcrNL8coo1H8FfoXa9aseaN+5KyK4ZBjIYThkTNcIsep36MzzkVs\nmbdqo95RsoXC+QrTvnx7kt9Jxn/LFg5t3653JCGEyPWk4BIGbcyyhUSE/sUP40Zimt9K7zjZxoha\nI1h6fCnOrq54zpihd5xsy9nZmYCAgAzrb8+ePVSoUCHD+hNCZB9ScAmDFRMTw7xRYxjk5kKD7umf\nOT03qVG8BsUsi/Fe/xYEh4RwOyxM70jZUnqX9TEyMuL8+fNPnjds2JBTp05lRDQhRDYjBZcwWO0/\naI9F2TzMX+uvd5RsyaOWB4FxfhSuUIH/eXnpHSfHSU5OTlU7GXslhAApuISBWrt1Lbv8drHph02Y\nGJnoHSdb6lKpC//c/If+fd9jTWAgydHRekfKMEeOHMHFxYUCBQrQvXt3unfvzoQJE1ixYgUNGzZ8\npu3TZ5n8/PyoVq0aVlZWODo6Mnny5Gfarly5EicnJ2xtbfniiy+eec3b25suXbrQu3dvrKys+OGH\nHwgNDaVu3bpYW1tTvHhxPDw8SExMBKBRo0YAvP322+TPn5+1a9cSGBiIg4PDkz7Dw8Pp1KkTRYoU\nwdbWFg8Pjwz/txJCGAYpuITBiX0YS98Bfen1WS/qla+nd5xsy8zEjMHVBxPldB4jS0tmLF6sd6QM\nkZCQQIcOHejbty+RkZF07doVX1/fVF3+s7S0ZNWqVURFReHn58eiRYvYuPHRzRgnTpxg+PDh/PTT\nT0RERHDnzh2uXLnyzPabNm2ia9euREVF0bNnT4yNjZk7dy537twhODiYgIAAFi5cCEBQUBAAR48e\nJTo6mq5duz7TV3JyMm3btqVkyZJcunSJq1ev0r1794z6ZxJCGBgpuITBaT+0PRb2FiwfK8vTpNeQ\n6kP45fhqxo76Hzv8M/jSrLc3aNqLD2/v1Ld/VdvXOHDgAElJSXh6emJsbEznzp2pWbNmqi7dubq6\nUrlyZQCqVKlC9+7d2b17NwDr1q3j3XffpUGDBuTJk4cpU6ZgZPTsr8h69erRrl07AMzNzXFxcaFW\nrVoYGRnh5OTE4MGDn/T3X0JCQrh27RozZ87EwsICMzMz6tevn5Z/CiFENiIFlzAoX/38PbvXB7Hp\nx00YGxnrHSfbsy9gT4vSLbCsHM+pkyc5duxYxnXu7Q1Kvfh4XcGV2ravERERgb29/TPfc3JyStW2\nBw8epHHjxhQpUoSCBQuyZMkS7ty586TfEiVKPGmbN29eChUq9Mz2T78OcObMGdq2bYudnR1WVlZ4\neXk96e+/hIeH4+Tk9EJRJ4TImeQnXRiMuzHRTBg7hY4tXKlfSf7Szygjao5gyV9LGD58OHPmzNE7\nTrrZ2dlx9erVZ7536dIlAPLly0dsbOyT71+/fv2Zdj179qRDhw5cuXKFe/fuMXTo0CdnxooXL054\nePiTtrGxsS8UT89fshw2bBiVKlXi3LlzREVFMW3aNFJSUlL1PhwcHLh8+XKqB9/nBklJSQQEBODl\n5cWsAQNYVaECOxo3JqJbN1I+/RR++QXk30tkU1JwCYPRqHdXCuY15+efZYLTjNTAsQEWJhaUa1EO\nX19fbt68qXekdKlXrx4mJib4+PiQmJiIr68voaGhaJrG22+/zfHjx/n777+Ji4vD+7kzaDExMVhb\nW5MnTx5CQkL4+eefn7zWuXNntmzZwr59+0hISGDixIn/WTzFxMSQP39+8ubNy6lTp1i0aNEzrxct\nWpSwV0zJUatWLezs7BgzZgyxsbHExcXligXH/fyCcHcfj5ubN+7u49m4ehOLf/uNQYMGPfn3MDIy\n4qqZGd8WK0bLa9cosWMHJnPnUmLMGLaMHKn3WxDijUjBJQzCnF+XcWJnMD9OGYexuYXecXIUTdMY\nUWsEK8NW0qVLF5YsWaJ3pHQxNTXF19eXFStWUKhQIdasWUOnTp1QSlG2bFkmTpxIs2bNKF++PA0b\nNnzmrNTChQuZOHEiBQoUYMqUKXTr1u3Ja5UrV2bBggX07NmT4sWLY2Nj88wdhS8blD9r1ix+/vln\nChQowODBg+nevfszbby9venbty/W1tasW7fumT6MjY3ZvHkz586dw9HREQcHhzdehim78PMLwtNz\nG9u3T+Xg3hGcvR5Mx8Hv8/H4CVSsWJHQ0FBCQ0OZMmUK3yxaRFBgIEknT3Lj+nUCr16l8ocf0n7d\nOuq0acOZM2f0fjtCpEm611LUNK0lMAcwBr5XSn313OvvA6MBDYgGhimljr6kH1lLMZe6//ABhSu/\nhfs75djku03vODlSbGIsjt84srL+Snq26ML148cxc3Z+7TbZaf2+/v37U6JECaZMmaJ3lEyRnY7F\n67i7j2f79qnY2ezkXv5B5MlrTfzt73F18cXfP3XH7vjdu6xcsIDv586lR48eTJw4kcKFC2dyciFe\nlKVrKWqaZgzMB1oClYAemqZVfK7ZeaCRUqoqMAX4Nj37FDlP5yEdKKil4LtaLiVmlrymeen/Tn/+\neBBAkrMz06ZO1TtShsoJxUhuEB9vQln7aVxP7kTRQrWIanqPuF7tOWcbSkxCTKr6qGxjw5cTJnDy\n5Ek0TaNixYpsnDULHo/jE8JQpfeSYi3gnFLqolIqEfgFaP90A6VUsFIq6vHTg0AJhHhs095NBKwP\nYP1vP2FiZq53nBxtWM1h/PDXCjp26cy3e/ZAUpLekTJMepfgEZkvKSmJiIv+REVMpnTR8VysfQKO\n9odf1xNjfRnnOc54BXhxI+ZGqvorXLgwPj4+bNu2jSFTprC6Rg24cCGT34UQby69BZc9EP7U8yuP\nv/cqA4Hf07lPkUPEJcTRs09Penj0oEHVBnrHyfFKWZeinkM93m5RnJvXrxOyYoXekTLM8uXL+fzz\nz/WOIV7hzp07NGvWDEubZPKWHMi5huvhfDMI8qK0xW8sbbmYA4MOEBkXSYUFFRiyeQhn7qRujFb1\n6tXZsXcv/VNSGN2uHTy1dqUQhiRdY7g0TesMtFRKffD4eS+gtlLqhfUpNE1rDCwA6iulIl/yupo0\nadKT525ubri5ub1xNmH4Wn7YktBdodw8dhNjY5lzKytsD9vOp9s/JcW/JEUiI9n5eDb0l8kp44Zy\ngux8LOLj42nWrBkuLi5MnzGdBgtduXkukdLH22NhnoKHR3PatGn0pP3NBzeZHzKfRYcW0cipEaPr\njaZ2idr/uZ+VoaH0a9OGwcWKscjXF8qUycy3JXKhwMBAAgMDnzyfPHlymsZwpbfgqgN4K6VaPn4+\nFkh5ycD5qoAvj4qzc6/oSwbN5yIb92+kY/OOBO0LosE7cnYrqyilqLigIt3tP2XqwNE8OH8es4IF\nX9o2O3/I5zTZ9VikpKTQp29f4h4+ZPUvq+n5W09SVAq/dvn1P9dIfZDwgGVHlvF18Nc4FXRiVL1R\ntC7bGiPt1Rdmfvv7b7q0bEmPokVZFRwMFnLHs8g8aR00n96CywQ4DTQFIoAQoIdS6uRTbRyBnUAv\npdSB1/QlBVcu8eBhDA6VyuLe0Y3Vs1frHSfXmXdwHkGX93B6djjeo0fTqWPHl7bLrh/yOVF2PRZt\nhg3jz337CAsOxjPAk4v3LuLX0w8zE7NU95GUksTa42uZsX8GCckJjKo3ip5VepLHOM9L2+84dYpW\nzZvTtmdPNnz11UvbCJERsrTgerzDVvz/tBBLlVLTNU0bAqCUWqJp2vdAR+Dy400SlVK1XtKPFFy5\nRP2OLfj7WBh3j58gj1nqf/GKjHE//j7Oc5yZZD2JbRu28fvvLx9WmV0/5HOi7Hgsxnz2GbOWLWP/\nzp0citvDsiPLCOwXiGUeyzfqTynFH+f/YOb+mZy4dQLP2p4Mrj4YK3OrF9oGnz9PK1dXls6ZQ+fO\nndP7VoR4qSwvuDKKFFy5w7er5jF0uBe//Pgt73XornecXGvE7yOw1Cz5vvf3/Pnnny9dizA7fsjn\nVNntWKz67jv6jBrFivnzadGhGVUWVWF3v91UKlwpQ/o/cu0IM/fPZFvYNgZVG4RnHU+K5y/+TJuQ\nkBDatm3LgQMHKFWqVIbsV4inScElDFZ09D2KudSgVrWK7FqzWe84udqp26dwW+FGl0tdsLG2eekd\nftntQ/6/uLm50bt3bwYOHPjCa5cvX6Zy5crcv3/fIKeXyE7H4tCuXdTp3JmPPvqIWRMn0nN9T5wL\nOvNF0y8yfF8X711kdvBsVh1dRccKHfm03qdULPz/U0HOmTOHn376ib1792ImZ9NFBsvSiU+FSItm\nPTpjbKSx/SdfvaPkehVsK1ClaBVKuJVg6dKlJOWgOble5XVzdTk6OhIdHW2QxVZ2Enf3Lq26daNJ\nhw7MmjiRHWE7CL4SzPhG4zNlf84FnfFp5cNZj7M4F3TG7Qc32q1ux97Le1FK4enpib29PaP79IFU\nLiouRGaRgktkiZ//+Ikju/bz87zZmJqa6h1HAB61PNgQuQF7a2uWDx+udxyRA0ybOJH6Njb4L11K\nXFIcw39m6DSrAAAgAElEQVQfzvxW88lrmjdT91sobyEmuE7goudFWpVpRb8N/WiwvAF7Lu9h+bJl\nbNq4kTkt2jyzaLaf36unRBEiM0jBJTLd/dj7DBgwkP6j+9C2xbt6xxGPtSnbhusx1ynargWTgoIg\nLk7vSKkWHh5Op06dKFKkCLa2tnh4eODt7U3v3r2ftLl48SJGRkakPHVm49y5c9SuXRsrKys6dOhA\nZGTkS9vevXuX/v37Y29vj42NDR1fcSen+H9Hjx5lyZo1LNy5EyNN48u9X1K1aFXalGuTZRksTC0Y\nVnMYp0ecZliNYfTd0Jc+2/rSa+QHjDq4lzOH67F7tzfbt0/F03ObFF0iS0nBJTJdyyEtsbG1YfHE\nxXpHEU8xNjJmeM3hmL99l+sRERxeskTvSKmSnJxM27ZtKVmyJJcuXSIiIoLu3bv/5+VApRQ//vgj\ny5cv59q1a5iYmDBy5MiXtu3duzdxcXGcOHGCmzdv8vHHH2fGW8kxkpOTGTRoEF988QXFixfnzJ0z\nzA+Zz9yWc3XJY2xkTK+qvTj14SmaODfhqzzLKO5Sixv5RmBi/GilubCwacybt0OXfCJ3koJLZKql\nfks5uOEg29dul/ExBmhgtYH8cXEjZV0b8dn69Wna1vvCBbTAwBce3q9Yz+5l7V/V9nVCQkK4du0a\nM2fOxMLCgjx58lC/fv3/HFSuaRp9+vShUqVK5M2blylTprBmzZoXtrt27Rr+/v4sXrwYKysrTExM\naNiwYZpz5iY+Pj5YWloycOBAlFIM9xuOV0MvShTQd+lcMxMz/lf3f9Q64MFl4xqYW5pT9J3/X+43\nLk5WuBBZ5/VT/QrxBvz8gvDx2U5MXALBx7+mfe8uvFX6Lb1jiZcolLcQnSp04r6VDRsnrCD5/HmM\nU3kLvXfJkniXLJnqfaW1/auEh4fj5OSEkVHa/150cHB48rWjoyOJiYncvn37hf5tbGywsnpxfifx\nosNnzjBl2jRCDhxA0zRWH1vN7djbeNR+YYU33eQzNoGAqZg4vUNE+FCMHXaRHN4Yc/NkvaOJXETO\ncIkM5ecXhKfnNrZvn8Kl+zsxs3Dg6LbSMlbCgI2oNYLg2F/Q8uXj22xwWdHBwYHLly+TnPzsh6Wl\npSWxsbFPnl+/fv2FbS9fvvzM16amptja2r7Q/927d4mKisrg5DmPiomhTefOvD1gAGXKlOFe3D0+\n2f4Ji9su/s+le7LSyJEtKF3ai1uXeqAa1CfZwoOSpcbi4dFc72giF5GCS2QoH5/thIVNo0z93lwL\nv4jJg/WcP/+FjJUwYNXsquFc0ImGHdz4IyJC7zj/qXbt2tjZ2TFmzBhiY2OJi4tj//79vPPOOwQF\nBREeHk5UVBTTp09/ZjulFKtWreLkyZPExsYyceJEunbt+sKlbjs7O1q1asXw4cO5d+8eiYmJBL1m\nke/c7NM+fbgbG8uGKVMA8Arwol35dtQpUUfnZM9q06YRc+e64+4+gWqaAzy4yLsd1DOLZguRWn5+\nQbi7p32qEym4RIaKjzfB6q0lnD/5OyWLfkhMZHVAxkoYOo9aHsSWPMPOLVuIiYnRO85rGRkZsXnz\nZs6dO4ejoyMODg6sWbOGZs2a0a1bN6pWrUrNmjV59913nymm/h3D1a9fP+zs7EhISMDHx+eZ1/+1\ncuVKTE1NqVChAkWLFn2mnXjkvK8vcwICmO3jg5WZGSFXQ/A95cv0ptP/e2MdtGnTCH//KRzesoSG\nQ1xY8eMiOYsp0uz/r+JMTfO2MtO8yFB12g/k1PlgrEyLcfnIziffd3efgL//FB2TiddJTE7Eea4z\nZbaWYWDvgfTp0ydbzW6e0xncsUhOpoqLCynlynF87VqSUpKo9V0tPq77Mb2q9tI73X86fO0wDTs2\npH+t/sz3ma93HJGNuLuPZ/uOiRQo3YT75/bJTPNCH9djrnMn+Sjx9x9w7e//n02+dOlxMlbCwJka\nmzKk+hDMapixYsUKveMIA3fgq684fvkyvgsWALAgZAEFzQvyfpX3dU6WOi52LjTo34AfV/3A38HB\nescR2Uh8vAn2NduSknw1zdvKGS6RIR4mPqTWF7UImxnGlDHe7NgbTVycMebmyXh4NJexEtnA9Zjr\nVJhbAeNvjDn852GcnZ0N66xKLmZoZ7jaly9PrRYt8Jo3j6v3r/L24rfZN2Af5W3L6x0t1UKvhjKz\neSOuRBZk79Wrb3TXq8h9XJp35thfQVgXnM6tcx/I4tUiayml6Lq6K3+M/YOvJ3790sWBRfbQc31P\nrq4Iw9XUnim//WZQH/K5mSEVXPv376dHjx6cPn0ac3Nz3lv7HhVsK/B54xcXQDd07y1sSuj4wzT3\nGMm3kyfrHUcYuCObN9O0Qwcc3Svw99bjgCxeLbKY104v9i3bR9OaTRkwYIDecUQ6eNTy4Hz5BL7a\nt0/vKMIAKaUYM2YM3t7emJubs/XsVg5fO8zYBmP1jvZG/vfuVN5pbMdSHx8u3LihdxxhwG5cuECH\nzp1xdTHBc9Ak3N0npLkPKbhEuiwKXcSKX1dgcsaE77/7XmaTz+bqlKiDbVmNlDx59I4iDJC/vz+3\nb9+md+/ePEx8yIitI1jQegEWphZ6R3sjdR3q8vDd4pQsVZJeH32kdxxhoBISEujcqxeVqtryzoyx\n9O/03hvdBCYFl3hjm05vYqKvN/G/xbN2zVqsra31jiTSSdM0PGuPpEStMnpHEQYmJSWFsWPHMm3a\nNExMTJi2ZxrV7arjXsZd72jpMqHJZCq9HU2wnx//XLyodxxhYJRSjBgxAvJpnOhpyugGn71xX1Jw\niTdy4MoBBqwbSszP5owY+AF16hjWRIfizXWr3I17b4XrHUMYmHGzZ3MT6NChAydvnWTJn0uY03KO\n3rHSrb5jfR64OlCyUR16T5qkdxxhYEbOns2+/fu42+ous1vOTtfZXCm4RJqdvXOW9r90wNzfmeK2\ntkz6PPsNlhWvZmFqwdCGXfWOIQxIXGgo38yYwRAvLwCG+Q1jYqOJFM9fXOdkGWOSmzcJDa5xYuNG\nLl26pHccYSBWHD7MwilTaDOhM8ULFadTxU7p6k/uUhRpcvPBTeotrYf90RqE+O8lfPdubEuX1juW\nyGCX7l3C2VqmhTAUet+lOKBJE9ZHRXHv0CFWHl2Jz0EfDg46iLFRzllBwm2FGwX2FqAIRfj+++/1\njiN0di8hgWJ169KlVXO2Wy1jV99dVC5S+Zk2j38u5S5FkfEeJDzg3dXv4vKwCXvWbmXLjBlSbOVQ\nTgWd9I7wWs7OzuzcufO/G4p0u7dtGz/+9RdfzZhBZFwko3eMZnHbxTmq2AKY6DqRE+VOsHHjRs6c\nOaN3HKGzdqNHk//WLfLWuEOPt3q8UGy9CSm4RKokpSTRfX13nPOUYsPXGxjTsydNe/bUO5bIpV53\nxicpKSmL0+Rsw2fOxKZiRYY2bcrYP8bSpVIXahSvoXesDNfYuTF2he1o0rMJk2QsV662+uBB9i1b\nxkKPvmw6t5nJjTNmjjYpuMR/Ukrxod+HxCXEEfVzJB+0bcsXCxfqHUvkUr179+by5cu8++675M+f\nn5kzZ2JkZMSyZctwcnKiWbNm7N69GwcHh2e2c3Z2JiAgAHj0//SXX35JmTJlsLW1pVu3bkRGRurx\ndgxa4t9/43/gAFMmTGB/+H42n9nMtCbT9I6VKTRNY2Kjifzt9DeBfn78LWdQcyWlFJ8MHMj7NWsy\np9BOpjaZSkHzghnSd7oLLk3TWmqadkrTtLOapr1wv6SmaRU0TQvWNC1O07RP0rs/kfWm751OSEQI\nZULLkJyczJwlS0Dm2xI6WblyJY6OjmzZsoXo6Gjee+89AIKCgjh16hT+/v4vPfuladqTeeJ8fHzY\ntGkTQUFBXLt2DWtraz788MMsfR/Zwep9+6hWvjwDmjdl6JahzHafjZW5ld6xMk2zUs2wtrJmYLkS\njPrgA73jCB0s/+gjip09i/vk94hPiqf/O/0zrO90FVyaphkD84GWQCWgh6ZpFZ9rdgfwAGalZ19C\nHz/+/SPf/vktXaK7sHvnbtauXYupqanesYQB+LeASe8jvf4trry9vbGwsMDc3Pw/t1myZAlTp06l\nePHimJqaMmnSJNatW0dKSkq68+QUKSkpzFi4kM+mT8fnoA/FLIvRrXI3vWNlKk3TmOQ6iV0dk9l5\n/z5LV6/WO5LIQlfPneOz+fNZOHc2n4VMwaeVT4aOVUzvGa5awDml1EWlVCLwC9D+6QZKqVtKqUNA\nYjr3JbLYjrAdjNoxijFFxzD/6/ls2bKFggUz5tSqyP6UUhnyyCjPX0J8nYsXL9KxY0esra2xtram\nUqVKmJiYcEOWd3ni999/J0+ePFSoVYHpe6ezoPWCXLGShHtpdxKLWNG6ST1GT0n7bOIie1JKMfR/\n/2P40KFsLHOFJiWbUM+hXobuI70Flz3w9AyJVx5/T2Rzf13/i/d932e47VhGDf2M9evXU6pUKb1j\nCQHw0g/+p7+XL18+YmNjnzxPTk7m1q1bT547Ojri7+9PZGTkk0dsbCx2dnaZGzwb+eqrr/jss8/w\n3ObJiFojKFuorN6RsoSmaUx0nci12heIun2bZatW6R1JZAFfX18uXLjAe+OH893h7/iy2ZcZvg+T\ndG6foRPDeHt7P/nazc0NNze3jOxepNKle5do+3NbPq0wkXEDP+fjbt2oVy9jK30h0qNo0aKEhYXR\npEmTl75erlw54uLi+P3332nevDlffPEF8fHxT14fOnQo48aN44cffsDR0ZFbt24RHBxMu3btsuot\nGLSdQUFERERg9pYZx3ceZ3Xn3HVprU3ZNkyynERTt7p4zZ3LgF699I4kMtGl+/f5+NNPWb50KZ/t\n+ozR9Ue/dFLfwMBAAgMD33xH6bwUUAfwf+r5WOCzV7SdBHzymr6U0N/d2Luq4vyK6vPtXyrzMmVU\nK3d3pVJS9I4ldGDIP5MbN25Ujo6OytraWs2aNUsZGRmp5OTkZ9qsWLFC2dnZqSJFiqhZs2apkiVL\nqoCAAKWUUikpKWr27NmqfPnyKn/+/Kp06dLKy8tLj7eSKll6LBISVPG6dVWv6VOV0zdOakfYjqzb\ntwHZcHKDquJTQxkVKqQ2HDigdxyRiap//LEq26yZ2nJ6iyo3r5yKT4pP1XaPfy5TXTOla6Z5TdNM\ngNNAUyACCAF6KKVOvqStNxCtlPr6FX2p9GQR6ReXFIf7KneqFKrOau/d2AIng4MxypNH72hCB3rP\nbi7+X1Yei99nzqTt9Ol8uHIQd+Ov8lOnn7Jkv4ZGKUW1JdUofNqFfBEP2PDrr3pHEplg+9mztKxZ\nk+C9u+m9qytzW86lVdlWqdo2rTPNp3tpH03TWgFzAGNgqVJquqZpQwCUUks0TSsGhAIFgBQgGqik\nlIp5rh8puHSUolLosb4HKckphH4TTnJMDOf27MHM2lrvaEInUnAZjiw7FkpRrnp1bN+qwNkqOzg2\n7BjFLItl/n4NlO9JX6bumEr4lHAOHDhAaVlZI0dRSlG8SROq2NrS7KOaBF0KYkvPLanePssLrowi\nBZe+Pt3+KQfDD+K825mLZ8+xdeNGLIsU0TuW0JEUXIYjq45F6MqV1P7wQ1ymV2dAra4Mrzk80/dp\nyFJUCm8vfpsqx6uQPzE/S5Ys0TuSyECzFy9m9GefcSR0F403tODAoAOUsSmT6u2l4BJpNvfAXBYf\nWkzdo3U5e/Is/v7+5MuXT+9YQmdScBmOrDoWNRo04L6VBQU7RBE8MDjHrZf4JtYcX8NX27/iwpQL\nHDt2DHt7uRE/J0iOi6NA9eoMdHPjXrNo7PPbM73Z9DT1IYtXizRZf2I9M/bNoNHpRhw7cgw/Pz8p\ntoTIhe4fPszJ0FAia1zMkYtTv6kulboQZxZHq+oVmN2mjd5xRAbZPGQIJa5cocsn77Hzwk68Gnll\n+j6l4MrF9lzawzC/YbS73o7g3cFs27aNAgUK6B1LCKGDFQEBlKpkS886rXGxc9E7jsEw0owY33A8\nVxvG8v0//3Dh0CG9I4l0ij9/nk9XrWLuN7P53x8f81Wzr7DMY5np+5VLirnUyVsnabTclaKB73D3\n6AmO/PknRYsW1TuWMCBySdFwZPaxSElJoVSZUtxrdY9LX1/K0eslvonklGQqL6xM3NoilLbIS4C/\nv96RRDrMmjSJ3Vu20O7boaz4ewV7++99o1UU0npJMb0Tn4ps6Fr0NVr/2AbzTaW5HHGOv/z9pdgS\nL5UblnLJzfz8gvDx2U5ExAWuxtzgvbd7SLH1EsZGxoxvNB6fez+wa2YwN48epUjVqnrHEm8gMjKS\nrxYuxO8PP97d8S5b39+aZb/n5AxXLhMdH43bwsZc/CEB48Rkjm/bRuESJfSOJbKRlAcPyFOyJP/7\nYiQzB43XO454Q35+QXh6biMsbBqYukKzP3G+MIT5s9rTpk0jveMZnKSUJCouqEj0zzbUK1wU302b\n9I4k3oCXlxc3btwgX+d8xCXFseTdN7/zVAbNi1dKTE6ky+JOnF58HytTMy4GB0uxJdLMKF8+qlYu\nz9pf/PSOItLBx2f7o2IrzwkwCoGYj7h46mvmzduhdzSDZGJkgldDL0o0K8im4GDiEhL0jiTSKPji\nRb5esIBuw7ux+p/VTGs6LUv3LwVXDufnF4S7+3gaNh1H8UGV2D/9AG85OHJ2/37yygB58YbGdelK\n+F9nuRh5We8o4g3Fxz8aUWJRaSyaW3kI/QSAuDi5O/FV3q/yPpGFz2Hu6MDYZcv0jiPSqP/EibzT\nvj1f/vMlE10nYpvXNkv3LwVXDvbvJYPtxzqx12YVt3+9hmlyNcZ/NAFjU1O944lsrMvgYeSxsODT\nBVn7F6LIOGZmSZhZXCXu3E6KxNeFuEerSpibJ+uczHCZGpsyrsE4HJrmJ2DJErmpJBsJ8PXlzMaN\n9OrfkJsPbjK0xtAszyAFVw4212cbYUUKQmVX+C0KHv5E5O0g5s//Q+9oIrszNaV+lUrs2n5IPnSy\nqZEjW+BYtR82pRy5cfBR4Vy69Dg8PJrrnMyw9X67N7HFLhIdc49du3bpHUekRnIyH8yaRcPmzZn1\n11R8WvpgYpT19wxKwZVDXYu+xt8lf4WbM2BnCUj+E2gPyCUDkTFmzJnHvcNHCTofpHcU8QbcW9bn\nwuXjlCthhWstH9zdJzB3bksZMP8f8hjnYVyjcVi6WTJr1iy944hU2DJ/Phf/+QeXXmWoZV+LxiUb\n65JD7lLMYRKTE1l1dBUzVkznvG8slnHluXt7M5D3SRt39wn4+0/RL6TIMUq9U4oybcqwfdp2vaOI\nNBo9cxrz5i3m0rEQiljZ6R0nW4lPiqfM7DI8/Pohu7dup7KLTBRrsBISaFC+PBb163C46nYODz6M\nU0GnDOla7lLMpR4mPmRByAJKzS3Foi+XcfqHmzSrVRer/LV4utiSSwYiIw3pP4TALYE8SHigdxSR\nRlvmz6NlVWcptt6AmYkZYxuPpVyVfMxu3hzkZIHBOjF5MmciIsjbLBbP2p4ZVmy9CTnDlc3FJMSw\n+NBiZgfP5u3Etwn96QoPoqL47vPP6dWvH35+Qcybt4O4OGPMzZPx8GgulwxEhrlz5w7FHIox/4/5\nDKk3RO84IpWCDgfRvJ4rYX+GUqJyDb3jZEtxSXGU/KIk0VPvcOCHH3irRw+9I4nnpaTQxcqKYp2b\n41ftCCeGn8DC1CLDuk/rGS4puLKpyIeRzAuZx/yQ+bg5ulHqXCm+mT4f1+bN2fDDDzLlg8gyNRrX\n4EHZB5z89qTeUUQqVe9UHSMTI0LXhOodJVvzOejDtLGrqZBowu49e/SOI55z+PBh2rRujfVnNkxr\nOY2OFTtmaP9ySTGHu/ngJmP+GEOZeWW4cO8Cy+st58rcK+zdupe/goPZ/ttvUmyJLOXZrgeR/pc4\nH3le7ygiFe7cv8ORbUeYNVYGfKfXBy4fkFI/mj1Hj3JHCi6DM3nyZOr1rI+9rT0dKnTQO44UXNnF\nlftX8NzqSYX5FYiOj8a/nT95t+elb/u+9OrVi6CgICpWqaJ3TJELtWnoyo27RswI+FbvKCIVRs0Z\nhW1ZW1yrueodJduzMLVgXOuB5C9bkhHTp+sdRzxl/6FDHAw5yG7b3cxtOdcg1oWVgsvAhd0NY/Dm\nwVRdVBUTIxM2tNvBwe+u0qJOc8zNzTl58iTDhw/HyEgOpdCHTY0alCxbls2rA0lRKXrHEa+hlOKX\nFWsYPnS43lFyjCE1hmDa0Ih1hw+TmJiodxzxWKcxYyjw7jv0culFpcKV9I4DSMFlsI7fPE4v317U\n/r42xSyLEdDzICErb+BWuzHx16+ze9Mmvv76awoXLqx3VCEYVK8ed47eYtcFmQjSkH2+dglx0WaM\nGzhO7yg5Rl7TvIzu3ANlU4BJP/2kdxwBrA4O5taRw0SVOoa3m7fecZ6QQfMG5vC1w0zbM429l/fy\nUe2PKJdcgwkz53AqKIhSZcqw/KOPaNi3r94xhXhG3Jkz5KtVi6bTWrH9w9V6xxGvUNytDk6aCcG7\n9uodJUeJSYjBbmhlSv1jw98hR/SOk7slJ2Nfty6JZcz4clR/BlQbkGm7kkHz2dTey3tp9VMr2q1u\nRz27eswoNgP/8f4Me+997B484M9Fizj3559SbAmDZF6uHG9VqMARv7NExUXpHUe8xOHwU1w79A+z\nPvhA7yg5jmUeS8b0GUzYpbMcOHBA7zi5mu/8+Vw/c4bidePo904/veM8Q85wZRE/vyB8fLYTH2+C\nmVkSI0e2oHXrhuw4v4Npe6Zx+d41OuRpS/zxOHzX+1KpUiWGDRtGhw4dMJWFpkU2sGrFCj7w+oQ5\nG79gSA2Zk8vQ1BnQiUshp7l27B8wgAHEOU10fDR2Xe1oaNSQrRu26h0nd1IKxxo1SCkG65YsoE6J\nOpm6O5mHywD5+QXh6bmNsLBHC8Ricw6b+kPIX+cGRucLEncqmetHj1PFwYFu779Ply5dKFeunL6h\nhUgjpRT2zvYU6lOIY1OO6R1HPCU+MZ78ZcsyunULpi78Xu84OdaErROY8d4MTh84jHPlynrHyXVO\nffstNTw9ab20K2t6/pjp+8vygkvTtJbAHMAY+F4p9dVL2vgArYBYoJ9S6oWL3Dm54Grh7sWOo12h\n1CoosJYikRbcv51E3NWr2Do60qJyZUZ16sQ7XbuCmZnecYV4YxMmTuCbwG8IXR9KxcIV9Y4jHpu6\nZDwzxiwg8swZjOVGm0wTFRdFveqFaRpTBJ9LV/SOk7soxfuFC3HK/gFb9l7ELn/mL1mVpQWXpmnG\nwGmgGXAVCAV6KKVOPtWmNTBCKdVa07TawFyl1Avn+TRNUy1aeDFyZItMX3rmZZf33mSficmJ3Iq9\nxfWY608ed25eIv70Ka4ej+BYzF1ORd3g7sVI1FUjiDEGk3ewK1iBInlLUcb2OusOLJDT+yLHOHPm\nDC51XRj20zBmtpypdxzxWMmahXm3YCl8dhzUO0qON2rVSBYNWMyZffsoXrOm3nFyjXP791OzUQM+\n2TyJ8a0mZck+s7rgqgtMUkq1fPx8DIBS6sun2iwGdimlfn38/BTgqpS68VxfChSlS3sxd657phVd\nL1zeg2f2maJSiHwY+UwRdePBDW5FXuX2vQjOP7xJxM0Ibty+QdTtKCzjLCl51Yb46Lzci00gKiaG\nh3fuYFyoEKXti+JarxZ7AsM5dewroDKPTgQ+4u4+AX//KZnyPoXQS1WXqlytcZUbi29gYmSid5xc\nb9/pfTR8pyHhF8KxL2avd5wc717cPZwbNKKBnQNbNvvpHSfXaN65OUceHiFicwR5jPNkyT7TWnCl\n97ehPRD+1PMrQO1UtCkB3HiuHc4N3yfRxIiPv/2KL39fho2poohtHhSPikKlFEopIu8lciXeBJRC\nQ4ECVAqFjBOxsTF+0k4pRUpKCrF347kTDVpKCndvR5FSNC8lizQlPjaGiMQYwowe0NF7Dnm+1IiN\njaXIQxtQNiQnJZOUmERiQgIJD+NQDx6gAVbWVtha21KtRDUcSjiQUsWCKEtLSpcpQ9XKlWlWtSol\nnlpe51GR9+tzRd44PDxapusfXwhDNKh3f2atmsK2c9toU66N3nFyvdFfj+at+m9JsZVFCpoXpPl7\ndfht5m/E376Nma2t3pFyrH+vVkVG3SP07wDG+EzKsmLrTaS34Ert6bHnK8CXbnfr1A40pTAihRjL\nfNzOb8I/NtrjDv6/C7NoRUrcv2eKNB7VlxqxJslcsDVG0zQ0tEf/NdIwilHEJZmBkTHxiUmkqATQ\njFHxZnDHAZIKUKToCQZ1bkmRgkWIeRjP5QcPyF+gAPmtrLApUgQnZ2cqFCuGc8GCaV4i4N+zdfPm\nTSAuzhhz82Q8PFpm+qVTIfTQqWRJPj4Vz6IDy6Xg0llUXBQHNhxg/ar1ekfJVZZ4fInf9zv539hR\nLPxuud5xcqRnrlaVagCly7D2qwQaFA/KtM/WwMBAAgMD37yDp88GpfUB1AH8n3o+FvjsuTaLge5P\nPT8FFH1JXwqUAqXc3cerzNKihdeT/Tz9yMx9CpGrJCcr2woVlF2PaurWg1t6p8nVRizwUAVKFFAp\nKSl6R8l1mnq+p8zKlJJ/+0zy5LPcPkRhaamw35/ln+WPSqjU10zpnfj0EFBW0zRnTdPyAN2ATc+1\n2QT0AdA0rQ5wTz03futpjy61NU9nrFcbObIFpUt7Zek+hchVjIx4r1YtYo8n8POxn/VOk2sppVjy\n826a9OlrEAv35jY/fj6X+JgYZmz6Ve8oOVJ8/OMLdCVHQ7V6cLUuAHFxxq/ZSl/puqSolErSNG0E\nsI1Ho8GXKqVOapo25PHrS5RSv2ua1lrTtHPAA6D/q/pzd5+Q6Zfa5PKeEJlv4siRLHR1ZcmenxhZ\ne6TecXKlZfvXkPhXGN/4yCUtPRQvUIyKjcvg/b9hbP3mVLruiBcvMjNLwqrUWqL+OoRmvu/JOCVz\n8+/WHNEAACAASURBVGRdc72OTHwqhMgUJV1cSHE2ZcOCxVSzq6Z3nFynfKfGJFy6x4U/ZW0/Pfj5\nBTH801+5fGkRmP8OkS0z/S783GTLlt0Mmf4xRnFmXDm8H3h0tWru3Kw7gSJrKQohDMLA9u25e/w2\ny/+SMyxZ7ULkRcIOnGJ0+/Z6R8m1fHy2c/nUArBpAraPlroKC5vGvHk7dE6WM5zPc4JrR89QzaoS\nrq7euLtPyNJi603IGS4hRKZITEzEzt6OlAEpXJt6DTMTWUUhq7Qd34c/Vuwk9uxZjCws9I6TK7m5\nebN7tzcYXQLTUuQpspKE8J64unoTGOitd7xsLTYxlmZv21IwuTi/nz6nWw45wyWEMAimpqb07N4T\nqzNWbD6zWe84uUZcUhxH1+7i3Vo1pNjSkZlZ0qMvUpyoRQlqpgwGLdmgxxhlF1O3fs6pUw+ZP+UL\nvaOkiRRcQohM06tXL+IOx7HsyDK9o+QaS/cs5c6FCL7z/lzvKLna03fER8XP5sjtJIr8X3v3HV/T\n/cdx/PXNtDcl9qjZGFW0KIlfBTFaVVSrQbX2rFq1R41WtUajVWpUSVrUSJoQ46IqihCCSIxYIbEj\n1UTG9/eHUFpU5rn35vN8PDzcM3LO+7bXyed+v9/zPS3ayB3x6XT+1nkWz/qKVuXLUqFjR6PjpIo8\nd0MIkWnq1atHHvs87Pp9F5HtInHK62R0JKs38+uZvNLWlQI1axodJVv75x3x16Irkssumldfq21w\nMss2zGcYsYdtGb93k8U9h1hauIQQmUYpRbdu3Sh1uhQ/BP9gdByrt/f8Xi7tuMTkj6V1yxy0bt0E\nf/8pmEwTmfpeV87vu8gUk2V1g5mTPef3sGnVJtxbuFOlShWj46SaDJoXQmSqs4GBVGvWjJJTyxE2\n9KhMwpmJ3Ca5ceSHI0SGR8p/ZzOjExLIX6UK+RsWZNf8NZQrUM7oSBYlWSfzomdDTow/wR+mnTg7\nOxsdSQbNCyHMS9kGDchRpgx2QfkIvBBodByrdfXOVUxrTAwdMFSKLTOk7O3p36ED0ftiGBEw0ug4\nFmfF4RVE7HWk8Ev1zaLYSgspuIQQmUspOjVpQtSRP2VOrkz06ab5JJ2xp1f3XkZHEU8wefp09J04\nNu89we/nfzc6jsWIvRvLcL9xxPgeYc6UKUbHSTMpuIQQmW7S6NHcPHWGjXs2cSfhjtFxrE5SchK+\nC36l/IsvUqBAAaPjiCewt7Nj+siRVNhvw2D/wSTrZKMjWYSZv80k/8HncSpTlg716hkdJ82k4BJC\nZLrnypenorMzhY+UYe3xtUbHsTrrj67ndPBZxnl0NzqK+A/9evbk4vGLJFxOYNmhZUbHMXtnb55l\nfuBCTm0+wJw337S4OxMfJgWXECJL9H/zTcKPnJM5uTLBnNnzcMiZE4/33zc6ivgPOXPmZMCAAZQ6\nUoox28YQEx9jdCSzNmLLCCqFvEiJ0qXpMNKyx75JwSWEyBIDPvqIoloTFBRExM0Io+NYjbBrYRzZ\nH8VbrVrKYHkLMXDgQAK37qFxrleYssNyxyRltl1nd7H75G7Or9nG6vbtwcHB6EjpIgWXECJL2Nna\n8n637pQ9W1a6UjKQ56Kx3Io4zxeT5Re3pShQoAB9ChYk109XWHJoCWHXwoyOZHaSdTJDNg3B/XAN\nXrW15eUxY4yOlG5ScAkhsoyHhwcXdl9gSdASGTCcAWLvxuLvtZ73natTtGhRo+OIVBgyeTLegfvo\nVGEAH236yOg4ZmfpoaXYJdqxYfMBxn/2mcW3boEUXEKILFSpUiWqV60O4bAjYofRcSze8qDlREQo\nhi6S6TYsTZEuXahdsyamFYcIuxaGX7if0ZHMxu3424zdNpZXLr9CYxcXnAcMMDpShpCCSwiRpbp1\n60a+4/lkTq500lozc/FMKleuTPXq1Y2OI1JLKRb37k3o1q30rTmOoZuGcjfprtGpzMK0XdNwLemK\n10Ivxo8fb3ScDCMFlxAiS3Xs2JHTgWH8tnuz3KGVDr+d+42rO64yeshoo6OINKrevTt1qlRhmacP\n5QuWZ/4f842OZLjTN07zXdB3XD9UgnL16lHTih7CLgWXECJL5c+fnyK1a1PkxPP8dPQno+NYrGnr\npmF7zZYOHToYHUWklY0Ni4YN47CfH31rj2XarmlExUYZncpQwwOG0716fzZ9u5gxEycaHSdDScEl\nhMhyvTt3JuTwOZYGyZxcaRF5O5KAXyNo+153HKxgMHF2VqdLF5xbtWLx4jV41PJg7LaxRkcyjCnC\nRNClIHZ7n6JS06a0rVvX6EgZSgouIUSWGz5wIHf/+ou7u25w4uoJo+NYnO/XTCd5z0X69ulndBSR\nAdZ8+im7ly9ncO3BbAzbyIHIA0ZHynJJyUkM8R9Cv6oj2Lt+A8s7djQ6UoZTWmujMwCglNLmkkUI\nkfmav/02pyJO0nlqc6a/Nt3oOBbjbtJdXFo1IPyK5srBQ0bHERnEw8ODKlWqUKxlMZYFL2NXj13Z\naiLbhQcWsuLwCq6tykHO+Hj27zD/u5iVUmitn/l/krRwCSEMMXPsWCJCQln1+0qSkpOMjmMxNgR5\nERx+g34ffGh0FJGBRo8ezZw5c3ir4lvcSbiDV4iX0ZGyzK24W4zfPp5hxXtxcvduVk2xzkl809zC\npZQqBHgDZYEIoJPW+uZj9vseaA1Ea62dn3I8aeESIpt5uVUrLhUN55tx82j1fCuj41iEzl3rs2Zz\nBHciI3GwszM6jshAHp06Ub5IEV4b2YV31r5DaP9QcjvkNjpWpvt488fcjLvJ7RG/UStvXj7Zt8/o\nSM8kK1u4RgEBWuvKwNaU5cdZArRMx3mEEFZqwqBBsD9B5uR6RocvB3PgUAzN2rSRYssKTW7cmPkL\nF1LRrgKNyzRm5u6ZRkfKdOHXwll6aCmdbFuyKzycwV7W27KXnhauUKCp1jpKKVUcMGmtqz5h33LA\nRmnhEkI8LCkpifIVynOj3Q3OfnaWQjkLGR3JrPX5uiMrh6zl+MVLlCxWzOg4IqMlJ9PO2ZnQ8uXZ\n+uMCan9bmwO9DlCuQDmjk2Wa171ep1HpRmyZtIX2zZvTd/hwoyM9s6xs4XpOa31/wpAo4Ll0HEsI\nkQ3Z2trS68NeFA8rzsojK42OY9Zuxt1k2Xof2r77thRb1srGhrkzZnBy504OhEYyuMFghgdYTgGS\nWltObyEkOgTnP505c+YMHwwZYnSkTPXUNmmlVABQ/DGbHnlst9ZaK6XS3Tw18aFJzlxcXHBxcUnv\nIYUQZq5Hjx7MnDWTRYGLGFDfOp6ZlhkW7VuE2q8YteNJozeENSjXti1tX3qJPsOGcWZ7ANW+rsb2\nM9txLe9qdLQMlZicyBD/Icxs9hmjuo1l6tSp2NvbGx3rqUwmEyaTKc0/n94uRRet9WWlVAlgu3Qp\nCiHSol3TplzXoXj+HEDN56znUR4ZJVkn49TNiVKnSrF/936j44hMdjM4mCKurnz1ww8UL/8Xk3dM\nJqh3EHY21jNuz3OfJ6uPrabCjXZ4L1rErcOHsbGxrIkTsrJLcQPQLeV1N2BdOo4lhMjG6r75JofO\n52LZ/sVGRzFLAacCuL3zNuNGjDM6isgCBWrVoufQoXwycSJvVGlPoZyF+O7Ad0bHyjA3/rrBpB2T\nGPHSRJZMncr8+fMtrthKi/ROC/ETUIaHpoVQSjkB32mtW6fstwpoChQGooHxWut/3ZIkLVxCZF/x\nCQnkdnLiFbeCbF0egoOtPK7mYbU/bcupBX9w82wktra2RscRWSAhIQGnqlX5fv58ytRxwm2FG8f7\nH7eKG0uG+A8hLjGOwNlhxDs6ctzPz+hIaZJlLVxa6+ta69e01pW11m735+DSWkfeL7ZSlrtorZ20\n1o5a69KPK7aEENmbo709rzZvTmR4bnzCfIyOY1Yirp8myvs0Ll3ek2IrG7G3t+fbzz9n/Cef4FzM\nmTervslE00SjY6Vb6NVQVhxewcuXnDkcFMS6OXOMjpRlrL8NTwhhEWaMGcOZ4yfx3jDP6ChmZfEX\nI4g6c475I2WwfHbTvn17HB0dWbVqFVOaTcErxIuQ6BCjY6XLsM3D+KTBcIbO/JLOXbpQpXJloyNl\nGSm4hBBmoUGNGhStXJlTW29y5sYZo+OYhbjEONZtPUnVVxtTtkgRo+OILKaUYtasWYwcORL7u3aM\nazKOIf5DsNThN/4n/Qm/Fk5hzwOUun6dZXPnGh0pS0nBJYQwG2OGDyc8LJKFBxYaHcUsrP75c44e\nP8Psydb5bDnx3xo3bkyrcuUY89pr9HmpD5diL7H+xHqjY6VaQlICQzcN5fOyA/nk59V8t3AhDmY+\nDURGS/Og+Ywmg+aFEFprqr5QleiG0UR9E5WtB89rrWnlWos/btlw/eAho+MIA90IDaVyzZpM+eYb\nKrmUobdPb472O0oOuxxGR3tmc/fOxSfMhxp7KhFz6hSLN20yOlK6ZeW0EEIIkaGUUowcNhK1V7E+\n1PK+xWekLacCCNt7lMUff2x0FGGwglWr0rNvXwZNnEi9og1wLubMl3u+NDrWM7t25xpTd06lf9n+\n/Lh6NTNWrDA6kiGkhUsIYVbi4uIoXqo4VYdVJXB0oNFxDKG15oWPXiDGJ4ZzYedQ6pm/RAsrpZOT\nKensTOUXXmDxgmk0WNSAw30P45TXyeho/2nArwNISkomePohunbtSr9+/YyOlCGkhUsIYdFy5MjB\nwP4DCV4XzMnrJ42OY4gtp7dwZuMZpk+YLsWWAEDZ2LDO05Odfn6c+P04ver2YtQW879z9Wj0UbyP\nenN0iw1RNjb06dPH6EiGkRYuIYTZiY6OpnyZUvSb8x6f985es89rrXEe7cyVH69w8cxF7Oys53Eu\nIv069O3L1qAgIkybqbGgOms6reHlUi8bHeuxtNa0WNGCwnfq4P3JIvbu20e9ChWMjpVhpIVLCGHx\nihUrRql69djgdZD4xHij42Spzac2E+ETwaQxk6TYEv+yau5ciInhu5WrmfG/GQzyG0SyTjY61mP5\nhPlwN/QsayYsY9isWVZVbKWFFFxCCLM0Y/hITgZH4BX4g9FRsozWmv4rZ5Ac7UiPHj2MjiPMkIO9\nPQHLl/PFmDG4lXDD1saW5cHLjY71L3eT7jLS7yMu+OemnLMzn8vnWboUhRDmq0D1apQqnZsi8W1w\ndExk0CA3WrduYnSsTOHru5Nxyz05FxpC/oJlmTt8pNW+V5F+o0aNIigoiEmLJtHh5w6EDggln2M+\no2M9MHvPbHZOWcWvhy5y9uhRShQsaHSkDCddikIIq+Dru5Nyuapz8ng0O450ZPPmqQwevAlf351G\nR8twvr47GTTYnyvJEVyPOM/5wO+t9r2KjDF16lSSk5NZP3cNLSq68enOT42O9MDFmItsXDCR3zcF\n4f3tt1ZZbKWFFFxCCLM0d+5mgg94k5yURNnaEwE4depT5s0LMDZYJpg7dzOnbRpx90wcZcq7khD/\nnNW+V5Ex7Ozs8PLywuubb3D1imXxwcWEXws3OhaBFwJp4lmfRO+/GNKlC+3btjU6ktmQgksIYZbi\n4+0AO+xzv8vZiwGQ7wIAcXG2xgbLBHHxtuSvP4UrJ05x69S0v9db4XsVGadIkSKs9fam3zofPP5s\nx0ebPzIsi6/vTmq815pXPZtxe2lu/ixZjpHLlhmWxxxJwSWEMEuOjokA3AmfBKeToEUzcIwhR44k\ng5NlvJjixyl8+DplKzTgZmz1B+ut8b2KjPVi69a8N2QI8xf6EXU0DP+T/lmeYd3GLbyzcjDH8p8m\ncU4/rpxw4IZui7//7izPYs6k4BJCmKVBg9yoWHEMkBMSJsHORHJ2r0Wf/i5GR8tQF2IucL6iicsh\np7A9PPbB+ooVP2HgwOYGJhOWYsGMGVSvVYsz6+0ZuLY3u87uyrJzR96OpPuO7sQklYct3eGvNcAm\nIiJmS5f4P8gkL0IIs3T/Dr1588bx55/J7Nt3lYp5KvJL0nJe182sYgb2W3G3cP/RnZrnX8CxRV60\n3kqpOBM5ciQxcGBLuUtRPLPf162jRO3a2G9/HrebbchzqxRVz7syqnenNH+OkpKTiP4zmsjbkVy8\nfZGLMRe5ePviI8vnbp2j8JV63DraGqJHQtHtEF0SkC7xf5JpIYQQFsHT05P1G9dzq8MtmldozpRm\nU4yOlC7xifG0+rEVZW3KsmHIBoKCgihbtqzRsYQFm7/0ZwaOGgyF3oIiReGVryhwugoL3xtHx3at\nHtn3dvztfxVSF2MuEhkb+WA5KjaKgjkLUjJvSUrmK4lTHidK5iv593JeJ0rnK81r7T8iaM86yLsR\nLjZ+cI4WLcbh72/Z/06fJrXTQkgLlxDCIvTs2ZNpY8cy/u3P+CxkBmXyl+HDuh8aHStNknUyPdb3\noECOAtgG2NKzZ08ptkS6bVwVDH/9AeGvw/E6cPAgN/83ga57O7Eh6Q0u3b70oLBK0kmPFE4l85ak\nUqFKNC3X9MFyibwlcLB1eOo5v/fz42DgLzQv04iAkL+LrXtd4i0z+y1bFCm4hBAWwdHRkZatWjFi\n1hz2mfxosrQJJfOVxP15d6OjpdqoLaM4d+scs1/8ijYD3QkLCzM6krAC8fF2EFMK2AF0gti+sP4n\nXkh0oNlbLz/SOpXfMX+6u+XnL1/OoMGDGdemDfW79MJm3jji4mylS/wJpOASQliM+V98wfJatVi7\n4Td+6fwL7Va1w+9dP+o61TU62jObu3cuG05s4Lcev1GtZXs6Dh5MgQIFjI4lrMD9O3shD7Ae6IUd\njSkc+yo96mTco3Xi4uIY0aMHG1avZkXnzryzbBnY2kqB9R/kLkUhhMXIUbw4fdzcmPLZZ9R3qs/C\ntgtp59WOMzfOGB3tmaw5toaZu2fi39WfRRsCuHXiBJ+9+67RsYSV+PvOXgB74HtaOEZxMvA7/Pfu\nzZBzHDt2jAZ163L5l1846OnJOytWgK0Mjn8WMmheCGFREi5fJt+LL/LxxIlM6dWL+X/MZ/4f89n9\n/m4K5ypsdLwn2nV2F2/+9Cabu26mRlFn8lerxoCaNfl8zRqjowkr4uu7k3nzAh507Q3q68reLycz\nJSiIkg0bMmvaNDq/+GKqj5ucnMyiRYsYM2YM06dPp2f79qjC5vvvLSukdtB8mgsupVQhwBsoC0QA\nnbTWN/+xT2lgOVAM0MBCrfXcJxxPCi4hxDMZ7eHBnMBAYkNDsbGxYUTACH4//ztbPLaQwy6H0fH+\n5fiV47gsc+GH9j/gVtGND774gp/mzuXmjh3YlCtndDxh7bTm8qpVfDB7Nn7h4RR5+WWmTp7MB/Xr\nP3Uc152kJFaGhPCttzchGzdSycaGn729qVq1ahaGN19ZWXB9BlzVWn+mlBoJFNRaj/rHPsWB4lrr\nQ0qpPMAB4A2t9fHHHE8KLiHEM0m6e5cX6tZl5qef0q5dO5J1Mu+ufZfE5ES83/LGRpnPaInI25E0\nXNyQya6T8ajlwdU//6R4xYp4tmxJr6VLjY4nspPERK4tWECfvXv5xdeXl11d6eTiQtGiRR/543fp\nEnO8vDjm54c6d466rVrRr0sX3nZ3x8Hh6XctZidZWXCFAk211lEphZVJa/3UslcptQ6Yp7Xe+pht\nUnAJIZ7Zr7/+ypAhQzh8+DA5cuQgPjEetxVu1C1Rl9ktZhsdD4CY+BiaLGlCpxqd+OTVTwCYMno0\na5ct42BwMBQtanBCkV3dunWLJUuXcurkSa5cuXLvz7FjRMfGEmtjQ10nJ/rHxPBGu3Y4zJsHdnKP\n3T9lZcF1Q2tdMOW1Aq7fX37C/uW4d69qDa117GO2S8ElhEiVt956i2rVqjFlyr3JFW/8dYNG3zei\nd93eDH55sKHZ7ibdxf1Hd54v9DyerT1RSnH9+nWqVqnCjmnTqPahZc4hJqzYokUQFgaJifD669C4\nsQyIf4oMLbiUUgFA8cdsGgMse7jAUkpd11oXesJx8gAmYKrWet0T9tETJkx4sOzi4oKLi8szvAUh\nRHYVGRlJrVq1MJlM1KhRA4CzN8/S6PtGzGk5hw7VOxiSS2uNxzoPYuJjWNtpLbY2935pde/enbx5\n8zJv3jxDcgkh0s5kMmEymR4sT5o0KUu7FF201peVUiWA7Y/rUlRK2QM+gJ/W+qunHE9auIQQqfbN\nN9/ww5Il7Pz9d2xTvo0fvHQQtxVurOu8jkZlGmV5ptFbRrM9Yjvbum0jl30uADZs2MDQoUMJDg4m\nT548WZ5JCJGxUtvClZ6RpRuAbimvuwH/arlK6WpcDBx7WrElhBBp1evdd7kUEcEbY8c+WFenRB1W\ntF9Bh586cOLqiSzN47nPkzXH1+Dzjs+DYuvatWv07duXJUuWSLElRDaVnoJrBtBcKRUGNEtZRinl\npJTyTdmnEdAVcFVKHUz5Iw9XEkJkGJu8eVnYqxe+np7sCA9/sL5FpRZM/9903Fe6ExUblSVZ1oWu\nY+rOqfh39adIriIAJGtNAw8P2rz1Fk2ayEzcQmRXMvGpEMLyJSXRsmFD9tvbE71rFzYPzS000TQR\n33BfTN1M5HbInWkR9pzfQzuve48aesnppQfre3z9Nd6ffkr0/v3kcXLKtPMLIbJWVnYpCiGEebC1\nZc28efwZHk7/r79+ZNOEphNwLuZM59WdSUxOfMIB0ifsWhjtvduz/I3ljxRbv50+zbJx41hev74U\nW0Jkc9LCJYSwGkt69uSDjRsJCQmhWrFiD9YnJCXQZlUbyhcoz4LWC546u3ZqXY69TMPFDRnbZCzv\n13n/oXMm8ZyLC/Vu3mSTyQTZ/DEoQlgbaeESQmRbPebPp+VrrzHs/fdJSkp6sN7e1p6fO/5M4IVA\npv82PcPOF3s3ljYr29CtVrdHii2Ad2bPJuHUKdZPmSLFlhBCWriEENYlISGBli1b4uzszFdfPXpz\n9P3H7ExtNpWuNbum7zxJCbTzakepvKVY2HbhI61mwcHBvNqoEV5t2uDu5ZWu8wghzJO0cAkhsjV7\ne3tWr17Npk2b+Pof47mc8jrh+44vwzYPY9uZbWk+h9aa3j69sVW2LGjzaBflyZMncXd35/sPP8R9\n4cI0n0MIYV2khUsIYZVOnz5No0aNWLx4Me7u7o9sM0WY6PRzJ7Z6bMX5OedUH3vC9gn4nfRje7ft\nj9z5GBkZSePGjRk1ahS9evVK93sQQpgvaeESQgigQoUKrFm9mm5vv82PW7c+ss2lnAtzWs6h9crW\nXIi5kKrjLjywkB+P/IjPOz6PFFs3btygRYsWfPjhh1JsCSH+RQouIYTVatiwIYP/9z88unRhe2jo\nI9u6OHdhQP0BuP/ozq24W890PJ8wHyaYJuDf1Z9iuf++C/L09eu4tGqFm5sbo0aNytD3IISwDtKl\nKISwbomJdHRxYf2lS/y0ejVv1KnzYJPWmoF+Awm9Gsqv7/6Kg63DEw+z7cw23l79Nj7v+FC/ZP0H\n66/cuUNFFxeeL1KEfT4+2NjI91ghsgPpUhRCiIfZ2fGTvz9dSpbkTVdXus2YQXLKlzulFHNaziGP\nQx4+2PABT/rSt/LISrqs6cLPHX9+pNjaEBJC+caNKXLjBoHTp0uxJYR4Irk6CCGsnsqTh2U7drB2\n8GC8v/ySZu+9x507dwCwtbFlZYeVhF0LY9z2cY/8nNaaz3d/zqgto9jmsY2m5ZoC8FdiIi2GD+eN\nRo1ooxSh3t7Y16qV5e9LCGE5pEtRCJGtXLtxg/79+hFy5AheXl688MILAFz58woNv29Iy7yvE7Yq\nB3HxNkRU9UNVuMZvfXZSKl8pAMLDw+ns7s7Z+HjW9ejBq+PHg62tkW9JCGEA6VIUQoinKFywIKtW\nruTjjz/G1dWVkSNHYjKZyGeThxFOk1hw/Bs2n6/DzueOcS4uD2pJe4J3nebKlSvMmjWLV155BQ93\nd6IDA3l10iQptoQQz0RauIQQ2VZ4eDhLly5ly+bNHD50CIfnq1AwT2WoehibqErcvlKYqydvkSNh\nF46OChcXF7744gsqVqxodHQhhMFS28IlBZcQQgDnjx/nnZ6jOH37Jlevn8YxV17y5yzHhcudqV/+\nAHv2fCWD4oUQD6S24LLLzDBCCGEpSlerRq68zkTumQrAXeB2yraCL56UYksIkS5yBRFCiBSDBrlR\nseKYR9ZVrPgJAwc2NyiREMJaSJeiEEI8xNd3J/PmBRAXZ0uOHEkMHNic1q2bGB1LCGFmZAyXEEII\nIUQmk2khhBBCCCHMjBRcQgghhBCZTAouIYQQQohMluaCSylVSCkVoJQKU0ptVkoVeMw+OZRSe5VS\nh5RSx5RS09MXVwghhBDC8qSnhWsUEKC1rgxsTVl+hNY6DnDVWtcGagKuSqnG6TinEJhMJqMjCAsi\nnxfxrOSzIjJTegqudsCylNfLgDcet5PW+k7KSwfAFriejnMKIRdFkSryeRHPSj4rIjOlp+B6Tmsd\nlfI6CnjucTsppWyUUodS9tmutT6WjnMKIYQQQlicpz7aRykVABR/zKZHpmLWWmul1GMn0dJaJwO1\nlVL5gU1KKRettSmNeYUQQgghLE6aJz5VSoUCLlrry0qpEtxrvar6Hz8zDvhLaz3rMdtk1lMhhBBC\nWIysenj1BqAbMDPl73X/3EEpVQRI1FrfVErlBJoDkx53sNSEFkIIIYSwJOlp4SoE/ASUASKATimF\nlRPwnda6tVKqJrCUe2PFbIAftNafZ0RwIYQQQghLYTbPUhRCCCGEsFaGzzSvlGqplApVSoUrpUYa\nnUeYN6VUhFLqsFLqoFLqD6PzCPOhlPpeKRWllDry0Lr/nKBZZE9P+LxMVEpdSLm+HFRKtTQyozAP\nSqnSSqntSqmjSqkQpdSglPWpur4YWnAppWyB+UBLoDrQRSlVzchMwuxp7t2sUUdrXd/oMMKsLOHe\nteRh/zlBs8i2Hvd50cDslOtLHa21vwG5hPlJAIZqrWsALwP9U2qVVF1fjG7hqg+c1FpHaK0TUehE\n1QAAAexJREFUAC/gdYMzCfMnN1iIf9Fa7wJu/GP1M03QLLKfJ3xeQK4v4h+01pe11odSXscCx4GS\npPL6YnTBVRI4/9DyhZR1QjyJBrYopfYrpT40Oowwe880QbMQDxmolApWSi2WLmjxT0qpckAdYC+p\nvL4YXXDJiH2RWo201nWAVtxr1n3V6EDCMuh7dwjJNUc8zQKgPFAbuAR8YWwcYU6UUnmANcBgrfXt\nh7c9y/XF6ILrIlD6oeXS3GvlEuKxtNaXUv6+AvzCvW5pIZ4kSilVHCBlguZog/MIM6a1jtYpgEXI\n9UWkUErZc6/Y+kFrfX/e0VRdX4wuuPYDzyulyimlHIDO3JtQVYh/UUrlUkrlTXmdG3ADjjz9p0Q2\nd3+CZnjCBM1C3JfyS/O+9sj1RQBKKQUsBo5prb96aFOqri+Gz8OllGoFfAXYAou11tMNDSTMllKq\nPPdateDeUxJ+lM+LuE8ptQpoChTh3niK8cB6HjNBs1EZhfl4zOdlAuDCve5EDZwBej80RkdkU0qp\nxsBO4DB/dxuOBv4gFdcXwwsuIYQQQghrZ3SXohBCCCGE1ZOCSwghhBAik0nBJYQQQgiRyaTgEkII\nIYTIZFJwCSGEEEJkMim4hBBCCCEymRRcQgghhBCZTAouIYQQQohM9n/XjwQf8XFARAAAAABJRU5E\nrkJggg==\n",
       "text": [
        "<matplotlib.figure.Figure at 0x10e942710>"
       ]
      }
     ],
     "prompt_number": 65
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## Sparse matrices\n",
      "\n",
      "Sparse matrices are often useful in numerical analysis dealing with large systems.\n",
      "Whenever the problem involves vectors or matrices with mostly zero values then one should take advantage of the sparsity.\n",
      "SciPy has an extensive library for sparse matrices, with basic linear algebra operations (such as equation solving, eigenvalue calculations, etc)."
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "To create a sparse matrix we have to choose the format that will be stored in:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from scipy.sparse import *"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 66
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# dense matrix\n",
      "M = array([[10,0,0,0], [0,.3,0,0], [0,2,1.5,0], [-3,0,0,100]])\n",
      "print M"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[[  10.     0.     0.     0. ]\n",
        " [   0.     0.3    0.     0. ]\n",
        " [   0.     2.     1.5    0. ]\n",
        " [  -3.     0.     0.   100. ]]\n"
       ]
      }
     ],
     "prompt_number": 67
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# convert from dense to sparse\n",
      "M_sparse = csr_matrix(M); # or csc_matrix(M)\n",
      "print M_sparse "
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "  (0, 0)\t10.0\n",
        "  (1, 1)\t0.3\n",
        "  (2, 1)\t2.0\n",
        "  (2, 2)\t1.5\n",
        "  (3, 0)\t-3.0\n",
        "  (3, 3)\t100.0\n"
       ]
      }
     ],
     "prompt_number": 75
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We can go back to the dense form:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": true,
     "input": [
      "print M_sparse\n",
      "print M_sparse.todense()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "  (0, 0)\t10.0\n",
        "  (1, 1)\t0.3\n",
        "  (2, 1)\t2.0\n",
        "  (2, 2)\t1.5\n",
        "  (3, 0)\t-3.0\n",
        "  (3, 3)\t100.0\n",
        "[[  10.     0.     0.     0. ]\n",
        " [   0.     0.3    0.     0. ]\n",
        " [   0.     2.     1.5    0. ]\n",
        " [  -3.     0.     0.   100. ]]\n"
       ]
      }
     ],
     "prompt_number": 69
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "The easiest way to create a sparse matrix is using the `lil_matrix` command. It creates an empty matrix and then we allogating non-zero entries by using matrix indexing. This way we avoid creating large dense matrices."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = lil_matrix((4,4)) # empty 4x4 sparse matrix using a linked list format\n",
      "A[0,0] = 1\n",
      "A[1,1] = 3\n",
      "A[2,2] = A[2,1] = 1\n",
      "A[3,3] = A[3,0] = 1\n",
      "A"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 70,
       "text": [
        "<4x4 sparse matrix of type '<type 'numpy.float64'>'\n",
        "\twith 6 stored elements in LInked List format>"
       ]
      }
     ],
     "prompt_number": 70
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We can also change the format based on the sparsity:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = csr_matrix(A); A # if the rows are sparse"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 71,
       "text": [
        "<4x4 sparse matrix of type '<type 'numpy.float64'>'\n",
        "\twith 6 stored elements in Compressed Sparse Row format>"
       ]
      }
     ],
     "prompt_number": 71
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = csc_matrix(A); A # if the columns are sparse"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 72,
       "text": [
        "<4x4 sparse matrix of type '<type 'numpy.float64'>'\n",
        "\twith 6 stored elements in Compressed Sparse Column format>"
       ]
      }
     ],
     "prompt_number": 72
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "---\n",
      "### Exercise 12\n",
      "\n",
      "Construct a $10000 \\times 10000$ `lil_matrix` $A$ and add values to it by setting all rows and columns $2000$ to $3000$ into random values. Also set random diagonal entries (you can do this without a loop!). Then use the `spy` function from `matplotlib.pyplot` to \"spy\" the matrix.\n",
      "\n",
      "Using a vector $b$ with $10000$ random values solve the system $A x = b$. \n",
      "To do so use the direcet solver `linsolve` from `scipy.sparse.linalg.dsolve`.\n",
      "Change to different sparse formats (CSR and CSC) and measure the time needed to solve the system. You can use the `time` package."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# type your solution here"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 73
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "##### Solution"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "#%load solutions/sparse_matrices.py"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 74
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## References:\n",
      "\n",
      "#### NumPy:\n",
      "\n",
      "- [NumPy tutorial](http://www.scipy.org/Tentative_NumPy_Tutorial) \n",
      "- [NumPy tutorial video](http://showmedo.com/videotutorials/series?name=i9KuJuNcG) \n",
      "- [NumPy: lock 'n load](http://mentat.za.net/numpy/intro/intro.html)\n",
      "- [NumPy for Matlab users](http://www.scipy.org/NumPy_for_Matlab_Users)\n",
      "\n",
      "\n",
      "#### SciPy\n",
      "- [SciPy getting started](http://www.scipy.org/getting-started.html)\n",
      "- [SciPy tutorial](https://docs.scipy.org/doc/scipy-0.14.0/reference/tutorial/index.html)\n",
      "\n",
      "#### More general\n",
      "\n",
      "- [General Python tutorial](http://docs.python.org/tutorial/)\n",
      "- [Dive Into Python](http://diveintopython.org/)\n",
      "- [Hans Petter Langtangen book](http://link.springer.com/book/10.1007/978-3-642-02475-7/page/1)"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Other packages\n",
      "Besides NumPy and Scipy, there are many other useful Python packages for scientific computing. Here is a short list:  \n",
      "\n",
      "- [matplotlib](http://www.matplotlib.org/) - plotting and visualization\n",
      "- [sympy](http://sympy.org/) - symbolic computation\n",
      "- [mpi4py](http://mpi4py.scipy.org/) - parallel computing\n",
      "- [petsc4py](http://code.google.com/p/petsc4py/), [pytrilinos](http://trilinos.sandia.gov/packages/pytrilinos/) - Python bindings for the \"big 2\" parallel scientific libraries\n",
      "- [pyCUDA](http://mathema.tician.de/software/pycuda), [pyOpenCL](http://mathema.tician.de/software/pyopencl) - GPGPU computing\n",
      "- [FENiCS](http://fenicsproject.org/), [FiPy](http://www.ctcms.nist.gov/fipy/), [PyClaw](http://numerics.kaust.edu.sa/pyclaw/) - solve complicated PDEs with very sophisticated numerical methods\n",
      "- [networkX](http://networkx.github.com/), [pygraphviz](http://networkx.lanl.gov/pygraphviz/) - graphs"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      " "
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      " "
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "*Copyright 2015, Maruan Al-Shedivat, Yiannis Hadjimichael, ACM Student Chapter.*"
     ]
    }
   ],
   "metadata": {}
  }
 ]
}