{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "# Finite Differences for Boundary Value Problems"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 1,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "import numpy as np\n",
        "import matplotlib.pyplot as pt\n",
        "\n",
        "import scipy.sparse as sps"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "We'll solve\n",
        "\n",
        "$u''+1000(1+x^2)u=0$ on $(-1,1)$\n",
        "\n",
        "with $u(-1)=3$ and $u(1)=-3$."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 2,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "#n = 9\n",
        "n = 200\n",
        "\n",
        "mesh = np.linspace(-1, 1, n)\n",
        "h = mesh[1] - mesh[0]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 3,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "A = sps.diags(\n",
        "    [1,-2,1],\n",
        "    offsets=[-1,0,1], \n",
        "    shape=(n, n))\n",
        "\n",
        "if n < 10:\n",
        "    print(A.todense())"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Create `second_deriv` as a matrix to apply the second derivative. Can only do that for the interior points!\n",
        "\n",
        "* change `shape` and offsets\n",
        "* Take `h` into account"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 4,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "second_deriv = sps.diags(\n",
        "    [1,-2,1],\n",
        "    offsets=np.array([-1,0,1])+1,\n",
        "    shape=(n-2, n))/h**2\n",
        "\n",
        "\n",
        "if n < 10:\n",
        "    print(second_deriv.todense())"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 5,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "factor = sps.diags(\n",
        "    [1000*(1 + mesh[1:]**2)],\n",
        "    offsets=[1],\n",
        "    shape=(n-2, n))\n",
        "\n",
        "if n < 10:\n",
        "    print(mesh[1:-1])\n",
        "    print()\n",
        "    print(factor.todense())"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 6,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "A_int = second_deriv+factor\n",
        "\n",
        "if n < 10:\n",
        "    print(A_int.todense())"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 7,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "A = sps.vstack([\n",
        "    sps.coo_matrix(([1], ([0],[0])), shape=(1, n)),\n",
        "    A_int,\n",
        "    sps.coo_matrix(([1], ([0],[n-1])), shape=(1, n)),\n",
        "    ])\n",
        "A = sps.csr_matrix(A)\n",
        "\n",
        "if n < 10:\n",
        "    print(A.todense())"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Next, assemble the right-hand side as `rhs`:\n",
        "\n",
        "Pay special attention to the boundary conditions. What entries of `rhs` do they correspond to?"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 8,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "rhs = np.zeros(n)\n",
        "rhs[0] = 3\n",
        "rhs[-1] = -3"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "To wrap up, solve and plot:"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 9,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "import scipy.sparse.linalg as sla\n",
        "\n",
        "sol = sla.spsolve(A, rhs)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 10,
      "metadata": {
        "collapsed": false
      },
      "outputs": [
        {
          "name": "stderr",
          "output_type": "stream",
          "text": [
            "/usr/local/lib/python3.5/dist-packages/IPython/core/formatters.py:92: DeprecationWarning: DisplayFormatter._ipython_display_formatter_default is deprecated: use @default decorator instead.\n",
            "  def _ipython_display_formatter_default(self):\n",
            "/usr/local/lib/python3.5/dist-packages/IPython/core/formatters.py:669: DeprecationWarning: PlainTextFormatter._singleton_printers_default is deprecated: use @default decorator instead.\n",
            "  def _singleton_printers_default(self):\n"
          ]
        },
        {
          "data": {
            "text/plain": [
              "[<matplotlib.lines.Line2D at 0x7f2c40aafa20>]"
            ]
          },
          "execution_count": 10,
          "metadata": {},
          "output_type": "execute_result"
        },
        {
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAEACAYAAABI5zaHAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztfWmQXcd13tcz2AbAACAAgiABbgJBQqQILjYpipbjiRaL\nUiTBikmX5KS02C67XJGjUuSdtgnIsezYtOwklisVW5KXisySGMkSHWuhlpFMU6JYIkgCFBeIoAAQ\nIAb7OoNtpvPjTHMu3tyl7+1zTt/3pr+qKQxm3rzTr2/fr7/+zum+xlqLhISEhITeQl/sBiQkJCQk\n8CORe0JCQkIPIpF7QkJCQg8ikXtCQkJCDyKRe0JCQkIPIpF7QkJCQg+CjdyNMX3GmMeMMV/ges+E\nhISEhGbgVO4fAPB9xvdLSEhISGgIFnI3xqwG8BYAf83xfgkJCQkJYeBS7n8G4NcApO2uCQkJCS1A\nMLkbY/4dgBFr7eMAzORXQkJCQkJEmNCzZYwxHwHwHwGcAzAAYBDAZ6217+54XVL1CQkJCQ1gra0t\nmoOVu7X2t621l1lrXwHgnQC+3knsmdemL6ave+65J3obeuUr9WXqzzZ/NUWqc09ISEjoQczifDNr\n7TcBfJPzPRMSEhIS6iMp9y7F0NBQ7Cb0DFJf8iL1ZzsQnFD1DmSM1YqV0Du4917giiuAO++ME/+B\nB4AXXwR+/ueBOXPitAEAJiaAviTFZiSMMbANEqqJ3BNai927gauuAjZsAO67L04b3vY24Ac/AM6e\nBTZvBgYH9dtw/Dhw443Av/4rsHKlfvyEuGhK7kkLJJRizx5g3744sX//94E77gC+85048QFg61ZS\n70uWAE89FacNf//3wPbtwEMPxYmf0J1I5J5QiNFR4PWvB/7wD/Vj79wJ3H8/8Fd/Rcr1pZf023D8\nOE1sV14JXHMN8Nxz+m2wFvjYx4A3vSmRe0I9JHJPKMSHPkTk8sQT+rEfewy4/XZg+XLg1a8GHnlE\nvw3f/z7wylcC/f3A2rXAtm36bfjmZO3Z7/wO2TIJCb5I5J6Qi6eeAj7/efp64gkieU388IfA5ZfT\n97fdFsea2boVeNWr6Purr46j3D/1KUrm3nIL8PTTwIkTuvE/+EHg2muB978fGB/XjZ0QhkTuXYDR\nUeBnfoYITwvPPEOEcs01wKxZ5L1rYseOdpH72rVxyH37dmrD3LnATTfpr2Aeegj49V8HPvMZYNcu\n3dgJYUjk3nJMTADvfjfw2c/qEtyOHVSCCADr1+tbM9n4t94KfO97wLlzum3YsuV8ct+2TX8Fs2sX\ncOml9P1rX6vru1tLn/mtbwXWraOJRgvf+Abwp3+qF68Xkci95fj4x0k1f/CDpKa1kFXON9wQh9xd\n/CVLgBUrqCRRE1nlvmQJsGCBbmLX2vPJ/fbbgYcf1ot/4ABgDLBsGfCKVwDPP68X+xOfoCqhhOZI\n5N4A990HbNyoE+vxx4F3vYsIVpPcs573DTcATz6pF7szPgCsWgXs3asXf/9+4NQpiuugnVQ9eJDs\nmIUL6f9XXqlrjfzgB/SZjQHWrNEj93PngC9+kcb76KhOzF5EIvcG+PjHgS9/WSfWzp2k3Nat01fu\nsWyZEyeAsTHgwgunfnbxxbqq+emnKZFoMltHtJOqu3YBl1029f+LL9ad4LZtI3IHiNy1bJmHH6bP\nff31tHFMGkeOACMj8nG0kci9Jg4coKTW1q3kh0vD3eCuzlojJnC+LfLKVwIvvECEqxX7ssvOJ9aV\nK3WJ7aWXgNWrz/+ZtnLPWjIAcMEFNPGdOqUTv5PctZT7Aw/QzuBbbwUefVQ+3kc/SrZnryGRe018\n7nPAm98MLF2qo2R27iSiGxykmDt3ysc8epS22y9bRv+fM4c8Vy1iy04sDtrkvm8fcNFF5/8shnLP\nkntfH7VJS2Vu20bHPwBTnrtGQvmBB4C3v52qtb77Xfl4jz0GfP3r+slyafQMuX/rW8Bdd8nHuf9+\nOsRq/Xp5H9qpNEeyWtaMI9dO5ax1DEHWEsrG1yT3kRFK4mZx5ZXUNi10kjuga81klfuyZUR+hw7J\nxjx0iFZNN99M5K6h3DdvplVprOMlpNAz5P7d7xLxStYBnzhBfuBb3qKTZHQ3tyNZbXLPYsUKPXLv\nTKYC+p77vn3Tyf3CCynRqoU8cl+5UqcfXBmkI3eXVJVere7dS0lsY2i8j4zITij79lHS9q67gK99\nTS7Ok08CP/mTcu+fB1Vyl6wTdps9/uAPZN4foKNfL7mESuI0koydN7cWuf/wh9OVsya5t8WW6ST3\n5cuJ3LWW7y6ZnoVWP+zbB8yeTVagg4bvvnfvlB3W308K/nvfk4v3+ON04uYb3iBL7o88QtaP5g5j\nVXL/kR+hRKQEtm+n8sRHH5VT1CMjUwNPw5ZxfrtDbOWupVrbSu4DA0R4WjdoTFsmq9odNGrd9+49\n/1jjK64gUSWFzZtp5+/rXkfWrtRGua1b6fgGDZvJQZXcjx+XS8pt3w5cdx0lYqR28WXJfe1aGojH\nj8vEAqYrt7VrdTbyxLZlOic1gCyRw4cp0auBPHJ37dCY5MbHafNaZ8WOli2zc+f01ZuGLZO9xwD5\n/nbKfcUKutekVuNbtlBpp+YmtGByN8bMNcY8YozZbIzZYoy5p+i1/f0y5D4+PjUYV66UqybIDrz+\nfqqD3rJFJhYwvc7ZqWdpWyAvoalF7tbSZ+wk1v7+KVtEA0Xkvnw5lcNKY2SELJG5c8//udYK5sCB\n8/cZADThSyeUO5W79DV3yh2QXSVs3Qr84i92Gblba08D+LfW2psA3AjgzcaYW/Neu2GDDLm/+CIN\ngnnzZEvFOisoLr9cdsnYqWAHBugQr5Mn5WIC9Jk6FeOFF+qQ+4kTZH3Mmzf9d1qqdWwMOHMGWLRo\n+u+0lHueJQPo2TIHDkxVaTksXy5fLZOn3KUm05MnqZ/XrZONNTJCds9P/zTw7W/r7VVhsWWstW6T\n8FwAswDkasv3vU+G3LdvJz8QkCf37MBbtoy2iEsh7wZftkxWObpyt84bW0u558V20Ewmrlhxfimo\ng9bqYc8eSt53QmuCO3hw+nWQHu/AdOUuOZm6zXKzZ8vG2rqVLJmLL6Yzip59lj9GHljI3RjTZ4zZ\nDGAvgAettblpg/XrZTzj7dvJDwRkSahzY8vy5XKDvfPQKI2YAKlWY2iVkIUWueeRioOWai2yZABZ\nJZlF0STnxIu0NXfwII21LJYulSd3Tc+9c6xJxXJ+O6B7+BuXcp+YtGVWA3i1MebavNetXk2DlttW\niKncpW70AweA+fOp7DILaeV+6ND55W8Og4OUzJQ+guDgwfz4gL5yz4OWLVN0HQYG6EvaHsmzZQYH\nya46fVourqbn3jmBSeVTOh/68sIL/DHyMIvzzay1x4wxwwDuAPD9zt9/+MMbsXAhPb7tne8cwtDQ\nEEvc55+nsygAXXJfvlzuYKOiZbm0ci9SzsZMJXQ7K1k04gN002ts/y8j9+XLddpQRO7A1AqmqJ84\nkHcdjKE2HTpEbeDG+Pj0RK7kZNo5gUkq9/e9j75ftqx6J+zw8DCGh4eD4waTuzFmOYCz1tqjxpgB\nAG8A8Ed5r924cSMef5w2DDDxOoDzlfuiRaQwR0dJ+XLBWl3lfvgwHRTViVjKHZiyZiTJvcpz/9a3\n5GI7tMWW6SxHdXArmOuuk4ufZ8sAU9aMBLkfPEietPPAAbqfz5yhYzjykuyh8TTIffv2qTN6fKyt\noaHzhe+mTZsaxeWwZS4G8A1jzOMAHgHwZWvtPxe9WKJWO+u5GyOj3k+coIObsjaJpIouIndp5V5G\n7hoVM2333LUSqmXXQSOpmmfLALJJ1U7xBND9LGWXdI41iWtrLd3L7louWyZvqTlwlEJusdbebK29\n0Vq73lpbegDAVVfxVsycPk2nGGaXcitW8JN73kFSM1W5S6KM3LX87rwTIbXbUHYdLrpI9jqMjpJF\n0pnvAWTJqdNvd5CaUPOUO/e9dfz41M5mQCcp7aB+cBj3mdhHjhAJZsvWJAZ/nqqQVDFlyr3Xyb0o\n/pIldL2l0RZbpqgfLriAxocUnCWTVwoqOeaz58pkITWhdpL7woVUj85ZNNB5H3eVcq+LVat4VfWR\nI3TTZyFhy+SpucWLSeVIbIkvU+4xEqqADrmXee6O1KTLAPNWaQ6LF1O115kzsm0oI3fpSa5sDEgq\nz5GRfOWuRe7G8Mc6dOj8+7inlfvixbwDU4vc85R7Xx9dOImLNZOVexGxzJ1LO3SlyzHLlHtfn7w1\nBlQr91jkLm3LFCl3Dc8d4LeAsn47QAniU6fkxQEQgdyd6uBSXzHJHZBLcMZS7lXkLu03lxELQNda\n0pKwtroN0hOsqw4ZHMz/vXQfHDiQXykDyCdUY3ruAP9E0nkfG0P/17Bm1Ml93jz6gFzPgYxN7lIq\nbiYr96L4gLzffOoUqfOysjvppKpTe3meN9DbtozWSZzumI3OsSZtywB6vnuUJzFxDk6XUM1Cqlqm\niNw1lfv8+TQwR0en/44DZeQqbUeMjwPHjuV/bgdpYjtyhKzDMkiTe9kEC+jYMmXKXYqY8u5lQKa/\njx2jKpY5c2RjddoygJ7vHo3cjx7lea/Dh3WU+/79049ABeSUdBG5S9b9AuUJzcWL+a5bHo4cIU+y\nv7/4NdLKPW8l2ImYqydAx5Yp89yliOno0fyJVYLci6wn7ms745Q7Z1K1yJbhtg+KbnqpwV6kYiRj\nAuXEsmgRbeaSOrK0yusG5FXr0aPV5B4z7wHET6hqk7uE5170GSWUe+d93PPKXZLcly6lzQOcGelj\nx4oHHveFcrvaishdSjmOjVHszhMhHfr7yRaSevpUFakB8qrVR7lLr2Cq+mHBAtq8J1VxUWbLuLNl\nJMpRiywxiWoZTXLvvJY9rdylyb2vj9+bPXo0/+ENEj706CiV/HU+hScbU2Lmd6RSlMgDeC21Tvgq\n95lO7sbI5h7KbJl582QeGHP6NOVc8oTFsmV0zTlXjEVjjXuVkGfLJOXuiaKbcdEiUtscsJbUal5p\nmgTRlql2QPasjSrlzL1PoTN+FblrJFSryF1yggP8VjCS1kzVdZAY884OyxMWEivGMuXOXQqZlHtD\naJD76Cip6OxpdQ4SROtD7lLKvYpcJVVrUu6E2PZUWZ07MGXNcKLIb3dYvJjvfgbKlbtknTvQ48qd\n8+YoSjxyknuRJQPEUe5Sg8OXVKSILbZiBfzJXbINvtdBog3uuOyi8Q7IKfcycue8n4FicueOM+Oq\nZTgHZl4pJMB7kYqSqYCMiq4id24V4+BDKpLEFluxAt2j3KUmuWPH6N7pK2GGXib3uXPJGuJ42tT4\nOFWXdX6unlbuXORurY4tU6bcuY9TAPzIXYJcYit3H2LVsGWqNjG1wXOXmuQcuZchhi2jRe4A5dY4\nYh05Qu/VuW8jKXcPnDpFs2zeVnFu5V404F1VC2f1QBW5cw90h4MHy+MCssq9bIXkIJ1Q9alznwnK\nvQyLFvGXw1ZNqprkzvX5iu7jpNw9UKb0uJV71cDjHOyxyP3oUT9ylyK2shWSQxsSqvPmUVke1/lI\nnYjpuVeNdUBm/Gkr97JJjEu551XKALRPYXxc/nTTRO4eqFIz3AMvFrn7qDZJ5eyj3AcH6aY4d06m\nDT7kbozcJDcx4X8dYtkyvUDux48Xf04usZaXTAXOf9C4JLq6WkaT3MsG3uBgbyh3nxs7tnJ3xCo1\nwfiQOyDnu584Qcqu7HwdIL4tI0HuZf0uodyLjlTmVO5F97H0ChRgIHdjzGpjzNeNMd83xmwxxvzn\nqr/RUu5chFtFOjNNuUuRu49yB2TLIX3JXWqS87kGwMyzZTj7+9w5KvksOmaD03MvstcGB2kilwSH\ncj8H4L9Ya68F8BoA/8kYs67sD+bPp84NLTcqI0FNW4ZrpneoIvcFC8jv5X68n69ylyAVa+nmLVJT\nWUhZEqdOkS1Sdpa7g1Q/1CH3WLYM93gHdG2Z48fpealFx2xwfb4iWwag+K0nd2vtXmvt45PfnwDw\nNIBVZX/jzsYInYl7NaFapR6NkalYiKncx8ZoB3Dn+dp5kFrSlm2B74RUP/iSey/aMprkXrUSl6yW\nAbqE3LMwxlwB4EYAj1S9lmNZWUaCnOpCO6FaNfgkYgJxlbuvJQPIEZuvJQPI2TI+eQdATrn7xO/2\nUsiic6IcuLijbExz5+nyMIvrjYwxCwHcD+ADkwp+GjZu3Pjy9319Qzh6dCgo5pEj5bWq3ZpQrRp8\nQFxyj0lqgCyx+U4wkp67Txuym+d8Vhqc8btduZclUwH63e7d4XHK7uMy5T48PIzh4eHg+CzkboyZ\nBSL2v7fWfr7odVlyf+ghHuW+Zk3+77R2qHLHspYu+sKF5a/jvsFOnyaSKDpm2MHlS86c8bNQfJGU\nO8HXlnEH2Y2OUg5GM363k3uVeFq0CHj66fA4ZfdxGbkPDQ1haGjo5f9v2rSpUXwuW+YTAL5vrf3v\nvn/AZcsUDYhuTaieOkU37ayKaVdzU0cWUjXedZS7lDVUh9ylqlV8r4NrA/cKxuc6LFhAk8r4OG/c\ntpA710q8jNw1bBmOUsgfA/AfALzOGLPZGPOYMeaOqr/juDmOHy8eEK6ihGOzi2ZC1Ue1u5gxyB2Q\nSSbWUe4Sni/QXcodkLPmqq5DXx9vQtDt9C2rUtJOqHLEaqrcuRBsy1hr/xVAxZaL6eAgd7fhIw/G\nTM2OVVvqq6CZUPXx2wF+cqlDKhLKuY5yl6rzbwu5X3KJ32tjVUy52HUm5DI48VSWO3CflSPH4OO5\nc/RrU8+dC1F2qAJ0sUJvjiqVy0EC4+N0KJjGYABmtnL3jS+1pG0DudeZ5CT6wTe+ZqkxQDt2BwZ4\nCNHHc5dW7t2yiakRFi4MP0lRg9zd6qDsfGuJDRZViEnuUp57HVsmtnJvg+cusZnIV41zxvY5Zhng\nu+5t8NwXLuwCz70pFizoDnLX3rF34oSfLaN51kYnJGyZ2KQG1E/qtsFz5yQI95Bqnx262sqdM6aG\ncq+qeutpW4bjw2mQu+9ZG8mWCY8fO6Hqs3nMIXadO8Bvy7gx4ONpdzO5+1S/OX+/KU6fpoq3vOcu\nuxg9S+6hyn1igv6+rMZXS7nHSKjGtmV6MaHq2/fAFLlzPoELiFstU3eC7VZyr7rOHI/aq7JXe9qW\nCVXuY2O0fCw7GlVLuS9cSHW/ExNhsQB/5c79HNVuKMFzkLJl6lhTc+ZMbSLibkOshGosS8h3xaRF\n7kB431bdxz2t3EPJ/eRJnV2cPgO+r492bmpk8h04qo2yiF2tUke5L1hAkzvnJhqgnnIHZKyZbiL3\nbs0zaRzvUUXuPe25h9oyWlv0fUmHS03G8tzr+M0S5F6HWLg30TjEJveJiXpt4B4DsSqWfCvEuFar\nbVDubvxy23pZdK1y1yL3OktG6Y0PnfFi2TJSyr3OhhgJa6bOBOfawNkPJ09SLXfVU5ik4sdU7r6C\nhmMy1cijVd3HLtkq9RxeICl3ljhcsWLEc4hN7nXiA/wVM+6BxXUO4YpJrkD8MRCD3HtFuQPy1kxS\n7gxxAN4tyz7xuH3nmOQ+MeHvuzpwE5vbrFZnazv3BFOX3LmvQyxbRpPcrdXJ1/mSu2TFTDRyHxiY\n2jTRBL7kLn3EQTYWl3L3ITnnO3MNjpjkfuIEJaR97QjXBu6cQ53JxbWhl8h9JtgyJ09WV9kBOspd\numImGrkbE2bN+M6Mobtgq2rpHTS2LHeCsxwyJrnXqZRxkEgot4Hc6+QdYtoyMcidgwx9y105lHtV\nnJ61ZYCwD1d2IqSDxhEHDpp+YDYmV7VG7BK8uqcLclsidZOpAP/qIbZyj2nL+Ix5DjL0vc6hfetj\nr/asLQPoKHfpIw4ctBOqLiYHuddNJg4M0NOYOM7KB+qTGtCbtkzdFYxbmXJsngPqbeKKUQrJRe5a\nyn3G2jJAuHKv6jwO5a5ty9QhGa7B4WL6JhON4VUdbbBl6hCbQ2zP3R2DGzrGHerudTh2jKdO21fQ\ncNisvvdXqpYJRAj5+mS8e125cw2OpsqZc/t5bGJt0obY1TJAnGcJAFPHL3DUafvGXbBAT7lrlGr3\ntC2jodxHR8PURZ1kTyi5j4/TzTJ/vt/ruZR7bFJpSqy9Zss0uQ4xDq1zGBwMtwXdmB8YqH4th5jx\nXaGFxvKtpe9Z5S5N7v39pC7GxprFAOrZMqEXysXytUe4Zv7YlkTdGncg7vELDrETqq4NMZS7i81l\ne5Y9DMfBKfcQseZ7nbWUe+vJ3RjzcWPMiDHmyTp/J51QBcJ9Om0/sM7NxTU4YqvW2PHb0oYmkyz3\n6Yx1x19o7DoT+5w5JNhCj+LVUO69ZMt8EsCb6v6RtHIHwiaQiQmydXxsEi5/v87NzUXuTZRzbHJP\ntsxUGzj6oerJQXngql7RjFlnB7i0KOwKW8Za+xCAw3X/Tku5N+1AnzPjHTiSPU2WxZpHHkjEbhq/\nDdUy3AnVmJPs2BhZmEVPDsqDZsFCNmYI6foUYgA85J42MbVYuWtXrtRVj7GVO9fAbINq7tY2cE0w\nsSYW7dWChigE/ASLxAF8WcySe+vp2Lhx48vfDw0NYcGCIbz0UrP30rhI2uQeY1kMxFfubUmoNiU3\na+sdOFbWhibXQbsE1yHGmA9dIWuIQmdxVRViFPXf8PAwhoeHmwXPIBq5A8Azz7RbuftWygAzU7n3\nmude1++ePZssO99SPp82xFLuTWJ3o6DxvacHBoAzZ6hUs86BdgCNBx+Lq+izDA0NYWho6OX/b9q0\nqV4DJsFpy5jJL2+EXKg63pmGcs8OhqaYyZ57XWKZN4+OPzhzJl4bgPgrmJjKncOa0y4i8P2cxlAh\nRRNhqGX9VIGrFPJTAB4GcLUxZqcx5n0+f9dUVZ89SyQ6Z071a0MSMHUGfOgpl8DMVu51icUY3oRm\nk4QqwNcPbjNPnYeFAPGVe2hs7WoZDatVaxdsFVhsGWvtzzb5u6adV2ezj5YtA0x9nrrLe4eZ6rk3\nVc3Omlm2LG4buBKadR8WAsQ5atph4UJgZEQ3rpYtAzTnDt/P1HRl4IuuPFtGK9GpneypS7JcFSux\nlXuT+K4NXDXesW2ZJhOsi6991LRDN5ZC1uUOSX6SVu5dWQpZ5wJplUICPEvGurZMt3vubSBW3wRY\nURtikSvAu3KI5bm31ZZpKtZ8KmWAKVv57Nn6MXwQndzbrNyb2jJNUXegz59PO2hDz/Nuopy5JpbT\np+lcEZ/8SSdiHLMs1YaQ1Us3e+7aq2MNW8a32CMkhg+i2zJJuU+h7g3mzvMeHW0e08WNpdyb2hGu\nDbHO1sm2IeYEwzXJxqxz1yoicNVVvmWroTlBH/QsuWvZMlqeu7ZyB/h2Ccby3NtArE0rZYC41SpA\n3OvAMblqVsvUPXU1RLn7krtkUjUquQ8MkN9Ztzbc19MCwhIwTWwZzVJIF1Pb9wTi2xGcbWjDBBOS\nUD1xIvyJSE2VezcdP6BlsyblDvJanW9cB1qeVjco99CYdR8Q4jB3Lv0bcvwq0B5ibVq+GnuCcbtk\nY1yHbtuh2sTfl1bu7oFCEohK7kCzDtQ6FiBGKaT2DVZ3qZoFB7G1hdxD2sC1Q7RpG2IoaBe3m0oh\ntWKdPOkvlnpWuQPNBojWsqft1TJAOMHFJtfQhGobyL0X2tAk/oIFdFRwSLWWti2jIdaS5z6JJh3Y\nVuUeEsvVe2urpyYTikOvKPeQhGpbyD1G3qWvj/Jm0ueeZxFaIKGRQ0ue+ySaeE5anadJ7mfONKv3\n5ii/jEnu3Z5Q5ayW6cZJNmT81Xk4Nke8Nir3nid3Sc89ZKbXtGWaKuhQ1RabXEOIlXOHbuwJphs9\ndyBs/Dlv2ufh2A4poeqPnif3uXPJE2yyxVcz2RNrE0u3Kkau+K4N3Votw9WGpuNAW9DMn08+f5Oj\ntZvYMkm5B0Ca3EOO4tWslmmqnLTPs8kiJVTjny3j2hAyBtzOzbrlsECYuGjq8zcpnwaa2TJpE1MA\nmnw4DbtkYoIUQp0BH0K0MTxPF3emK/deSKiGtMEp2iblsDGsSI0zqULiJOU+CemEqotRtwNHRynR\n01Y/0IHjxo5JbLHjA2HEOn8+bSA6dy68DU0n2VBrLvQadMtqta4tkxKqgZC2ZYBmg0F7U0cMzzMk\nLtAO5c6x9T6kDcbwlSLGVO4xJpamnzmEdNuo3FNCNQMN5V43BhCu3GPYMrGVcwixzppFpaMcp2I2\nTagC8fshdHIJrVjqFkGjlUNLnvskmhJvHS+8yUVqombcZ2miJEMPjmqK2EfuhsR3bYhJrK4NIUnV\nps9PzcaPpdx73ZYZHa13P09M1DurqfW2jDHmDmPMM8aY54wxv1Hnb7USqnVjNBl4c+aQR3/mTL2/\nc/FilEJ2s3LvlTaEJDQBnnLYGMpdm9zr2jL9/XRPj435/83YGDBvnn+urtXkbozpA/AXAN4E4DoA\n7zLGrPP9+7qek7XttWVcLI0lo0O3e+4hk4trA0cZYJ1dknltiLXXgCN+TM+9zbYMUJ87NLjJFxzK\n/VYA26y1O6y1ZwHcB2CD7x/X/XCnT5PXOmuW/99oJVSbxgLilULGVO4hz0/laoMj1qaqGQg/goBj\ngptpyr1pKWSTPJokuTet2fcBB7mvArAr8/8XJ3/mBemZsUkMoDtUBRBXNYbGDnl+KlcbQpOpXG2I\nuXqJ6blrCpq6tgxQfyXeJuVeQ/8WIk/z5KYgNm7c+PL3Q0NDGBoaUiH3pscKd4Nyzz7Nqr+//t/H\nVO6hdgRXG0KItQ1t6FbPven117Q+pZV7Hv8NDw9jeHjY/00KwEHuLwK4LPP/1QD25L0wS+4OdROq\nTZX77t31/qbJEg7QV+59fVMDpIkCjance4FYXRtCqmViW1MnTgAXXdTsb2N57rt2Vb8uL14T7pBW\n7q4ix1mDTvg6bNq0yf8NM+CwZR4FcJUx5nJjzBwA7wTwBd8/rptQbUrump57k2VW6A5F7aVxNm7T\nTUShXjMQTmwhRw9wtSHEFsnGb3odYtXYa1qf587R4YHz5tWPJSk++/vpUYmnTtVrlw+Cyd1aOw7g\n/QC+AuBn8bJfAAAgAElEQVQpAPdZa5/2/XstW6ZJtYy2565tjzR9QIiD20RUp1Qsi15S7jHbEFKC\nC4RXy3QDubv7uW7iXFq5A3JJVQ5bBtbaLwG4psnftjmh2g2lkEDzG8xVHs2e3Syui338eLMTBbmI\ndWSk+d9ztCG0WoajDe46uAeXa8WPZcto3s9a/LRsWb2/q0L0HarOc/ddUmolVLulFDI0ZsyEZuz4\nrg3dXi0T2oaZUC0TshJvo/j0QXRynz2bfKfTp/1er9V5MWwZ7Rsstucdm9Q42xAzoeraEOM6DAzQ\nvdvk4Rma1TJN7y8NW6ZnyR2ol1TVVO5a1TITE/T5m54t0nRpHFs5x55cgN5IqLo2aHvfgO7DcBy0\n72dp5S51eFhryN33w/XiJqYmZ8eHxgT4SCUp9/htCPG+Qyf5JrHHx+s/DCcbT8uW0VLuEgnVGUPu\nbfXoAJ5NLNo+v0Mi93YkVGOuoJqsGkZH6z8c20E7h1ZXudedsHrelpEkd83suubAc2h6Y7dBuaeE\nalxyDy2HBfTHvKYtkzz3QNTxnJp6WmNj5G37QisBA8S5uVzcpNx5LJGQzVxc/dBkDJw6RUUNHOWw\ndRCyWpgzh/q6Tl1/qpaJBOmEal8fedp1d8JqJnti2DKxlTtHQjXmuSoOoU+E4rgOTfshVjI3NIlb\nl3RDxFpKqAZA2papG8Mle5qc8R3DlpnJyn3hQiLVOquyLDiqZYD4/dA0PufKpW7ckDFfd4Xcdlsm\nJVTRnNzrDEDtZE+sg6NiK3cOYunro2sVcjLhTCZ3rlVD2wVN222ZkCObizBjyL1OjBjJnpmq3EOJ\nxbUhhNxDE6pA84qZ0OenOsRW7k08d80x3+ZNTE0PG6xCK8hdOqEK1LtITWMA3VUK2QvKPaQN7pGN\nXBNMkza4sRbyJCgg/hhou6Bp8yamkMPXytAKcpdOqAL1LlLIwNOszHFoqlxjK3eOhGpIG0ZH6aCt\nJg85yWtDkyMIOFcv3eS5d4stk0ohA9FLtkxfH50ZXSdBEutJPDGVO8fzU0PbwJVMDWlD7D6I6bmH\njnmNahlXdnn2bPVrm64Ee165tymhGmLL1I0FzEzPneP5qaFt4CLWNrQhZvwmsbX3djS1ZeqcneP2\nDMyqeZB6Uu6TaLK9t24MbbKNtYkp5oFVXHaEa8NMJ/eYYyCGLdPELmkazzeW9jHhVWgFufsmVM+e\npeqCJg8kqDMYYiR7Qm7wgQHarXfuXL2/6/b6ao42cFTKAM2rZWLnHWaK5x5aAefDT1obpXzRCnL3\nTaiGVBbU8ei6zZZxO/a0bzBg6trV3UTERWpA908wMROqsXbHdku1DJCUexB8Z64QQtBW7nVm4pi7\nBLk2EWmXf2bRloRq02oZjja41azvQ28443dDKWSILSOt3FtJ7saYO40xW40x48aYm5u+jwa5a5VC\nuliayr1JzPFxIoIm+YtONE2oxSb3XmtDrNVbjOMH6sQ8e5YsyyZ2LiCv3N3qt+nBc0UIVe5bALwD\nwDdD3sSX3EMGhEZFTjZW3YSqtvftlqmhm2eaxAZSQrUtbejWVWMTm7XpWJdW7v39VDU2Nlb/b8sQ\nRO7W2mettdsABFGE9MwI1BuAM0G5c8R06Gbl3isJVaBZP8Ta66C5Kzt0F3IdZ0FDfPqiFZ679MwI\ntL8UUls9xVbOsUkNaM8E0+2T7Ny5ZHv4bPThitvG6rdQX5/bd68stzfGPAjgouyPAFgAd1trH6gT\nbOPGjS9/PzQ0hKGhIQD+CZlQcm/jJqYzZ6jSJHQzT92kFie5NqmWaAux9lIbYnnu2fPVlyzx+xtN\n5R5SKQPoKPfs5xkeHsbw8HCzN8qgktyttW8MjjKJLLlnkZ35y54I000J1RdfrBeL4+CoOuTSBsW4\neHG8+EBvVcu4NsTKfThy8iF3l+Bs8ryEzng+4LBlNMSn46es8AWATZs2NXpPTlumMT351mlrJVQ1\nbRkuBa29cSqLlFANqzOP7blrb6Jy1z5E0Gjn0DSVOxdCSyF/yhizC8BtAP7JGPPFpu/lYytoJlS1\nqmW4EpspoVr/7zjbMG9efd+Zuw11++HcOb5y2DrjTzvHpFX91lMJVWvtP1prL7XWDlhrL7bWvrnp\ne/nYClqdF1pFoT3Qgfo3dmxy5axUaUO1jDHNKmY4VzBNNvZwlcPWHfOh/V6nNlxTuYc8B6JVyp0T\nPh+OY4eqz2DQTvbMROV+7BgfsQ4MTPm4dcA5wQHdt4LhvAZ1bZnQz9zfT7k6n2NLQsd6KoUMhI8t\nE6Jy5syhrfJnzpS/ztrwRFss5V4noRdbuR87xpdQbXq2DmdCFWiWVO0Vcm/zalWr+q2nPHdOSNsy\ngJ93durU1I6xpoih3OtaAr2k3Ju2IbZyd89P5fC8gfoVUzOF3LtBufc8uUsmVAH/1UHogK9Tdsk1\n0Bctqq/cY5L70aO85F6X2Nz+gqbnjeShyREQCxfSijJG/G62ZQD/Ma9V/ZZsmQJIV8u4GFUDUDuT\nz6nc65A7py3SjcrdTeIcycSmbZDog7q7lDmPX/Adf9pWJEe1TFLuAfBRXqGDwmcAcgy8OqWQXAq6\nri0TS7UBlPzkKsFr2gZuSwaofw2OHuWbYIG4k0sMcvftbw07NzROTyt3LVumagByJNnmzqUlf1Xy\nFuDbRNJEucci97aoZm5yj63cY3rusci9V2yZnlbuGraMz0zPsVStU/PMZcvUrdSISe7cpNakDRLK\nve414M47xJzgYtkyGtUy8+eTUBsfL36N28A2b16zGD2v3KWrZXxuPk5VcfRo9eu44rlab5/VAsBr\nCTQhFU47okkbpMg9Zj/U9dyTcveDK7Utu7bu/JqQM+N7VrlXfTi3SaXpzAj4K/duGngOdXdIxtxE\n1BblHrsNbVDu3Uzuvisljnus6vNpHxPug9aQe5Xq4Dg9UVu5+ww8TgXpG9Nt1OLcel/H720LucdO\nqHIr97lzyTrwXb1xk7vPShXQT6hyFC1Ik/uMtmW4Zt+qwcDlQ2or9zoxx8aAWbPCz5DPog65SpB7\nkzr/2LYMt3I3Jt51WLzYv/+17zGOSTQp9wD4KvfQGD7KnWPAt1m5Syln38HJTWoA3by+yhGQq5aJ\nldTOtsH3OsS0ZTji+kxkExM8FWk+5K5RS18HrSH3qplLS7nHsGW0lXtsW0Qi/pIl9ci9LcpdIrFc\nR0FzHj9w8iSRaRU07zFHuqG7gDVsGd+DDX3RKnIvuzE4BkQbPXfOG1xzmdqJ2OS+eDFw5Ij/69tA\n7rEnOc74fX3+m300E6pcq0Rpcp89m75OnWr+Hp1oDblr2DKayt3Hg3Q7NUOWc1nUmVBiK3fuyaWu\nLdOGapnY/cA9ufiMP2t5N+755NA4+lia3AF+a6Y15F61LNHy3LkGvE/1gIvFtVNzJtsybfDcmxw/\nEGsFw3G0dSd8xt/YGCnUWZVPb+aJx3k/S5M7d1K1NeTuzls/fTr/992m3H0G3tGj/k+L943pQy69\naMvU9dyPHAEuuIC3DXV8Z0DmOvj2gyNZzoopH0GjfX49l+2ZlHsgyqwZjqVc2zx3iYOjusWWkSqF\n9E1IHTnCO7ECJE4GBvxv0JjKPVY5Kje5nzhRPpl2k3Jv8rCXMoQ+IPuPjTFPG2MeN8b8X2NMUDeW\nLUu4NiK0SbkfOcJL7jPZlpk1i3Yv+y5rufveIXY/+Cr3XiD3/n665mWTabckVIH6q88qhCr3rwC4\nzlp7I4BtAH4r5M3KKmY4Os9tky97Qr1mnTu3cu+WahkJxQrUuzm4LTEH335w5wBxHnsM+OcepPIe\nmuQOVAs2zYRqaGHEkiX1Kr6qEETu1tqvWmvdoug7AFaHvF+VLRNK7lU7+M6doxtuYCAsDtBuco+9\niUiCWFwbfG6O8XFSe9wJVaD+6onz2GNgZtkyPjG1lDuHWKhbFFAFTs/95wB8MeQNymwZjaNxnfXD\nccO1mdwlbuwLLgAOH/Z7rcTKAainWgcH+R5vl8UFF/iRq8QGJmBm2TJA9UpJS7kfPhyeoOdW7pUF\nScaYBwFclP0RAAvgbmvtA5OvuRvAWWvtp8rea+PGjS9/PzQ0hKGhofN+X2bLaBw2pH0UQC/ZMr4D\nc3ycKjW4avuz8CV3iWSqw5IlfpOc5OrF96hpCXLfvbs6rqZy5+rnqmQnx5hy99Dw8DCGh4fD3gwe\n5G6tfWPZ740x7wHwFgCvq3qvLLnnocyW4dqmX6XcuQbeggVU1nn2LJWc5eHoUWDlSp54QFzl7ktq\nro+57QjXhtjkXke5S+UdYtoyTz9d/hoJ5V5ll2gody5y37VruvDdtGlTo/cLrZa5A8CvA3i7tbag\nQt0fZbYMx7IHqFbunMfgVi0Ze8lz9yU1KcUK+PvNUpUyQD3l3gZrihO+FWKaezs4lfvx48WltpzK\nnQuhruP/BLAQwIPGmMeMMX8Z8mZltszhw8DSpSHvTiib6bV37HGTjO8mmpi2jDS5t0G5x7RlYiv3\nKnLnEmkOPsqd43O6Utu8ssuJCZ7+5C6FDNoEbK1dy9UQoNyWOXRIR7lrZ/I5Sba/n0rrTpwoH2gS\nN/bCheSll9lQUrEd2kLuO3ZUv04qoepjBwJ0HVat4o0dg9y1SiFdrGPHptvDx49Tv4ceqVD38Lsq\ntGqHatH25fFxIiyOi6TluQP65O4TU+JMEYBsKB/lIeU1A+3w3GMnVN3jFrXOXMkiFrlrKPeyWFzj\nqW22DCuWLiWF3glnX3CUrml57i6WNrlX+fyjo/Q4tjJV1xQ+g1NKsQL+ykdqAxMQvxQS8LsOsTYx\nadsyEsq9E1znFPU0uS9bBhw8OP3nXJYMkJS7pC3io1oPHaLrLIE22DKxlTvg1w8SKyifg8M0bRl3\npDbXLuCie+vw4aTcK1Gk3DkHRNlg4JxEXCzpXW11Y0qSio9qPXSIJzGehzq2jJRq9k2oSttTMa6D\nz0FemsrdrcSlj9TmEgvz508dS8GBVpH7smXF5M41EMsGw8GDvKqyTMmcOUMXkuOogyyqlKM0qfgo\ndylyr1MKKancu6FqSOI69PeXn4ppLT+5l31W7j6WJndjeI8gaBW5L10qb8uUKXduci/zIJ0lw72Z\np2j149DLyr0Ntozrg6qjh2P3g5Q9ViZoxsYobzZvHl+8IisX4Bcy0uQO8FozrSL3JUuo88bHz/85\np3Ivsy0klHsVuXOjbLAD8rbITCf32bMpYV11pjv3WMui6jqMjdE9xn0iJVBuS3GrdgBYvrx4vHPv\n50jkHoD+furAzg/HOSguuKBY2c4EcpcmlZi2zOAgVQN1ioNOSJI70P7EshNLEkdAlJGtBLkvWwYc\nOJC/UkrKvWXI8905bZnly2kw5EHTc5ck9zJbRpLcY9syfX20wUR7C3wnqpKq1sZNLB88KBfbkW0e\nJMh9YIBEYd5KSVO5c32unib3PN+d05ZxyqIzoz8xwRsHKF8lSJF7Ud7C4cAB6gMJxLZlgGrVOjEh\ncyJiFlX9cOwY+c6czy/NoiqxLLlq0FbuLmbehCKh3PP6lasUEuA9gqCV5C6p3OfMoa3CnRfpyBFS\nfRxPZXe48MJiFSOlHttsyzjFKnGDO1SRu9sq3t8v14Yq5S55DYDqPpCcYLWVu4uZN+a5iwcuvBDY\nv3/6zznvZc4jCFpH7nkXintQ5M30Ejdc0WAA4tkyBw7Es2VOnqTJde5cmfiuDWXEKm3JuDbEUs5A\n9cpBktxjKfe8mNzK/aKLgJGR6T9Pnrsn8pQ7t12SR7oS5L58OX2WvE0dsWyZgwdlbZkyYpW2ZIDy\nnAqgQ+5V/SDpeQMzU7nnxeQe64ncA5Gn3LmX8lrkPns2WT15N7okuR8+XLxLMGYJnga5l62WAD3l\nHtOWKcv1uPgzQbmPjBAhc8Y5fJietZxFIndPaCn3PFtGQtEuX55PNvv3y8SbPZs85aKKEY2EatEG\nnpmk3MtuUGlyX7GifIKTTqi2RbmPjFBfcGHWLBq/2VjuxFou+6enyb1TubuzqTmfuZlHuFI3XFFS\ndWSE9xF7WRRZM2fPUh241Lkq8+ZROeLYWP7v26DcDx2Kr9yl+2FwkK51jOtQltDXTqhyK3dgujXj\nkrZcD1ufUdUyrsyIc8OFli1TFAsA9u7lH3gORYP94EG6uSQ2rziUJRPbQO7798v1u0NVQlVauRtT\nrt6lE6rayr0opga5c5ZBAjOsWobbkgGKbRlNcpdU7kUVM5LJVIeyZWUbbJl9+3iX6nnwSahKkjtA\n427fvvzfSR8BcfIkrRw6oancT5wge7DzqUmhWLHifHLntvk6bZ8QhD4g+8PGmCeMMZuNMV8yxgTT\nVZ5ylziPIqZyP3eOPqMU0RbZMhqkUkZsbVDu3D5sHnxsGenrsGJFHHLv6ytO6GomVPftI5XNvUq9\n6KLz+5Wb3C+5hN6z6mwiH4Qq9z+21t5grb0JwP8DcE9ogzqJSWLTSx4BSNV/500kBw7Q5+TcMJVF\nkS0jWePuUJRABuRLAF382MrdpxxVuh/KyF16ki8af9wWRjZe5zWXsGSA6bbMvn28fdnXB7ziFcAP\nfsDwXiF/bK3NPs56AYCSY/r90Lms6wVbpjOWpN8OxLVlVq7MrwUGdJV7UcWOBrmvWEHXvOgAM40V\nVBG5nzrFX6DQibwJdmyMrgn38wtcvM7JRGqF1knuL7xAZMyJtWtbQO4AYIz5r8aYnQB+FsDvBTeo\n7/yl/Z49/EQYO6Eq6bcDxcpRQ7mvXEmTVx40yN2d2VJ0Zr8Guc+eTWO4SL1rkXveCkryREiHPCXt\nLBmJuAsWkNWZrQ7SUu7btwNr1vDGWLsW2LYt/H0qyd0Y86Ax5snM15bJf98GANba37HWXgbg/wD4\nlfAmna88t2wBXvUqjnedwoIFpKpGR6d+pknuGso9lud+8cXASy/l/06D3IFia8ZauhbS5A4UT3Ln\nzlGyT7ocs0i5ayW1O8ff/v10L0jAmOljXpPcJZS7I/d7723+PpWur7X2jZ7v9Q8g331j0Qs2bpz6\n1dDQEIaGhnJft2IFsHs3sG4dkfsHPuDZAk8YM2WXXHYZkby1Mg8viKHcy2yZa6+ViwvQ5/ryl/N/\np0Xurs87b7ojR+gaS55t47ByJU1y69ef//PDh8l65KqLLkJscu+cXHftAi69VDbmwYPA6tX0/5ER\n4g9udFbLPP88P7mfODGMr351GL/3e8Cf/Enz9wlK6RljrrLWOndoA4Cny16fJfcy3HIL8N3vAv/m\n3wDPPQdcd11IK/PhCOCyy6YUrcSSMesBu/ffu5ey4lKYybYMUKzcNSwZh4svzu8HjdUTUFwKqRF/\n2bLpeZedO+lek4yZveYjI8BP/AR/HGd3TUzQKuyll/g/1113DeHee4ewYQPwmc8AzzyzqdH7hOqH\nP5q0aB4H8AYALBr79tuBhx+mpcnq1TKKOlvVsXu3nJKeP59UWra0SUO5F9kyGgnVPFvm2DHqB4mE\nWifK9hZokXvRJKdRBgm0U7lLknunFSQ1kc+dO3Ve1M6dwKpVlGPhxKpVtMr8x38E3ujrm+QgtFrm\nTmvtemvtjdbaDdbaAre1Hm6/Hfj2t2X8dodsFcsTTwA33CATx8XKkk2sahlN5d5ZreKW5ZKJPIei\nIx80lXsRuWsr987rsHevfB/kiYudO2VtmRUrzu9vKc8dmPLdJfx2gETQmjXAX/91RHKXwiWX0Oz4\n2c8C118vEyM7GDZvBm68USYOMJ3cpZX74sV0Jk/nRggNYhkYoK/OXarSnmsWRbX2bmOLBooSyxo1\n7gBdg7lzp1cN7dgBXH65bOwYyv3aa4Gnnpr6fzeTO0BJ1f37w6ylVpI7QOr9c5+TI/ebbwa+8x36\nfvNm4KabZOIA+sq9rw+44goafA5jY0T2GsSSZ81I39xZFNkybVDuL744lfSTRp41s2MHjQ1JXHop\nxclCWrmvXw88+SR9f+oUFUlIPfFrzRrge9+TSaY6XHUVcNttYadNtprcz56Vs2Ve/3rgG98AzpwB\ntm6VtWWyq4SzZ+nUN2kFvWYNDT6H55+nm1ry8XIOeclEbeXeVlvmhz+UV84OeeSuEX/VKsqxuNMN\nx8dpsl+1Si7m9dfTfTwxMXWdpSzAX/kV4M//HHj6aTlyf8c7gA99KOw9Wk3uc+bQ8kQCl1xCA+DT\nn6bvJR+YfP31wOOP0/f79hH5SJPsmjXnK/dt2+T6shN5xKZJ7m1R7nm2jIYt4tBJ7taSgpaO39cH\nXHMN8Oyz9P+9e2nFKFmCumQJxdi+nQokJFfGN91EK4V//mf+DUwOr3kNEXwIWkvuN94IPPig3Pkr\nAKn3e++V9dsB4NZbgUcfpe/37pX12x06lbs2uefZMprkHlu5L1lCeY/sRjlAl9w7yyH376fqLe6T\nEvOwbh3wzDP0vZYl56yZL31Jpgwyi7vvpsnyyitl44SgteTe10d17pJ43euoUkbSbwfI39+yhSyg\nJ54Arr5aNh5Ay8WY5N6p3KU91874e/ZMf9SgZimkMdPP2ZmYIKLTVO7Z+JoTS5bcta79DTcQud9/\nP3DXXbKxfuzHaC+ORg6rKVpL7hoYGqKbUJrcFy6kGX7LFuCBB4C3vlU2HhDXlun03K2lRKIWuQ8O\nks22Z8/5P9eslgGmT3IjI9Q2iX0beehMqmv6/bGU+6c/TX7/q18tH++WW+RjhGBGk/vSpcA995C/\nJY1bbwX+5V+Ar30NeMtb5ONdeSUpJncy4bZtlIHXQCepHTw4tflDC52HL7nSUOkzXbLotKc0lTNA\nZPfEE3Hix1Du69dTkvPOO+WPd+gGzPguuOcenRv+lluAj36UVgkam1jmzaPE7a5dRGoHD+raIllS\n0yyDdLj6ajq6wmHbNlrNaN70nSsYjTLELF71KiJYd3y2Zvy1a2nVcPas3vW/6iqq77/zTvlY3YAZ\nT+5auOUWGuRve5teTJdUdfW4GmWQAJHanj1TuyM1k6kOncp961aZM4rK0LmC0VbuCxZQv7uqFU1b\nZt48Kn187jnaXKRB7rNm0c7222+Xj9UNSOSuhPXryWt9+9v1YjrfXdNvB2jFMG8ekQkQh9w7lfvW\nrXJ7JoqwZg3ZBA6a5OqQ3dyjPbmsWwe885009m6+WSfmDTfoHHHRDUjkroQ5c4hoNSplHFzFjDa5\nGwO89rXAQw/R/9ui3LXJ/bbbgEcemfq/NrkCRHbOd9e2ha69llaL//APeqvGhCkkcleEZqUGQFbQ\nX/4l7abTJHfgfHLfsUOf3NesoUegnTtH/4+l3EdHaVMNoE+uwBS5j4yQTaaZUP7d36UigsFBvZgJ\nU0jk3sN405vImvibv6HlsSYcuZ88CXzlK8CP/7hu/IEBmkx37KA27N6tVy3kYMyUerc2jnJ3FTMf\n+ADw3vfqWhaLFiVijwnB/Z8JbcDKlcAdd+jHXb+e7JiPfYw2fGgrVoAssG3b6Pjja66R3e1chNtu\nowPq1q+nctDFi3XjX3YZHaT11FPAJz+pGzshLhK5J4hg1iwitnvuoYcOxMDatbRyGRzUt2QcbrsN\n+PCHyRZ5//v14xsD/OqvAhs26DwoJaE9SOSeIIbXvpZ875AHDoTguuuAv/s78vt/9EfjtOGWW+hc\noWef5XmifRPcfXecuAlxkcg9QQy//Mt01EKs3YK/8Au0M/UjH6G2xMDixbSCeM97ZE8eTUjohLGd\nz+GSCmSM1YqVkJBF9uHkMfDSS3SIVyoHTGgCYwystbVHMIumMsb8qjFmwhjT4jPSEmYqYm9qufji\nROwJ+ggmd2PMagBvALCj6rUJfBgeHo7dhJ5B6ktepP5sBziU+58B+DWG90mogXQD8SH1JS9Sf7YD\nQeRujHkbgF3W2i1M7UlISEhIYEBltYwx5kEA2Y3zBoAF8DsAfhvAGzt+l5CQkJAQGY2rZYwxrwLw\nVQCjIFJfDWA3gFuttftyXp9KZRISEhIaoEm1DFsppDHmBQA3W2sPs7xhQkJCQkJjcG4vsUi2TEJC\nQkIroLaJKSEhISFBD2Ibw40xdxpjthpjxo0xhc9hMcbcYYx5xhjznDHmN6Ta0+0wxlxgjPmKMeZZ\nY8yXjTG55wtO9vdjxpjNxphIR3a1E1VjzRgzxxhznzFmmzHm28YY5Se/dhc8+vM9xph9k+PxMWPM\nz8VoZzfAGPNxY8yIMebJktf8j8mx+bgx5saq95Q89WMLgHcA+GbRC4wxfQD+AsCbAFwH4F3GmHWC\nbepm/CaAr1prrwHwdQC/VfC6k9bam621N1lrf0qvee2G51j7eQCHrLVrAfw5gD/WbWX3oMa9e9/k\neLzZWvsJ1UZ2Fz4J6stcGGPeDGDN5Nj8JQD/q+oNxcjdWvustXYbyn34WwFss9busNaeBXAfgA1S\nbepybADwt5Pf/y2AIuJOeY98+Iy1bB/fD+D1iu3rNvjeu2k8esBa+xCAsmKUDQD+bvK1jwBYbIwp\nfbZb7CcxrQKwK/P/Fyd/ljAdK6y1IwBgrd0L4MKC1801xnzXGPOwMSZNlFPwGWsvv8ZaOw7gSDov\nqRC+9+6/n7QRPj15VElCM3T2925UcGXQkb8lG5zuttY+4PMWOT+bsRneig1jvrjMWrvXGHMlgK8b\nY5601r7A2c4uhc9Y63yNyXlNAsGnP78A4FPW2rPGmF8CrYrSaqgZanNlELlba0Mfw/AigGzSajWA\nPYHv2bUo68/JZMtF1toRY8xKANM2ik2+x97Jf18wxgwDuAlAIne/sbYLwKUA9hhj+gEsSvs2ClHZ\nnx1991cA/ptCu3oVL4LGpkMlV2rZMkW+26MArjLGXG6MmQPgnaDZPmE6vgDgvZPfvwfA5ztfYIxZ\nMtmPMMYsB3A7gO9rNbDl8BlrD4D6FgDuAiWuE/JR2Z+TIsRhA9JYrIJBMVd+AcC7AcAYcxuAI86m\nLYS1VuQLlPDbBWAMwEsAvjj584sB/FPmdXcAeBbANgC/KdWebv8CsBR03MOzAB4EsGTy5z8C4H9P\nfl2FZuoAAACmSURBVP8aAE8C2AzgCQDvjd3uNn3ljTUAmwC8dfL7uQA+Pfn77wC4Inab2/zl0Z8f\nAbB1cjx+DcDVsdvc1i8AnwIp8dMAdgJ4H6gq5hczr/kLAD+YvLdvrnrPtIkpISEhoQcRu1omISEh\nIUEAidwTEhISehCJ3BMSEhJ6EIncExISEnoQidwTEhISehCJ3BMSEhJ6EIncExISEnoQidwTEhIS\nehD/HwvSsjVYZXShAAAAAElFTkSuQmCC\n",
            "text/plain": [
              "<matplotlib.figure.Figure at 0x7f2c4a6b44a8>"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "pt.plot(mesh, sol)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 136,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": []
    }
  ],
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.5.1+"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}