solving heat eqn by discritizing and then as an ODE

PDE.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "using BenchmarkTools, DifferentialEquations, Plots\n",
    "pyplot();\n",
    "using PDEOperators"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "heateqsystem (generic function with 1 method)"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "function heateqsystem(t, u, du)\n",
    "    A(t, u, du)\n",
    "#     scale!(du, length(u)^2)\n",
    "end\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "test (generic function with 1 method)"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "function test(N)\n",
    "    h = 1 / (N-1)\n",
    "    x = collect(0 : h : 1)\n",
    "    u0 = -(x - 0.5).^2 + 1/12\n",
    "    tspan = (0., 2.)\n",
    "    prob1 = ODEProblem(heateqsystem, u0, tspan)\n",
    "    sol1 = solve(prob1, dense=false, saveat=[0.03, 1.0],force_dtmin=true)\n",
    "    return x, sol1, prob1\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "0.010101010101010102"
      ],
      "text/plain": [
       "0.010101010101010102"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "A= LinearOperator{Float64}(2,2,100,:Neumann,:Neumann;bndry_fn=(1,-1))\n",
    "x, sol, prob = test(100);\n",
    "gap = 1/(100-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "BenchmarkTools.Trial: \n",
       "  memory estimate:  138.14 KiB\n",
       "  allocs estimate:  2171\n",
       "  --------------\n",
       "  minimum time:     1.914 ms (0.00% GC)\n",
       "  median time:      2.105 ms (0.00% GC)\n",
       "  mean time:        147.116 ms (0.00% GC)\n",
       "  maximum time:     1.968 s (0.00% GC)\n",
       "  --------------\n",
       "  samples:          38\n",
       "  evals/sample:     1"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "@benchmark solve(prob, dense=false, saveat=[0.03, 1.0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 59,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<img src=\"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAYAAAByNR6YAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAIABJREFUeJzs3Xt8HXW97//3zKy1kjRZK73Qa5K2QCnQAuVWbSiCILQCCiiNbBEVbAW2dLtPe7bbs4/HDfhQ5IG7OZ7tbyNCEXCrXAqoKNoqCAI2QgsUaYFCuTXphd7X/TYzn98fk6wmbS6T1cuaT/p+Ph55dDVN0u/01STffL+zZgwRERARERHRQWNWegBEREREQw0nWEREREQHWaiSf7nruti8eTOi0SgMw6jkUIiIiIgGJCJIJpOYMGECTLPvdaqKTrA2b96MpqamSg6BiIiIaNDa29vR2NjY559XdIIVjUYBABs2bMDo0aMrORQqQy6Xwy233IKbbroJ1dXVlR4O+SQi2JrI46bWH6Pl2n/ELieCHTnBjhywPSfYngV25IGdOcHOPJAo9P5xhoWA4RGgPgIMrzIQC3uPo2HvcTQCRMPe72tDQG0YqAsDw0IGhlne+9eEvF/DJlSsYrsiyNrwXhwgUwTSjiBdBNJFIGUDaRtIFgTJIpAoAoluj+MFYE9esKcI7MkDbi9PMQqZwMgqYFQVMLrGwFHdfj2qBhhTbWCEVcADP/khvvu//geG1/JzTxt+7dRt+/btmDJlSmkO05eKTrC6vqBGo1HEYrFKDoXKEIlEUFVVhVgsxi8SAZEuCtrTQEdasCkNbM70/HVLRrA1CxTdYcC0/42HXvDer8YCxtYAY2oMjB4GnDwKGF1t4KhqYFS1gVFVKD0eWeVNrCJW8CdEh8Lwg/RxRASpIrC7AOzMATvz3iS3a2K7PetNeLdlgfVxwfYcsCPXNSmrBU66Bf/9O6/FuGFAwzADE4YBDbV7f22oBZpqDYytASzzyOwVRPzaqVsulwMw8A+FRiUv05BIJFBfX494PM4JFtEAXBF8mAU+SAreT3m/fpACNqYE7WlvYrU73/N9RlZh7zfbYcC4YQbGDwPG1XT+Osz75lsX5jdfDRzXm3xtzeydLG/NCLZkgc1pweYMsKnz14K79/1CBkqTraY6YGKdgUl1wOSogUmdj2v5f4DIF79zl4quYHXJZrOcYCmUzWbxT//0T/jRj36EmpqaSg9nSEgWBO8mgXcTgneTPR+/n+z5TbM+Akzq/GY5e6yJprq930AbO1cxakK9f9Ps0S7MdlpYpoExNUAUWfznN/v+3BPxJmIdaaC9awKe8lY2N6aAtg9dtKcBp9uP10dVA8dEDRwT6/y12+PGWq6AHUxD+WtnoVDABx98AMdxKj2UA2JZFiZNmoRIJLLfn2WzWV8fIxATrP7OwqfgMk0TjY2N7DdI6aLg7QSwfo9gQ6LzJQ5sSHgrEl1qQyh9g7ukycTkaM8Vh+FV5X/DYzvdBupnGN727lHVwKmjev9/YrveSlfXiuj7ScF7nZP6v2510ZEGuuZfEdP7v3hczMCUmIEp9d7j44d7ky9TwflzQTJUP/86Ojpw1VVXIZPJVHooB8WwYcPwwAMPoKGhocfr/XbjFiHRISDiffN6Y4/gjd2CN+PehGp9XNCR3vt2o6qAKfUGpsSAKTEDx9UbODYKHBvzvkFqOPGbhqa84612bUgI3kkINiSAt+PeDwTvJgC78ztHjQUcVw8cX+9NuE7sfDl+uPeEBjoyuK6LG264AVu3bsV3vvMd9eeW5XI5fPvb38b48eNx55139phUqdoizGQynGAplMlk8OUvfxn3338/hg0bVunhVISId+7Tut2CtbsE63YL3tjjTaySRe9tIiYwtR44friBLx1nYmq9gePrgan1BkZWV+YbENvpdjj6VVkGjqsHjqvf//+o7Xrn/3X90PBW5w8Qz693saVz8cKAt4V94ggD04YbmD7CwEkjgWnDjSP+fK+h+Pm3Y8cOvPzyy/je976HU089tdLDOSgWLlyIb33rW9i5c2ePKx34XaELxATLsqxKD4HKYFkWmpubj5h+u3KC13YL/r5T8PddgrW7vYlV10SqNgRMG2Fg2nDgs5PN0jeWo6PBO3/lSGs31FS6X8g0cGzMW2m9eJ8/ixe8VduuHzTe2CP41fsulry2922OjgInjTBw8kgDp4z0fp1a733cI0Gl+x0Ke/bsAYB+rwulTdex7N69u8cEy2+3QEywqqqqKj0EKkNVVRUWL15c6WEcdK5450St2SlYs0uwZqfgtV17t/YiJjBtBHDyCAOXTzJx0kjvp/OJdXrORRmq7Y4UQe5XHzEwa6yBWWN7vj5d9CZbXT+YrN0luP9tF5s6P6+qLGDacOCUkQZOHbX35UDONQyqIPcrl+t6z8AZSpPGrmPpOrYufucsgZhgpdNpbhEqlE6n8dnPfhaPPfYYamtrKz2cshRdwbrdwEvbBS/v9CZTr+4UpG3vzxtqgRkjDXzxOBOndP60fVw9EFb+k/ZQaHck09ivNmzgzNEGztznmtI7c94PMH/vfHl1J/DQuy5ynU9Cm1znnah/2lEGzuh8GTeMn39Unvfffx/XXHMNXnnlFRx99NFYs2bNoD9GOp0e+I0QkAlWOByu9BCoDOFwGC0tLWr62Z2TqVXbBas7J1R/3yXIO975IicMB04bZeAzk02cOsrAjJEGRtfo/kLeF23tqKeh1G9UtYGPTzDw8Ql7X2e7gvV7UFpBfmWH4P+tdbGr8zpvE4YBp3dOts48ysDM0QbGKpp0DaV+2sRiMXz3u99FPB7Ht771rbI+ht9ugZhg9XadCQq+SCSCBQsWVHoYvRIRvJMAXtguWNX58soOQdYBTAM4cThwxlEGrp5i4oyjvK2II+nE2yC3o4EN9X4h08D0kcD0kQa+MMV7nYh3Yv1LOwQv7xC8tEPwX6+72OFdVBsT64CZow3M7JxwzRxtIBoJ5uf0UO8XBP/xH/+Bt956C3fddRcA7xyxKVOm4K233sLZZ5+NZ555puyP7XfOEogJFrcIdUqn07jgggvw5JNPVnyZe09e8OJ2wd+2CV7ofNnZ+dPusTHvC+8Vk018ZIyB046wyVRvgtSOBu9I7GcYRulacFcc7b2ua9K1arv3+b9qu+C7a1ykit6q9PQRwKwxBj46xsSsMQZOHB6MJ5wcKf3eTQj29HEv0wMxPAIcE+u/44IFCzB16lTcfvvtGD58OO69915cdtllGDly5AH//aq2CLmCpVMkEsHixYsPez8RwdtxYOU2wcoPXfx1q+B17wksGFEFfHS0gYXTvS+oM0cbGFWhSyEEWaXa0cHBfp7uk66WY7zXOa7gzT3e6vXftrl4YZvgp285cAWIhb0J11ljvZePjjEQq8Aq15HQb0dOcNzDdq83ND9QlgFsvTqEo/r52j58+HDMmzcPP/3pT7Fo0SL8+Mc/xkMPPXRQ/n5VK1jch9ap6zyCQ63geNsBz28VPLdVsPJDb3XKAHDSSOBj40z86wwDzWO8E9B5cc6BHa52dGiwX9+sbtuLXzneuzhkquitbv1tm/f14z/Xubj5Ze9ryMkjgdljTXxsnIGzxxloqjv0Xz+OhH5HVRt4+3OhQ7aC1d/kqsvXv/51XHrppTjxxBMxevRonHbaaQfl71d1DlYqleIWoUKpVAof/ehH8cILL6Curu7gfdyi90Xw2S3ehOrF7YKcAwwLAc1jDNw43cTszp8+6wN6jkXQHap2dHiw3+DUhQ2cN8HAeZ0n0ot4F0dd+aHgrx+6+PNmFz9+w/uziXXwJltjDZw73sQJww/+D21HSr+BtvEOtRNOOAHHHHMMrrvuOtx+++0H7eOmUilfbxeICZb2S+ofqaqrq9Ha2nrA/eIFb3XqL1u8l5d2CBwBRld7X+hunen9dDljlKH+8ghBcbDaUWWw34ExDO9WPscPN3Bt5yrX9qz3dej5D70f7B58x4UjLkZXA+eMN3DuOG/CddLIA7/eHfsdPl/96lexcOFCzJs3D4B3FfapU6cin88jHo+jsbERX/ziF/H973/f98f0263PCdbbb7+NL3/5y9ixYwfq6+tx3333Yfr06fu93T333IPbbrsNruvi/PPPxx133DHoLb9QKBDzPBqkUCiEuXPnDvr9kgXvi9jTm72Xl3cKXPGeen3ueG9Z/5xxh+YnR/KU246Cgf0OvtE1Bj5ztIHPdJ5A37WS3vWD3/98wUXRdTGqyvs6dd4EA+eNNzFtxOC/TrHf4fP000/ja1/7WmleMmzYMHR0dBzQx/Q7Z+nzra6//npcd911uOaaa/DII4/gmmuuwapVq3q8zXvvvYdvf/vbePnllzF27FhcdtlluOuuu3DjjTcOarDJZJJbhAolk0mceOKJeOONNxCNRvt8u7zjfaF6cpPgz5u9cyEcAcYPA86fYOCGEy2cO9679QYnVIeH33YUTOx36NWFDcxpNDCn884vWds7h6vrB8PFf/MmXGNqgPPGG/hEg4kLJhg42se2GPsdeps3b8b555+PkSNHYsWKFQf1YyeTSV9v1+sEa9u2bVi9ejX++Mc/AgCuuOIKLFy4EBs2bMCUKVNKb/fII4/g0ksvxbhx4wAAN9xwA2699dZBT7BqamoG9fYUDDU1NVi2bNl+/VwRrNkJPLnJxVObvOX2rAMcVe19IfrycSbOm2BiKk9Ir5i+2pEO7Hf41YR6nseV7lzh+vNm7+WG571nKh4TBT7RYOCCCSbObzB6PRmb/Q69CRMm4M033zwkH9tvN7O3V7a3t2P8+PGlZTDDMDBx4kRs3Lixx9tt3LgRkyZNKv1+8uTJ+71Nd/l8HolEoscLABSL3t1yc7kccjnvqnHZbBb5vHcho0wmU3qcTqdRKBRKj7veN5VKwba9+5skk8nS40QiAcdxSo9d14WIIJFIQETgum5pHI7jlB7btl2apdq2XTqprVgslq6BUSgUSo/z+XzpDtv5fB7ZbHbIH1MoFMKpp54Kx3GwKS2487UsPvenAsb+3MYZv7Jxy8sOTAP49ik2Vl0KfHh1CD+dlcOCqYLjhxs9ji8oxzQUO/V2TI7joLm5Gfl8fsgc01Ds1NcxiQiam5uRy+WGzDFp6xRyC7iw0cS3T8rj2Ysc7PxiCA+eU8TcBu9criv/7GDMf9s481c2/uWvWTzVXkTBkdJYmpubS5+LQTmmg9HJcRyISI/HIrLfYwC+Hruu6+tx1/0C+3rsOI6vx/uOvetYu3cyzV6nTvuTXqxevVqmTp3a43UzZ86Up556qsfrFi5cKLfeemvp9+vWrZOmpqbePqSIiNx0000CYL+XK6+8UkREFi1aJIsWLRIRkfnz58tNN90kIiLz5s2TJUuWiIjInDlz5O677xYRkVmzZsnDDz8sIiLTpk2T5cuXi4hIQ0ODrFy5UkREotGorF27VsT7l5L29naJx+MCQOLxuLS3t0vXP8PatWslGo2KiMjKlSuloaFBRESWL18u06ZNExGRhx9+WGbNmiUiInfffbfMmTNHRESWLFki8+bNKx3n/Pnzh/Qx/dPib8hj65NiXrlEjvr/2gV3FQQ/yUnTj9vlWy/acuZVi+XHd92j6piGYqe+jul73/ueRKNROf/884fMMQ3FTn0d03333SfRaFROOOGEIXNMQ63T6ed/Sr728za56qmiWD/cKrirILU/LUjVv/xBbnxkndROPkXdMfXX6cEHH5QzzjhDHnzwQclkMiIismrVKsnn82LbtqxatUps25Z8Pi+rVq0SEZFMJiMvvfSSiIgkk0lZs2aNiIjs2bNHXnvtNRER2blzp7z++usiIrJt2zZZv369iIhs2bJFNmzYICIimzZtkvfee09ERDZu3CgbN24UEZH33ntPNm3aJCIiGzZskC1btoiIyPr162Xbtm0iIvL666/Lzp07RUTktddekz179oiIyJo1a+Rvf/ubnHHGGVJbW9uj04svvlhq1B+j8x162LZtG6ZMmYJdu3YhFApBRDB+/Hg8//zzPbYIf/CDH+Cdd97BnXfeCQD4/e9/j1tvvRXPP/98r5O5fD5fmhV3zQSbmpqwadMmTJgwoTSTrq6uRjabhWmaqKqqQiaTgWVZqKqqQjqdRjgcRiQSQTqdRiQSQTgcRiqVQnV1NUKhEJLJJGpqahAKhZBIJFBbWwvLspBIJFBXVwfDMJBMJhGNRiEipctEOI5Tuqq8bdvIZrOIRqOwbRu5XA51dXUoFosoFAqora1FoVBAsVhEbW0t8vk8HMfBsGHDkM/n4bouampqhtQxvbXHxTPbI/j9Rgd/2QpkHQNjwgVcNDGEiyeF0TwyhzHDdB3TUOzk55gA7xzKpqYmVFdXD4ljGoqd+jomy7LwzjvvoLGxEbW1tUPimIZip65jSiRTWJ+pwp+3mPj9B0W07TBRdA0cXevi4okmLppo4fS6FMYO13NM+3ZyHAcXX3wxTj/9dFx33XWIRCJwHKe02uO6bo/HlmWVVun6e9y1qjfQYwAwTbPPx47jwDCMAR+bpgnDMJDL5XDPPffgpZdewiOPPIKxY8eW/u8Vi0UcddRRiMfj/Z4/3usECwA+/vGP45prrimd5H7bbbdh9erVPd7m3Xffxdlnn93jJPc5c+Zg4cKFff6F3SUSCdTX1w84SKqMvOM9e+aJjYLft7vYkADCpnfphIuaDHyy0cT0Mp5BQ0R0JEsVvfO2lrcL/tDu4v0UUGUB544zcMlEA5dMNHFsha8hVY6//e1vWLx4cWlLUbtIJILW1lbMmjWrx+v9zl36nGCtX78e11xzDXbu3IlYLIZ7770XJ598MhYsWIBLL70Ul156KQDg7rvvxm233QbAm5Tdeeedvi/T0DXI9vZ2NDY2+nofOrQ2pwVPtAue2OjiyU2CtA001gKXNJm4eKKB8ycYqOu8jx8nyHqxnW7sp1v3ftFoFG/Fgd+3u/h9u/dDbdEFjq8HLplo4pImAx8br+cagKlUCps3by6tHmllmiYmTJjQ64VgOzo60NTUVP4E63Do+k+2e/duDB8+vFLDOKJJ5zP+frvRxW8/EKzeITAN4KwxBi6eaOCSJhMnj+x9lcp1XWzevBkTJkzwf9IfBQLb6cZ+uvXXL1nwLmnzROeEa0sGqI8AFzUZ+PREExc1GRhRpWOyNVTt2bMHI0aMGHCCFYgrfHKL6fAqOIJntgh+84Hgtx+4aE97N0H9ZJOBfz7JwkVN/m6QbBgGYrEY+ynEdrqxn2799YtGui54akJE8PKOzh+AN7r4wtMOLMM7TeOySQYum2T6uu4WHVx+P+8CsYLFLcJDL1Hw9vt//YGLJzYKEkVgUh1w2SQTn55k4JxxBiLW4D5RuU2hF9vpxn66lduvIyX43UYXj28UPLVJUHCBU0YCl082cfkkE6eO4oJFOXbkBP/6goP/mm2hJjTwv5+qLcI9e/agvr6+UsMYsrZnvVWqx953S5+Mp43yJlWXTzZxSh9bf36JSOnZLvyk1oXtdGM/3Q5Gv2RBsKLD+6H5dxsF8YJ3o+rPTjbx2ckGzhprwFJy3lalLX3TxVefc7DmsyHMGDXwv1k8Hsfw4cN1bBFWcI435HSkBL9638Vj7wue3er9u35snIHbP+r9hDMpevA+4aTzgnhdT5cmPdhON/bT7WD0i0YMzDvGwLxjTBRd7+T4X70vePhdFz9cC4ypAS6f5E22zpsw+B2KI8mzW7wT8pNFATDwv5PfOUsgJlipVIonuR+A95OCR99zsew9wQvbBGETuKDBwE/OtnDpJANjag7NJ1YqlfK1TErBw3a6sZ9uB7tf2DRwQYOBCxqAH51l4sVtgsfe974v3PUmMDwCXDbJwLyjTVzYaKCKk60euhYjkkV/b9915fyBBGKLkF8kBu/dhOCR91w88p538+QqC/hko/cJ9OlJBuoj/AQiIjqSiQj+vgudP4C7eHOP94SmSzsnW3MbDVT7OOdoKPsgKZj8oHebn4fOt/C5Ywd+Vq7fuUsgVrC67vFD/duY8pZ/H3rHu5xCtQVc3GRg0UkWPjXRQPQwT6ocx8Gbb76JE044AZZlHda/mw4M2+nGfrodrn6GYWDGKGDGKAu3nGFi3W7gkfdcLHvXxc83OIiFvZWtK481cWHDkbmN+NzWvWtMflew/M5ZAjHBSqfTGDFiRKWHEUhbMoJl77p46F3vzu1VFnBRo4H/eYo3qeq66GclpNNpNDc3o6OjgyuQyrCdbuynWyX6GYaBk0YCJ420cPMZFl7f3fkD+7su/nuDgxFVwGcmGfiHY02cN8FA6Ag5Qf7ZrS6mjwDeTwKJor8Nva6baA+EW4QBtDsvePQ9wQPvuHhmi8AygLmNBq48xsSlkwzEuP1HREQHgYjgtV3AQ52TrXcS3gnynzvaxOenGGgeYwzpJ1Kc8HAR508w8av3XfzjNBP/fvrAK4p+5y6BuASwbduVHkLFZWzBQ++4uOyPNsb+3MZ1zzkwANx1toUPrw7ht3NDuPo4M1CTK9u20dbWxn4KsZ1u7KdbkPoZhoFTRhn43kwLb38uhFWXW7h6ionH3ncx+3EHxzxo499edPDarqH3bP8PM4L1ceCc8QZiEf9bhH67BWKClc1mKz2EinBcwZObXFzzjDep+oc/O/gwC/zgoyY2fSGEJy8JYf4JZmBvi5DNZtHS0nLE9tOM7XRjP92C2s8wDJw52sSSWRY2fj6Epy+xMKfRxF1vujjlURszHi3iB6866EgNjclW90sZRcMGkj7vUe23G7cIDzMRwau7gJ+/7eKX77jYkgGOiwFXH2fiqmNNTKkP5mSKiIiOTAXHu6jpzze4ePwDQd4Bzptg4OopJq44Wu9pK//0Vwd/6HCx4cowzvudjfHDgF+eP/Cp6dwiDJgtGcF//N3BjMdsnPaYjZ+97WLe0SZeuMzC+s+F8O+nW+omV7ZtY8WKFUdEv6GG7XRjP9209YtYBj49ycRDnwhh69Uh3HOOd57S/GcdjP25jav+bGN5uwvb1bWy9exWF+eO877vRsNDdIswl8tVegiHRMYWPLDBxUV/sNH4Sxv/Z7WLE4cb+N1cC5u+EMJ/nmXhI2NMtScQ5nI5LF68eMj2G8rYTjf2001zv/qIgWuPN/HUJSFs/HwIN59uYs1OwUXLHUx8wMY3XtBxvtaunHdy/znjvWnQYCZYfrtxi/AgExH8bZvg3re861UlisDZYw18aaqJlqMNDA/o+VRERETlEBG8vAP4WeepLztywOlHAddONfH5Y02Mqg7e973ffuDi0j86ePfKEI6OGfjH5x28uN3FS58JD/i+qrYIi0Wf08YA25wW3LbGwYnLbJz1uIMVHYJ/PsnE258L4blLQ/jqCeaQm1wVi0UsW7ZsSPQ70rCdbuyn21DrZxgGzhht4P+dZWHzF0L49YUWJtYaWNTmYsIvbLQ8aeOJjcHaQnx2q6CxFpgc9X4fDcP3Se5+uwViglUo+DyqgCk4gsfec3HJchtND9i45WUXZ4428OTFFt77hxC+c6a+86oGo1AooLW1VW2/Ixnb6cZ+ug3lfmHTwGWTTfxqTgibvhDCbR8xsT4u+NQKbwvxf69ysCFe+YnWs1sE54zbe42vaBhI+Jzv+u3GLcIyvLFbcM96Fz9728X2HPCR0Qa+crx3BVzeA5CIiGgvEcErO4F71rv4xQYX8QJw7ngD84/3noU47DDfDzFVFAy/38Z/zTZx/YneCfs/fM3Bt1a7SF87xLYINcziM7bg/rdczH7cxrRHbNz7lourppj4+xUhvHB5CNefaB1xk6tCoYClS5eq6Ec9sZ1u7KfbkdbPMAycfpSB/5ptYcsXQvj5eRZMAF96xsGEX9i48a8O1uw8fGs9bR8KHAHOGbd3ChQNG8jY3vUpB+K3WyAmWEHeh351p2DhX73/BNf8xUFtyLvj9uYvhPDDZgsnjzyyJlXdDbXzCI4kbKcb++l2JPerCRn4whQTf/5UCBuuDOHGad5tak57zMZHfm1j6ZsuUj7vCViuZ7cKjqoGThi+93WxiPdryscVGPx24xZhLzK24MF3BD95w8WL2wXjaoCvHG9i/vEmjokduRMqIiKig812BU9sFNz1pos/tAtqw8BVx5q44UQTpx118L/nnvNbG6OrgUcv3HtR0eXtLi5a7qD98yE01vX/d6raIszn85UeAgDgzT2C/9HmoOEXNhY862BkFfDYBRY2XhXC92ZanFztI5/Po7W1NTD9yD+20439dGO/nkKdJ8Y/8ckQ3v98CItPNvFEu4vTf2Vj1m9s3P+Wi6x9cNaCcrbghW2Cc8b3/H4e7Tz1ys+J7n67BWKC5ThOxf5u2xUse9fFeb+zceIyG7/Y4OL6E01suDKEP1wUwmeONhE2ObHqjeM4aGtrq2g/Kg/b6cZ+urFf3ybWGbjlDAvv/0MIv7rQQn0YuOYvDhp+aWNxm4O39hzYROvF7YKC2/P8K8A7BwsAkj62J/12O2K3CHfnBUvfdPGjdS7a097FQP9xmveMhiqLEyoiIqIgeCfhnbLz0/UuduaBCxsMfG2aiU9NNBAa5ALId1928IO/u9j1pRCsbu/7flJw9IM2/nSxhQsa+l974hZhH97a45203tR565pPNBhY81nvYqBXTTE5uRqEfD6Pm2++mcvcCrGdbuynG/sNzrExA7d/1ELHVSH87OMWkkXgM39ycMyDNr73ioMPM/7WiTK2YHmH4OxxRo/JFbB3i9DPxUb9dhv4ttGHgeu6h/zvWL9H8K3VDh59TzC6GviXU7wT6MYN44SqXK7roqOj47D0o4OL7XRjP93YrzzVIQNfPM7AF48z8fIOwR2vO/jeKy5ufsnFhY0GWo42cdkkAyP3uTXP+j2CO99wcd9b3jW4/vs8a7+PPZhzsPx2G/JbhFsygltecrF0vYuGWuDfT7PwhSkGqg/zhc2IiIjo4NqdF/xyg4uH3xU8t1VgGcAFDQZajjFRHwHufMPFk5u8yzLMP97E9SeYOLqPJ6xV3VPEklkmFk7ffwLWnaotwkNxR/FEQfDt1Q6mPGRj2Xsubv+IifUtIcw/weTk6iDRfEf4Ix3b6cZ+urHfwTOiysCN0y1/PrMwAAAgAElEQVT85dPerXl+2Gwi6wALnnUw70kHGRv4749baP98CLd9xOpzcgV03o/QxwqW326B2CI82H73gYuvPOsgWQT++SQT/2vG0LvRMhEREe01fpg32bpxOrA1I4gXgOOH+//eH4v4m2D5NaS2CEUEra+5+MYLLj410cAds60BLxhGRERENOPRIs4ZZ+JHs4fQFmE2mz3gj1FwBF99zsG/vODimzNM/HoOJ1eHWjabxYIFCw5KPzq82E439tON/YIpGjaQ8HEdLL/dArFFaJoHNs/bmRNc8aSDlR8K7jvXwpenBmLeOOSZponGxsYD7keHH9vpxn66sV8w+T0Hy2839VuEb+4RfGqFjXgB+NWFFs4ex/+wRERENDife9LG7gLwp4v7X3tStUWYyWTKer+NKcHsx21UmcALl4U4uTrMMpkMWlpayu5HlcN2urGfbuwXTH5PcvfbLRBbhJbV/wllvbFdwVV/dlAXBp77dGi/C4vRoWdZFpqbm8vqR5XFdrqxn27sF0zRsIFkYeCLiPrtpnaL8P+scnDbqy7+8ikLs7lyRURERAfg31c7uPctF+1Xhft9O1VbhOl0elBv/9QmF7eucfGdM0xOrioonU5j7ty5g+5Hlcd2urGfbuwXTH5PcvfbLRCzk3C4/9lid9uygqufdnD+BAPfnBGI4R+xwuEwWlpaBtWPgoHtdGM/3dgvmLrOwRpoY89vN1VbhK4ILlnu4KUdglevCGE8b9RMREREB8EvN7j4wtMO0teGMKyfW+oNyS3C1tdcLO8Q/OzjFidXAZBOp9Hc3MxlboXYTjf20439ginauTCVKPT/dqq2CCORyIBv8+I2F//2ootvnGLik02BGPYRLxKJYPHixb76UbCwnW7spxv7BVPXBGug87D8dlOzRfjZP9nYkBCsvjyEiMXVKyIiIjp4XtouOPPXNl76TAinHzVEtghTqVS/fy4ieG6r4PJJJidXAZJKpTB9+vQB+1HwsJ1u7Kcb+wVTrHNhKjnA/Qj9dgvEBKu6urrfP39zD7AjB3xsHCdXQVJdXY3W1tYB+1HwsJ1u7Kcb+wVTaYtwgHOw/HYLxJXcQ6H+h/H8hwLTAGaN4QQrSEKhEObOnVvpYVAZ2E439tON/YKpdJL7AOdgDTRn6RKIFaxkMtnvnz+3xcVpowxEI5xgBUkymURjY+OA/Sh42E439tON/YJpWAgwjYG3CP12C8QEq6ampt8/f26rcHswgGpqarBs2bIB+1HwsJ1u7Kcb+wWTYRioCw38LEK/3QK/RdiREryf4vlXQRQKhdDc3FzpYVAZ2E439tON/YKr62ru/VG1RZhIJPr8s+e2ekt1Z3OCFTiJRAKxWKzffhRMbKcb++nGfsEVDQ98krvfboGYYNXW1vb5Z89tFRxfD4yp4QQraGpra9HW1tZvPwomttON/XRjv+CKhg0kBjgHy2+3QGwRWpbV5589t9Xl9mBAWZaF6dOnV3oYVAa20439dGO/4IqGB94i7G/O0l0gVrD6Wm7bnRes3Q2cPS4Qw6R9JBIJGIbBZW6F2E439tON/YLLzzlYqrYI6+rqen39XzvPv+IKVjDV1dWhvb29z34UXGynG/vpxn7B5eccLL/dAjHBMozeJ1DPbRVMGAYcHT3MAyJfDMNALBbrsx8FF9vpxn66sV9w+TkHy2+3QEyw+rpoV9f1r/ifMJiSySTq6+t5sTyF2E439tON/YLLzzlYqi40Go3uv0SVtQWrd/ACo0EWjUYRj8d77UfBxna6sZ9u7BdcfiZYfrsFYoIlsv9y3AvbBEUX+BhPcA8sEUEikei1HwUb2+nGfrqxX3D5Ocndb7dAzF5SqdR+r3tuq6A+AkwfUYEBkS+pVApNTU299qNgYzvd2E839guuaNhA3gEKTt+TKL/dAnEdrFgstt/rnt8qmD3WgGVyizCoYrEYfwJTiu10Yz/d2C+4omHv12QRGNXH5a56m7P0JhArWI7j9Pi97QpWbuP5V0HnOA7WrVu3Xz8KPrbTjf10Y7/g6j7B6ovfboGYYKXT6R6/f3UnkCry+ldBl06n0dzcvF8/Cj620439dGO/4PIzwfLbLZBbhM9tdVFlAWeO5gQryHizUr3YTjf20439gisW8eYdyaIA6H0OomqL0LbtHr9/bqvgo6MNVFmcYAWZbdtoa2vbrx8FH9vpxn66sV9wlVaw+rmau99ugZhgZbPZ0mMRwXNbBWdzezDwstksWlpaevQjHdhON/bTjf2Cq2uClehni9Bvt0BsEXa/aNfbcWB7judfaRCNRtHR0VHpYVAZ2E439tON/YKrzsc5WKouNNp9uW193Hvq6oxRnGAFnW3bWLFiBZe5FWI73dhPN/YLrpBpYFio6xys3qnaIszlcqXH2zsfHlVdocGQb7lcDosXL+7Rj3RgO93YTzf2C7ZouP9zsPx2C8QWYV1dXenxjpxgeAQI8wKjgVdXV4d169ZVehhUBrbTjf10Y79gi4b7Pwer+5ylP4FYwSoW9x7JjhwwmqtXKhSLRSxbtqxHP9KB7XRjP93YL9gGuuGz326BmGAVCnvX4nbkBEdVc/VKg0KhgNbW1h79SAe20439dGO/YIuGjX7PwfLbzZAK3hApkUigvr4e8Xi8dOGuT6+wYQB4fG4gdi+JiIjoCPLpFTZMA/jNnN7nIb3NXXoTwBUsnuCuRaFQwNKlS/lTmEJspxv76cZ+wTbQSe5+uwVigtXzHCzBaG4RqsDzCPRiO93YTzf2C7Zo2Oj3JHe/3QKxD1dbW1t6vJ0rWGrU1tZixYoVlR4GlYHtdGM/3dgv2LyT3Ps+e6r7nKU/gVjByufzAICiK4gXwJPclcjn82htbS31Iz3YTjf20439gi0W6f9ZhH67BWKC5TgOAGAnLzKqiuM4aGtrK/UjPdhON/bTjf2CbaDLNPjtFqhnEa7dJTj5URttl1qYNTYQcz8iIiI6gix908VXn3PgLAjBNPbfUVP1LMKu5bbtOW+uxy1CHfL5PG6++WYucyvEdrqxn27sF2zRzhs+p/pYxVK1Rei6LgDvEg0Atwi1cF0XHR0dpX6kB9vpxn66sV+wdU2w+tom9NstUFuEP37dwddXuijMD8HoZVmOiIiI6FB6fquLj/3WwRstIZwwXPkWYdedqbsuMsrJlQ68I7xebKcb++nGfsEWDXtzkL4u1eC3WyAmWF14DSwiIiKqpK4twsQBXmg/UFuEV/3ZxpYM8PSnAnH9UyIiIjrCbM8Kxvzcxq8utHD55P3XoVRtEWazWQDeFuFormCpkc1msWDBglI/0oPtdGM/3dgv2AY6yd1vt0BMsEzTG8b2nPASDYqYponGxsZSP9KD7XRjP93YL9iqLCBs9n0Olt9ugdoibPplEddONfGdM61KDYmIiIiOcKN+VsS/nmLim6fuPx9RtUWYyWQgIqVnEZIOmUwGLS0tyGQylR4KDRLb6cZ+urFf8EXDQKKPLUK/3QIxwbIsCxkbyDnAaG4RqmFZFpqbm2FZXHHUhu10Yz/d2C/4+rsfod9ugdki3GVEcfSDNv54kYULGwMx7yMiIqIj0OzHbUytB+49d/+rGqjaIkyn09jB+xCqk06nMXfuXKTT6UoPhQaJ7XRjP93YL/j6W8Hy2y0QE6xwOFy6DyEv06BHOBxGS0sLwuFwpYdCg8R2urGfbuwXfNEwkOzjQqN+uwVmi/A3H9bhS884yFwbQk2Iq1hERERUGV/5i4039gBtlw2RLcLaEDi5UiSdTqO5uZnL3AqxnW7spxv7BV80bPR5HSxVW4SRSISXaFAoEolg8eLFiEQilR4KDRLb6cZ+urFf8MUifZ+D5bdbYLYIv/FqLV7eIVj1Gd6HkIiIiCrn9lcd3Paqi11f2v98K1VbhKlUqvM2OZUeCQ1GKpXC9OnTkUqlKj0UGiS20439dGO/4IuGgUQB6G0Nym+3QEywqquruUWoUHV1NVpbW1FdzXDasJ1u7Kcb+wVfNGzAEe8C6Pvy2y0QE6xQKIQdvNGzOqFQCHPnzkUoxG1dbdhON/bTjf2CL9q5M9jbeVh+uwVigpVMJrE9x2tgaZNMJtHY2IhkMlnpodAgsZ1u7Kcb+wVfrPM89t4mWH67BWKCFamqxq48r+KuTU1NDZYtW4aamppKD4UGie10Yz/d2C/4+lvB8tstEOuTaTcEV3gOljahUAjNzc2VHgaVge10Yz/d2C/4omFvwSdREAA9F39UbRF+sNM7I59bhLokEgnEYjEkEolKD4UGie10Yz/d2C/4+lvB8tstEBOstOEtt3GLUJfa2lq0tbWhtra20kOhQWI73dhPN/YLvv7OwfLbLRBbhHuK3jyPW4S6WJaF6dOnV3oYVAa20439dGO/4KsNeRuDvU2wLMvy9TECsYLVsScLA8DIqkqPhAYjkUjAMAwucyvEdrqxn27sF3yGYaAuDCQL+19oVNcWIaowsgqwTG4RalJXV4f29nbU1dVVeig0SGynG/vpxn46RMNAopcVLL/dAjHB2pHn9qBGhmEgFovBMDgx1obtdGM/3dhPh2i49y1Cv90CMcHaliryBHeFkskk6uvrebE8hdhON/bTjf10iEUMJIv7bxGqutBo3AlzBUuhaDSKeDyOaDRa6aHQILGdbuynG/vp0NcKlt9ugZhg7cwJr4GlkIggkUj0erdxCja20439dGM/HaJhIFHY//V+uwVigrU963KLUKFUKoWmpiakUqlKD4UGie10Yz/d2E+Hvlaw/HYLxHWwdhdNbhEqFIvF+BOYUmynG/vpxn461IUNpG13v9fHYjFf7x+IFaxUERjNFSx1HMfBunXr4DhOpYdCg8R2urGfbuynQ5UFFHpJ5LdbICZYAC/ToFE6nUZzczPS6XSlh0KDxHa6sZ9u7KdDxAQK+y9g+e4WiC1CgBMsjXizUr3YTjf20439dIiYva9gqdoiBHijZ41s20ZbWxts2670UGiQ2E439tON/XSIWL2vYPntFpgJFi/ToE82m0VLSwuy2Wylh0KDxHa6sZ9u7KdDX1uEfrsFYoswbAJ14UqPggYrGo2io6Oj0sOgMrCdbuynG/vp0NcWoaoLjY6q8n9vHwoO27axYsUKLnMrxHa6sZ9u7KdDxALy2rcIR0R6OQIKvFwuh8WLFyOXy1V6KDRIbKcb++nGfjpETAMFZ/8rt/vtFogtwjHDrEoPgcpQV1eHdevWVXoYVAa20439dGM/HSIWIAAcAULdNtnq6up8vX+vK1hvv/02zjrrLEydOhUzZ87s8z/CM888g5qaGpx66qmll3JO2hsR5hVtNSoWi1i2bBmKxV7uJUCBxna6sZ9u7KdDpHOGtO+J7n679TrBuv7663Hdddfhrbfewje/+U1cc801fX6A448/HmvWrCm91NTU+PqLuxse5tVsNSoUCmhtbUWh0MvdMCnQ2E439tON/XQoTbD2maL47bbfBGvbtm1YvXo1rr76agDAFVdcgfb2dmzYsOHARtqPcXwKoUq1tbVoa2tDbW1tpYdCg8R2urGfbuynQ1Xn2Uv7rmD57bbfBKu9vR3jx49HKOSdnmUYBiZOnIiNGzf2+gHeeecdnH766Zg5cybuuOOOfv+yfD6PRCLR4wUAasU7YSyXy5VOHstms8jn8wCATCZTepxOp0uzx3Q6XVqqS6VSpTP7k8lk6XEikSjdNyiRSMB1XYgIEokERASu65bG4ThO6bFt20gmk6XHXXfPLhaLpcvkFwqF0uN8Po9MJlN63LVVOpSPqVAo4I477ih9zKFwTEOxU2/HlEwmsXTpUuzevXvIHNNQ7NTXMaXTaSxduhS7du0aMsc0FDv1dUyZTAZLly7Fjh07hswxDcVOFrzj2BlP9jgm309OkH2sXr1apk6d2uN1M2fOlKeeemrfN5V4PC579uwREZH29nY5+eST5aGHHtrv7brcdNNNAu+csR4vZy28XUREFi1aJIsWLRIRkfnz58tNN90kIiLz5s2TJUuWiIjInDlz5O677xYRkVmzZsnDDz8sIiLTpk2T5cuXi4hIQ0ODrFy5UkREotGorF27VsR7GoC0t7dLPB4XABKPx6W9vV26/hnWrl0r0WhURERWrlwpDQ0NIiKyfPlymTZtmoiIPPzwwzJr1iwREbn77rtlzpw5IiKyZMkSmTdvXuk458+fP+SPKZVKSUNDg/zbv/3bkDmmodipt2P6/ve/L3PmzJFPfOITQ+aYhmKnvo7pZz/7mcyZM0dOOOGEIXNMQ7FTX8f01FNPyZw5c4bUMQ3FTj/642uCuwpSO+mkHsf00ksvlY6nPxARuf/++2XGjBkyY8YMue222yQajUqxWBQREdd1ZezYsfL222/3+4FERG699VZZuHBhn3+ey+UkHo+XXrr+oZet+VBERLLZrGSzWRERyWQyksvlREQknU6XHqdSKcnn86XHhUJBRESSyWRpzIlEovQ4Ho+Lbdulx47jiOu6Eo/HxXVdcRyn9I9k23bpcbFYlEQiUXqcTCZFRKRQKEgqlRIRkXw+X3qcy+UknU6XHmcyGR4Tj4nHxGPiMfGYeExKj+mZjoLgroKsak/0OKbdu3f7n2Dt69xzz5V7771XRESWLVsmZ5xxRq/vvHnzZnEcpzSYs846S+65555+/8Luuma1T7+1zff7UHDkcjlZsmRJ6T8t6cF2urGfbuynw4vbHMFdBXl1h9vj9du2bfM1wer1WYQ/+clP8JOf/ARTp07Fbbfdhnvvvbf0ZwsWLMDjjz8OAHj00Udx8sknY8aMGZg1axYuvPBCXHvttf72JrsZHuKzCDVyHAdtbW2lvWnSg+10Yz/d2E+HiOld/Krg9ryUlN9uhsg+lyg9jBKJBOrr67Ft1x6MHlFfqWEQERER9fDGbsG0R2w8/2kLs8ftXY/qmrvE43HEYrE+3z8Qt8qBzWuBaJTP53HzzTeXnplBerCdbuynG/vpEOnjMg1+uwViguW6vBehRq7roqOjg/0UYjvd2E839tOhrwuN+u0WiC3CgZbZiIiIiA6nDzOCcb+w8fgcC5+epHSLkHcU14l3hNeL7XRjP93YT4e+tgj9dgvEBIuIiIgoSPraIvSLW4RERERE+7BdQfgeG/eda+HLU5VuEXbdQ4h0yWazWLBgAfspxHa6sZ9u7KeDZQAG9t8i9NstEBMs0wzEMGiQTNNEY2Mj+ynEdrqxn27sp4NhGIhYQMHpudHntxu3CImIiIh6EbuviFvOMLHoZKv0OlVbhJlMptJDoDJkMhm0tLSwn0Jspxv76cZ+ekTM/U9y99stEBMsy7IGfiMKHMuy0NzczH4KsZ1u7Kcb++kRsfY/B8tvN24REhEREfVi8gNFXD3FxHdnKt0iTKfTlR4ClSGdTmPu3LnspxDb6cZ+urGfHhFz/xUsv90CMcEKh8OVHgKVIRwOo6Wlhf0UYjvd2E839tOjty1Cv924RUhERETUi9MfK2LWGBN3nM0tQjqM0uk0mpub2U8httON/XRjPz0iloGC23MdStUWYSQSqfQQqAyRSASLFy9mP4XYTjf204399OjtMg1+u3GLkIiIiKgXF/7exqgq4MFPhEqvU7VFmEqlKj0EKkMqlcL06dPZTyG20439dGM/PXp7FqHfboGYYFVXV1d6CFSG6upqtLa2sp9CbKcb++nGfnr0tkXotxu3CImIiIh6ceVTNnbngT9erHSLMJlMVnoIVIZkMonGxkb2U4jtdGM/3dhPj962CP12C8QEq6amptJDoDLU1NRg2bJl7KcQ2+nGfrqxnx69bRH67cYtQiIiIqJe/OPzDl7c7uKlz+y9eruqLcJEIlHpIVAZEokEYrEY+ynEdrqxn27sp0dvK1h+uwViglVbW1vpIVAZamtr0dbWxn4KsZ1u7Kcb++nR270I/XYLDfwmh55lWQO/EQWOZVmYPn16pYdBZWA73dhPN/bTo7cVLL9zlkCsYHGZVKdEIgHDMNhPIbbTjf10Yz89ensWoaotwrq6ukoPgcpQV1eH9vZ29lOI7XRjP93YT4/etgj9dgvEBMswjEoPgcpgGAZisRj7KcR2urGfbuynR1UvW4R+uwVigsWLremUTCZRX1/PfgqxnW7spxv76dHbCpaqC41Go9FKD4HKEI1GEY/H2U8httON/XRjPz0iprHfBMtvt0BMsCp4rVM6ACKCRCLBfgqxnW7spxv76RGxAFcAx93bym+3QEywUqlUpYdAZUilUmhqamI/hdhON/bTjf30iHTOkrqvYvntFojrYPE2OTrFYjH+BKYU2+nGfrqxnx5dE6y8A9R0zpj8zlkCsYLlOM7Ab0SB4zgO1q1bx34KsZ1u7Kcb++kR6bymaPcVLL/dAjHBSqfTlR4ClSGdTqO5uZn9FGI73dhPN/bTo7RF2G1O5bcbtwipbLxZqV5spxv76cZ+evR2DpaqLULbtis9BCqDbdtoa2tjP4XYTjf204399Ohti9Bvt0BMsLLZbKWHQGXIZrNoaWlhP4XYTjf204399Ohti9Bvt0BsEfJiazpFo1F0dHRUehhUBrbTjf10Yz89IpZ3W5yCKwC8x6ouNMplUp1s28aKFSvYTyG20439dGM/Pap6OQdL1RZhLper9BCoDLlcDosXL2Y/hdhON/bTjf30KJ2D1W2L0G+3QGwR1tXVVXoIVIa6ujqsW7eu0sOgMrCdbuynG/vp0duzCP3OWQKxglUsFis9BCpDsVjEsmXL2E8httON/XRjPz16O8ndb7dATLAKhUKlh0BlKBQKaG1tZT+F2E439tON/fTo7TINfrsZUsEbIiUSCdTX1yMej/Nio0RERBQo6aKg7j4bvzzPwueneGtSfucuXMGishUKBSxdupT9FGI73dhPN/bT40BWsAIxweI+tE48j0AvttON/XRjPz1C3qWvekyw/HbjFiERERFRH6ruKaJ1lokbp3vLWaq2CPP5fKWHQGXI5/NobW1lP4XYTjf20439dIlYPVew/HYLxATLcZyB34gCx3EctLW1sZ9CbKcb++nGfrpEzJ6XafDbjVuERERERH0Y//MivjbNxLdP5xYhHSb5fB4333wz+ynEdrqxn27sp0uV5i1C13UHfiMKHNd10dHRwX4KsZ1u7Kcb++my7xah327cIiQiIiLqw0mPFHFhg4n/26xwi5B3FNeJd4TXi+10Yz/d2E+XiNlzi9Bvt0BMsIiIiIiCKGIZyDuD3+zjFiERERFRH875rY3JdcDPzgsBULZFmM1mKz0EKkM2m8WCBQvYTyG20439dGM/XfbdIvTbLRATLNMMxDBokEzTRGNjI/spxHa6sZ9u7KfLvldy99uNW4REREREfbj8jzZsF/jdJxVuEWYymUoPgcqQyWTQ0tLCfgqxnW7spxv76bLvFqHfboGYYFmWVekhUBksy0JzczP7KcR2urGfbuyny75bhH67cYuQiIiIqA9f+YuNN/cAKy9TuEWYTqcrPQQqQzqdxty5c9lPIbbTjf10Yz9dqiyjxwqW326BmGCFw+FKD4HKEA6H0dLSwn4KsZ1u7Kcb++ninYO1d7PPbzduERIRERH14RsvOPjtBy7e/Jw3seIWIR1y6XQazc3N7KcQ2+nGfrqxny77PotQ1RZhJBKp9BCoDJFIBIsXL2Y/hdhON/bTjf10iZhA3un2e5/duEVIRERE1Ifvr3HQ+pqL7V9UuEWYSqUqPQQqQyqVwvTp09lPIbbTjf10Yz9dIiZQ6LaC5bdbICZY1dXVlR4ClaG6uhqtra3spxDb6cZ+urGfLvueg+W3G7cIiYiIiPrwkzccfO2vLpwFCrcIk8lkpYdAZUgmk2hsbGQ/hdhON/bTjf10iZgGXAGczmth+e0WiAlWTU1NpYdAZaipqcGyZcvYTyG20439dGM/XSKdtx7s2ib0241bhERERER9WPaui8895WDPl0Oojxi6tggTiUSlh0BlSCQSiMVi7KcQ2+nGfrqxny5VXStYnc8k9NstEBOs2traSg+BylBbW4u2tjb2U4jtdGM/3dhPl0jnTKlri9Bvt9AhGs+gWJZV6SFQGSzLwvTp0ys9DCoD2+nGfrqxny6lCVbnCpbfOUsgVrC4TKpTIpGAYRjspxDb6cZ+urGfLl0nuec7V7BUbRHW1dVVeghUhrq6OrS3t7OfQmynG/vpxn667LuC5bdbICZYhmFUeghUBsMwEIvF2E8httON/XRjP10iltep0HkdLL/dAjHB4sXWdEomk6ivr2c/hdhON/bTjf102fckd1UXGo1Go5UeApUhGo0iHo+zn0Jspxv76cZ+uuy7Rei3WyAmWBW81ikdABFBIpFgP4XYTjf20439dNn3Su5+uwVigpVKpSo9BCpDKpVCU1MT+ynEdrqxn27sp8u+K1h+uwXiOli8TY5OsViMP4EpxXa6sZ9u7KfLvudg+Z2zBGIFy3GcSg+ByuA4DtatW8d+CrGdbuynG/vpsu8Wod9ugZhgpdPpSg+BypBOp9Hc3Mx+CrGdbuynG/vpsu+9CP124xYhlY03K9WL7XRjP93YT5dQ52WvVG4R2rZd6SFQGWzbRltbG/spxHa6sZ9u7KeLYRiImHsvNOq3WyAmWNlsttJDoDJks1m0tLSwn0Jspxv76cZ++kQsIN+5Rei3WyC2CHmxNZ2i0Sg6OjoqPQwqA9vpxn66sZ8+EVPphUa5TKqTbdtYsWIF+ynEdrqxn27sp4+3Reg9VrVFmMvlKj0EKkMul8PixYvZTyG20439dGM/fSLW3gmW326B2CKsq6ur9BCoDHV1dVi3bl2lh0FlYDvd2E839tOn+xah3zlLIFawisVipYdAZSgWi1i2bBn7KcR2urGfbuynT/ctQr/dAjHBKhQKlR4ClaFQKKC1tZX9FGI73dhPN/bTp/sWod9uhlTwhkiJRAL19fWIx+O82CgREREF0kd/beOUkcDd54R8z124gkVlKxQKWKacKz0AABZHSURBVLp0KfspxHa6sZ9u7KdP9y1Cv90CMcHiPrROPI9AL7bTjf10Yz99qqzBn4PFLUIiIiKifly83EaNBTx6obItwnw+X+khUBny+TxaW1vZTyG20439dGM/fSLm3lvl+O0WiAmW4ziVHgKVwXEctLW1sZ9CbKcb++nGfvp0PwfLbzduERIRERH14+qnbXSkgWc+xS1COgzy+Txuvvlm9lOI7XRjP93YT5/uV3JXtUXoum6lh0BlcF0XHR0d7KcQ2+nGfrqxnz4R0yhtEfrtxi1CIiIion58faWDZ7a4+PsVYV1bhLyjuE68I7xebKcb++nGfvp03yL02y0QEywiIiKioOp+L0K/uEVIRERE1I+bX3KwdL2LjquUbRFms9lKD4HKkM1msWDBAvZTiO10Yz/d2E+f7luEfrsFYoJlmoEYBg2SaZpobGxkP4XYTjf204399Ol+L0K/3bhFSERERNSPH6118M0XXWS+omyLMJPJVHoIVIZMJoOWlhb2U4jtdGM/3dhPn4gF5DtXsPx2C8QEy7KsSg+BymBZFpqbm9lPIbbTjf10Yz99IqYBVwDHFd/duEVIRERE1I9fbHBx9dMOMteGUMwk9WwRptPpSg+BypBOpzF37lz2U4jtdGM/3dhPn0jnbKng+p+zBGKCFQ6HKz0EKkM4HEZLSwv7KcR2urGfbuynT2mC5fifs3CLkIiIiKgff2h3cfFyBx1XhRB1uEVIh1g6nUZzczP7KcR2urGfbuynT/cVrAPaIvz617+OyZMnwzAMrFmzpt8PcM899+C4447Dsccei69+9asoFouDGzWASCQy6PehyotEIli8eDH7KcR2urGfbuynT/dzsPx263WL8Nlnn8UxxxyDs88+G7/+9a9x6qmn9vrO7733HmbPno2XX34ZY8eOxWWXXYa5c+fixhtv9PWXc4uQiIiIgu6FbS5m/cbB368IYVLoALYIzznnHDQ2Ng74Fz7yyCO49NJLMW7cOBiGgRtuuAEPPPDAoAeeSqUG/T5UealUCtOnT2c/hdhON/bTjf30iZgGAG+L0G+3AzoHa+PGjZg0aVLp95MnT8bGjRv7fPt8Po9EItHjBQC6FtFyuRxyuRwA72aK+XwegHfV1K7H6XQahUKh9LhrSzKVSsG2bQBAMpksPU4kEnAcp/TYdV2ICBKJBEQEruuWxuE4TumxbdtIJpOlx13/oMVisbT/WigUSo/z+Xzp6q75fL50M8ihfEzV1dX4/ve/D8MwhswxDcVOvR2TYRhobW2F4zhD5piGYqe+jsmyLLS2tsK27SFzTEOxU1/HFAqF0NraikKhMGSOaSh26n5MIXjHtDuZ9r+1K/2YNGmSvPLKK33++cKFC+XWW28t/X7dunXS1NTU59vfdNNNAmC/ly9+8YsiIrJo0SJZtGiRiIjMnz9fbrrpJhERmTdvnixZskRERObMmSN33323iIjMmjVLHn74YRERmTZtmixfvlxERBoaGmTlypUiIhKNRmXt2rXSuRUq7e3tEo/HBYDE43Fpb2+Xrn+GtWvXSjQaFRGRlStXSkNDg4iILF++XKZNmyYiIg8//LDMmjVLRETuvvtumTNnjoiILFmyRObNm1c6zvnz5/OYeEw8Jh4Tj4nHxGMaAsf0+1XrBXcVBFM/Jq+//nrpePoDEZH7779fZsyYITNmzJCf/vSnpT8caIJ1++23y/XXX1/6/RNPPCGzZ8/u8+1zuZzE4/HSS9c/9IYNG0REJJvNSjabFRGRTCYjuVxORETS6XTpcSqVknw+X3pcKBRERCSZTEqxWBQRkUQiUXocj8fFtu3SY8dxxHVdicfj4rquOI5T+keybbv0uFgsSiKRKD1OJpMiIlIoFCSVSomISD6fLz3O5XKSTqdLjzOZzJA/pkQiIRMmTJDt27cPmWMaip16O6bt27dLQ0ODbNmyZcgc01Ds1Ncx7dy5UxoaGmTz5s1D5piGYqe+jmnXrl3S0NAgHR0dQ+aYhmKn7sf07p6i4K6CPLY+KRs3bvQ1wer3OliTJ0/u9yT3d999F2effXaPk9znzJmDhQsX+lo96zrJfefOnRg5cqSv96HgsG0bq1atwsyZMxEKhSo9HBoEttON/XRjP322ZgTjf2Hjt3MsnBXdg1GjRpV3kvv111+PxsZGdHR0YO7cuZgyZUrpzxYsWIDHH38cAHDMMcfglltuwezZszFlyhSMHj0a119//aAHzv9gOoVCITQ3N7OfQmynG/vpxn76dL9Mg99ugbiSe3t7u69nLVKwJBKJ0kScl9nQhe10Yz/d2E+fVFEQvc/GA+dbOLtqM5qamnRcyb22trbSQ6Ay1NbWoq2tjf0UYjvd2E839tOn+5Xc/XYLxPqkZVmVHgKVwbIsTJ8+vdLDoDKwnW7spxv76RPutkXod84SiBWsrutgkC6JRAKGYbCfQmynG/vpxn76GIaBsAkUXPHdLRATrLq6ukoPgcpQV1eH9vZ29lOI7XRjP93YT6eI6W0R+u0WiAlW15XASRfDMBCLxdhPIbbTjf10Yz+dIpa3Rei3WyAmWF2Xxyddkknvhpfspw/b6cZ+urGfTl0rWH67BWKCFY1GKz0EKkM0GkU8Hmc/hdhON/bTjf10qupcwfLbLRATrApeiosOgHS7gSfpwna6sZ9u7KdTxATyjv85SyAmWF13zSZdUqkUmpqa2E8httON/XRjP50ipreC5bdbIK6DxSvZ6hSLxfgTmFJspxv76cZ+OnWd5O53zhKIFSzHcSo9BCqD4zhYt24d+ynEdrqxn27sp1PENFBwxHe3QEyw0ul0pYdAZUin02hubmY/hdhON/bTjf106toi9NuNW4RUtlgsxisRK8V2urGfbuynk8otQtu2Kz0EKoNt22hra2M/hdhON/bTjf106roOlt9ugZhgZbPZSg+BypDNZtHS0sJ+CrGdbuynG/vp1LWC5bdbILYIebE1naLRKDo6Oio9DCoD2+nGfrqxn04RE0gUlF1olMukOtm2jRUrVrCfQmynG/vpxn46dZ3krmqLMJfLVXoIVIZcLofFixezn0Jspxv76cZ+OnVtEfrtZkgFr3aWSCRQX1+PeDzOZxISERFRYN3wnIOXdgie+kTG19wlECtYxWKx0kOgMhSLRSxbtoz9FGI73dhPN/bTKWIBeVd8dwvEBKtQKFR6CFSGQqGA1tZW9lOI7XRjP93YT6euyzT47cYtQiIiIqIB/O9VDh58x8Wai7N6tgg5i9epUChg6dKl7KcQ2+nGfrqxn05dzyL02y0QEyzuQ+vE8wj0Yjvd2E839tOpa4vQbzduERIREREN4D/+7uC7r7jY+BlFW4T5fL7SQ6Ay5PN5tLa2sp9CbKcb++nGfjp1rWD57RaICZbjOJUeApXBcRy0tbWxn0Jspxv76cZ+OnWdg+W3G7cIiYiIiAbw0/Uu5j/rYFdLBiNHDOcWIR06+XweN998M/spxHa6sZ9u7KdTpHPGlMopehah67qVHgKVwXVddHR0sJ9CbKcb++nGfjp1TbDyRX/duEVIRERENIDfvO/i8j85eOfSDI4dp2SLkHcU14l3hNeL7XRjP93YT6eI5f2ayCp6FiERERFRkHVtEdo+d3a5RUhEREQ0gOe2uDjndw5Wz83gzElKtgiz2ez/397dxUZVp3Ec/5057SkN7YxLVrT2hdqKoRQsEWwsaNJqCL2gMdGUlcRookUxJmxsAskmJFblCmNv8AJEl4S4XljxLWJiRAzBpCYaMArGBIrtFCzirpuemem083KevdBWkL6cnp3Of57d3yfhouEYn+QL5u955pwxPQIFkEwm0dXVxX4KsZ1u7Kcb++n0+4rQ32q3IA5YoVBBjEHzFAqFUFVVxX4KsZ1u7Kcb++nkhCwAQNbn0YkrQiIiIqI5nP1FsOpIBp+0JbBx+Z90rAjHxsZMj0ABjI2NobOzk/0UYjvd2E839tNpakU4pugpQtu2TY9AAdi2jZaWFvZTiO10Yz/d2E+nyacIuSIkIiIiypGRMcEt/8jgrfUJbFmlZEWYSCRMj0ABJBIJbNq0if0UYjvd2E839tNp8g6WqqcIi4uLTY9AARQXF6Ozs5P9FGI73dhPN/bTafKA5VlFvq7nipCIiIhoDqmsoOTvGexfm8D2tVwR0gJKJBJoaWlhP4XYTjf20439dCr+7cQU0/RdhI7jmB6BAnAcB93d3eynENvpxn66sZ9OlmWhOAR4Ia4IiYiIiHKm7FAau1fE8bf1S3SsCOPxuOkRKIB4PI7Gxkb2U4jtdGM/3dhPL8cGYsmUr2sL4oC1aNEi0yNQAIsWLUJvby/7KcR2urGfbuynlxMChCtCIiIiotypeTONrbfEsbdNyYowFouZHoECiMViqKqqYj+F2E439tON/fRybCA+oWhFWFpaanoECqC0tBR9fX3spxDb6cZ+urGfXg6fIiQiIiLKrTVH0mguj+PgJiUrQtd1TY9AAbiui3A4zH4KsZ1u7Kcb++nl2Bbi44pWhIsXLzY9AgWwePFi9Pf3s59CbKcb++nGfnrNZ0Xo76oFZtu26REoANu20djYaHoMCoDtdGM/3dhPrxIbyHiWr2sL4g4Wb5Pq5LouLMtiP4XYTjf204399HJCQCKV8XVtQRywysrKTI9AAZSVlWF4eJj9FGI73dhPN/bTy7EBD/62bgVxwLIsf7fbqLBYloVwOMx+CrGdbuynG/vp5YSAjChaEfJlazrFYjFEIhH2U4jtdGM/3dhPLycEJNOKVoTl5eWmR6AAysvLMTo6yn4KsZ1u7Kcb++nl2EBW04rQ4LtO6b8gInBdl/0UYjvd2E839tPLCVlIe/6uLYgDVjweNz0CBRCPx1FdXc1+CrGdbuynG/vp5YSAcZ8nrIJ4Dxa/JkencDjM/wNTiu10Yz/d2E8vxwYylr97UwVxByubzZoegQLIZrM4e/Ys+ynEdrqxn27sp5cTAtI+sxXEASuRSJgegQJIJBJoaWlhP4XYTjf204399HJsYCLLFSEtMH5ZqV5spxv76cZ+ejkhICWKVoSZjL93SlBhyWQy6O/vZz+F2E439tON/fQqsaHrKcJkMml6BAogmUyis7OT/RRiO93YTzf208sJAamsvwcUCmJFyJet6VReXo6LFy+aHoMCYDvd2E839tPLCQEpT9FX5fA2qU6ZTAYff/wx+ynEdrqxn27sp5djW/B8vmGjIA5Y4+PjpkegAMbHx9Hd3c1+CrGdbuynG/vp5czj1GSJwbedua6LSCSC0dFRPklIREREBe2tAQ9/Ofpv4K9/nvPsUhB3sNLptOkRKIB0Oo2+vj72U4jtdGM/3dhPL8ff9zwDKJADViqVMj0CBZBKpdDb28t+CrGdbuynG/vpxRUhERERUY4du+Rh4xFFK0Ke4nVKpVJ47bXX2E8httON/XRjP73mcwerIA5Y3EPrxM8R6MV2urGfbuynF1eERERERDl26p+CtW/+omdFODExYXoECmBiYgK9vb3spxDb6cZ+urGfXiXaVoTZbNb0CBRANptFf38/+ynEdrqxn27sp9d8XtPAFSERERGRD0MxQe0hrghpgU1MTKCnp4f9FGI73dhPN/bTS92LRj3PMz0CBeB5Hi5evMh+CrGdbuynG/vpxacIiYiIiHIslhKE9ytaEfIbxXXiN8LrxXa6sZ9u7KeXuhUhERERUaEr5oqQiIiIKPeK9v0L2R1zrwiL8jjTdSbPdleuXDE5BgWUTCaxc+dOvPTSSygtLTU9Ds0D2+nGfrqxn24lKRdj+P0MMxOjd7AuXLiA+vp6U/96IiIiokAGBgZQV1c34+8bvYO1ZMkSAEA0GkUkEjE5CgXgui6qq6sxPDzMFa8ybKcb++nGfrqNjo6ipqZm6gwzE6MHrFDo10+LRSIR/iFTLBwOs59SbKcb++nGfrpNnmFm/P08zUFERET0f4MHLCIiIqIcs3t6enqMDmDbaG1tRVGR0W0lBcR+erGdbuynG/vp5qef0acIiYiIiP4XcUVIRERElGM8YBERERHlGA9YRERERDmWlwPWuXPnsH79etx+++246667cPbs2Wmve/3117F8+XLU19dj27ZtSKfT+RiPZuGn3fHjx9Hc3IyVK1eisbERu3btgud5BqalP/L7dw/49Wsf7rvvPtxwww15nJBm47fft99+i9bWVjQ0NKChoQHvvPNOniel6fjp53keuru7sXLlStxxxx1oa2vD+fPnDUxLV9uxYwdqa2thWRa+/vrrGa+b9dwiedDW1iaHDh0SEZG+vj5Zt27ddddcuHBBKioqZGRkRDzPk46ODnnllVfyMR7Nwk+7U6dOycDAgIiIJJNJ2bBhw9Q/Q2b56Tfp5Zdflq6uLolEInmajubip18ikZBbb71VTp48KSIimUxGrly5ks8xaQZ++r377rvS3NwsqVRKRERefPFF6ezszOeYNI0TJ07I8PCwLFu2TE6fPj3tNXOdWxb8gPXTTz9JeXm5pNNpERHxPE9uuukmOXfu3DXX7d27V5566qmpn48ePSobNmxY6PFoFn7b/dEzzzwjzz33XB4mpNnMp9+ZM2fk3nvvlfPnz/OAVSD89jt48KBs3brVxIg0C7/93nvvPWlqahLXdcXzPNm5c6c8++yzJkamacx2wJrr3LLgK8Lh4WFUVFRMvSvCsizU1NQgGo1ec100GsWyZcumfq6trb3uGsovv+2udvnyZbz99tvYvHlzvsakGfjtl06nsW3bNhw4cAC2bZsYlabht993332HkpISbN68GWvWrMGjjz6Kn3/+2cTIdBW//To6OtDa2oqbb74ZFRUV+PTTT/HCCy+YGJnmaa5zCz/kTjnjui46Ojqwa9curFu3zvQ45NPzzz+PBx98EA0NDaZHoQAymQyOHTuGAwcO4PTp06isrMTTTz9teizy6auvvsKZM2dw6dIl/Pjjj7j//vuxfft202NRDiz4Aau6uhojIyPIZDIAfv0gbTQaRU1NzTXX1dTUYGhoaOrnwcHB666h/PLbDgBisRja29vxwAMPoLu7O9+j0jT89jtx4gT27duH2tpa3HPPPXBdF7W1tbwLYth8/tvZ1taGyspKWJaFRx55BF988YWJkekqfvsdPnx46uGSUCiExx57DJ999pmJkWme5jq3LPgBa+nSpbjzzjvxxhtvAACOHDmCqqoq3Hbbbddc99BDD+GDDz7A5cuXISLYv38/Hn744YUej2bht108Hkd7ezva29uxe/duE6PSNPz2O3nyJIaGhjA4OIjPP/8c4XAYg4ODuPHGG02MTb/x22/Lli348ssv4bouAOCjjz5CU1NT3uela/ntV1dXh+PHjyOVSgEAPvzwQ6xatSrv89L8zXluWbiPhv3u+++/l7vvvluWL18ua9eulW+++UZERJ544gl5//33p6579dVXpa6uTurq6uTxxx+feqqCzPHTbs+ePVJUVCRNTU1Tv/bs2WNybPqN3797k3744Qd+yL2A+O13+PBhaWxslNWrV0t7e7tEo1FTI9NV/PQbHx+Xrq4uWbFihaxevVo2btw49VQ2mfPkk09KZWWl2LYtS5culfr6ehGZ37mF30VIRERElGP8kDsRERFRjv0HwNYMspqDdvYAAAAASUVORK5CYII=\" />"
      ]
     },
     "execution_count": 59,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x, sol, prob = test(100); plot(x, sol(0.3))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "4-element Array{Float64,1}:\n",
       " -1.83333 \n",
       "  3.0     \n",
       " -1.5     \n",
       "  0.333333"
      ]
     },
     "execution_count": 37,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "T = Float64\n",
    "grid_step = one(T)\n",
    "aorder = A.boundary_length-1\n",
    "first_order_coeffs = PDEOperators.calculate_weights(1, (0)*grid_step, collect(zero(T) : grid_step : aorder* grid_step))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 60,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<img src=\"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAYAAAByNR6YAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAIABJREFUeJzs3Xl4VdW9//HP2fskJJCEClKRhKEo8yAIkQTrULVgHRESeusPW1qstgoK0VbrBLbqVSvRqvWq6LW09rYmqNhKC1qqOCRAGGUeREqY5zNPe6/v7w+ao1aGnUPCPuvk83oenmcLiGvlzWm/2SvZxyMiAiIiIiJqMobbCyAiIiLKNBywiIiIiJqY163/sFIKO3fuRH5+Pjwej1vLICIiInJMRBAIBNCpUycYxrHvU7k2YO3cuROdO3d26z9PRERElLL6+noUFRUd89ddG7Dy8/MBAJs3b0aHDh3cWgalKBqN4sEHH8TUqVORk5Pj9nKokdhPb+ynN/bThxLBa58KHlhqwx8H7hxo4L/OOIj+vc5OzjHH4tqA1XAsmJ+fj4KCAreWQSnKzs5Gq1atUFBQwP+B0BD76Y399MZ+eli2XzDxYxu1ewVju3vwxDATnfM82Ls3AQAn/PImj1uPafD7/Wjbti18Ph8HLCIiIkoL+6OCe+sUZqxX6Hca8PRwE9/q9PnXWjmdX1y7g9UgEokcdYHxeBz/+te/YNu2C6tqOqZpomvXrsjOznZ7KU0qEolg0qRJeOaZZ5Cbm+v2cqiR2E9v7Kc39ktPthK8sF7hviUKSoCnSg3c0teA1/jynapIJOLoz3N9wDraV+Bv374d119/PcLhsAsranqtW7fGn/70JxQWFrq9lCZjGAaKioqO+x0UlL7YT2/spzf2Sz8f7VaY+LGNTw4CP+rlwSPFJr6ee/QjQKfd0u6IUCmFn/zkJ9i9ezd++ctfan8+HY1Gcf/99+PMM8/E888/zxcUERFRmtgREty12MYfNwvO6+DBs+cbKO5w/P+f1uaIMBwOf2mB+/fvx7Jly/Dwww9j0KBBLq6s6UycOBH33nsvDhw4kDHfMRkOh/GDH/wAM2fOROvWrd1eDjUS++mN/fTGfu6L2YKnVin8arlCay/w8oUmxvf0wHDwXE6np2vJMe22225Dt27d4PF4sGLFiuRv2LRpE4YPH46ePXuiuLgYa9ascfRrTpmm+aV/Pnz4MAAc99kSumnYy6FDh1xeSdMxTROlpaVf6Ud6YD+9sZ/e2M9df69XGDDLwr1LFH7c28DGsV78qJfhaLgCvjq3HEtywCorK8NHH32Erl27fuk33HzzzbjpppuwceNG3HXXXRg/fryjX3OqVatWX/pnpVSjNqCDhr007C0TtGrVChUVFV/pR3pgP72xn97Yzx2f+gXXzLNwxVwbnfM8WDnaiydLTXytVePeTcZpt+SAdeGFF37lrtHevXuxZMkSjBs3DgAwZswY1NfXY/Pmzcf9tcYIhUKN+v3pZOvWrbj44ovRtm3bjDnOdCoUCmHkyJFa92vJ2E9v7Kc39ju1QgnBfXU2+s2ysOKAoPpSE/+4wkS/dqm9TZ/Tbsf9Sq76+nqceeaZ8HqPfKmWx+NBly5dsG3btuP+2tHEYjH4/f4v/QAAy7IAHPli8Gg0CuDI+/w03O2xbdvRdcPX6p/oWkS+ct3w3zzRtVLqS9dt2rTBQw89hFdffTW5T6VUcl1fvBYRxONxAEfOb2OxGIAjoRp+PhQKIZE48gCzYDCY/NgEAoHktd/vT67B7/dDKQURgd/vT37cGj62tm1/6eMcCASS18FgEACQSCSSf1ni8XjyOhaLJc+ZY7FY8ttSGzplZWXh2muvTe4vE/YEHPn224Z9ZPKesrKycM011yRfH5mwp0zsdKw9ZWVl4eqrr04+6DAT9pSJnY61p6ysLFx11VXJb3rKhD2lY6doNIaqTxV6VyXwxCc2fj7QwNIr47imyILH40l5Tw1zzwnJf+jatassX75cRESWLFkiPXv2/NKvFxcXy/z584/7a0czdepUAfCVHzfccIOIiEyZMkWmTJki69atk27dusmCBQtERGTz5s2ya9cuERHZsGGD7N27V0RE1q5dKwcOHBARkVWrVsnhw4dFRGTFihUSCARERGTp0qUSDodFRKSurk5isZhYliV1dXViWZbEYjGpq6sTEZFwOCxLly4VEZFAICArVqwQEZHDhw/LHXfcIT/+8Y/lwIEDsnbtWjl06JC0a9dOFi1aJCIir7/+uvTp00dERHbs2CGfffaZiIhs27ZNtm3bltzTrbfeKiIiZWVlMn36dBERGTFihMyYMUNEREpKSqSqqkpERPr27Stz584VEZHCwkKpqakREZH8/HxZvXq1iIgAkPr6evH5fAJAfD6f1NfXS0PW1atXS35+voiI1NTUSGFhoYiIzJ07V/r27SsiIlVVVVJSUiIiIjNmzJARI0aIiMj06dOlrKws2W7ChAlf6iQiMmHCBJk6dSr3xD1xT9wT98Q9pdWeRk+8T7o9s1HwYlzOfKhG7qmc0WR7atiHz+eT4znugLVnzx7Jz8+XRCIhIiJKKTnjjDNk06ZNx/21o4lGo+Lz+ZI/Gha4ZcsWERGJRCISiURk3bp1cu6558qaNWtERMSyLNl0yJKl+5Qs3p2Quj1Hv16y13Z0vWSvnbzefNgWy7KS6z/W9f79+6VDhw5y4MABsSxLKisr5Yc//GHy98yfP1/OOeccERGxbVts2/7SdcOeVq5cKSIioVBIotGoiIgEg0GJxWLJ63g8LiJHhryGj63f709e+3y+5H/X5/OJbduilBKfzydKKbFtOxndsqzkdSKREL/fn7xuGELj8bgEg0EREYnFYsnraDQqoVAoed0wqDZ0CgaDct555yWH3EzYk8iRQbthH5m8p2AwKMXFxXLw4MGM2VMmdjrWnhr6NXximQl7ysROx9pTMBiUoUOHJteZCXtKl067fBG5vcYSc0ZMzv5TTP62zW7yPe3YscPRgPWV52B169YNs2fPTn5N0cUXX4zx48dj/PjxmDVrFh599FEsWbLkhL92Ig3Pkdi/fz/at2+f/Pn169dj3LhxePXVV9G7d2/sjwrOeNWCaoandZkeYPc4L07POfE57C233IKzzz4bU6ZMQa9evfDaa69h8ODBAID3338fkydP/tJ3X37Rf+4pEyQSCcyePRujRo1CVlaW28uhRmI/vbGf3tiv6SkR/G6j4O7FNsIWcP9gA5MHGGhlpvZ1Vsdz4MABnH766c6fg3XzzTdjzpw52L17N0aOHIn8/Hxs3rwZL7zwAsaPH49HHnkEBQUFeOWVV5L/8vF+zakT/eU6PceDTWO9OBxv9B99Ql/LhqPhCjjyGItrrrkGffr0QYcOHZLDVUuVlZWF8vJyt5dBKWI/vbGf3tivadXtU5j4scLifYLrz/Lg8WEmCts0/WDVwOlQnBywXnjhhaP+hl69eqG2trbRv+ZUMBg84Zs9dy9ovg+UU71790b37t1x00034fHHH3d7Oa4LBoMYNmwYFi1ahLy8PLeXQ43EfnpjP72xX9PYGxHcU2fjfzcIBrYDPrjKxAVnNv+7pTR8Ef6JuP6+LTq9Fc6Pf/xjWJaFsrIyAEe+26KoqAjl5eVYu3YtioqK8Itf/MLlVZ4aOTk5qKys1KoffY799MZ+emO/k2MpwdOrbfSssvDGVsGz5xtYep33lAxXgPO5xfW3ynH87Y5p4L333sMtt9ySvD3YunVrbN++3eVVucPr9WLkyJFuL4NSxH56Yz+9sV/q3tupMKnGxtpDwE29DTxUbDj+Up+m4nRucf0OVsOzL9LZzp070bt3byxbtgyTJ092ezlpIRAIoKioSIt+9FXspzf20xv7NV59UPDd+RYumWOjbbYHS67z4vkLzFM+XAHO5xbXbx/l5ua6vYQT6tSpE9avX+/2MtJKbm4uqqurtehHX8V+emM/vbGfc1FLMH2VwiMrFAqygN9fbGLc2Z7kQ3bd4LSb6wOWTkeE9Dmv14vS0lK3l0EpYj+9sZ/e2O/ERARvbxNMrrWxLQhMHmDg/sEGCrLd/6Y3bY4IGx6V36DhrQMaHl+fCRr20rC3TOD3+1FQUPCVfqQH9tMb++mN/Y5v42HBlfNsXPOOjbMKPPhkjBe/HmamxXAFfHVuORbXbx+1adPmS//cqVMnZGdnY8aMGfjxj3+s/UPYEokEZsyYgezsbHTq1Mnt5TSZNm3aoLa29iv9SA/spzf20xv7HV0wIXhouULlKoXC1sAbl5kY1c3d48CjcdrtK09yP1UanuR+tCehLly4EBUVFck3mdRddnY2KisrUVJS4vZSiIiI0oqI4E+fCn62yMbBGHD3OQZ+fo6BXG96DVYNjje/fJHrd7AabpV+UUlJCd555x3s3LkTSimXVtY0DMNAp06dMu5hck7/glF6Yj+9sZ/e2O9zKw8IJtXY+HC3YHQ3D6aXmOiWn56DVQNtjgiPNXjk5eWhZ8+ep3g15FReXh7q6+szbnBsKdhPb+ynN/YDDkYFDyxV+J91Cj3bAu98x8S3i/T4OmWn3VwfsNLtbJWc8Xg8KCgoYD9NsZ/e2E9vLbmfrQQvbzjyFjdxBTx+noFJ/QxkN8ObMjcXp91cHxf5oDU9BQIBtG3blv00xX56Yz+9tdR+tXsUznvLws0f2biqiwcbx3pxx0BTq+EK0OhBo/n5+W4vgVKQn58Pn8/HfppiP72xn95aWr/dYcFdi238fpNgyOke1FxjoPQM1+/vpMxpN9d36NI3MdJJEhH4/X720xT76Y399NZS+iWUYPonR96Uec42wQvfNLHoWlPr4QpwPre4vstgMOj2EigFwWAQnTt3Zj9NsZ/e2E9vLaHfu9sVBr5u4eeLFX7Qw8Cm73pxUx8DpqHXceDROO3m+hFhS/8WVV0VFBRk/GdfmYz99MZ+esvkflsDgjsW2nhjq+DCjh68domJge31H6q+yOnc4vodLNu23V4CpcC2baxZs4b9NMV+emM/vWViv4glmLbURp9qC4v2Cf50iYn3r8q84QpwPre4PmCFQiG3l0ApCIVCKC0tZT9NsZ/e2E9vmdRPRPDGZwp9qi389wqFyf0NrC/34r/OMjL2MRROu6XlW+UQERFRelt/WHBbjY13dwiu6OzBU6UmerTNzKHqi5zOL67fwbIsy+0lUAosy0JtbS37aYr99MZ+etO9nz8uuHOhjQGzLGwJCP46wsScy70tYrgCnM8trg9YkUjE7SVQCiKRCMrLy9lPU+ynN/bTm679lAh+v1GhV5WF/1mnMG2IgdVjvLiqq+ujxCnltBuPCImIiOi4lu0XTPzYRu1ewXe7e/DrYSY657WMO1b/iUeE1Kwsy8K8efPYT1Pspzf205tO/fZHBTd/aGPomxYCCcF7V5r486XeFjtcARodEUajUbeXQCmIRqOoqKhgP02xn97YT2869LOU4LdrjjyFvWqLwm9KDSwf7cXFnVwfG1zntBuPCImIiCjpg10Kk2psrDoITOjlwSPFJjrkttw7Vv9JmyPCRCLh9hIoBYlEAtXV1eynKfbTG/vpLV377QgJrv+nhYvetpHr9WDRKBMzLvRyuPoPTru5PmDF43G3l0ApiMfjqKysZD9NsZ/e2E9v6dYvZgseXWGjV5WF+TsF/3uhiZprTBR3cH1ESEtOu/GIkIiIqIX62zaFybU2tgSASf0MTD3XwNda8Y7V8WhzRJguEzw1Tjwex0svvcR+mmI/vbGf3tKh36d+wTXzLFw5z0bnPA9WjvbiyVKTw5UDTru5PmCl2xk0OZOuX0NAzrCf3thPb272CyUE99bZ6FttYcUBQfWlJv5xhYl+7ThYOeW0G48IiYiIMpyIoHqL4I5FNvZFgZ8NNPCLQQZaezlYNZY2R4SxWMztJVAKYrEYKisr2U9T7Kc39tPbqe63+qDg0jk2vvtPG0NO92BtmRe/GmpyuEqR026uD1i2bbu9BEqBbduora1lP02xn97YT2+nqt/hmGByrY1Bb1jYERb8/XITs0d40b2Ag9XJcNqNR4REREQZRInglQ2CX9TZiNjA/YMNTO5vINvkYNUUeERIzSoWi2HatGnspyn20xv76a05+y3eq1Dylo0bP7QxosiDDeVe/Pwck8NVE3LazdvM6zghpZTbS6AUKKWwfft29tMU++mN/fTWHP32hI/csXplo2BQe+DDq018s6Pr91AyktNuPCIkIiLSVEIJfrtGYepSBdMAHh5q4KbeBkyDd6yaizZHhOn8buJ0bDq8GzwdG/vpjf301lT9/rlDYfAbFioWKnzvLAMbx3rx074mh6tm5rSb60eERERE5Ny2oODOhTaqPxMMP8ODJdeZOPd0DlXphkeEREREGohagic+UXhkhULbbODxYSbGne2Bx8Ph6lTS5ogwEom4vQRKQSQSwY033sh+mmI/vbGf3hrbT0Twl38p9Jtl4cFlCrf2M7BhrBc39DA4XLnAaTfXjwgNw/UZj1JgGAaKiorYT1Pspzf201tj+m08LLi91sbc7YIRhR7MudxE769xqHKT09cdjwiJiIjSTCAueGi5wpOrFQpbA0+Wmri2K48D04E2R4ThcNjtJVAKwuEwysvL2U9T7Kc39tPb8fqJCP64WaFXtYWn1yjcN9jA2nIvRnXjcWC6cPq6c/2I0DRNt5dAKTBNE6WlpeynKfbTG/vp7Vj9VhwQTPrYxkd7BGO+4cH0YSa65nOoSjdOX3c8IiQiInLRwajg/qUKz69T6NUWeGa4iUsLXT9gomPQ5ogwFAq5vQRKQSgUwsiRI9lPU+ynN/bTW0M/fyCE59fa6FFl4dVNCk8MM7ByjJfDVZpz+rpz/YgwKyvL7SVQCrKyslBeXs5+mmI/vbGf3rKysjC47BZc9E4WVhxUGN/Tg/8uNtGxNY8DdeD0dccjQiIiolNkV1hw1yIbf9gsGHq6B88MN1ByBu9Y6YRHhNSsQqEQSktL2U9T7Kc39tNP3BY88YmNXlUW/lav0P2jx/DPb8c4XGlImyPC7Oxst5dAKcjOzkZFRQX7aYr99MZ+enlnu8JtNTY2+YFb+xq4b6DCgtxeyGnFfjpy+rrjESEREVEz2BoQVCy08eZWwUVnevDMcBMD2vHrrHSnzRFhMBh0ewmUgmAwiH79+rGfpthPb+yX3iKWYNpSG32qLdTtE/z5EhPvXfn5cMV+enPazfUjwpycHLeXQCnIyclBZWUl+2mK/fTGfulJRPDm1iN3rXaFgTsHGrhnkIE2WV++a8V+enPajUeEREREJ2ndIcFttTb+sUNwZWcPnio1cXZbHgdmIm2OCAOBgNtLoBQEAgEUFRWxn6bYT2/slz78ccEdC20MfN3C1oDg7ZEm3r7ce9zhiv305rSb60eEubm5bi+BUpCbm4vq6mr20xT76Y393KdE8IdNgrsW2wgkgF8ONVAxwEAr88R3rdhPb0678YiQiIioEZbuE0yssbFwr+C73T14YpiJojweB7YU2hwR+v1+t5dAKfD7/SgoKGA/TbGf3tjPHfujgps+tFA820LIErx3pYk/X+pt9HDFfnpz2s31I8I2bdq4vQRKQZs2bVBbW8t+mmI/vbHfqWUpwfPrFO5fogAAvyk18NO+BrxGanet2E9vTrvxiJCIiOgYFuxSmFRjY/VB4Ee9jrwpc4dcHge2ZDwipGbl9/vh8XjYT1Pspzf2a37bg4Lv/dPCxW/baO31YNEoEy9d6G2S4Yr99KbNEWFeXp7bS6AU5OXlob6+nv00xX56Y7/mE7MFlasUHl6u0CYLeOUiE9/v4YHhabq7VuynN6fdXB+wPE34l5ZOHY/Hg4KCAvbTFPvpjf2ax5xtCpNrbXwWACb1MzBtiIG22U3/MWY/vTnt5voRIR+0pqdAIIC2bduyn6bYT2/s17Q2+wRXz7Nw1TwbXfM8+GSMF0+Wms0yXAHspzun3RwNWH/7299w7rnnYtCgQejfvz9mzpwJANi7dy8uv/xy9OjRA/3798cHH3zQ6IXm5+c3+t8h9+Xn58Pn87GfpthPb+zXNEIJwb11NvrNsvDJQcGsy0y8e4WJvqc1750l9tOb024nPCIUEYwbNw7vv/8+Bg4ciK1bt6J3794YPXo07r77bpSUlGDu3Lmoq6vDddddh88++wxZWVmOF+rSNzHSSRIR+P1+5OXl8Ta3hthPb+x3ckQEVVsEdy6ysS8K3D3IwF3nGGjtPTUfS/bTm9O5xdEdLI/Hg8OHDwM48tXz7du3R6tWrVBVVYWf/OQnAIDi4mJ06tQJCxYsaNRCg8Fgo34/pYdgMIjOnTuzn6bYT2/sl7pVBwWXzLHxX/+0MeR0D9aWefHgEPOUDVcA++nOabcTDlgejwevvfYaRo8eja5du+Kb3/wmZs6ciUAggEQigY4dOyZ/b7du3bBt27aj/jmxWAx+v/9LPwAk73ZFo1FEo1EAQCQSQSwWAwCEw+HkdSgUQjweT14nEonkZi3LAnDkbLTh2u/3w7bt5LVSKvmZg4hAKZVch23byWvLspJnrJZlJT+YiUQCoVAIABCPx5PXsVgM4XA4eR2JRDJ+TwUFBQiHw2jVqlXG7CkTOx1rTwUFBQgGg8jJycmYPWVip2PtqaCgAIFAAK1bt86YPTV3p32hBH66IILBb1jYERLMvjiB2SO8KMpJnPI9FRQUwOfzJR9YyU567cnx0a6cQCKRkIsuukgWLFggIiKLFy+Wjh07yu7duyU7O/tLv7e8vFxefvnlo/45U6dOFQBf+TFu3DgREZkyZYpMmTJFREQmTJggU6dOFRGRsrIymT59uoiIjBgxQmbMmCEiIiUlJVJVVSUiIn379pW5c+eKiEhhYaHU1NSIiEh+fr6sXr1a/v0wVamvrxefzycAxOfzSX19vTR8CFavXi35+fkiIlJTUyOFhYUiIjJ37lzp27eviIhUVVVJSUmJiIjMmDFDRowYISIi06dPl7KysuQ+J0yYkPF7sixLRo8eLffff3/G7CkTOx1rT5ZlyfDhw+WFF17ImD1lYqdj7cmyLDnrrLNkzpw5GbOn5urUqbBIZqyzpe1LQTGePSSPr7Dk/16rdnVPlmVJmzZtZOXKleyk4Z62bt2a3NPxnHDAqqurkx49enzp54YOHSrvvPOOtG7dWnbt2pX8+eLiYnn33XeP+udEo1Hx+XzJHw0f6E2bNomISCQSkUgkIiIi4XBYotGoiIiEQqHkdTAYlFgslryOx+MiIhIIBCSRSIiIiN/vT177fD6xLCt5bdu2KKXE5/OJUkps205+gCzLSl4nEgnx+/3J60AgICIi8XhcgsGgiIjEYrHkdTQalVAolLwOh8MZvyefzyd5eXmyd+/ejNlTJnY61p4a+u3bty9j9pSJnY61p4Z+Bw4cyJg9NUenj3bGZfCsqODFuFw/Py4b96THnhr6HTx4kJ003NO2bdscDVgnfKucPXv24Oyzz8bixYvRp08fbN68Geeddx5WrFiBBx54AN26dcO0adNQV1eHUaNGYevWrY6+yJ1vlUNERM1hT1hwd52N320UDGoPPDPcxDc7uv5UIsoQTueXE34X4RlnnIEXX3wRY8eOhWEYUErh2WefRZcuXfDYY4/hhhtuQI8ePZCdnY1XX321Ud9BCCB5zkl6sSwLdXV1KC4uhtfr+vNqqZHYT2/sd3QJJfjtGoWpSxW8BvDc+QZu6m3ATPFNmZsL++nN6dzi+ps9b9++HYWFhW4sgU5CIBBAnz59sG7dOj7LRUPspzf2+6p/7jjypszrDgM39zHw0FAD7XPSa7BqwH5627FjB4qKik54B8v1AYtHhERElKptQcEdC23M+kxw/hkePDPcxODT03OwoszgdH5x/VCaR4R6siwL8+bNYz9NsZ/e2A+IWoJfLbPRu8rCx3sEf7jYxIdX6zFcsZ/enHZzfcBqeGYF6SUajaKiooL9NMV+emvJ/UQEb21V6DvLwq+WK0zsZ2BDuRfjehjaPBW9JffLBE678YiQiIi0sOGw4PZaG/O2Cy4v8uCpUhO9vqbHUEWZQ5sjwoanqJJeEokEqqur2U9T7Ke3ltYvEBf8fJGNAa9b2OgTvDXCxN8u13e4amn9Mo3Tbq4PWA2PrCe9xONxVFZWsp+m2E9vLaWfiODVTQq9qi08u0bh/sEG1pZ5cU1XfY4Dj6al9MtUTrvxiJCIiNLO8v2CSTU2Pt4jKPuGB08MM9E1X9+hijKHNkeEnOD1FI/H8dJLL7GfpthPb5nc70BUcMtHNobOtnAoJvjHFSaqL/Nm1HCVyf1aAqfdXB+weAatJ34Ngd7YT2+Z2M9WgufX2uhZZeGPmxWmDzOwYowXlxa6/n9TTS4T+7UkTrvxiJCIiFz18e4jT2FffgAY39ODR4tNnNE6c+5YUWbR5ogwFou5vQRKQSwWQ2VlJftpiv30lin9doUFN7xn4Zt/tWF6PKi9xsQrF3kzfrjKlH4tldNurg9Ytm27vQRKgW3bqK2tZT9NsZ/edO8XtwW/XnnkOHDudsGMC0wsGmWi5AzX/y/plNC9X0vntBuPCImI6JSZV69we62NTX7g1r4GHhxi4LRWmX3HijILjwipWcViMUybNo39NMV+etOx32d+wah3LFw+10bH1h6sGO3F08PNFjlc6diPPue0m7eZ13FCSim3l0ApUEph+/bt7Kcp9tObTv3CluDRFQqPf6LQIQf48yUmxnb3aP2g0JOlUz/6KqfdeERIRERNTkTwxlZBxUIbu8PAnQMN3DPIQJusljtYUWbQ5oiQ7yauJ74bvN7YT2/p3m/tIcG3/2aj7B82BrbzYE2ZFw8Xmxyu/i3d+9HxOe3m+hEhERFlBl9c8OBShWfWKHTLB94eaeLKLq5/Hk/kCh4REhHRSVEi+P0mwd2LbQQSwH2DDVQMMNDK5B0ryjzaHBFGIhG3l0ApiEQiuPHGG9lPU+ynt3Tqt2Sfwvl/sfHDBTa+1cmDDeVe/GKQyeHqONKpHzWe026uHxEahuszHqXAMAwUFRWxn6bYT2/p0G9fRHDvEhsvrRf0bwe8f5WJi87k3ycn0qEfpc5pNx4REhGRY5YSPL9O4f4lR75V/Vek+WNSAAAgAElEQVRDDfykjwGvwTtW1DJoc0QYDofdXgKlIBwOo7y8nP00xX56c6vfgl0K575p4bYahbJveLBxrBcT+5kcrhqJrz+9Oe3m+hGhaZpuL4FSYJomSktL2U9T7Ke3U91ve1Dws8U2/vypYNjXPVg8ysDQDq5/fq4tvv705rQbjwiJiOioYragcpXCw8sV2mQBj51n4vs9PDBa8FPYibQ5IgyFQm4vgVIQCoUwcuRI9tMU++ntVPSbs02h/ywLDyxRuKm3gY1jvRjf0+Bw1QT4+tOb026uHxFmZWW5vQRKQVZWFsrLy9lPU+ynt+bst9knmFxrY0694NJOHrw1wkTf0zhUNSW+/vTmtBuPCImICKGE4OEVCtM/UejYGqgsMTG6W8t+U2aio+ERITWrUCiE0tJS9tMU++mtKfuJCF77VKF3tYXKVQp3DzKwrtyLMd8wOFw1E77+9KbNEWF2drbbS6AUZGdno6Kigv00xX56a6p+nxwQ3FZrY8EuwbVdPagsMdG9gENVc+PrT29Ou/GIkIiohTkUEzywROG5dQpnFwBPl5oY2dn1Aw0iLWhzRBgMBt1eAqUgGAyiX79+7Kcp9tNbqv1sJZixXqFnlYXfbVJ4tNjAqjFeDlenGF9/enPazfUjwpycHLeXQCnIyclBZWUl+2mK/fSWSr+FexQm1igs3S8Yd7YHj51nolMbHge6ga8/vTntxiNCIqIMtjssuHuxjZmbBIPaA88ON3F+R96xIkqVNkeEgUDA7SVQCgKBAIqKithPU+ynNyf9EkpQ+YmNXlUW/rpN8D/nG1gyysvhKg3w9ac3p91cPyLMzc11ewmUgtzcXFRXV7OfpthPbyfqN3+Hwm01Ntb7gJt7G/jVUAPtc3gcmC74+tOb0248IiQiyhD/CgjuWGTj9c8E55/hwTPDTQw+nYMVUVPS5ojQ7/e7vQRKgd/vR0FBAftpiv309p/9Ipbgl8ts9Km2ULNH8Oq3THx4NYerdMXXn96cdnP9iLBNmzZuL4FS0KZNG9TW1rKfpthPbw39WrdujdlbFabU2tgRBib3N3D/YAP52Rys0hlff3pz2o1HhEREGtpwWHBbjY13dghGFnnwm1ITvb7GwYqoufGIkJqV3++Hx+NhP02xn74CccHtH4bR+09hbDhsY/a3Tfz9cg5XOuHrT2/aHBHm5eW5vQRKQV5eHurr69lPU+ynHxHBHzcLfr7YxuFYFu7sE8SDw3PQOtv1z5Opkfj605vTbq4PWHy3dj15PB4UFBSwn6bYTy/L9wsm1dj4eI+g7BsePDHMwGliIjeLw5WO+PrTm9Nurr86+aA1PQUCAbRt25b9NMV+ejgQFfz0IxtD3rRwOC6Yf4WJ6su8OE2C7Kcxvv70ps2DRvPz891eAqUgPz8fPp+P/TTFfunNVoIX1yvct0TBUsCTpQZu6WsgyzjymTP76Y399Oa0m+t3sFz6JkY6SSICv9/Pfppiv/T10W6FobMt3PKxwqiuHmz6rhe39zeTwxXAfrpjP7057eb6gBUMBt1eAqUgGAyic+fO7Kcp9ks/O0OCce9ZuOCvNrIMDxZea+Lli7z4eu5Xv96D/fTGfnpz2s31I0I+A0tPBQUF/OxLY+yXPuK24KnVCr9arpBrAi9dYOKHvTwwjvOFtOynN/bTm9O5xfU7WLZtu70ESoFt21izZg37aYr90sPceoUBr1u4p07hRz0NbBzrxYTexnGHK4D9dMd+enPazfUBKxQKub0ESkEoFEJpaSn7aYr93LXFL7j2HQvfmWujU2sPlo/24jfDTXytlbNv/2Y/vbGf3px241vlEBGdImFL8OgKhcc/UeiQAzwxzMTY7h4+D4lII9q8VY5lWW4vgVJgWRZqa2vZT1Psd2qJCGZtUehTbeGxlQp3DjCwvtyL755lpDRcsZ/e2E9vTru5PmBFIhG3l0ApiEQiKC8vZz9Nsd+ps+ag4LK/2Sifb2NgOw/WlHnxULGJNlmp37ViP72xn96cduMRIRFRM/DFBdOWKjyzRuEb+cBvSk1c0cX1z2mJ6CTxiJCalWVZmDdvHvtpiv2ajxLBKxsUelZZmLFe4aGhBlaXeZt0uGI/vbGf3rQ5IoxGo24vgVIQjUZRUVHBfppiv+axZJ/C8L/Y+NEHNi7t5MGGsV7cPchEK7Npv4id/fTGfnpz2o1HhEREJ2lfRHBPnY2XNwgGtAOeGW7iwjNd//yViJqBNkeEiUTC7SVQChKJBKqrq9lPU+zXNCwleGa1jZ5VFl7fKnhmuIGl13mbfbhiP72xn96cdnN9wIrH424vgVIQj8dRWVnJfppiv5O3YJfCuW9auL1WYWz3I09hv7WfCa/R/M+0Yj+9sZ/enHbjESERUSNsDwruXGTjtS2Ckq978OxwE0M68EGhRC2FNkeEnOD1FI/H8dJLL7Gfptiv8WK24L9X2OhVbeH9XYKZF5n4+Bp3hiv20xv76c1pN9cHLJ5B64lfQ6A39mucOdsU+s+y8MAShZ/2OXIc+P2eJ35T5ubCfnpjP7057cYjQiKiY9jsE0yutTGnXnBZoQdPl5rocxqPA4laMm2OCGOxmNtLoBTEYjFUVlayn6bY7/iCiSOPXeg3y8LqQ4LXLzPxznfSZ7hiP72xn96cdnM0YMViMUycOBE9evTAgAEDMG7cOADApk2bMHz4cPTs2RPFxcVYs2ZNoxdq23aj/x1yn23bqK2tZT9Nsd/RiQj+/KlC72oLT65S+MUgA2vLvRj9jdTelLm5sJ/e2E9vTrs5OiKcMmUKLMvC008/DY/Hg927d6Njx4645JJL8P3vfx/jx4/HrFmz8Nhjj6Gurs7Rf5hHhESUTj45ILit1saCXYLrunlQWWKiW376DFVElB6a7IgwFArh5ZdfxsMPP5z8DK5jx47Yu3cvlixZkrybNWbMGNTX12Pz5s2NWihvkeopFoth2rRp7Kcp9vvcoZhg0sc2Br9pYXdYMO87Jt74tjethyv20xv76a3Jjgg//fRTtGvXDo888giGDh2KCy64APPnz0d9fT3OPPNMeL1eAIDH40GXLl2wbdu2Yy7I7/d/6QdwZIADjry3T8P7+0QikeQGwuFw8joUCiW/PTIUCiW/kj8YDCbffDEQCCSv/X5/8lae3++HUgoiAr/fDxGBUiq5Dtu2k9eWZSEQCCSvg8EggCPfOdCw3ng8nryOxWIIh8PJ60gkkvF7Ukph69atyZ/PhD1lYqdj7Ukphc8++yz5ezJhT43tZCvB08tD6Fll4XebFB4cmMCK6wx8u9CT9ntq6Nfw85ncKRP3pJTCli1bkr8nE/aUiZ2OtSfHR7tyAkuXLhUAMnPmTBERWbZsmbRv314WLlwoPXv2/NLvLS4ulvnz5x/1z5k6daoA+MqPG264QUREpkyZIlOmTBERkQkTJsjUqVNFRKSsrEymT58uIiIjRoyQGTNmiIhISUmJVFVViYhI3759Ze7cuSIiUlhYKDU1NSIikp+fL6tXr5Z/H4NKfX29+Hw+ASA+n0/q6+ul4UOwevVqyc/PFxGRmpoaKSwsFBGRuXPnSt++fUVEpKqqSkpKSkREZMaMGTJixAgREZk+fbqUlZUl9zlhwgTuiXvintJ4Ty/+Y6UMeSMheDEuV7+5X3aGlPZ7ysRO3BP3lI57atiHz+eT4znhgLVv3z4xDEMsy0r+3NChQ+W1116T/Px8SSQSIiKilJIzzjhDNm3adNQ/JxqNis/nS/5oWODWrVtFRCQSiUgkEhERkXA4LNFoVEREQqFQ8joYDEosFktex+NxEREJBALJdfj9/uS1z+dLrtvn84lt26KUEp/PJ0opsW07+QGyLCt5nUgkxO/3J68DgYCIiMTjcQkGgyIiEovFktfRaFRCoVDyOhwOZ/yeIpGITJo0SQ4fPpwxe8rETsfaUyQSkYkTJyb/u5mwJyeddoWUXP9ORPBiXAa/Hpd3Pg1ouadIJCK33npr8r+baZ0yfU+RSERuueWW5H83E/aUiZ2Otafdu3c3zYAlIvLtb39b5syZIyIiW7Zskfbt28v27dvloosukldeeUVERKqrq2XIkCFO/rjkQgHInj17HP87lD4ikYhMmTIl+ZeZ9NLS+sVtJdNXWlLwSlzaz4zL82stsWzl9rJS1tL6ZRr209uePXscDViOvotwy5YtmDBhAvbv3w/DMPDAAw9gzJgx2LBhA8aPH48DBw6goKAAr7zyCgYMGODoaJLfRUhEp8L8HQq31dhY7wN+0sfAr4YYaJeTvl/ATkTprUkfNNq9e3e89957WLVqFVauXIkxY8YAAHr16oXa2lps3LgRS5YscTxcfVHDF6eRXiKRCG688Ub201RL6PevgGDMuxYu+5uN9jkeLLvOi9+eb2bEcNUS+mUy9tOb027eZl7HCRmG6w+TpxQYhoGioiL201Qm94tYgsdXKjy6UqFdK+CP3zLxvbM8afWg0JOVyf1aAvbTm9NufC9CIsoIIoLZWwUVC23sCAMVAwzcO8hAfnbmDFZE5D5t3ouw4VkWpJdwOIzy8nL201Sm9Vt/WHD5322M/oeN3l/zYPUYLx49z8zY4SrT+rU07Kc3p91cPyI0TdPtJVAKTNNEaWkp+2kqU/r544JfLVd4apVClzzgrREmru6SWceBR5Mp/Voq9tOb0248IiQi7YgIXt0s+PkiG744cM8gA3cONJDjzezBiojcp80RYcOj70kvoVAII0eOZD9N6dxv2X7BN/9q4/vv27igowfrx3px37lmixqudO5H7Kc7p91cPyLMyspyewmUgqysLJSXl7OfpnTsdyAquLdO4cX1Cn1PA+ZfYeKSQtc/R3SFjv3oc+ynN6fdeERIRGnNVoIX1ivct0TBVsCDQwzc2s9AltFy7lgRUfrgESE1q1AohNLSUvbTlC79PtylMORNC7d+rHBdVw82fdeLyQPMFj9c6dKPjo799KbNEWF2drbbS6AUZGdno6Kigv00le79doYEP1tk4/8+FRR38GDhtQaGfd31zwfTRrr3o+NjP7057cYjQiJKG3Fb8NRqhV8tV8g1gf8uNvHDXh4YGf7YBSLShzZHhMFg0O0lUAqCwSD69evHfppKx35z6xUGvG7hnjqFCb0MbBzrxYTeBoero0jHfuQc++nNaTfXjwhzcnLcXgKlICcnB5WVleynqXTqt8UvmLLQxl/+JfjWmR68fpmJ/u04VB1POvWjxmM/vTntxiNCInJF2BL89wqFX3+i0CEHmD7MRHn3zH8KOxHpTZsjwkAg4PYSKAWBQABFRUXspyk3+4kIqrco9K6y8PhKhTsHGFhf7sXYswwOVw7x9ac39tOb026uHxHm5ua6vQRKQW5uLqqrq9lPU271W3NQcFutjX/uFFzdxYMnS02cVcChqrH4+tMb++nNaTceERJRszscE0xbpvDsGoXu+cBvhpv4TmfXb6ATETWaNkeEfr/f7SVQCvx+PwoKCthPU6eqnxLB/25Q6FVt4aX1Cg8PNbCqzMvh6iTx9ac39tOb026uHxG2adPG7SVQCtq0aYPa2lr209Sp6Fe3T2HixwqL9wmuP8uDx4eZKGzD48CmwNef3thPb0678YiQiJrU3ojgnjobL28QnNMOeGa4iQvO5B0rIsoMPCKkZuX3++HxeNhPU83Rz1KCp1fb6Fll4Y2tgt+eb2DJdV4OV82Arz+9sZ/etDkizMvLc3sJlIK8vDzU19ezn6aaut/7OxUm1dhYcwi4qbeBh4oNnJ7D48Dmwtef3thPb067uT5g8bk3evJ4PCgoKGA/TTVVv/qg4M5FNqq2CEq/7kHdKBNDOvDvRHPj609v7Kc3p91cv3fPB63pKRAIoG3btuynqZPtF7UEDy+30bvawoJdgpkXmfjoGg5Xpwpff3pjP71p86DR/Px8t5dAKcjPz4fP52M/TZ1Mv7f/pTB5oY1/BYDb+xt44FwDBdkcrE4lvv70xn56c9rN9TtYLn0TI50kEYHf72c/TaXSb5NPcOVcC1e/Y6N7vgefjPHiiRKTw5UL+PrTG/vpzWk31wesYDDo9hIoBcFgEJ07d2Y/TTWmXzAh+MViG/1nWVhzSPDGZSbmfcdEn9M4WLmFrz+9sZ/enHZz/YiQz8DSU0FBAT/70piTfiKCP38q+NliGweiwD2DDPz8HAO5Xg5WbuPrT2/spzenc4vrd7Bs23Z7CZQC27axZs0a9tPUifp9ckBw8ds2rn/PxrAOHqwr92LqEJPDVZrg609v7Kc3p91cH7BCoZDbS6AUhEIhlJaWsp+mjtXvYFQw8WMbg9+0sDcqeOc7Jl7/thfd8jlYpRO+/vTGfnpz2o1vlUNEsJXg5Q1H3uImroCp5xqY1M9AtsnBiojoi7R5qxzLstxeAqXAsizU1tayn6a+2K92j8Kwt2zc/JGNK7t4sGGsF3cMNDlcpTG+/vTGfnpz2s31ASsSibi9BEpBJBJBeXk5+2kqEolg9Pif4gfvWRj+FxsCwcfXmJh5sRdntuZgle74+tMb++nNaTceERK1MAkleHq1woPLFLIN4OFiAzf2MmAaHKyIiE6ER4TUrCzLwrx589hPM+9uVxj4uoWfL1a4uNV2rB0D3NzH5HClGb7+9MZ+etPmiDAajbq9BEpBNBpFRUUF+2lia0Aw+l0LI/5uo0OOBx+NjOPTh65Ca4m5vTRKAV9/emM/vTntxiNCogwWsQSPrVR4bKVCu1bAr4eZ+N5ZHsfvBk9ERF+mzRFhIpFwewmUgkQigerqavZLUyKCNz5T6FNt4ZEVCrf3N7BhrBfXn23A4/Gwn+bYT2/spzen3VwfsOLxuNtLoBTE43FUVlayXxpaf1gw8u82xvzDRt/TPFg9xotHzzORl/X5XSv20xv76Y399Oa0G48IiTKEPy745TKF36xW6JIHPFVq4qouPA4kImpK2hwRcoLXUzwex0svvcR+aUCJ4PcbFXpVWXhurcK0IQbWlHlxdVfjmMMV++mN/fTGfnpz2s31AYtn0Hri1xCkh2X7Bd/8i40fLLBx4ZlHnsJ+72ATOSd4U2b20xv76Y399Oa0G48IiTS0Pyq4t05hxnqFfqcBTw838a1Orn++RESU8bQ5IozF+BweHcViMVRWVrLfKWYpwW/X2OhZZeG1LQpPlRpYPtrb6OGK/fTGfnpjP7057eb6gGXbtttLoBTYto3a2lr2O4U+3KUw9E0Lk2oURnfzYONYL27rb8KbwlPY2U9v7Kc39tOb0248IiRKcztCgp8vsvF/nwqKO3jw7HAD533d9c+NiIhaJB4RUrOKxWKYNm0a+zWjmC14bIWNXlUW3t0hePlCEwuvNZtkuGI/vbGf3thPb067eZt5HSeklHJ7CZQCpRS2b9/Ofs3k7/UKt9fY2BIAJvYzMO1cA19r1XTPs2I/vbGf3thPb0678YiQKI186hdMqbXx122Cb53pwdPDTfRvxweFEhGlC22OCPlu4nriu8E3rVBCcF+djX6zLKw4IKi61MT8K5tvuGI/vbGf3thPb067uX5ESNSSiQiqtwjuXGRjbxT42UADd59joE0W71oREemMR4RELll9UHBbjY33dgmu6erBkyUmuhdwsCIiSmfaHBFGIhG3l0ApiEQiuPHGG9kvBYdjgsm1Nga9YWFHWPD3y028NcJ7Socr9tMb++mN/fTmtJvrR4SG4fqMRykwDANFRUXs1whKBL/bKLh7sY2wBTw81MCUAQayzVN/14r99MZ+emM/vTntxiNColNg8V6FSTUKi/cJrj/Lg8eHmShsw+NAIiLdaHNEGA6H3V4CpSAcDqO8vJz9TmBvRDBhgYVhb9mI2YIPrjLxx0u8rg9X7Kc39tMb++nNaTfXjwhN03R7CZQC0zRRWlrKfsdgKcFzaxUeWKpgeIDfnm/g5t4GzBTeN7A5sJ/e2E9v7Kc3p914REjUxN7fqTCpxsaaQ8BNvQ08VGzg9Jz0GKyIiOjkaHNEGAqF3F4CpSAUCmHkyJHs9wX1QcF351v41hwbBdkeLLnOi+cvMNNyuGI/vbGf3thPb067uX5EmJWV5fYSKAVZWVkoLy9nPwBRSzB9lcIjKxQKsoDfX2xi3NkeeDzpN1g1YD+9sZ/e2E9vTrvxiJAoRSKCt7cdeabVtiAweYCB+wcbKMhO38GKiIhODo8IqVmFQiGUlpa22H6bfIIr59m45h0bZxd4sKrMi18PM7UZrlp6P92xn97YT2/aHBFmZ2e7vQRKQXZ2NioqKlpcv2BC8NByhcpVCoWtgTe/beLarul9HHg0LbVfpmA/vbGf3px24xEhkQMigj99KvjZIhsHY8AvBhn42UADuV69BisiIjo52hwRBoNBt5dAKQgGg+jXr1+L6LfygOCit238v/dslJ7hwfpyLx4419R6uGpJ/TIR++mN/fTmtFujBqxXXnkFHo8Hs2fPBgDs3bsXl19+OXr06IH+/fvjgw8+aPRCc3JyGv3vkPtycnJQWVmZ0f0ORgUTP7Zx7psW9kcF715hYtZlXnTN13ewatAS+mUy9tMb++nNaTfHR4Rbt27F9ddfDxHBXXfdhVGjRuFHP/oRunTpgmnTpqGurg7XXXcdPvvsM0ffwsgjQkpXthK8vEFwT52NuAKmnWtgUn8DWWnyFHYiInJPkx4RKqVw44034plnnkGrVq2SP19VVYWf/OQnAIDi4mJ06tQJCxYsaNRCA4FAo34/pYdAIICioqKM61e7R+G8tyzc/JGNq7p4sHGsFxUDzYwbrjK1X0vBfnpjP7057eZowKqsrMT555+PIUOGJH/uwIEDSCQS6NixY/LnunXrhm3bth31z4jFYvD7/V/68UXRaBTRaBQAEIlEEIvFABx5U8WG61AohHg8nrxOJBIAjpyHWpYF4MjGG679fj9s205eK6UgIvD7/RARKKWS67BtO3ltWVbyA2hZVvK8NZFIJL89Mx6PJ69jsVjyzR9jsRgikUjG7yk3NxevvvoqDMPIiD1tPRTFD963MPwvNiCCmmtMPHdeHKeZn+9Dtz0d7+9ebm4ufv/73yffUysT9pSJnY61p9zcXMycOTN5WpAJe8rETsfaU25uLn73u98lvxstE/aUiZ2OtSfHR7tyAqtWrZKSkhKJx+MiInLRRRfJm2++Kfv375fs7Owv/d7y8nJ5+eWXj/rnTJ06VQB85ccNN9wgIiJTpkyRKVOmiIjIhAkTZOrUqSIiUlZWJtOnTxcRkREjRsiMGTNERKSkpESqqqpERKRv374yd+5cEREpLCyUmpoaERHJz8+X1atXy7+PQaW+vl58Pp8AEJ/PJ/X19dLwIVi9erXk5+eLiEhNTY0UFhaKiMjcuXOlb9++IiJSVVUlJSUlIiIyY8YMGTFihIiITJ8+XcrKypL7nDBhAvekyZ5ilpILf/WGZP+PX9rPjMv5d/+v3D91mtZ7ysRO3BP3xD1xT+myp4Z9+Hw+OZ4TDljPPfecdOzYUbp27Spdu3aVVq1aSYcOHeS5556T1q1by65du5K/t7i4WN59992j/jnRaFR8Pl/yR8MCN23aJCIikUhEIpGIiIiEw2GJRqMiIhIKhZLXwWBQYrFY8rph6AsEApJIJERExO/3J699Pp9YlpW8tm1blFLi8/lEKSW2bSc/QJZlJa8TiYT4/f7kdSAQEBGReDwuwWBQRERisVjyOhqNSigUSl6Hw+GM35PP55O8vDzZu3evtnv666cR6V0VF2NGTG5eEJUDEZVxnUSO/nevod++ffsyZk+Z2OlYe2rod+DAgYzZUyZ2OtaeGvodPHgwY/aUiZ2Otadt27Y5GrAa/Rysiy++GJMnT8aoUaMwfvx4dOvWLflF7qNGjcLWrVsb9UXuBw8exGmnndaYJVAasG0b69evR+/evZPHTLrYGhBULLTx5lbBBR09eGa4iXPaZ9bXWJ2Izv2I/XTHfno7dOgQ2rVrd8Ivcj+pJ7k/9thjuOGGG9CjRw9kZ2fj1VdfbfSbV/Ivl55M00S/fv3cXkajRCzBYysVHlup0K4V8H/fMvFfZ+n3FPamoGM/+hz76Y399OZ0bnH9Se719fUoKipyYwl0EnR6zIaI4M2tR+5a7QwDdwwwcO9gA3lZLW+waqBTP/oq9tMb++lt+/bt6Ny5c/PewWoKeXl5bi+BUpCXl4f6+vq077fukOC2Whv/2CH4TmcP3r3CRI+2LXewaqBLPzo69tMb++nNaTfXB6yWeDyTCTweDwoKCtK2nz8ueHCZwtOrFbrkAX8ZYeKqLi3zOPBo0r0fHR/76Y399Oa0m+vvRcgHrekpEAigbdu2addPiWDmRoWeVRaeX6cwbYiBNWVeXN3V4P+YfUG69iNn2E9v7Kc3p91cv4OVn5/v9hIoBfn5+fD5fGnVb+k+waQaG7V7BWO7e/DEMBOd8zhUHU069iPn2E9v7Kc3p91cv4Pl0tfY00mSLzyd1237o4KbP7RRPNtCICF470oTr13q5XB1HOnUjxqP/fTGfnpz2s31Aavhkfikl2AwiM6dO7vaz1KC366x0eM1C1VbFH5TamD5aC8u7uT6X+u0lw79KHXspzf205vTbq4fEfJbVPVUUFDg6mdfH+xSmFRjY9VBYEIvDx4pNtEhl3esnHK7H50c9tMb++nN6dzi+qf6DW+iSHqxbRtr1qw55f12hATX/9PCRW/byPV6sGiUiRkXejlcNZJb/ahpsJ/e2E9vTru5PmA1vDs26SUUCqG0tPSU9YvZgkdX2OhVZWH+TsH/Xmii5hoTxR1c/yuspVPdj5oW++mN/fTmtJvrT3Lnk2zpRP62TWFyrY0tAWBSPwPThhhom807VkREdOo5nV9c//Tfsiy3l0ApsCwLtbW1zdrvU7/g6nkWrpxno3OeBytHe/Fkqcnhqgmcin7UfNhPb+ynN6fdXB+wIpGI20ugFEQiEZSXlzdLv1BCcF+djb7VFj45KKi+1C0KKmoAAB2cSURBVMQ/rjDRrx0Hq6bSnP2o+bGf3thPb0678YiQ0oaIoGqL4M5FNvZFgZ8PNHD3IAOtvRysiIgoPfCIkJqVZVmYN29ek/VbdVBwyRwb//VPG0NO92BtmRe/HGpyuGomTd2PTi320xv76U2bI8JoNOr2EigF0WgUFRUVJ93vcExwe42NwW9Y2BkW/P1yE7NHeNG9gINVc2qqfuQO9tMb++nNaTceEZIrlAhe2SD4RZ2NiA08MNjA7f0NZJscrIiIKH1pc0SYSCTcXgKlIJFIoLq6OqV+i/cqlLxl48YPbYwo8mBDuRc/O8fkcHUKnUw/ch/76Y399Oa0m+sDVjwed3sJlIJ4PI7KyspG9dsTFvxogYVhb9lIKMGHV5t49VtedGrDwepUS6UfpQ/20xv76c1pNx4RUrNLKMFzaxWmLlUwPcBDQw3c1NuAaXCwIiIivWhzRMgJXk/xeBwvvfTSCfu9t1Nh8BsWptQqfO8sAxvHevHTviaHK5c57Ufpif30xn56c9rN9QGLZ9B6OtHXEGwLCsb+w8Ilc2y0zfZg6XVe/M83TbTP4WCVDvg1IHpjP72xn96cduMRITWpqCV44hOFR1YotM0GHh9mYtzZHng8HKyIiEh/2hwRxmIxt5dAKYjFYqisrEz2ExH85V8K/WZZeHCZwq39DGwY68UNPQwOV2noP/uRXthPb+ynN6fdXB+wbNt2ewmUAtu2UVtbC9u2sfGw4Iq5Nq59x8bZBR6sKvPi18NMFPBNmdPWF/uRfthPb+ynN6fdeERIKQsmBA8tV6hcpVDYGniy1MS1XXkcSEREmYtHhNRsRAS/WxdDp5cP4Terbdw32MDaci9GdeNxoC5isRimTZvG15+m2E9v7Kc3p928zbyOE1JKub0EaoSVBwSTamx8uNtAN98mzP3RQPTqkO32sqiRlFLYvn07X3+aYj+9sZ/enHbjESE5cjAquH+pwvPrFHq1BZ4ebuKyQtdvgBIREZ1S2hwR8t3E05utBC+ss9GzysIfNik8MczAyjFefLN9nO8Gr7FoNMp+GmM/vbGf3px2c/2IkNJXzR6FSTU2lu0HftDDg0fPM9Gx9ZGvseL3vhARER0bjwjpK3aFBXctsvGHzYIhp3vwzHADpWe4frOTiIjIddocEUYiEbeXQP8WtwXTP7HRq8rC37cLXrzAxKJrzaMOV5FIBDfeeCP7aYr99MZ+emM/vTnt5voRoWG4PuMRgHe2K9xea2OjD7ilj4FfDjVwWqtjP3LBMAwUFRWxn6bYT2/spzf205vTbjwibOE+8wsqFtqY/S/BhR09eGa4iYHt+SwrIiKio9HmiDAcDru9hBYpbAmmLrXRd5aFuv2CP11i4v2rnA9X4XAY5eXl7Kcp9tMb++mN/fTmtJvrR4Smabq9hBZFRPDG1iN3rXaHgTsGGrhnkIG8rMbdtTJNE6WlpeynKfbTG/vpjf305rQbjwhbkLWHBLfV2Ji/U3BFZw+eKjXRoy2PA4mIiJzS5ogwFAq5vYSM54sLKmptnPO6ha1BwV9HmJhzufekhqtQKISRI0eyn6bYT2/spzf205vTbq4fEWZlZbm9hIylRPCHTYK7FtsIJIAHhxioGGAgx3vyd62ysrJQXl7OfppiP72xn97YT29Ou/GIMEMt3SeYWGNj4V7Bd7t78OthJjrn8TiQiIjoZPCIsIXaFxHc9KGF4tkWQpbgvStN/PlSb5MPV6FQCKWlpeynKfbTG/vpjf30ps0RYXZ2tttLyAiWEjy/TuH+JQoA8JtSAz/ta8BrNM9dq+zsbFRUVLCfpthPb+ynN/bTm9NuPCLMAB/sUphYY2P1QWBCLw8eKTbRIZfHgURERE1NmyPCYDDo9hK0tT0o+N4/LVz0to3WXg8WjTIx40LvKRmugsEg+vXrx36aYj+9sZ/e2E9vTru5fkSYk5Pj9hK0E/v/7d17dJT1ncfxz8wkgZibV265EINQIGBUQIPoWeIFoyVChXDsBXQFlG7XC2nV065daNW9qKTHI7uromsPxe0uwQsIFhRRhBIRFIpgq4DkJuEi6DyTydyeeb77R5opSEKeQJJnvpPP6xzOmSRPyO/xDZ4vz2/mmajgN59YeGyHhbRk4KW/82DWUBfcrp67atW3b19UVlayn1Lspxv76cZ+utntxi1CZdbUWXigOooaH3BvoRsLxriRlcLtQCIiop6gZovQ5/M5vQQV9nkFZetMTF4XxeB0F/40LQmV4z2ODVc+nw85OTnspxT76cZ+urGfbna7Ob5FmJqa6vQS4po/IviXnRae2mVhwDnAihs8uC3fBVcPbge2JTU1FVVVVeynFPvpxn66sZ9udrtxizBOiQiWfyH42dYojgaBh4vceLjIjXO64C7sREREdGbUbBEahuH0EuLOJ8cF162J4vYNUYy90IU/T0/Cr8Z44mq4MgwDmZmZ7KcU++nGfrqxn252uzm+RZiWlub0EuLGNyHBP39k4T8/tXBJJrDuZg8m5Tg+A7cpLS0N1dXV7KcU++nGfrqxn252u3GLMA5YInjpM8HPt0URiAILrnDjvkI3Ujzxc8WKiIiIuEWoxtYjFopXRjFnUxQ35bjw+Ywk/OxST9wPV4ZhwOVy9fp+WrGfbuynG/vppmaLMD093eklOOJwc8sVq5c+F1x2AbC5zIMJAxyfd21LT09HfX19r+2nHfvpxn66sZ9udrs5PmA5fbuBnhaxBP+xx8KCjywkuYH/muDG3OFueLrpTZm7i8vlQmZmZq/rlyjYTzf20439dLPbzfFLJr3pRmsbvrRw+asmfrrVwg8vcePzGUmYN9KjbrgCWrplZWX1qn6JhP10Yz/d2E83NTcazcjIcHoJ3a6uSfDTD6JYcUAwob8L26d6cPmF+oaqE2VkZMDr9faKfomI/XRjP93YTze73Ry/guXQixh7RNAUPPpxFMOXm/jjYcGyEg82lekfroCWboZhJHS/RMZ+urGfbuynm91ujg9YTU1NTi+hy4kIVtVaKFxh4tEdFu4tdOOz8iT88BJ3wuy5NzU1ITc3NyH79Qbspxv76cZ+utnt5vgWYaLdA+vzbwT3V0extkFwU44Lb5Z68J1zE2OoOlFmZib/9aUY++nGfrqxn2525xbHr2BFo1Gnl9AlfGHBw1ujGPWKic+8gpWTPPhDgg5XQEu3PXv2JEy/3ob9dGM/3dhPN7vdHB+w/H6/00s4KyKCl/dZ+E6ViWf2WPjl5W58Oj0Jtw5OnO3Atvj9fowfP159v96K/XRjP93YTze73fhWOWdh5zHBvX+MYvNhwfSLXXjqKg8GZyTuUEVERNTbqXmrHNM0nV5Cpx0PCv5hcxRjXjNxPCRYf4sHVTck9arhyjRNVFdXq+xH7Kcd++nGfrrZ7dbhgBUMBjF16lQMGzYMRUVFuPHGG7Fv3z4AwJEjR1BaWoqhQ4di1KhReP/99zu90EAg0OnvcUrUEjz35yiGLjfx8j4Li65yY+e0JFyf7fic2uMCgQDKy8tV9aO/YT/d2E839tPNbrcOtwiDwSA2bNiAm2++GS6XC4sXL8aKFSvw3nvv4a677kJeXh4WLlyIbdu24Xvf+x4OHDiA5OTkDn+wti3CPx6ycO+WKHYcA/5+mAv/Os6D/uf0nitWRERE1IVbhH379sUtt9wSe8J2cXExampqAADLly/HvHnzAADjxo3DoEGDsHHjxk4tNN4vkTY2C2a+a+KaN6LwuFyovtWD//67pF4/XJmmiXXr1sV9P2ob++nGfrqxn25dtkX4bU8//TSmTJmCY8eOIRKJYMCAAbGv5efno66urs3vC4VCMAzjpF8A8M033wBouVIWDAYBtFx+C4VCAIDm5ubYY7/fj3A4HHsciUQAtNz0q/WEfT5f7LFhGLGXUxqGAcuyTrqDrmVZsXVEo9HYY9M0cewbH57aFcWw/zPxh3oLS671YPN3BaPTWy4NhsPh2CsJQqEQmpubY49bLx/G2zm1vn+SaZqxG6VFIpHYeXTmnILBIB544AF4vd6EOadE7NTeObX2a/25iXBOidipvXNq7df6cxPhnBKxU3vnFAwGcf/998d+biKcUyJ2au+cbG/tSic8/vjjUlxcLH6/X7766itJSUk56evl5eXy4osvtvm9CxYsEACn/Jo5c6aIiMyfP1/mz58vIiKzZ8+WBQsWiIjI9OnTZdGiRSIiMmnSJFmyZImIiBQXF8vy5ctFRGTkyJGydu1aERHJzs6WLVu2iIhIRkaG7N69W/66DSr19fXi9XoFgHi9Xqmvr5fW/wS7d++WjIwMERH5zbrdkvT4p+JZEpZb/2e/fOfyq0REZPny5VJcXCwiIkuWLJFJkyaJiMiiRYtk+vTpsfOcPXt23J3Tli1bJDs7W0RE1q5dKyNHjuQ58Zx4TjwnnhPPiefUyXNqPQ+v1yunY3vAevLJJ2XMmDHy9ddfxz53zjnnSGNjY+zjcePGydtvv93m9weDQfF6vbFfrQusr68XEZFAICCBQEBERJqbmyUYDIqIiN/vjz1uamqSUCgUexwOh0VExOfzSSQSERERwzBij71er5imGXscjUbFsizxer1iWZZEo9HYfyDTNOVPXxoydV1E8HxYrnktKLuOWRKJRMTn84mISDgclqamJhERCYVCscfBYFD8fn/scXNzc9ycU+vjSCQihmHEHp/tOYXDYVm2bFns90mEc0rETu2dUzgclqVLl8Z+n0Q4p0Ts1N45tfZr/X0S4ZwSsVN75xQOh+W3v/1t7PdMhHNKxE7tndPRo0dtDVi27oNVWVmJl19+GevXr8d5550X+/ydd96J/Pz82JPcp06dipqamk49yf3gwYMYOHBgh8d3p2ZT8G87LTyxy8JFfYGnrvJgRoEroW8Uerb8fj9uuOEGrF+/HmlpaU4vhzqJ/XRjP93YT7fGxkYMGjSowye5dzhgNTQ0IDc3FwUFBcjIyAAA9OnTB1u3bsXhw4cxc+ZMHDhwACkpKVi8eDFKSkpsLTAeXkUoIni1RlDxQRSHmoGfXerGLy5zIy2ZgxURERGdyu780uGbPefk5LT7ppT9+/fHW2+9dearBGJPNutpn34tuG9LFO8cFHw314V3bvHgkiwOVnaFw2EsXboUs2bNQkpKitPLoU5iP93YTzf2083u3OL4HTJbn8HfU7xhQUV1FEWvmKhtEqy+yYPVpUkcrjopEomgqqqqx/tR12A/3dhPN/bTzW63XvNehJYIfrdX8PCHUfgiwCOXu1Ex2o0+Hg5WREREZI+a9yJsvTdFd/roqGDCqiju3BhFySAXPitPws8v83C4OguhUAiVlZU90o+6Hvvpxn66sZ9udrs5PmC13sCrOxwNCO7eZGLc6yb8puC9yR78/rok5KRzsDpb0WgU1dXV3dqPug/76cZ+urGfbna7JeQWoWkJnv2zhV9utwAAj451Y94IN5LcHKyIiIjozPXaLcL3Gy1c8ZqJ+7ZYKC9w4fMZSfjHQg+Hqy4WCoWwcOFCXuJWiv10Yz/d2E83u906vE1Dd7Msq0t+n4YmwYMfRvG/+wXF/Vz4cKobYy9yfH5MWJZloaGhocv6Uc9iP93YTzf2081uN/VbhKGoYNEuC4/vtJCeDPz7lR7MGuqCm3dhJyIioi6mZouw9Z2vz8T2oxYKV5hY8JGFe4a78fmMJNw5zM3hqgcEg0FUVFScVT9yDvvpxn66sZ9udrs5vkV4Nv5pm4UUN/CnaUkYeR6HKiIiIooParcIvWHBRb8zsegqN+4d5emGFRIRERGdTM0WYSAQOKPvW1sviFjArYMdP4VeKRAIYM6cOWfcj5zFfrqxn27sp5vdbo5PJ273mS1hZa2Fyy8ABmdwa9AJbrcbOTk5Z9yPnMV+urGfbuynm91uKrcIw9GW7cGK0W4sGMPtQSIiIuoZarYIm5ubO/09GxsFRgSYku/48nut5uZmlJeXn1E/ch776cZ+urGfbna7OT6heDydvwK1slYwOB0oOr8bFkS2eDwejB8//oz6kfPYTzf20439dLPbTd0WoYgg7/cmbst34+mr+YeTiIiIeo6aLUK/39+p4z/+CmjwA1MG88ntTvL7/bjppps63Y/iA/vpxn66sZ9udrs5PmAlJyd36viVtRbOTQGuHcgBy0nJyckoLy/vdD+KD+ynG/vpxn662e2mbouw6JUIRp/vwrIS1TehJyIiIoUScovwgCHYdRyYwpuLOs7v92P8+PG8xK0U++nGfrqxn25qtghTUlJsH7uytuW9B0tzuD3otJSUFFRUVHSqH8UP9tON/XRjP93sdlO1RViy2kRqEvBmKbcHiYiIqOep2SJsamqyddzxoGDTIeGrB+NEU1MTCgsLbfej+MJ+urGfbuynm91ujg9Yffv2tXXcmnpBVICyPMeXTGjpVllZabsfxRf20439dGM/3ex2U7NFOO1tE182Ax9M4fYgEREROUPNFqHP5+vwmKApWNfA7cF44vP5kJOTY6sfxR/20439dGM/3ex2c3zASk1N7fCYDQcFfpO3Z4gnqampqKqqstWP4g/76cZ+urGfbna7qdginLcpivUHLeydkQSXi1exiIiIyBlqtggNwzjt10UEq+stlOW5OVzFEcMwkJmZ2WE/ik/spxv76cZ+utnt5viAlZaWdtqv7zgGfOkHyvI4XMWTtLQ0VFdXd9iP4hP76cZ+urGfbna7Of6SPI/Hc9qvr6q1kMU3d447Ho8HhYWFTi+DzhD76cZ+urGfbh3NLa0cv4LV0aW2N+os3JzrQrKbA1Y8MQwDLpeLl7iVYj/d2E839tNNzRZhenp6u19raBJ8/BVvLhqP0tPTUV9ff9p+FL/YTzf20439dLPbzfHJ5XRPXF9dZ8HjAm7O5dWreONyuZCZmckXHijFfrqxn27sp5vdbo4PWKe7YdcbdYJrB7hwXh/+IYw3Pp8PWVlZvFGeUuynG/vpxn66qbnRaEZGRpuf90cE7xwUvnowTmVkZMDr9bbbj+Ib++nGfrqxn252uzk+YLV3n9O3vxSEokAZ794el0QEhmG024/iG/vpxn66sZ9udrs5Pr00NTW1+fk3ai0MPxcYmsUrWPGoqakJubm57faj+MZ+urGfbuynm91ujt8Hq63bzFsiWF0vuGOo4/MftSMzM5P/+lKM/XRjP93YT7eO3t6vleMTTDQaPeVzHx4RHAkAtw7m1at4FY1GsWfPnjb7UfxjP93YTzf2081uN8cHLL/ff8rn3qgTXNAHGN+PA1a88vv9GD9+fJv9KP6xn27spxv76Wa3W1xuEb5RZ+GWXBc8vHt73OIblerGfrqxn27sp5uaLULTNE/6uMYn+OQ4cCtfPRjXTNNEdXX1Kf1IB/bTjf10Yz/d7HZzfIoJBAInffxGrYVkNzAph1ev4lkgEEB5efkp/UgH9tON/XRjP93sdnOJQy9lMAwDWVlZ8Hq9J11um/Rmy2T41i2O714SERERnaS9+eXbHL+CdeKlNiMseK+Rd2/XwDRNrFu3jpe4lWI/3dhPN/bTTc0WYTAYjD1e1yCIWHz+lQbBYBAVFRUn9SM92E839tON/XSz2y2utghnvWti53HBrmnJTiyJiIiI6LTUbBFGIhEAgGkJ1tQLbs1zfElkQyQSQVVVVawf6cJ+urGfbuynm91ujk8z4XAYAFB9WHA8BJTx7u0qhMNhVFZWxvqRLuynG/vpxn662e0WN1uED22NYuleCwd/mAS3i0MWERERxR81W4Stk+CqWguT81wcrpQIh8N44YUX+C8wpdhPN/bTjf10s9vN8QErEolgr1fwmRco4/Ov1OBzCHRjP93YTzf2081ut7jYInyhJg2/2G7h2MwkpCXzChYRERHFJzVbhKFQCKvqBDcMcnG4UiQUCqGyshKhUMjppdAZYD/d2E839tPNbjfHB6xjzVFsPiR89aAy0WgU1dXViEajTi+FzgD76cZ+urGfbna7Ob5FuOTjrzF3exoafpCE7DQOWURERBS/1GwRrqk1MeZCF4crZUKhEBYuXMhL3Eqxn27spxv76aZmi3BDI/jmzgpZloWGhgZYluX0UugMsJ9u7Kcb++lmt5vjW4R4+it89IPzccWFHLKIiIgovqnZIhyQKrj8AqdXQZ3Fd4PXjf10Yz/d2E83u90cH7Buyha4ePd2IiIiSiCObxFW7f4a0wvPdWIJRERERJ1id4swqQfXdJLWuW540lcwDMcvpFEnBQIBPPjgg3jyySeRmprq9HKok9hPN/bTjf10O3LkCIC/zTHtcewK1hdffIEhQ4Y48aOJiIiIzsr+/ftRUFDQ7tcdu4J1/vnnAwDq6upaXk1IqhiGgdzcXNTX15/2EinFJ/bTjf10Yz/dvF4v8vLyYnNMexwbsNzulm3BrKws/gFTLDMzk/0UYz/d2E839tOtdY5p9+s9tA4iIiKiXoMDFhEREVEX8yxcuHChYz/c48HEiRORlOTYTiWdBfbTjf10Yz/d2E83O/0cexUhERERUaLiFiERERFRF+OARURERNTFOGARERERdbFuH7D27t2Lq6++GsOGDcO4ceOwZ8+eNo978cUXMXToUAwZMgRz585FJBLp7qWRDXb6bdiwAVdeeSVGjhyJwsJCPPTQQ7Asy4HV0rfZ/fsHtLztw3XXXYdzz+V7g8YLu/0++eQTTJw4ESNGjMCIESPw6quv9vBKqS12+lmWhYqKCowcORKXXnopSkpKsG/fPgdWSye67777kJ+fD5fLhZ07d7Z73GlnF+lmJSUl8tJLL4mISFVVlYwdO/aUY7744gsZOHCgNDY2imVZUlZWJosXL+7upZENdvp9/PHHsn//fhERCQQCMmHChNj3kLPs9Gu1aNEimTNnjmRlZfXQ6qgjdvr5/X65+OKLZdOmTSIiYpqmHDlypCeXSe2w0++1116TK6+8UsLhsIiIPProo1JeXt6Ty6Q2bNy4Uerr62Xw4MGyY8eONo/paHbp1gHr8OHDkpGRIZFIRERELMuS/v37y969e0867oknnpB77rkn9vGaNWtkwoQJ3bk0ssFuv2/7yU9+IgsWLOiBFdLpdKbf7t275dprr5V9+/ZxwIoTdvstWbJEvv/97zuxRDoNu/1ef/11KSoqEsMwxLIsefDBB2X+/PlOLJnacLoBq6PZpVu3COvr6zFw4MDYfSJcLhfy8vJQV1d30nF1dXUYPHhw7OP8/PxTjqGeZ7ffiQ4dOoQVK1Zg8uTJPbVMaofdfpFIBHPnzsVzzz0Hj8fjxFKpDXb7ffrpp+jTpw8mT56Myy67DLNmzcLRo0edWDKdwG6/srIyTJw4EQMGDMDAgQPxzjvv4Ne//rUTS6ZO6mh24ZPcqcsYhoGysjI89NBDGDt2rNPLIZt+9atf4bbbbsOIESOcXgqdAdM0sX79ejz33HPYsWMHsrOz8eMf/9jpZZFN27dvx+7du/Hll1/i4MGDuP766zFv3jynl0VdoFsHrNzcXDQ2NsI0TQAtT6Ktq6tDXl7eScfl5eWhtrY29nFNTc0px1DPs9sPAHw+H0pLSzFlyhRUVFT09FKpDXb7bdy4Ec888wzy8/NxzTXXwDAM5Ofn8yqIwzrz/8+SkhJkZ2fD5XLhRz/6ET744AMnlkwnsNtv6dKlsReXuN1u3HHHHXj33XedWDJ1UkezS7cOWP369cMVV1yBZcuWAQBeeeUV5OTk4JJLLjnpuGnTpmHVqlU4dOgQRATPPvssbr/99u5cGtlgt19TUxNKS0tRWlqKRx55xImlUhvs9tu0aRNqa2tRU1ODzZs3IzMzEzU1NbjoooucWDb9ld1+M2bMwLZt22AYBgDgzTffRFFRUY+vl05mt19BQQE2bNiAcDgMAFi9ejVGjRrV4+ulzutwdum+p4a1+Mtf/iLFxcUydOhQGTNmjOzatUtERGbPni0rV66MHff8889LQUGBFBQUyF133RV7RQU5y06/xx57TJKSkqSoqCj267HHHnNy2fRXdv/+tTpw4ACf5B5H7PZbunSpFBYWyujRo6W0tFTq6uqcWjKdwE6/YDAoc+bMkeHDh8vo0aPlxhtvjL0qm5xz9913S3Z2tng8HunXr58MGTJERDo3u/C9CImIiIi6GJ/kTkRERNTF/h9R9kb7CDY1DQAAAABJRU5ErkJggg==\" />"
      ]
     },
     "execution_count": 60,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sum(first_order_coeffs.*sol(0.01)[1:4])/(1/100) # the first derivative at different times\n",
    "derivs = []\n",
    "for i in 0:1/99:1\n",
    "    push!(derivs, sum(first_order_coeffs.*sol(i)[1:4])/(1/100))\n",
    "end\n",
    "plot(x,derivs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [],
   "source": [
    "# pyplot()\n",
    "# p = plot(1)\n",
    "# @gif for pt = 0:gap:1\n",
    "#     push!(p, 1, sol(pt))\n",
    "# end"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## KdV equation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = LinearOperator{Float64}(1,4,100,:Neumann,:Neumann;bndry_fn=(2,1));\n",
    "C = LinearOperator{Float64}(3,4,100,:Neumann,:Neumann;bndry_fn=(2,1));"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "KdV (generic function with 1 method)"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "function KdV(t, u, du)\n",
    "    du3 = zeros(du)\n",
    "    C(t,u,du3)\n",
    "    A(t, u, du)\n",
    "    @. du = -1*du*u\n",
    "    du .= du.+du3\n",
    "end\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "test_kdv (generic function with 1 method)"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "function test_kdv(N)\n",
    "    h = 1 / (N-1)\n",
    "    x = collect(0 : h : 1)\n",
    "    u0 = -(8(sech.(2(x+8))).^2 + 2(sech.((x+1))).^2)\n",
    "    tspan = (0., 10.)\n",
    "    prob = ODEProblem(KdV, u0, tspan)\n",
    "    sol = solve(prob, dense=false, saveat=[0.03, 10.0],force_dtmin=true)\n",
    "    return x, sol, prob\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\u001b[1m\u001b[33mWARNING: \u001b[39m\u001b[22m\u001b[33mInstability detected. Aborting\u001b[39m\n"
     ]
    }
   ],
   "source": [
    "x, sol, prob = test_kdv(100); #plot(x, sol(0.30))"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Julia 0.6.0",
   "language": "julia",
   "name": "julia-0.6"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "0.6.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}