iceSheet_v0.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "a000c010",
   "metadata": {},
   "source": [
    "# Glacier and ice sheet flow\n",
    "\n",
    "<img src=\"Figures/MargerieGlacierPic.jpg\" width=350 height=300 />\n",
    "<center><i>Margerie glacier (Alaska)</i></center>\n",
    "\n",
    "# Part 0: Information\n",
    "This is the notebook associated with the lecture on Finite-Difference Methods for modeling ice sheet flow, which was presented by Vincent Verjans on November 21st 2022. \n",
    "\n",
    "The first part of the notebook is a summary of the derivations that we made in class. The second part of the notebook is the model code.\n",
    "\n",
    "Feel free to use the code and play around with it. For example:\n",
    "- Investigate issues of numerical stability. How do the time step and the spatial resolution affect the stability of the model? How can you relate the Courant–Friedrichs–Lewy condition to this problem?\n",
    "- Vary the values of the ice rheology parameter (bigA in the code). Does it change the time needed to reach a steady-state? Does it affect the ice sheet profile?\n",
    "- Change the accumulation rate (bdot in the code). Try to implement a space-dependent accumulation rate.\n",
    "- Etc.\n",
    "\n",
    "Any questions: vverjans3@gatech.edu\n",
    "\n",
    "# Part 1: Derivation of the ice flow equation and Finite-Difference scheme\n",
    "\n",
    "We consider the problem of the flow of ice. This first part of the notebook gives the derivation of the finite-difference scheme for simulating glacier flow. The details of the derivation are explained during the lecture, and on the whiteboard. For those of you who are interested in the derivation of the momentum equations, I recommend:\n",
    "- \"Fundamentals of glacier dynamics\" by C. J. van der Veen\n",
    "\n",
    "The essence of glacier flow physics lies in the notion of conservation of momentum: stresses acting on ice (gravity, friction with the bed, etc.) balance each other and cause ice deformation. The momentum conservation results in a complex system of equations: the \"Stokes equations\". But the system can be simplified to the \"Shallow Ice Approximation\" (SIA) if:\n",
    "1. We use the hydrostatic approximation in the vertical (vertical stress varies much more along the vertical than horizontal direction)\n",
    "2. We consider an ice mass along a flowline (2D) and neglect lateral forces\n",
    "3. The ice is frozen to its bed (i.e., all the ice flow is due to deformation under gravitational forces)\n",
    "\n",
    "The shear stress in the vertical reduces to:\n",
    "$$\n",
    "\\tau_{xz}(z) = -\\rho g (h-z) \\frac{\\partial h}{\\partial x}\n",
    "$$\n",
    "where $g$ is gravity (9.81 m s<sup>-2</sup>), $\\rho$ is the ice density ($\\sim$ 917 kg m<sup>-3</sup>), $h$ is the ice thickness, $x$ is the horizontal coordinate, $z$ is the vertical coordinate. This is literally the force that shears the ice!\n",
    "<div>\n",
    "<img src=\"Figures/glacierSIAsketch.png\" width=\"500\"/>\n",
    "</div>\n",
    "\n",
    "The constitutive relation relates stress to deformation. In other words, it quantifies how much ice deforms under a given shear stres. Deformation is characterized by strain rates ($\\epsilon$ $[$m s<sup>-1</sup>$]$)). Glacier ice is a visco-plastic material:\n",
    "<div>\n",
    "<img src=\"Figures/graphGlenGlowLaw.png\" width=\"300\"/>\n",
    "</div>\n",
    "\n",
    "This is known as Glen's flow law:\n",
    "$$\n",
    "\\epsilon_{ij} = A \\tau^{n-1} \\tau_{ij}\n",
    "$$\n",
    "with $n\\approx 3$ and $A$ depends on ice properties (temperature, fabrics, impurities ...). Note: here, we assume $A$ constant and set $A=3.5\\times 10^{-24}$s<sup>-1</sup> Pa<sup>-3</sup>.\n",
    "\n",
    "By definition:\n",
    "$$\n",
    "\\epsilon_{xz} = \\frac{1}{2}\\left(\\frac{\\partial u}{\\partial z}+\\frac{\\partial w}{\\partial x}\\right)\n",
    "$$\n",
    "where $u$ is the horizontal velocity and $w$ is the vertical velocity. Because of hydrostatic equilibrium, $\\frac{\\partial w}{\\partial x} << \\frac{\\partial u}{\\partial z}$. And $\\epsilon_{xz}$ simplifies to:\n",
    "$$\n",
    "\\epsilon_{xz} \\approx \\frac{1}{2}\\frac{\\partial u}{\\partial z}\n",
    "$$\n",
    "\n",
    "Along a 2D flowline (in x-z dimensions), we have:\n",
    "$$\n",
    "\\left\\{\n",
    "    \\begin{array}\\\\\n",
    "        \\epsilon_{xz} = A \\tau_{xz}^{n} \\\\ \n",
    "        \\epsilon_{xz} = \\frac{1}{2}\\frac{\\partial u}{\\partial z} \n",
    "    \\end{array}\n",
    "\\right.\n",
    "$$\n",
    "\n",
    "And thus:\n",
    "$$\n",
    "\\frac{\\partial u}{\\partial z} = -2A\\left[\\rho g (h-z) \\frac{\\partial h}{\\partial x}\\right]^{n}\n",
    "$$\n",
    "\n",
    "Integrating and setting $u\\left(z=0\\right)=0$ (ice frozen to the bed):\n",
    "$$\n",
    "\\int_{0}^{z} \\frac{\\partial u}{\\partial z}\\,dz = \\int_{0}^{z} -2A\\left[\\rho g (h-z) \\frac{\\partial h}{\\partial x}\\right]^{n}\\,dz \n",
    "$$\n",
    "$$\n",
    "\\longrightarrow u(z) = \\frac{2A\\left(\\rho g\\right)^{n}}{n+1} \\left(\\frac{\\partial h}{\\partial x}\\right)^{n} \\left[\\left(h-z\\right)^{n+1}-h^{n+1}\\right]\n",
    "$$\n",
    "\n",
    "Now, we want to find the depth-averaged velocity $\\overline{u}$. Thus, we integrate a second time and we divide by the ice thickness ($h$).\n",
    "$$\n",
    "\\frac{1}{h}\\int_{0}^{z} u(z)\\,dz = \\frac{1}{h}\\int_{0}^{z} \\frac{2A\\left(\\rho g\\right)^{n}}{n+1} \\left(\\frac{\\partial h}{\\partial x}\\right)^{n} \\left[\\left(h-z\\right)^{n+1}-h^{n+1}\\right]\\,dz \n",
    "$$\n",
    "$$\n",
    "\\longrightarrow \\overline{u} = \\frac{-2A(\\rho g)^{n}}{n+2}\\left(\\frac{\\partial h}{\\partial x}\\right)^{n} h^{n+1}\n",
    "$$\n",
    "\n",
    "Now, we use the continuity equation. This equation ensures mass conservation. It stipulates that any change in ice column height must be balanced by a corresponding influx/outflux of ice and/or by a source term. The continuity equation is given by:\n",
    "$$\n",
    "\\frac{\\partial h}{\\partial t} = \\frac{-\\partial}{\\partial x}\\left(\\overline{u}h\\right)+\\dot{b}\n",
    "$$\n",
    "\n",
    "where $\\dot{b}$ is the snowfall rate ($[$m s<sup>-1</sup>$]$). It thus represents a source term in our equation. This has the form of an advection equation. But notice that $\\overline{u}$ depends on $h$. Thus, using the formulation derived for $\\overline{u}$, the continuity equation can be re-written as a diffusion equation:\n",
    "$$\n",
    "\\frac{\\partial h}{\\partial t} = \\frac{\\partial}{\\partial x}\\left(D \\frac{\\partial h}{\\partial x}\\right)+\\dot{b}\n",
    "$$\n",
    "\n",
    "And the diffusion coefficient is:\n",
    "$$\n",
    "D = \\frac{2A(\\rho g)^{n}}{n+2}\\left|\\frac{\\partial h}{\\partial x}\\right|^{n-1} h^{n+2}\n",
    "$$\n",
    "\n",
    "We evaluate our variable of interest on a discretized grid: ... $h_{i-2},\\: h_{i-1},\\: h_{i},\\: h_{i+1},\\: h_{i+2}$ ... (where $i$ is a space index)\n",
    "\n",
    "But we evaluate the diffusion coefficients on a staggered grid: ... $D_{i-\\frac{3}{2}},\\: D_{i-\\frac{1}{2}}, D_{i+\\frac{1}{2}}, D_{i+\\frac{3}{2}}$ ...\n",
    "\n",
    "We write our system as a finite difference equation in a semi-implicit scheme:\n",
    "$$\n",
    "\\frac{h_{i,t+1}-h_{i,t}}{\\Delta t} = \\frac{1}{\\Delta x}\\left(D_{i+\\frac{1}{2}}\\frac{h_{i+1,t+1}-h_{i,t+1}}{\\Delta x}-D_{i-\\frac{1}{2}}\\frac{h_{i,t+1}-h_{i-1,t+1}}{\\Delta x}\\right) + \\dot{b}_{i}\n",
    "$$\n",
    "Notice that this is semi-implicit because $D$ depends on $h_t$.\n",
    "\n",
    "Setting all the $t+1$ terms on the left-hand-side:\n",
    "$$\n",
    "h_{i,t+1} - \\frac{\\Delta t}{(\\Delta x)^2}\\left[D_{i+\\frac{1}{2}}\\left(h_{i+1,t+1}-h_{i,t+1}\\right)-D_{i-\\frac{1}{2}}\\left(h_{i,t+1}-h_{i-1,t+1}\\right)\\right] = h_{i,t}+\\dot{b}_{i}\\Delta t\n",
    "$$\n",
    "\n",
    "We re-write this equation as:\n",
    "$$\n",
    "-\\alpha_{i} h_{i-1,t+1} + \\beta_{i} h_{i,t+1} - \\gamma_{i} h_{i+1,t+1} = \\theta_{i}\n",
    "$$\n",
    "where:\n",
    "$$\n",
    "\\alpha_i = \\frac{\\Delta t}{(\\Delta x)^2} D_{i-\\frac{1}{2}}\n",
    "$$\n",
    "$$\n",
    "\\beta_i =  1+\\frac{\\Delta t}{(\\Delta x)^2} \\left(D_{i-\\frac{1}{2}}+D_{i+\\frac{1}{2}}\\right) \n",
    "$$\n",
    "$$\n",
    "\\gamma_i = \\frac{\\Delta t}{(\\Delta x)^2} D_{i+\\frac{1}{2}} \n",
    "$$\n",
    "$$\n",
    "\\theta_i = h_{i,t}+\\dot{b}_{i}\\Delta t\n",
    "$$\n",
    "\n",
    "And the full system can be written as a tridiagonal matrix system:\n",
    "$$\n",
    "\\begin{bmatrix}\n",
    "   \\beta_1 & -\\gamma_1 & 0 & 0 & 0 & ... & ... & ... & ...& 0 & 0 & 0 \\\\\n",
    "   -\\alpha_2 & \\beta_2 & -\\gamma_2 & 0 & 0 & ... & ... & ... & ...& 0 & 0 & 0 \\\\\n",
    "   0 & -\\alpha_3 & \\beta_3 & -\\gamma_3 & 0 & ... & ... & ... & ...& 0 & 0 & 0 \\\\\n",
    "   . & . & . & . & ... & ... & ... & ... & ...& . & . & . \\\\\n",
    "   0 & 0 & 0 & 0 & 0 & ... & ... & ... & ...& -\\alpha_{n-1} & \\beta_{n-1} & -\\gamma_{n-1} \\\\\n",
    "   0 & 0 & 0 & 0 & 0 & ... & ... & ... & ...& 0 & -\\alpha_{n} & \\beta_{n} \\\\\n",
    "\\end{bmatrix}\n",
    "\\begin{bmatrix}\n",
    "    h_1 \\\\ h_2 \\\\ h_3 \\\\ . \\\\ h_{n-1} \\\\ h_{n}\n",
    "\\end{bmatrix}\n",
    "= \\begin{bmatrix}\n",
    "    \\theta_1 \\\\ \\theta_2 \\\\ \\theta_3 \\\\ . \\\\ \\theta_{n-1} \\\\ \\theta_{n} \n",
    "\\end{bmatrix}\n",
    "$$\n",
    "We need to solve this sytem at each time step $t$ for the vector of unknown $h$.\n",
    "\n",
    "First, we invert the coefficient matrix and solve for the $h_{t}$ vector. Then, we update the $\\alpha, \\beta, \\gamma$ coefficients, because they depend on $D$, and thus on $h_{t}$. With the updated $h_{t}$ and $\\alpha, \\beta, \\gamma$ coefficients, we solve for $h_{t+1}$ at the next time step. We proceed iteratively in this way.\n",
    "\n",
    "What boundary conditions should we use? We can assume that the most upstream grid cell corresponds to an ice divide. An ice divide is the summit of an ice sheet, where the surface slope is 0. Thus, we enforce a Neumann boundary condition:\n",
    "$$\n",
    "h_1-h_2=0\n",
    "$$\n",
    "In the matrix equation, this translates to:\n",
    "$$\n",
    "\\beta_1 = 1 \n",
    "$$\n",
    "$$\n",
    "\\gamma_1 = 1 \n",
    "$$\n",
    "$$\n",
    "\\theta_1 = 0\n",
    "$$\n",
    "At the end of the domain, we assume a strong ice removal mechanism, such that ice thickness is equal to 0. For example, the ice sheet can be assumed to reach the ocean at the end of the domain. Thus, we enforce a Dirichlet boundary condition:\n",
    "$$\n",
    "h_n=0\n",
    "$$\n",
    "In the matrix equation, this translates to:\n",
    "$$\n",
    "\\beta_n = 1 \n",
    "$$\n",
    "$$\n",
    "\\alpha_n = 0 \n",
    "$$\n",
    "$$\n",
    "\\theta_n = 0\n",
    "$$\n",
    "Now, we are ready to code our finite-difference ice sheet model!"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "32bc3bbb",
   "metadata": {},
   "source": [
    "# Part 2: Example Code"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "846876ef",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAsgAAAGoCAYAAABbtxOxAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAC5AElEQVR4nOydd3yV1fnAv+fe3Ow9CSEhhL1lqYAgTqiTusWB2IrV1lptq2hV1LYWRx2/1j0AxVVUcKKICMreICvMQAgJ2etm3dx7fn+cm0kWJGSQ5+vn+N685x3nzQXyvU+e8xyltUYQBEEQBEEQBIOlrQcgCIIgCIIgCO0JEWRBEARBEARBqIYIsiAIgiAIgiBUQwRZEARBEARBEKohgiwIgiAIgiAI1fBo6wG0F8LDw3V8fHxbD0MQBEEQBEFoJTZu3JiptY6ovV8E2U18fDwbNmxo62EIgiAIgiAIrYRS6lBd+yXFQhAEQRAEQRCqIYIsCIIgCIIgCNUQQRYEQRAEQRCEaoggC4IgCIIgCEI1RJAFQRAEQRAEoRpSxUIQBEEQTnPy8/NJT0/H4XC09VAEodWw2WxERkYSGBh4wueKIAuCIAjCaUx+fj7Hjh0jJiYGHx8flFJtPSRBOOVorSkuLiYlJQXghCVZUiwEQRAE4TQmPT2dmJgYfH19RY6FToNSCl9fX2JiYkhPTz/h80WQBUEQBOE0xuFw4OPj09bDEIQ2wcfH56RSi0SQBUEQBOE0RyLHQmflZP/siyALgiAIgiAIQjVEkAVBEARBEAShGiLIgiAIgiC0C5xOJ3/961+JiIggICCAq6++mszMzLYeluDmtttuw2az4e/vX9leeeWVyv45c+ZgsVhq9N94442V/UlJSSil8PPzq+zv1q1bWzxKo7S6ICulblBK/ayUyldKldfR31MptUApledua5RStmr9I5VS65RSRUqp/Uqpm2udH6mU+kwpVaCUylBKPa2Ukg8CgiAIgtDOmTVrFp9//jlr167lyJEjANxyyy0tdn2pA918pk6dSmFhYWW7++67a/QnJCTU6P/www+Pu0ZiYmJlf8X73BK05PvbFnWQc4BXAB/gjeodSqkI4Gf3/tuAQmAY4HT3BwGLgOeAccB4YIFSar/WerX7Mu8DBUA3IAz4FsgGnj6VD3XCrFwJP/988ud7ekJkJERFmW1kJEREgIeUthYEQRA6Jm+88QaPPfYYCQkJADzzzDP06tWLpKQk4uPjaxybk5ND165dWbVqFcOGDavcP378eC666CIeffRRJkyYwBlnnEFSUhJLly7l4YcfZsaMGa35SEIzePDBB9m9ezeff/555b6lS5cyefJkUlNTWb9+PRdeeCGzZ89m5syZZGRkUFBQ0CL3bnWb0lp/B6CUmlBH9/3AYa3149X2baj2+iqgGHhGa62B75VSC4DpwGqlVA/gQqCX1joPyFNKPQ08QnsT5KVL4bHHWv66YWFVwlxdnqt/3aULxMSAl1fL318QBEFo13z7J0jb0jr36nIGTHqxacfm5eVx+PBhRowYUbmvZ8+eBAYGsm3btuMEOSQkhGuvvZa33nqLl19+GYA9e/awevXqGlHLd955h4ULF7JgwQKKi4ub+UStwxrnBrJ1TqvcK1SFcLZ1ZJOP//TTT/nss88IDw/nyiuvZObMmfj7+1f2Jycn06VLF2w2G2PHjuVf//oXPXr0qHGNs846i7KyMgYOHMjjjz/OhAkT6rzX9OnT6d+/P6mpqURHRwPw1ltvMWXKFPz8/ACTlrNo0SI2b96MzWar8zonQ3sLN54H7FVKfY6JEB8BntZav+/uHwpscstxBZuAW6r152mt99fqj1dKBWqt86vfTCk1HSPXxMXFtfjDNMhDD8Ff/3ry5xcXQ0YGpKfDsWNmW9Eqvt661Wxz6vlLFhkJ3bqZFhtb9bp68/Y++TEKgiAIQhPJzzc/ooOCgmrsDw4OruyrzfTp07n88sv597//jbe3N2+//TaTJk0iJiam8phrrrmG888/HwBfX99TNPrOwT333MPTTz9NREQEu3btYtq0adxxxx2VH0jGjx/PL7/8Qq9evUhPT2fGjBlcdNFFbN26FT8/P8LDw1m9ejXDhw/H4XDwzjvv8Ktf/Yq1a9cyZMiQ4+7Xs2dPxo8fz9y5c5kxYwY5OTksWLCAFStW1Dhu1qxZx/25aS7tTZDDgVHA9cDVGGH+Uil1SGu9AggA8mqdkwtUrB9YXz/uY2r8DdNav4E7zWPkyJGa1sTDo3npEN7eEBICffo0fmxZWU2ZTk2FI0eq2sGDJt2jLpEOD68pzHFx0KOHafHxRrKlvqYgCEKHoakR3dYmICAAMJHk6uTm5ta7TPA555xDTEwMn3zyCTfccANz587ljTdqZG8eF3nuCJxIRLc1qR7dHzhwIC+88AITJkxgzpw5eHl5VabGAHTp0oU333yToKAg1qxZwwUXXIC/vz9nn302AJ6entxzzz188cUXzJ8/v05BBrjzzjsrU2PmzZtH//79a4zDYrEQGxvb4s/a3gS5AFittf7E/fX3SqlvgSuAFe7++FrnBFMlvgVA7Y8QwdX62g/b18H2NSd/vtUGQaGmBVbb+gUeL6yenialoton6jqx2yElxUhzcnJNiT5yBFavhqysmuf4+BhRri7N1bchISLQgiAIQqMEBwcTFxfHpk2bOOOMMwA4cOAA+fn59coTmCjy22+/jb+/P1arlUsvvbRGv8Ui8/RPFRXf25q/2K9CKYVSqt7+ims01D958mTuueceli9fzttvv82dd95Z5z1amvYmyFuAXnXsr/jObQV+XatvmHt/RX+QUipBa32gWn+SOye5/bBnC3wx++TPdx5XAMRg9YDAkJrSXLGteB0UDmFR5rXFWnWun5+JSDcUlS4shKQk0w4erNoePAirVkFubs3jAwOrZDkhAXr1gt69zTYuDqzW424hCIIgdE6mT5/O008/zXnnnUdYWBgPPvggEydObDAKfOutt/LQQw/xxBNPMG3aNKzyc+WU8dFHHzFp0iSCg4PZu3cvf/7zn7niiivwdqdjfv311wwdOpSYmBhycnJ46KGHCA8Pr4war1mzBn9/f/r160d5eTnvvvsuy5cv56mnnqr3njabjdtuu4377ruPvXv3MmXKlFZ51lYXZKWUFbABnu6vK5JcS4HXgZ+VUpOBL4BzgYupmmC3AHhGKfVX4CVMnvJVwEUAWuuDSqkl7mNux1SxeNB93fbFVdNNO1nKHVCYB3lZkJcN+dk1twXu7bFksy0pOv4aVg8IiYDQSAiNMtIcWq2FRZl+j2pJ7/7+MGiQaXWRm1tTnCu2+/bB4sUmd7oCm61KmitahTx37y4VOQRBEDoZFXmmo0aNorS0lIsuuoh58+Y1eE5wcDDXXHMN7733Hp999lkrjbRz8tprr3H33XdTWlpKZGQkv/71r3n88ccr+5ctW8Ydd9xBXl4egYGBjB07lu+//75yEt/Bgwd59NFHSU1NxdvbmwEDBvDll1/WSJmoizvuuINnnnmGqVOntniucX2ohsLap+SGSt0G1BU67aG1TlJKXQv8A1Om7SDwhNZ6frXzRwEvA4OBVOAxrfW8av2RwGsYaS4F3gFmaK1dDY1r5MiResOGDQ0d0rEpLYb8HCPLuRmQnQ7ZxyDrWM1taa0ZvkqZyHOFPIdHQ0RXiIgx28iYutM66sLlMvnP+/bB3r1mW73Z7VXHeniYqHOFOPfpA/36Qd++Jhda0jYEQRCaxK5du+jfv39bD+OU8vjjj7Nq1SoWL17c1kMRTgF2u52oqCgWL17MmDFjTvj8hv4OKKU2aq2PS/puizJvc4A5DfTPB+Y30L8eOLOB/nRMVFmojpcPRPgYqa0PraGooJY0u0U6Ow3SkmHHOii21zzPx6+mNNcn0BZLVS70uecef++0tCpZri7QP/9sUjsqqEgFqRDmim2fPiAzlAVBEDoVx44d48033zxucp5weqC15sUXX6R///4nJccni/wOW6hCKSOzfoEQ17vuY7QGez5kHHW3lKpt+lHYub5+gY7sBl3iICoWusSa12FdTKqHUhAdbdq4ccffMy0Ndu+GxMSq7erV8NFHpr+CuLjjxXnAAFP7WaLOgiAIpxX3338/r7/+Orfccstxk/OEjk96ejoJCQlERkYyf369sdNTQqunWLRXTvsUi9aiPoFOT4FjR0xOtKO06nirh4kyR8W6xTnOLdHdjFDbPBu+X3GxiTZXF+eKbfWoc3AwDBxoZLn6NjpaxFkQhNOazpBiIQgN0SFSLITTHKXAP8i0HnX8YXS5ICfdiHJaMqQdrnqduLlm9Fkpk/McFQvR3SE6Hrp2h+geENnVVODw8YEhQ0yrjtZw9KiR5V27YMcO2LkTPvsM3nyz6rjgYCPLtcW5a1cRZ0EQBEHopIggC62LxWLSKsK6wIBRNfu0hoKcmtKcdthsV31rItMVeNhMpDk6HrpWbz2MnCtVle98wQU175GRUSXMFduFC+Gtt6qOCwoysjxokJHvwYNNCw09Zd8aQRAEQRDaByLIQvuhomJGYCj0OaNmX4U8H02C1ENw9KB5nXIANi2vWRc6ILimOMckQLeeJpXDYjWr/0VGwnnn1bxHbXHesQM+/bRmxDkmpkqWK8S5Xz/w8joF3xBBEARBENoCEWShY1BdnvsNr9nnLDc5zkeT3ALtbltWwLKFVcfZPE2EuVtClTTH9DQTBj1sEBEBEyaYVoHWpjTdL7/Atm1m+8svsHSpWcIbTEm6vn2PF+e4OEnTEARBEIQOiAiy0PGxerhzlLvDiFrl44oKIOUgHNlvos1H9sOebbByUc3zu8QZYe6WYKS5W4KJQnt6mXzkrl1h4sSqcxwOMzmwujivWWOqalQQHAxnnFGz9e9vlv4WBEEQBKHdIoLcVhRmgD3j5M+3eIBXIHgHgM1XIpX14RsAvYeYVp2SIneKxn444hbnw3tg3Q9QsaaMshhxjusNsb2qtl3izCqAFZP7rr++6rr5+bB9u5HmrVthyxZ44w0ocq9kaLOZSYDVpXnoUCPTgiAIgiC0C0SQ24oDP8HW/7XMtZTViLJXteYdaAS68nVA1THeQWC1NX7d0xlvX0gYYFp1ykoh7ZBbmvfB4X1ucV5SVW/Z5mUizLG9a8pzSCQEBsKYMaZV4HSaaPOWLVXtm29gzpyqY+Lja0rziBEm31k++AiC0In46KOPePnll9m6dStFRUWUl5fX6N+6dSszZsxgy5YtpKWl8fPPP3POOedU9hcXF3PrrbeyZcsW9u/fz5NPPskjjzzS2o9xWrNkyRIeeeQRtm/fjre3N9dddx2vvPJKZf/+/fv5y1/+wtKlSwHo378/P//8Mzab8Q673c5DDz3E/PnzKSgoIDY2lg8++IBhw4a1yfPUhwhyW9FjHET0PfnzXQ4oLYCSAijNr9qWFkBustmWFgL11Ln29AffEPBxtzpfB3c+kfb0grg+plWntNhIc/JeOLwXkvfBL6vhpy+qjvELrBZp7g3d+0L3PkbG+/Uz7YYbqo5PS6spzVu2wOefV4l4RIQR5Yo2fLjkNQuCcFoTEhLC3XffTXFxMdOnTz+u39PTk6uuuoonn3ySM888flFdpRRjxozh7rvv5qGHHmrx8TkcjkrR64wsW7aMa665hrfeeovLL78crTU7d+6s7M/IyGDcuHFMnz6dOXPm4O/vz+bNm7FarYBZFW/y5Mn4+Piwbt06YmNjOXDgAH5+fi0yvpZ8f0SQ24qASNNOJS4XOOxQkl9Tpotz3S0HinIgL8V8rZ3HX8MroEqWfUPAN8w0v3B3CwObz6l9jvaAlw/0HGhadQpya0rz4b3w89dQ7F6kRClTx7lCluP7mddhXczqfpMmmVaB3W7SMzZuNG3TJvj+exOFBggPN6JcXZrj40WaBUE4LZjonuuxbNmyOvv79+/f4KIn3t7e3HfffZWvG+PBBx9k9+7dfP7555X7li5dyuTJk0lNTWX9+vVceOGFzJ49m5kzZ5KRkUFBQcEJPNHpxUMPPcTvfvc7rrnmmsp9w4dXTZx//vnniYuL4/HHH6/cN3Jk1RocixcvZuXKlRw5coRQd9nUhISEeu/Xlu+PCPLpjMVSlXLRGNplBLo4p0qca7/OO2K2tVdf9PQzouwbXlOcK177hJjyaqcjAcGmnnP1ms5aQ2aqSc1I2g2H9sChRFj7fdUxfoFuaXaLc/e+Jvrs5wejR5tWQXGxkeZNm6rE+dlnoeJXj6GhNaV51Cjo3l2kWRCEOnltjpMDh1pnFd2E7orf3dZ+//2fPn06/fv3JzU1lejoaADeeustpkyZUhnVdDqdLFq0iM2bN7da9Dg181tKytJa5V7enl2IDp/U6HF2u51169YxceJEhg8fzuHDhxk0aBDPPfdcpQT/+OOP9O7dmyuvvJKff/6Zbt268eCDD3LTTTdV9vfs2ZOnn36a2bNnExAQwPXXX88TTzxR5/e2Ld8fEWTBoCzgE2Qa8fUf53IaSbZngj3LbIvcW3smZCRCmb3mOUqBT6iR5YBI8K/VfEKNzJ8uKAURXU0bMaFqf7HdRJgPJZqWlAhLPzXpG2CqaXSNN7Lcoz/0GADxfY1Mn3WWaRWUlJjJgBXCvHEjPP+8qa4BJtJ85plGlita5Cn+jYUgCEIHo2fPnowfP565c+cyY8YMcnJyWLBgAStWrKhx3KxZswgKCmqjUbYPcnJycLlcvPnmmyxatIh+/frx3HPPcckll7Bnzx6Cg4PJzMxk/fr1fPzxx3z66af8+OOPXH755XTv3p1zzjmHzMxMtm/fzsSJEzl8+DDJyclccskl+Pn58be//e24e7bl+yOCLJwYFmtVZLg+HCVQlFlNorNMxQ57JqTthKIV1MiNrrxmhTRH1BRo78DTIxrq4wd9zzCtApfTrBRYIc2HEmHneljxddUxUbFuYXZLc49+ph70yJGmVVBWZsrNrV8P69aZ7bffmlQbMFHlUaOqxHnECAhowm8XBEE4rWjPEd224M477+Thhx9mxowZzJs3j/79+zNixIjKfovFQmxsbKuOqSkR3dYmwP3zYtq0aQwZYipDPfTQQzz77LOsWrWKSy65hICAAEaPHl2ZgnHRRRcxadIkvvjiC8455xwCAgKwWq3885//xMvLi969e/P73/+eDz74oE5BhrZ7f0SQhZbH5g1B3UyrC2e5keXC9FotA5LXmzzp6nh4QUAUBHSp1txf+4aa6HdHxWKtWvFvdLU6y3lZcHA3HNwFB3eatmZxVX9YlBHm+AGQ4JbnkMiqNIvf/c4cV1hoUjOqS/Mnn5g+pUxd5gppPvNMs8iJ1GkWBKETMXnyZO655x6WL1/O22+/zZ133lmjXymFOh2CNM0kKCiI+Pj4Or8XFfvOOOMM9u3b12B/XTT0/W2r90cEWWh9rB4Q2MW0unCUHC/PBccg9wgc2QSuamV/rDbwd8tyYG15Du+4qRtBYXDGWNMqKMwzaRmV0rwLNi6vygkPCjWinDAQeg4y29BIGD/etAoyM40oV0jzokUwd67p8/Y2gn322VWtWz0fdARBEFoYp9OJw+GgzL1SaUlJCQBeXl4opdBaU1paWnl8WVkZJSUl2Gy2ykoJpaWlaK1xuVyUl5dTUlKC1WqtNz/VZrNx2223cd9997F3716mTJlyip+y43L33Xfz0ksvceONN9KnTx+ef/55vL29GeMubXrnnXcybtw4Fi5cyBVXXMHy5ctZvHgxDz74IABXXXUVM2bMYObMmTz55JOkpKTw6quvcvvtt9d7z7Z6f5SuPeGqkzJy5Ei9YcOGth6G0Bgul8l5Lkg10lyQBvlpZlt4DJyOqmMtHiZFIzAaArtW23Y9fdI2Sorc0ryzSpyT91ctdhISCT0HQMIgI809B5qJhdXRGpKTYe1asxrgmjUmp7nih1BMTE1hHjECfDpB5RJBOE3YtWtXg5Uf2hNz5sxh2rRpx+0/ePAg8fHxJCUl0aNHj+P6Z8+ezW233QZAfHw8hw4dqtE/depU5lSvPV+L/fv307t3b6ZOncrs2bMr9y9btowLL7zwuHrMnRWtNTNnzuSNN96gpKSEYcOG8cILL9SIDM+fP59HHnmEI0eO0KNHD2bOnMm1115b2b9t2zZ+//vfs2nTJsLDw5k2bRqPPvpo5Qecumju+9PQ3wGl1Eat9cjj9osgG0SQTwO0y1TcKHALc8ExI9J5R83X1SPPnn41hbnydRewdvAUg9Jik55xYAfs3w77d0BqUlV/ZExVlLnnIBN19vWveY2yMrMSYIUwr1kDBw6YPg8Ps/rf6NFV0pyQcHp84BCE05COJMhthd1uJyoqisWLF1dGQ4X2Q3PfHxHkZiCCfJrjcpmJgvlH3S3V3Y5CUXa1A5WZMBjkFuegGHfrZqLOHZWiAjiw08jy/u3mdUaK6VMKouONLPceDL2GmNJzHrV+HZmeXjPKvG6dyXEGs6jJmDEwdqzZjhhh0jUEQWhzRJAbRmvNU089xcKFC1m/fn1bD0eoRUu8PycjyJKDLHQOLBZ3bnIUxNRaztJRUk2aj1a1Y7vBWZXrhldAlSwHdTOvg7uZOs/tPXrqGwCDzjKtgvxstzRvN9Hmbavg5y9Nn83TRJZ7DzHC3GuwKVt3+eWmgVm8ZOdOWL3atJUrzUqAYCb6jRxZJcxjxxqJFgRBaEekp6eTkJBAZGQk8+fPb+vhCLVoy/dHIshuJIIsHId2mehy7hGz2mDeEXdLqVnr2eZTTZyrbf0jOlaFjYoFTvZtg72/wL5fjEA73B8SgkKrZLn3EBNxrp2akZ4Oq1aZtnIlbNhg0jUAeveuKcz9+nXcSZSC0IGQCLLQ2ZEUi2Yggiw0Ga2hJK9KlvNSjETnu5fsrsDqZSLMwXEQHFvVfILbf8S5gnKHWdxk7zYjzPu2wdEk06cUxCRUCXOfM6Bbz5rSW1JiJvytXFklzZmZpi8kpEqWx483EWcvr9Z+QkE47RFBFjo7IsjNQARZaBFKC6uizbnJVa0kr+oYr4CawlzRPP3abtwnQmGeyWXeu9Utzb9AQa7p8w1wy/JQI8y9BteMMmsNe/fWFOZdu0yft7dZLXDcOCPMo0eDv3/tuwuCcIKIIAudHRHkZiCCLJxSivMgL7mmNOcmg6O46hjfsCpZDomDkO4QGGPqRrdntIZjyZC4Bfa4W/I+s19ZIK53lTD3OQOiutWMoGdmwooV8NNP8PPPZmETlwusVhg+vEqYzzkHwsLa5BEFoSMjgix0dkSQm4EIstDqaG1qOucmQ+7hauKcAi53PWeL1eQzh3Q3Ldi99WnZNedbnKICE1lO3GqEee82KHZXvAgKrZLlPkNNyTnPaqkVBQVm0l+FMK9dW1WTeeDAKmEeN04WMRGEJiCCLHR2RJCbgQiy0G5wOU1FjZxDNVtxTtUxPsEmt7lCnEO6m9J0lnYabXY54ch+d5TZLc1ph02fh81M+Os3DPqNMNLsX+0DQEmJmexXIcwrVxqJBlN/ecIE0847T4RZEOpABFno7IggNwMRZKHdU5IPOYeNLOe6pTn3SNUCKBYPU0EjNN7dehhxtrXTVe/ysowsJ26B3RtNxQxnuUm/iO0FfYdD/+HQbziEVVuWvLwctm0zwrx8uWk57g8PPXtWyfKECWYVQEHo5IggC50dEeRmIIIsdEhc5TWjzdlJppXmuw9QZnXAkHgjzKE9jDy3x0VPSoth33bYvcm0PVvMUtpgajD3Gw59h0H/EdC1R1W1DJfLCPOyZfDjj0acc3NNX69eVbI8YQJ07drqjyUIbY0IstDZEUFuBiLIwmmD1iYdI/ugW5jdW3tG1TG+YTUjzaE9zL72VH7OWQ6H9lQJ8+5NJuoMEBAMfc8wKRkDRkF836rJjE7n8cKc564i0qdPzZSMLl2Ou60gnG50JEF+8MEH+eqrr0hOTsbf359LL72Up59+mtDQUADeffddXnvtNXbt2oXVamXUqFE888wzDB48GIA1a9bw97//nQ0bNlBSUkKvXr149NFHmTx5chs+1enDwIEDOXToUOXXTqeTkpISNm7cyPDhw/nd737HvHnzapxjt9v597//zf3338/777/PnXfeWaO/uLiYSy+9lC+++OKUjVsEuRmIIAunPaWF1YT5IOQkQd5RwP1vgFeAEeWwnhCWYFp7kuaKahm7NkLiZiPMqe5/qH38THR5wCgYMNKsAlixVLbTCVu3GlletswIc747wj5wIFxwAVx4IZx7LgS2w8i6IDSTjiTIDz/8MNdeey2DBg0iNzeXW2+9FU9PTz53r9L58ssv06dPH8aMGYOHhwdPPvkk77zzDgcOHMDHx4dvvvmGzMxMfvWrXxEWFsYXX3zBjTfeyE8//cSoUaOaPT6Hw4HNZmv2dU4X/va3v7Fw4UJ27NhRZ//333/PJZdcwqFDh+hax2/w8vLy6Nq1K7Nnz+a6665r9njqe39EkJuBCLLQKXGUmAoa2Qch6yBkHzCVNLTL9HsH1SHNoW075urkZMDODbBrA+xYD0cPmv1ePiYlo/9II8w9B9YU5s2bYelSWLLElJgrLjZl5UaNqhLm0aNl4RLhtKAjCXJtvv76a6ZMmUJexW+BamG32/H392fTpk0MGzaszmNGjRrFlClTuO+++47re/DBB9m9e3elgAMsXbqUyZMnk5qayvr167nwwguZPXs2M2fOJCMjg4KKScKdnPLycmJjY3nooYf44x//WOcx11xzDS6Xi88++6zO/v/+97/84x//IDk5uU6xban352QEuZ1OeRcEoVWweUNEH9MqKC9z5zMfgKz9kHUAUreaCC6AT4hbmhOqxNknuE2GT0gEjP2VaQC5mSbCvHM97NwIH71k9nv5mOoYFcI8dIhZue+BB0wJudWrjSz/8APMmgX//Cf4+JjayxXCfMYZRqIFoYPz0Ppifslxtsq9BodY+deok58o/MMPPzBkyJAG+319fenVq1ed/WlpaezYsaPea0yfPp3+/fuTmppKdHQ0AG+99RZTpkzBz88s3uR0Olm0aBGbN29uvejx+jnm3+HWIKQ7jLrthE9buHAheXl53HrrrXX2p6Wl8cUXX/DVV1/Ve43XX3+d22+/vd7va1u+P60uyEqpG4DfA0MBX611nWNQSj0NPADcorWeV23/SOAVYBCQCsys1R8JvAZcBJQA7wAPaV0REhMEoUE8PCGit2kVOErMP9ZZByB7v4k2p2ymMj3DNwzCe0J4LwjrZcTZ5t36Yw8Oh9ETTQPIz3YL8wYjzf/7r9lv8zI5zIPOgoFnwrhzTF7yP/5h8pV/+qlKmGfMMOeEhJi85QsugIsuMhMA20v6iSCchnz66ae8+eabLF++vM7+PXv28Nvf/pZ///vfBAQEHNdvt9u5+uqrueKKK7jgggvqvEbPnj0ZP348c+fOZcaMGeTk5LBgwQJWrFhR47hZs2YRFNTO68+3Mq+//jrXX389wcHBdfa//fbbxMbGctFFF9XZv3LlSnbu3Nlg7nFbvj9tEUHOwQiuD/BGXQcopc4EfoUR4Or7g4BFwHPAOGA8sEAptV9rvdp92PtAAdANCAO+BbKBp1v8SQShs2Dzhsi+plXgKHGnZlREmvfB4XWmTylTci6sl5Hm8F5mhUBLK0dgA0PhrItMA7Mk9u6NsGMD7FgHH/2f2e/jZ6LLg86EgWfBpZfC5ZebvtRUk47xww+mVfyqsEcPuPhi084/H+r5ISEI7Y3mRHRbi/nz53PnnXfyxRdfMHz48OP6d+7cyUUXXcRf/vIXfve73x3XX1BQwKWXXkpkZCTvvvtug/e68847efjhh5kxYwbz5s2jf//+jBgxorLfYrEQGxvb/Ic6EU4iotua7N+/nx9++IHVq1fX2e9yuXjzzTe56667UPUEEl5//XUuvvhievTo0eC92ur9aXVB1lp/B6CUmlBXv1LKC3gbmA58WKv7KqAYeEab5OnvlVIL3MeuVkr1AC4Eemmt84A8dyT6EUSQBaFlsXlDVH/TKijJN7Kcuc+0Ixtg/4+mz+pZlZYR3ttEnP0iWjcKGxAMoy4wDUxVjJ3rYfs6I8yb3JGqgBAYOKoqwjxlCtx0k0kz2bcPvv8evvsO3n8fXn/dpF6cdZaR5YkTTfqGh2SwCcLJMHv2bP785z/z5ZdfMnbs2OP6N23axKRJk3j00Ue55557juvPzs5m0qRJJCQkMG/ePDwa+bs4efJk7rnnHpYvX87bb799XJUFpVS9ktdZef311xk6dChnnXVWnf3ffvstqamp3H777XX2Z2dnM3/+fD766KNG79VW7097/Bf8cWCp1np1HQ88FNika84s3ATcUq0/T2u9v1Z/vFIqUGudX20/SqnpGLkmLi6u5Z5AEDor3oEQM8w0MEJZmF4lzFn7IHEx7PrafXyQEeaIPialI6xX66ZmBIXB6EmmAWSmwY61Rpi3r4E1i83+sCgTWa6IMN99t2kOB6xZA4sXG2F+4gl4/HETTb7gAiPLF18M3bu33jMJQgfm//7v/3jiiSf47rvv6qw6sXLlSi677DKefvpppk+fflx/WloaF110EcOHD+edd97B2oR5Azabjdtuu4377ruPvXv3MmXKlBZ5ltOVsrIy5syZw9///vd6j3n99de56qqriIiIqLN/7ty5hIeHc9lllzV6v7Z6f9qVILvzi68FzqjnkACg9lTWXCCwkX7cx9QQZK31G7jTPEaOHCnlPAShpVEKAqJM6+GOBLnKzYqAFdKcuRdSNrmPt5gltCN6G2kO72POba3oTXgXOPdK07Q2y2FvX2uiy5t/hp/cuXLR8TDkbBg8GkacCePGwd//DllZJg3ju++MNH/6qTm+T5+q6PKECeDv3zrPIwgdjHvvvRcPDw/OO++8GvsLCwsBeOSRR8jLy+P+++/n/vvvr+xftGgR48aN4/XXX2f79u0cOHCATyv+/mHKxz388MP13veOO+7gmWeeYerUqZJr3AifffYZxcXF3HTTTXX2p6Sk8PXXX/PDDz/Ue4033niD3/72t036AANt8/60WZk3d4rFkopJekopT2Aj8JdqaRhJwCMVk/CUUi8C8VrrydWucx9mIt9wpdRkYI7WOrha/zBMFDnYnXZRJ1LmTRDakNJCI8oZeyFjj4k0O4pNn1egEeZwtzS31QRAlwuS9xph/mWNmfhXWmzyqnsPMbI8dDT0HGQWLdEadu+uii4vW2bKyXl6wvjxcMkl8KtfQd++MtlPOKV05DJvrYXdbicqKorFixczZsyYth6OUIvmvj8dqg5yHYIcDxzATKirIAQoAr7QWt+klJoGPK617l7tOu8CTq31NHcO8gGgp9b6gLv/duBRrXWDWeAiyILQjnC5IO+IkeXMPUac84+aPmWBkDgTXY7sC5H9wC+89cfoKIM9W+GX1bB1FRzcaaTYN8DkLw8eDUPHQFSsEeDSUlNzedEi03buNNfp0cOI8iWXmCoZvr6t/yzCaY0IcsNorXnqqadYuHAh69evb+vhCLVoifenQwiyUsoK2DAVKL4BKn7X6ACiax2+GngG+EBrnaWUCgb2uve9hKlk8TlwUUUVC6XU95g0i9upqmIxW2s9q6FxiSALQjuntMCkZGTscYvzPigvMX2+YVWyHNHXpGlYLK07voJcd3R5NWxbDRluoY+IMekYQ8aYSX/+7l8PHjpkRPmbb0xaRlGRWZhkwoQqYe7du767CUKTEUGun/T0dBISEoiMjGT+/Pk1qiMIbU9LvT8dRZBvA2bX0dVDa51U69gkqqVYuPeNAl4GBmPKwD3WQB3kUkwd5BmN1UEWQRaEDobLZVYBTN9tWkYiFLl/AWXzcU/8c0tzeC/waMVV8Sryl7etNsK8fR0UF5pIcs/BcMY5cMZYs8KfxWqiyz/9VCXMiYnmOj17VqViTJhgFi8RhBNEBFno7HQIQW6viCALQgdHa7BnVslyeqJZNhsNympW/6seZfZpxYk4znLYtx22rYItK2H/L2a8AcEwZDQMPcekYwS7U0UOHKhKxVi61OQue3ubyhiXXw6XXQYxMa03fqFDI4IsdHZEkJuBCLIgnIZUTP5LTzTinLUPnA7TFxgNkQOqajm3Zh5zQa6JLm9ZAVtXQJ478t2jvzu6fI6Z+Gf1gJISWL7cRJa//BIOHjTHDh9eJcvDh7d+SonQYRBBFjo7IsjNQARZEDoBTodZ/S99NxzbZbaOItPnFwFR1YTZv5XKy7lccCjRyPKWFWbin8tpJvsNPtvI8tAxENbFRJ137oSvvjKyvHq1OT862ojy5ZebKLNM9BOqIYIsdHZEkJuBCLIgdEIq8piP7YL0nWZbWmD6fEMhsn+VMAfGtI4wFxWYMnJbVsKWnyE73eyP7QXDx8Pwc6uiy5mZVZHl776DggJJxRCOQwRZ6OyIIDcDEWRBENAa8lLg2E5I32WEuTjH9HkFGlGO7A9dBkBwrCk5d6rHk7wPtq40C5Xs3mTymf2DYOhYGHGu2foHQVmZmej35ZfHp2JcdhlMngxnnCE1lzshIshCZ0cEuRm0tiDn5G8ip2BTC11NVW4rf/QpVU+fqra7nj5UZf9x+6vdU6m6+tznNHb/WmNQta+j6tlf7bWq4x7VRlft+sp9HQWqntdUrOVe/bUFi7JhsXhiUZ5mW/218kQp6ylZA15oJ2gNBcfcwrwTju0Ge4bp8wqALgMhaiB0GWRymk/1n4WiAlNzefNPsHkF5GebKhh9hprI8ohzISbBHFtXKkb37kaUJ0+Gc84Bj3a1mKpwihBBFjo7IsjNoLUFObfwF/IKtjXzKhrt3lJ9qxvoA9D19bn363r217p+vWM4Zfevef2671/fGE7Vn3N1nDRXiLTV4ouH1QerxRer1RcPq2+t1z6oUx2BFFqewnQjzGk7TCvKMvt9QmoKc0DkqR2Hy2kqY2z+CTb9BEm7zf7IGBjmTsUYMBI8vSA93cjywoVmZb/SUggLM2kYkyebZbClhNxpS0cS5L/97W988MEHZGVl4e3tzfjx43n++eeJi4sDYP78+TzxxBOkpKQAMHDgQP75z39y7rnnApCYmMjUqVPZu3cvDoeDbt268ac//Ynp06e32TOdTqSlpXHvvfeydOlSysvLGTZsGC+88AJDhw4FwOl08vDDD/PBBx+Qm5tLfHw8M2fO5Jprrqm8xrvvvsvTTz/N4cOHCQ0NZdq0acycOfOUBptEkJuBpFic/lT9WXcZddZGnHWFQOvar124tAOXq8y0ite6rOa2zn2lOF3FOJ1FuLSj3jFZLN54uKXZavWpfO1h9cfDIwCbNaDytYlYS7S6XVERYU7bDsfcwlziXtHeL8IIc4U0+4Wd2rFkpZk0jE0/mRzmshLw8jET/Sqiy8HhUFho8pUXLjTSnJtrJvVNnGhk+bLLIDT01I5VaFU6kiDv3r2b6OhogoKCKCoq4pFHHmHNmjWsWrUKgOTkZDw8PIiOjsblcvHJJ58wbdo0UlJSCA4OJi8vr3JxCavVyrZt27jwwguZN28eF198cbPH53A4sNlszb5OR+Wqq66ioKCAjz/+GH9/fx555BE+/PBDDh8+jFKK//u//2PWrFn8+OOP9OnTh88//5zrrruObdu20a9fP7Zu3cqIESP47LPPuPzyy0lMTOS8887jySef5I477mj2+Op7f05GkOX3a0KnoUourTUzTE4xLpcDp6uYcmcRTlcRTmdRrdfFOF1FlJcXUOo8RrnLjtbO465jUTY8rAF4ePibrdUfm4dboK0BlUJtsXiJSLcWSkFgF9P6XFiVw1whzMkbYP8yc2xAtJHl6EEmwuwV0LJjCesCF15rWlkJ7FgPm5YbYd7woxlr7yEw8jw4+zy4+mpwOEwJuYULTVuwAKxWOPfcqlSM2NiWHacgNEC/fv0qX2utsVgsJFYsnAPEVvvzqLXGarVSVFREcnIywcHBBAUFERRUVeNcKZMyl5iYWKcgP/jgg+zevZvPP/+8ct/SpUuZPHkyqamprF+/ngsvvJDZs2czc+ZMMjIyKCgoaOnH7jDs27ePP/zhD4S6P0T/5je/4dlnnyUrK4vw8HD27dvHhAkT6Nu3LwCTJ08mLCyM7du3069fP/bv309kZCRXXHEFYN7vCRMmsHXr1jrv15bvjwiyIJxiLBYbFosNm0dgk47XWuNyleJwFlDuLKC8vJByZwGO8gLKnYWUOwspKUulvLygzui0RXli8wiq1gJrfO3hEYhFWVv6MQUwEhrczbR+k0C7IOewOx1jOySthL1LAGUWLokebFpkX7B6ttw4PL1h2DjTbtdweI+R5A0/wgcvmhYdb2R51Pnwfy/Bf/4DGzZUyfIf/2jaiBFwzTWm9erVcmMU2ox/ZWaSWFbWKvfq6+nJQ+EnVmP8gw8+4K677iI/Px8PDw+ef/75Gv2HDx9myJAhFBQU4HK5uP766xk8eHCNY4YMGUJiYiJlZWUMGjSIG2+8sc57TZ8+nf79+5Oamkp0dDQAb731FlOmTMHPzw8waQOLFi1i8+bNrRc9nvM0HNrdOvfq3g9ue7BJh/71r39l3rx5XHXVVfj7+/PGG29wzjnnEO5+j++44w5uuukmdu7cSd++fVmwYAHl5eWMHz8egIkTJ9K1a1cWLFjAlVdeyc6dO/npp5+YPbuuBZbb9v0RQRaEdoZSCqvVG6vVG4ho8Finq7RSoCsk2lGeh6M8H0d5HsWlR3G6io47z8MaUKdAe9pCsHmEYLW0oKx1ZpQFQuNNG3CpyRvO2g+pv0DqNtj5Fez4HKw2Ux0jejBED4GQuJarkKEUdO9r2tW/M6kYG5bBhqXwzXvw5WwICoXhE2DUeTDzUfjnP2HPHiPKn30GDz1k2tChVbJcLdInCC3JlClTmDJlCmlpabz99tvHyW9cXBy5ubnY7Xbmz59PaWnpcdfYtm0bDoeDn376iZ9++qlSpmrTs2dPxo8fz9y5c5kxYwY5OTksWLCAFStW1Dhu1qxZNSLTnZWxY8cyd+5cIiIisFqtxMbGsmjRosr+hIQExo0bx6BBg7BYLHh5efHee+8RGWnmZPj5+TFt2jSmTp1KUVERTqeTBx98sN70l7Z8fyQH2Y3kIAunKy6Xo1KYazRn1T6ty2ucY7X64ekR4hbmYDxtoZVfe1gDJIWjpXCUmAl/qduMNOcdMfu9At2pGIOh65BTt8pfUYGphrHhR5O/XFxo8paHjDaR5WHjITAEDh82ovzJJ7BypTl34MAqWR44UMrHtWM6Ug5ybSryiSsmdNXFgAEDeOGFF5g4cWKd/XfffTdBQUH861//qrN//vz5PPzww+zdu5f//Oc/zJ49m02bTJWpZcuWccEFF1BeXt7p/91zuVz06tWLSZMmMWvWLLy9vXn33Xd56KGH2L59O1FRUZUTJD/88ENiY2NZs2YNkydPrswBnz17Ng888ABfffUVo0aN4vDhw9xwww1cfPHFPPnkk3XetyXeH8lBFgThOCwWG16eYXh51j1JTGuN01WEw5FLWXkOZY4cyspzcDhyKSpJxlG+nepVQBRWbLbgagJttl62MGy2ECxK/llpMjZv6DbcNICibJOKkboNUrdDkpmYREC0iS53HWrymG0tVHHCNwDG/sq0cofJW96w1Ajz+qUmit1vGJx5Idx8A/zpT5CSYnKVP/kEnnwSnnjCRJMrZHnIEJFlocUoLy/Hbrdz9OjRegW5vLycvXv31ivIFf31MXnyZO655x6WL1/O22+/zZ133lmjvyKPubOTnZ3NwYMHueeeewgMNCmDv/3tb3nwwQdZs2YNV155JRs3buT3v/893bt3B2DMmDGMGzeORYsWcfHFF7Nx40bOP/98zjrrLADi4+O56aabeOedd+oV5LZ6f+QnmSB0cpRSeFj98LD64cPxq65p7cRRnlcpzmWOHBzubVFJMi5d/debCptHMF62MBN1toW5X4dh8wiSHzKN4RsKCeNN09pElFN/Me3Actiz2NQ9jugLXc8wwhzSvWWE1MNmlrQeOgZu/xsc3GlEed0PMPdp03oNNrJ8zZXwhz9AWlqVLD/1FPzjHyZPuUKWhw8XWRaajMvl4pVXXuG6664jMjKSI0eOcM899xAfH185ee/dd99lzJgxJCQkYLfbeeGFFzh8+DDnn38+AN999x3BwcEMGzYMpRTffPMN8+bN4z//+U+997XZbNx2223cd9997N27lylTprTK83Y0wsPD6dOnD6+88gqzZs3Cy8uLd999l4KCgso0mLFjx/L+++9zxRVXEBMTw9q1a1m2bBkvvvhiZf+9997Lxo0bGTFiBMnJycybN4/hw4fXe9+2en9EkAVBaBClrG7ZPT56Y6LPJZQ5sihzZFPmyKLUkUWZI4uikkM1JhEq5YGnh1uaPavLczgeVt/WfKSOgVJmtb7gWOh/iVlBL2M3HN1q2uYPTPMJhuihRpa7DmmZ6hhKQcJA0677AxxNgrVLYN338MELpnXvC2ddCJdfDHfdBRkZJmf5k0/g2Wdh1izo0QOuvx5uuEEiy0KT+Oabb3jyySex2+0EBwczYcIElixZgod7UZs9e/bw6KOPkpmZia+vL0OGDOHrr79mwIABAOTn53P//fdz6NAhPDw86NGjB//+97/5zW9+0+B977jjDp555hmmTp0qucYNsHDhQv7617/SvXt3HA4HvXr1Yv78+SQkmAWKnn32WR544AHOPPNM8vPziYqK4s9//jO33HILADfeeCNHjhzhhhtuIC0tDX9/fy655JLjJmLWpi3eH8lBdiM5yILQsmitKXcW1pDmMkc2pY4sHI5sNK7KY60WX7w8I/CyhZutZwRetgg8rP4Sda6PomyTipGyxWzL7ICC8J5uWT4DwnqBpYUXo0lPgXVLjDDv2WL2de1hZPnMC6FHf8jOhi++gI8/hiVLwOk0aRg33miE2V0CSmgdOnIOcmtht9uJiopi8eLFjBkzpq2HI9Siue+PLBTSDESQBaH10NqFozyXUkcWpWWZlDoyKrcuV0nlcRaLF162CLw8w91bI86SrlELl8tUxzi6xUSXM/cBGjz93LnLZ5jmG9Ky981Oh/U/GFneucGUtYuMMaJ81oXQawhkZZkJfh9+CD/9ZFJHhg0zUeXrrzfLXwunFBHkhtFa89RTT7Fw4ULWr1/f1sMRatES748IcjMQQRaEtsdEne1uYc6oIc5Op73yOKVs7mhzJN7u5uUZJRHnCkoLTN7y0S1wdBsU55j9oQnQbRjEDIewhJYrJQeQn2Nyltd+b1byc5ZDaCScfTGMnmQWKTl6FObPh48+grVrzXmjRxtZvvZacNc5FVoWEeT6qaiSERkZyfz58xkxYkRbD0moRku9PyLIzUAEWRDaN+XOIkodmUacK+U5nXJnYeUxVosPXp5RbmmOwsszEi/PyM5d11lryDkEKZsgZTNk7AU0eAdBzBlGlqOHgGcL5oHb880qfmsWw5aVpkJGRFc4eyKMnggJAyApyaRgfPQRbN1q8pMnTDCyfPXVEHaKl+buRIggC50dEeRmIIIsCB2TcmcRpWXplJQdo6TsGKVl6ZSWpdeYIOjpEVIpzhVbT1soqiUjqB2FknyThpGyyWzL7KCsENUfYtzR5cDolptQZ883keVV31ZFlqNijSiPngTd+8Du3UaWP/zQLFDi4QEXXww33wxXXgm+MomzOYggC50dEeRmIIIsCKcPWmsc5bk1pLmkLJ0yRxYVNZ2VsuHtGYW3Vxe8Pbvg4xWNly0Si6UTFfdxOSEj0USWj2yqWqgkIMqIcsxwI87WFlq+tTDPlI1b/S1sX2fu3zXeiPKYSRCTYKLJH35oWnIy+PubiPLNN8N554FVlkk/UUSQhc6OCHIzEEEWhNMfl8tBqSPTiHNpGiVlprlcFbWcFV6eEXh7RuPtGYWPVzTenl3cy353AgrTq2Q5bQe4HODhbapixI4yEWYv/5a5V362mdy36lvYtcGkgsT2qoosd4kzk/rmzTN5y/n50LUrTJkCt9xiysYJTUIEWejsiCA3AxFkQeicVEabS1MpLktzi3Nqjdxmm0cw3p7R+Hh1cUecu2LzaCFRbK+Ul5pV/Y5sNK0410zqi+wHsSOh2ygIiGyZe+VkmHzl1d9B4mazL74/nHOpWeXPJwC++greew8WLYLyciPIN99shDnm+AVuhCpEkIXOjghyMxBBFgShOuXlhW5hTjWR5tI0ysqzK/s9rAH4eMXg4xWNj1cM3l7Rp++CJ9oFWQcgeQMkr69KxQiOc8vyiJaripGVZmR55SLYv93kQg8808jyWRdCUanJV543D9asMf3nn29k+eqrIaAFFko5zRBBFjo7IsjNQARZEITGcLpKKSlNo7j0KCVlqRSXplDmqJJmm0cIPl5dK5u3VzRWi1cbjvgUUZBmZPnIRkjfZdIjfELcsjwSugxsmbzlo0mw4mvTjiWDzROGnwvjLoUzxsGhw0aU582D/fvBx8dM6rv1VrjoIjPZTxBBFjo9IsjNQARZEISTweksobjsqJHmUrN1lOdV9nvawmtKs2f06TURsLTA5Cwf2WCqYpSXgs3H5C13c0eXm1tCTmvY94sR5VXfmvxlv0BTY/mcS6HvMFi3zojyRx+Zlfy6djWiPG0a9OnTMs/aQRFBFjo7JyPInbDGkSAIQsthtXrj75NARPA5xEZdR5+4P9G3+1+I6zKFiJAJeNlCsRcfIC3rWw4efYfdSf/iQMqbpGYuIq/wF8ocOXToQIVXAPQ8F879M1z3Fpw/A+LHQPpuWPlfmP9b+OFfsHepKTF3MihlFhqZ9hC8ugRmvALDxhlhfmIa/PFXcHAtPPgnSE01K/cNHw7PPmuWtR47Ft56y0z0EzoELpeLMWPGoJTiyJEjdR7z4IMPopRi3rx5dfYvWrQIpRS//e1vT+VQOxVZWVlMnTqVLl26EBQUxJQpU8jJyansdzqd/P3vf6dHjx74+/szbtw4tm3bVuMan332GUOGDMHf35++ffsyf/781n6MJnEahTEEQRDaBx5WPwJ8exPg2xuoWCGwgOLSoxSXHKGo9Ag5BZvJzl8HgNXqh69XN3y8YvD1jsXbq2vHXNzE6umupTzM5C1n7oPD6+DwWljzOqx9E6IGQNyZpiqGb+iJ38PDZuR42DgoKTI1ln/+Cr6cA5+/DXF9YNxl8O47UFJuosqzZ8Mdd8C998I115io8vjxYJEYUXvlhRdewLeB+tfr1q1j0aJFRNez+mJeXh733nsvY8eObdFxORwObLYWKnvYAbn11lvx9vZm7969OBwObrzxRm655Ra++uorAJ5//nnmzZvHDz/8QLdu3Xj88ceZOHEie/bsISAggDVr1nDzzTfz1VdfMWHCBL7++muuvvpq4uLiOOuss5o9vhZ9f7TW0rRmxIgRWhAEobVwuZy6qCRVZ+Wt10eOLdB7Dv9Hb9//uLs9ofclv6pT0r/U2fmbdUlphna5XG095JPH5dI664DWmz7UeuGftH73Oq3fvV7rRY9oveNLrQuONf8euZlaL3pf679N0fq6QVpfP0Trf92l9apvtS4p1nr1aq2nT9c6MFBr0LpHD62feELrpKTm37uds3PnzrYewgmRmJioExIS9ObNmzWgk5OTa/SXlJToQYMG6VWrVunu3bvr995777hr3H777XrWrFl66tSp+je/+U2993rggQf0FVdcUWPfDz/8oAMCAnRhYaH+8ccftdVq1e+++67u0aOH9vf3b5mH7IAUFhZqpZTesmVL5b5ly5ZpQCe5/x6NGjVKv/jii5X9ZWVl2maz6blz52qttf7rX/+qJ0+eXOO6EyZM0NOmTavzni31/jT0dwDYoOvwQokgC4IgtAFKWfDx6oKPVxcINOlv5c4iiktT3FHmFPLs28kp2AiA1eKNj3csvt5x+HrH4ePZtePkMisFoT1MG3YD5B4xUeXD62Dje6aF9oC4s0wL6nri9wgKg0lTTDuaBD99YdqLfzH5ymN/BX+9B55/HhYuhHfegZkz4fHH4YILTFT51782E/1Ocz50HiBZ21vlXrHKjxutCU0+3uVycfvtt/Pss88SHBxc5zGPP/44559/PqNHj66z/7vvvmPLli28/vrrjaZXTJ8+nf79+5OamloZjX7rrbeYMmUKfn5+gEkbWLRoEZs3b2696PGf/gRbtrTOvc44A158sdHDqstjBS6XC4CtW7fSvXt3XC7XcSljWmu2bNnCrbfeWme/y+ViSz3P2pbvTwf511UQBOH0x8Pqe1xqRpkjk6LSIxSVJFNckkx60V4AFFa8vbri6x2Ln3ccPt6xHafMXHA304ZcbSpiHF5n2paPTAvqBt3Phu6jzXEnStd4uOGPcN3vYftaWPY5/LgQFn8M3XrCuVfCJx9Bnh3mzoU5c+CmmyAoyNRVnj7dSIPQ6rz00kt06dKFq666iqSkpOP6N2zYwPz58+sVqvz8fO666y4+/fRTPJpQxaRnz56MHz+euXPnMmPGDHJycliwYAErVqyocdysWbMICgo6mUc6bfD392fChAk8/vjjzJkzB4fDwVNPPQWY7zvA5Zdfzssvv8yll15KbGwsM2fOxOl0VvZfdtllTJo0ie+//57zzjuPL7/8kpUrVxIfH1/nPdvy/RFBFgRBaKcoZVb28/KMICRgGGCizEUlyRSVHKaoJJnsvDVk5a0CwMsWgW+1KLPNIxilVFs+QuMEdIGBV5hmz4Jktyxv+xS2fQLBsUaU48dAYN35pvViscKQMaYVFcDqxbD8c3j/efjgRRg6Bi6+Eh58AFatNrnKs2fDq6/CmWcaUb7+erPc9WnEiUR0W5N9+/bx73//m/oqSpWVlTFt2jRefvll/Ot5T/7yl79w/fXXM2zYsCbf98477+Thhx9mxowZzJs3j/79+zNixIjKfovFQmxs7Ik9THNpQkS3LZg3bx73338//fv3x9vbmz//+c8sWbKE8PBwAGbMmIHdbufiiy/Gbrczbdo0+vfvX9k/YcIEXnvtNe6//35SUlIYN24cN9xwA/v27av3nm31/kiZNzdS5k0QhI6Iy+WguPRopTAXlR6uXDrbw+pfKct+3t3x8oxq/8JcQXEuHFoDh1abihgAofHQfYwR5uas4lc9BSPrGPgFwJhLYMKVENLVTOx74w3YudMsPHLzzR06qtxRyrzNmTOHO++8kwD3Yi8ul4ucnBxCQkL4xz/+wSWXXEJCQgKhoVWTO3NycvD19eWKK67g/fffJz4+nry8PKxWKwCFhWZFzC5dutQZkQYzsSs2NpaPP/6Ye++9lzvvvJO77roLgGXLlnHhhRdSXl5+Cp+84/L1119zzTXXkJqaWmdKTGZmJrGxsSxYsIBJkybVeY1Ro0ZxwQUXMGvWrDr7W+L96RB1kJVSNwC/B4YCvlprj2p9twK/A/oDTmA98IDW+pdqx4wEXgEGAanATK31vGr9kcBrwEVACfAO8JDW2tXQuESQBUE4HdBaU+pIrxLmksOVdZmtFh98vbvj59MdP+/4jiPM9iy3LK8ylTEAwnpB/Ggjy35hJ3ddlxO2rzNR5XU/QFmJScE4/2pTCeOXnfD66/C//0FpaYeNKncUQS4qKiI7u2rhnSNHjjB69GjWr19Pv3798PHxITU1tcY5o0eP5oEHHmDKlCmEhYWRlpZWQ5buv/9+PDw8eO655+jatf7c9hkzZrB48WISExM5evRo5a/rRZBrkpiYSEREBMHBwWzcuJEpU6Zw8803M3PmTADS0tIoKSkhPj6e5ORk7rrrLoqKivjhhx9QSlFeXs62bdsYOnQohYWFPPfcc7z55pts3bqVqKioeu/b3PfnZAS5LVIscjCC6wO8UasvAJgJrALKgceAxUqpBK11sVIqCFgEPAeMA8YDC5RS+7XWq93XeB8oALoBYcC3QDbw9Cl9KkEQhHaAUgpvzyi8PaMIDRwFQFl5HkXFSdhLkigqPkRBkYnIWi3ebmGOx9c7Hu/2Ksx+YTDgUtMK001UOWl11QS/iL5VsuwT3PTrWqwwZLRpRQWw+jtYugDefQY+eMEsbf3AH+CFF6qiyr/9Ldx3X4ePKrdHfH19a5R2q5CeLl26VKZUdOtWMyfdarUSEhJCWFhY5bG1r+nh4dGgHAPccccdPPPMM0ydOrXT5xo3xE8//cRjjz1GXl4eMTEx/OEPf+Dee++t7D9y5Ag33ngjKSkpBAQEcO211zJr1qzKf1ecTifTp08nMTERpRTnnXceK1asaFCOoW3enzZLsVBKTQCWVI8g13GMH1AIDNdab1ZKTQOeALq7S3OglHoPKNdaT1NK9QAOAL201vvd/b8BHtFa92hoPBJBFgShs+Aoz8NefMgtzEmUlZtC/xXCXCHNRpjbca3g/NQqWc49DCiI6g/xY80kP6+TjPIe3gM/fAo/fwn2AugSB+dfBeOvgJ17akaVR42CO+9s11HljhJBbkvsdjtRUVEsXryYMWPGtPVwhFo09/3pECkWlTdumiBfAXwIdNFaFyilXgTitdaTqx1zH3CL1nq4UmoyMEdrHVytfxiwCQjSWtdYRkkpNR2YDhAXFzfi0KFDLfJsgiAIHQlHeT72ahHmsnLza26LxRs/7+74+fTAzycBL1t4+4wwA+QdMaKctAryj5rocMxw6HEOdBtuFjE5UcpKYO33RpZ3bQSrB4yYABdcDd36wvsf1MxVvuUW+P3vYcCAFn+85iCC3DBaa5566ikWLlzI+vXr23o4Qi1a4v3pKCkWTUIp1Qd4C/iz1rrAvTsAyKt1aC4Q2Eg/7mNqCLLW+g3caR4jR46U2YqCIHRKbB6BBAcMIThgCFAlzEUlh7AXH6SgKBEwk/78fBLw90nAz6cHNo/Ahi7bugR1g6HXwpBrICcJDvxsZDl5Pdh8TX3lhHPMSn5NjYp7esO4y007ehCWfmZKxq1bAuHRcN6vYfkSSDwAr71mlrN+5RU4/3y45x647DJoQqkxoe1IT08nISGByMjIdrvkcWemLd+fdhlBVkoNAL4HXtJaP1Nt/4ucfAQ5WGtdW54rkRQLQRCEuilz5GIvPkBh8UHsxQdwuooA8LSF4++OLvv5xGO1eLfxSGvhcsGx7XBghSkf5ygG3zBTMi5hHIR0P/FrOspgw1ITVf5ljZHtYeeYiX2x/WH2HFMmLjkZ4uLgrrtM3rK7zFVbIBFkobNzWkSQlVLDMRPr/q61/k+t7q3Ar2vtG+beX9Ef5J7Ud6Baf1JDciwIgiDUj6ctGE/bcEICh5sqGWXHKCw+gL34IDkFW8jOXw8ofLxi8PPpgb9PAj7e3bCoNv4RY7FA9BDTyn8DRzbCwRWw6xvY+aWpsdxjHPQYC35NFFibJ4yeZNqxZPhxASxbCM/dC6GRcOG1sHEtrFgD//0vPPSQWa3vhhvgD3+Akcf9HBYEoR3SFmXerIANU4HiG6BiVkMpMAb4CnjQnf5Q+9xgYC/wDPASppLF58BFFVUslFLfY9IsbqeqisVsrXXdBfbcSARZEAThxHFpJ8UlR9wR5gMUl6YAGqU88PPujr9vL/x9euFpC2s/+csl+WZy38EVkLHH7IsaYPKVu58Nnn4ndj1nOWz6yazUt22VyVU+80KYeAM4vUzaxbvvgt0OZ51l0i+uuQa8vFr+2epg165d9O3bF4ulHU+4FIRThMvlIjExsf1P0lNK3QbMrqOrh3v/uUBRrb5faa1/dp8/CngZGIypg/xYA3WQSzF1kGdIHWRBEIRTj9NVYipkFB+gsHg/ZY4sAGweQfj79MLftxd+Pj2wWlpHDhul4JgR5YMr3JP7bBA3CnpOgC6DTRT6REg9BN9/bKLK9gKI6w0X3wBDxsHH8+Hll2HPHoiMNGXi7rwTapUua2kOHz6MUoqoqChsNlv7+aAiCKcQrTUOh4Njx46htSYuLq7O49qNILdXRJAFQRBanjJHDoXF+yks2oe9+CAuXQZY8PWOw9+3JwE+vdrHgiVaQ/YB2L8cDq6EskKTr5wwHnqOh8CG6+geR2kxrFwE330ESbvAx8+UibvwOth9wKRffPWVEfBf/xruvRfGjoVT8H1wuVxkZmaSl5cnC14InQoPDw+CgoIIDw+v9zcoIsiNIIIsCIJwajHpGMkUFu2jsHgfJWXHAFMdoyq6nICH1adtB+p0mHzl/cvg6BYjzxF9TVS5+9ng6dvIBaqhNezdatIvVn8H5Q4YdBZcfD2ExcMbb5rqFzk5pqbyn/8MV18t1S8EoZUQQW4EEWRBEITWxVFeUC26vB+nq4SKyX7+vr0I8O2Dt2eXto0uF2XDwZ9h3zKTgmH1gu5nGlk+kZJxAHlZZlLf9/+DzFQzqe+Ca2D0JfDFN2bFvr17TfWLP/7RVL+QVd0E4ZQigtwIIsiCIAhth9YuiktTKqPLxaVHAfCwBhDg24cAvz74effAYrG11QAhc69JwUhaaUrG+UVAz3Mh4VwIiGz6tVxOM6nvu4+qJvWdfRFMuhkSD8G//w3Ll5vFR377W5N+0f0kStIJgtAoIsiNIIIsCILQfigvL6SgeC8F9r3Yi/fj0mUo5YG/TwL+vn0I8O2DzSOgjQZXZuoq718Oqb8A2kSTe51vFiTxOIFV+1IPweKPTGS52A59h8Glt4AlCF58CT7+2Bx39dUm/eLMM0/JIwlCZ0UEuRFEkAVBENonLl1OUfEhCooSKSjag6PclLX39owmwK+POxUjum1SMeyZcOAnk69ccMyUiEsYD70vhOATqE5RVGgk+dv3IT0FImLgV1Og95nw1jtmSeu8PDOR7/774corwWo9ZY8lCJ0FEeRGEEEWBEFo/2itKXVkUGBPpKBoL8WlyUBFKkZvAnz74OeT0PqpGNoFaTth7xITXXY5IbIf9L4A4s5uelTZ5YT1P8I378HuTab6xXm/hnFXwleL4cUX4eBBSEiAP/0Jpk0Df//GrioIQj2IIDeCCLIgCELHo9xpp7BoLwVFeyks2lctFaMXgX79CfDtg9Xayktgl+Sb9Iu9P0BBarWo8gVm9b6msn8HfP0urFlsls0+8wKYOAV2JcHzz8Pq1RAcDL//vclTjog4VU8kCKctIsiNIIIsCILQsXFpJ0XFSRQUJZJv3025swCw4OfTg0C/fgT49m3dvGWt4dgOI8qH15rocERfI8rdRzc9qpyVZib0LZkP9nzoOcjkKbsC4KWX4LPPwNsb7rgD/vIXiD0BCReETo4IciOIIAuCIJw+aK0pLk0h376LAvtuysqzAfDxiiXQrx+Bfv3wtIW23oBK8uGAO6qcXy2q3Ot8CKl7ha/jr1EEP30J38yD1CQIi4KJN0K3IfDfV2Cee1HZW26BBx+Evn1P2eMIwumCCHIjtLYg/+LaxS+uHS10NYWq43Xdexs/pv7/Vz+27v/Xd3Rd90M1fdz1jaNp42549HV/lxobT91jaNZ4VMPjqfp/vd/RRp+mrhE0/q42ftX6vhsne2bl/1VD4zl+/E1/L5v2favvbqrafxYUuLfH/dfWq8MJQFXecoUsl5SlAeDlGWVk2bc/Xp6RrfN+aQ3HdlaLKpdDRB/oc7FZhMTahNxplwu2rDDpF9vXmjzlC6+DIRPgrTnw5ptQWgrXXAMPPQTDhp3qpxKEDosIciO0tiAnu1JI1iktci2Ndm9P7v9VfwLqfl19b3099R1z/P8xPyAaGXdjY6r/7g3vrb6n4WOaMraaZ53891E4nakt0zZseOGFF554KfcWTzwrX1f0Vb32wENkuwUpc+SQb99NQdFuikoOA+DpEUKAXz8C/Qbg4xXTOt/vknxTAWPvEhNV9g4y6Re9LwS/sKZdIykRvngHVn1rqlqcewWMvhw++swsZ52fD5MmwcMPw7hxp/Z5BKEDIoLcCJJiIbQlWjem2Q39v/rRTfmgUfPrxj7wNG1MLfNhrPK1bnw8DX3gaPzDV8Pf28aurJv6n9a4au0rw0EppZTpMkrNK0opxYmL+rBgwQtPfPHFT/niV7n1w0/54o8vvvhiOZFV3QTA1FvOL0qkwL4be/EBNC5sHkEE+g0gyH9Q65SP0y5I3Q6J35klrpWC2FHQbxJE9q/49VLDHEuGL+fCsgVmOeuzLoLzr4VFy8yEvowMUyLu4YfhV79q2jUFoRMggtwIIsiCILQl5bq8UpZLKXMLtPtrXUYJpRRRhF3bsVNEGY4a5ysUPni75dkPP3zxd2/9lB+B+OOlvNro6ToGTmcJBUWJ5BXuoLB4P+DC5hFCkP8AAv0G4e0ZdepluSAd9iyGfUuhzA7BcdB3IvQ4B2xNqMaRmwmL3jeT+ooLYchouPgmWLMNnnsODh+GoUNN6sU110gtZaHTI4LcCCLIgiB0JBzaQWE1Ya7YVt/nxFnjHC+8CCSAIBVAoAogkAACVSBBBGBTbbSEczul3FlMgX03efYd2IsPABpPW1i1yPIJLC19UgMohYMrIfFbyDkENl/odR70vRgCujR+flEBfD/f1FPOzXRXvpgKe1LhmWdg927o3Rv+9je46Sbw8Di1zyMI7RQR5EYQQRYE4XRCa00pZdixU6jt5FNAvi4gnwLydD5FFNc43gfvKmGulGez9VCdW57KnXby7bvJL9yBvSQJ0HjZIgj0H0iQ30C8PMNP3c21hoxE2P0tHF5n0jFizjBR5a5DobG0mrJS+OkL+GK2ScPoGg+X3QbZTpj1NGzaZER55ky44QaJKAudDhHkRhBBFgShM+HQ5RS4pTmPfCPPuoA8CiihpMaxgQQQrIIIJZhgFUyICiaIgE6Z8+woLyTfvpN8+47KCX5enlEE+Q0gyH8wnraQU3fzomxT/WLvEijONZHkvhdDz/PA07fhc11OWPs9LHwHknZBaCRccguU+MI/n4KtW6FfPyPK110Hls733gqdExHkRhBBFgRBMJTpMnek2USbc8klR+eRT0HllEULFoIIJEQFE6KCCMGIsz9+nabihqM8n3z7TvIKd1BcegQwdZaDAwYT6DcQD2sj0nqyOMtNibjEbyFjD9h8TD3l/peAXyPRbK1h22r4/G3YsQ4CQ+HSW6HYx4jyjh0wcCA8/jhcdZWIsnDaI4LcCCLIgiAIDVOuneSRR47OI0fnkkMuOToXO0WVx3jg4ZblIEJUMGEqlDBCTvsc57LyPPIKfyGvYBuljgwUFvx9exPsPwR/3z5YLKcoTSVrP+z8Gg6tNl93Hw0DLoOwhMbP3bMFPnkNtq6sEuVCT3jqX7BrFwwZYkR58mSpeiGctoggN4IIsiAIwslRpsvIJY9snUuuzqsU5xJKK48JJtDIsgolXIUSRuhpKc1aa0rKjpFXuJW8wu2UOwuxWLwI9BtAsP9QfL3jTk2E3Z4JuxeZFAxHMUQNMKIcM6zxPOXELfDJq7BtlVuUp0Kugn/Ngj17zEIjjz8Ol18uoiycdoggN4IIsiAIQstSpIvJ1FlkkU2mziZLZ9eYHBhIgJFlFUY4oYSpEDyVZxuOuGXR2oW9+CC5hdsosO/CpR3YPIII8h9MsP8QvDwjWv6mZUWw7wfYtQiKsiCwqxHlhHFgbeR7e5wo3wrZGp6aBfv3w8iR8MQTUkdZOK0QQW4EEWRBEIRTT7EuNrJcTZqrp2gEElAZZY5U4YQRhofq+JUVXK4y8u27ySvcRqG7bJy3ZzTBAUMI9BuEzcO/hW9YDofWwM6vIPsgeAWaCX19J4J3YMPn1iXK6Q7419OQlARnnWVE+eKLRZSFDo8IciOIIAuCILQNxbqELJ1NJkaYs3Q2hdgBMxkwlBAiVbi7ReCHb4eeCGgqYWwnt2AbJWWpgMLfpxchgcNMvnJLfiDQGo7tMHnKKZvAaoOEc6H/pRDUteFza4vyJbdAWpkpD3f4MEyYYGoqjxrVcuMVhFZGBLkRRJAFQRDaD8W6mHSdSbrOJENnkkFW5cInvvgQqcKJUBHuKHNoh40yl5RlkFe4ldyCbZQ7C7BafAkOGEpwwDC8WzoFI++IEeUDP4PLAd1GwMArIbJvw+dVF+WgUJh0KxzKMznKGRlw/fXwz39Cz54tO15BaAVEkBtBBFkQBKH94tIussmplOZ0nUkhhYCJMocRQqSKIEKFE6Ui8FOnqMTaKUJrF4XF+8jJ30xB0R7AhY9XN0IChhHoPxCrpQWXCS/ONctZJy6G0gKIGgiDfw1dBjWcMrF7M3z6qikTFxQKE2+GzUnwwovgcMDvfgePPgoRpyC3WhBOESLIjSCCLAiC0LEo0sVkVApzBplkV0aZA/Cni4oiWkXSRUXhr/zaeLRNp9xpJ7dgK7kFmyl1ZKKUjSC/gQQHDsPXK7bl0kscJWbRkZ1fQXEOhPeCQVdBt+GNi/L//mvqKEfGwAU3waIV8Pbb4OsLDzwA990Hfh3ney50XkSQG0EEWRAEoWNTEWU+pjNI1cdI0+mUUQaAP341hbkDLGiitaa49Ag5BZvJL9yBS5fhaQsjJGAYQf5DW25in7MM9i+H7Z+DPQNC4mDQryHu7PoXCtEatq6CD1+EpN0Q3w/O/jW89wksWADR0WYi37Rp4NG5lyoX2jciyI0ggiwIgnB6obUmh1zSdDqp+hjHdHplbWY/fGsIcwD+7VqYna4y8u07yC3Y4l7iWhHg24fggGEE+PZGtcSy365yOLgSti+E/KMQGA2DJkOPc6C+hU5cLli1CD7+D6SnwMAzoe8E+L/XYNUqs3z1rFlwxRVS8UJol4ggN4IIsiAIwumN1ppc8iqjy2k6nRJKADPxr4uKJFp1IUZ1wV+1cNm1FqS0LJOcgs3kFW6l3GnHwxpASOAIQgKGY/MIaP4NXC5IXgu/LICcQ+AXAQMvh17n1V9LudwB3/8PPnsD8rPhzAshpD888zwkJsLYsfDsszB6dPPHJwgtiAhyI4ggC4IgdC601uSRT5o+RqpbmIvdC5kEEkCMiiZGRROtotrlqn9aOyko2ktO/gYKi/cDFgL9+hEaOApf7+7Nj4hrDSmb4ZfPIHMv+ASbRUd6XwQ277rPKbbDV3NNKyuFc6+E4kB45jlIS4Nf/xr+9S/o20jlDEFoJUSQG0EEWRAEoXNjIsz5HNWppOhU0vQxynGiUEQSTowlmq4qmnBCsbRESkMLUurIIid/I7kFm3G6SvCyhRMSOJLggKFYLfXIbFPRGtJ2wPYFkLYdPP2h/yXQ71fgWU+1kLwsE03+/n/gYYPzr4UDefDCS1BcDPfcAzNnQnBw88YmCM1EBLkRRJAFQRCE6ji1k3SdSYpbmLPIBsATT7qqKLq6I8wB7Sgdw+VykG/fQXb+BopLU1DKRrD/YEICR+Hj1aX5N8jYY1IvUjYZUR54hVmdr76Icloy/O8/sHIR+AfBedfD8m2m4kV4uIkmT5tW/2RAQTjFiCA3ggiyIAiC0BAluoSjOo0UncZRnVq5RHYgAXRVXeimutJVdcFDtY+qDcWlR8nOX09e4Xa0LsfHqxuhgaMI9BuApb5Jd00l6wBs+RiObgHvIFP1os8F9ecoH9wFH74EW1dCeDSMugLe+ghWroSRI+E//4Gzz27emAThJGg3gqyUugH4PTAU8NVae9TqvxWYCUQDvwB3a603VusfCbwCDAJSgZla63nV+iOB14CLgBLgHeAhrbWroXGJIAuCIAhNpSJ/OUWnclSnkqrTKaccK1a6qihiVTdiVUy7WLDE6Swmt3Ar2fkbKHNkYbX4EBwwjNDAkXjaQpp38fTdRpSP7QTfMBh8FfSaUH/Vi+1r4b3nTGm4vsMgdDA8+yIcPQpTp5qKF11aINItCE2kPQnyRCAU8AHeqC7ISqlzgO+AXwPLgXuBPwO9tdb5SqkgYB/wHPAiMB5YAFyktV7tvsb3QAEwDQgDvgXe1lo/3dC4RJAFQRCEk8WpnaTpdJJ1Cod1SuUqf6GEEKdiiLXEEE5Ym5aS01pjLzlIdt4GCop2A5oA376EBZ3dvEl9Wpvc5C0fQeY+8I+CoddA/Dl1p064nPDjAvjw/6AwF8ZcBklF8PKr4OUFjz0Gf/wjeNYTjRaEFqTdCHLljZWaACypJchzAYvW+hb31wo4BDyqtZ6rlJoGPAF01+6BK6XeA8q11tOUUj2AA0AvrfV+d/9vgEe01j0aGo8IsiAIgtASVESXD+sjJLtSSCcTjcYHb7qpGGJVDDGqS5tWxnCU55Odv4Gc/A04XcV4e0YTFnQ2gf4DsSjryV20ourFlo8hJwmCYmDodRB3JtQ1qdGeD5++Bt9+CF4+MHoyfLEcvv7GVLl46SWYOLE5jykIjdJRBHkLMEdr/WK1fZ8D+7XW9yulXgTitdaTq/XfB9yitR6ulJrsPj+4Wv8wYBMQpLXOrzWG6cB0gLi4uBGHDh1q0WcUBEEQhBJdSoo+ymGdQoo+ShkOLFiIVlHEuoW5rSb6uVwOcgu3kZW3hjJHJh7WAEKDziQkYAQeVp+Tu6h2weF1sPV/kJcCIfFwxnUQU88S1ikH4N1nYMtK6NoDEsbBf96EvXvNAiMvvAAJCc16TkGoj44iyPuBf2itZ1fbNxdwaK1/q5R6G/DQWk+t1j8N+JvWupdS6hb3+d2r9VdElWO11kfqG49EkAVBEIRTjUu7OKYzTHRZp5BPAWBSMeItsXRXcYSooFYfl9aawuJ9ZOWtwV58wFS/CBhKWNDZeNnCTu6iLhckrYCtn0DhMQjvDWdcD9GD6z5+009GlFMPwZCxUBwML70MDgf85S/w0EPg53fSzygIdVGfILePqbZVFAC1/2UIBvZX64+voz+/Wn9d51f0CYIgCEKbYVEmchxNFGcxgjydz2GdwiHXYTa5trGJbQQTSHcVR7wljlCCWyVvWSlFgG9vAnx7U1J6jKy8NeTmbyYnf8PJ5ylbLJAwHuLHwL5lZsGRJf+AqIEw4iYI61nz+OHjYchoWPQBfPqqWWjkqftg9T745z9hzhx47jm4/npZtlo45bS3CPJc95hudX+tgCTgsWo5yI/XihC/Czhr5SD31FofcPffjslhlhxkQRAEod1i10Uc0skk6WSO6XQ0mgD8iXfLcjihrTrJz1FeSE7+erLzN+B0FeHt2cWdpzzo5PKUnWWw5wcjyqX50GMcDLsB/MKPPzY3Ez76DyxbAIEhMHgSzPkENm2CSZPg1VchPr7ZzygI7SbFQillBWyYChTfABWJV6XAWEzViSuBn4E/An+lqopFMLAXeAZ4CRgHfM7xVSzygNupqmIxW2s9q6FxiSALgiAI7YViXcIhncwhncxRnYZG44evW5ZjiSSi1WS5Ik85O28NpY5MPKz+hAaeSWjgKKzWk1ilr6wIti+EXd+AAvpfBoOuBFsdOc8HdsCcpyFxM8T3A8+e8OLLJn3jH/8w1S6sJzmpUBBoX4J8GzC7jq4eWuskdx3kx6mqg3xXrTrIo4CXgcGYOsiPNVAHuRRTB3mG1EEWBEEQOiKlupTD+ghJOpkUnYoLFz740F11o4eKI0pFtsrS1yZPeT9ZeauxFx/AojwJDRpFWNBoPKwnkRtcmAGbP4SklWaxkaHXQq/zwVJLeLWGVYtg3r8hOx2GXQjLd8G338GoUfDmmzB0aMs8pNDpaDeC3F4RQRYEQRDaO2XaQbJblo/oozhx4oM3CSqenpZ4wlopDaO4NI3M3J/Jt+9EKQ9CAoYTHjwWm0fgiV8scx9sfM8sOhLUDUbcDF3POD7PuKTIlIX76l2zbHXsGHj5HcjOhr/+FR59FHxOsvKG0GkRQW4EEWRBEAShI+HQ5RzRKezXSRzRR3HhIohAelriSVDxBKqAUz6G0rJMMnNXklu4DQUEBQwlPPgcvGyhJ3YhrSF5HWx6HwqOQfQQI8oh3Y8/Nmk3vP64Sb/oMxL2FcNHH0OvXvDGG3DeeS3xaEInQQS5EUSQBUEQhI5KqS7loD7MAVcSaaQDEEE4PS3x9FDd8VEnkSt8ApQ5csnKW0VOwSa0dhHkN5DwkHF4e0ae2IWc5bDnO9j2KTiKoOcEGHo9+NZaEtvlhO8+gg9fMl/3PQ9mfwr798NvfgPPPgshzVxGW+gUiCA3ggiyIAiCcDpQqO0c0EnsdyWRQy4KRYyKpqeKJ07FYlOnrsKro7yQrLzV5ORvwKXLCPDtS0TwOHy8Y07sQqWFptpF4rdg8YCBV5jJfLZaop+ZCu88BRuXQbfeYA+Ft+ZAeDj85z9wzTVSEk5oEBHkRhBBFgRBEE43snUO+11JHNBJ2CnCAytxKpaeKp4YFX3KJveVO4vJzl9Ldt5anK4S/HwSiAged+K1lPPTYPP7ZmU+nxAYdqNZaa/6uLWGdUtg9r8gNwsGToCFP8OmzXD55fDKK9CtW4s/o3B6IILcCCLIgiAIwumK1ppjpLPflcRBfZgyyvDGi56qB30sPQlRwafkvk5XKTn5G8jKW025046PVywRIePw9+l1YqJ8bJeZyJe1HyL6wlm/OT4/uagAPngJvv8YQqPAuy+8Phs8POBf/4K77jKLlwhCNUSQG0EEWRAEQegMOLWTI/oo+/VBDusUXLiIIJy+lp70UN2xKVuL39PlcpBbsIXMvJU4yvPw8epGVOj5+Pk0uIZXTbQL9i83E/nK7NDvV6Y0XO36yYlb4I0n4Mg+6H0mrNwPy5bD2LEwdy707Fnn5YXOiQhyI4ggC4IgCJ2NYl3Cfn2QPa595JKPBx4kqO70sfQkgvAWLxmntZOcgi1k5PxEuTMfP58eRIacj6/3CaRAlBaY+sl7l4JPMIycCt3PrplrXO6AL2bDZ6+DhydEjYA33ofycnjpJZg2TXKTBUAEuVFEkAVBEITOitaaDDJJdO3noD5EOeUEE0QfS096qh4tXgXD5Sonp2ADGTk/43QVEeDbh8iQ8/H2imr6RTL2wLq3ITvJlIU7cxoEdq15zNEkeOvvsGMdRPeFDamwchX8+temJFx4HctcC50KEeRGEEEWBEEQBHBoBwf0Ifa49pNBJhYsxKlu9FE96aq6tOjEPqerjOy8tWTmrcLlKiHQbxCRIRPw8gxr2gVcTtizGLZ8DE4HDLwSBk02UeMKtIbln8PcZ0wZOd8BMPtDCA2F2bNh0qQWex6h4yGC3AgiyIIgCIJQkxydyx7Xfvbpg5RSih++9FY96WNJwF/5t9h9nM5iMvNWkZ23FpcuJzjgDCJCzsXTI6hpFyjONZP4Dq4A/0g483aIGVbzmMw0ePUR2L4WogfBd1tg9274wx/gmWdkFb5OighyI4ggC4IgCELdOLWTw/oIe/R+UnQqALEqhgGqL11VlxbLVS4vLyQjdwU5+ebncUjgCCKCx+Hh0UQZT9sOa9+G/KMQOwpG3QZ+1dIoXC5Y9D58+CLYfKEkCj76FPr3h3nzYPjwFnkOoeMggtwIIsiCIAiC0DiFupBE134S9T5KKCGIQPpb+tBbJbRYBYyy8jwycpaTW7AFpTwICzqTsKCxeFibEOV1lsOur8xqfCgYejX0uxSs1RZISd4H/30YknZBl6EwfylkZsKTT8Jf/wpWa4s8h9D+EUFuBBFkQRAEQWg6Tu3koD7ETtceMsnChge9VU/6W/oQpAJb5B6ljiwycpaTV/gLFosX4UFjCAsajcXSBBEvzIANcyF5PQR1M7WTowZU9Zc74JNXYeHbEBAGyRb47gcYNw7efRfi41vkGYT2jQhyI7S2IBdoB3YcLXAlhap8VfWVqtZbvZ969tc+r75zj79+4/ds6TJBgiAIQvsiXWeyy5XIQX0YFy5iVDQDVF+6qa4t8jOgpOwY6dk/UlCUiM0aSFTYhQT6DWratY9shHWzwZ4BfSfCsCk1l6xO3AIvPwTHjkDoYHj/K7P/5Zfh5pulHNxpjghyI7S2IH/pSmah63Cr3a89oGq9borMH7+vKcfUf15TPgA05byqZzh+X33n1nUeqq79jZ+nqo2o9nkneq2TO6/6CBoaQ/3n1TkGVd+fjYbHWXMcTblnUz78NX5eUz5YNuW8uj94untUw+ed6PvQ6P2Ou96J/T092fOqn9PgMSIr7ZoiXUyi3sdu116KKSYAf3f6RU+8lGfjF2gEe/Eh0rK+o6QsFR+vbnQJm9i0GsrlpabSxa5vwD8CxtxVM5pcUgTvPgM/fArB3WBzBqzfCNddB6++aipeCKclIsiN0NqCfETbSdFFzbqGdv9fV77G/VpX6zdf6yacV/3/dZ1bdY26z63zPF2zn6ae1+A46z6v5tjrHyMncV7NY+o6r6nXq/yGNGn8VftO7H51nVd738meV/2cho6pf19dfx4F4cQ40Q8qde87/oNgxdcWVLVm+izu/fW99lUeBONJkLKZLZ4EK0+C8SQAG5ZOJPdO7eSQTmanK5F0MvHAg16qB/0tfQlRTaxMUQ9aa3ILt5Ke/QPlzkKC/IcQFXoBNo8mpHUc2wWrX4WCY2YlvjNuqBlN3rgcXn8MCvLB1gs+XAhRUWYFvgsuaNa4hfaJCHIjSA6yILQtWtct6RVfw/EfoOo+pmEpr6u/qec19IH0ZM+r/YGn8hjdtPOa+gG40ftVO49Gzm3qeXV9QGyov9HzdP3PV/cYmvKB+/h7uQCX++y6Xrvcr2vs05oiysmljELKqY0FCMTIcrDyJAibe2v2RSkfIvE+LSU6U2ezy5XIAZ2EExfRqguDVX9iVHSzfiPgdJWSmbuSrLxVgCI8eCzhQWMbz092lMCWj2D3IgiIgjF3Q2S/qv78bLNU9fqlEBwPP+6Gffvh73+Hhx4CS8vVgRbanhMWZKXUMyd5rxe11kdP8tw2QwRZEARBaAnKtYs8HORSRq4uI8+9zaXm69oi7YmFaHzopvyIUb7E4Es35UcQttMitaRYl7BH72OXaw9FFBNKCEMsA4hXcc1afKTMkcux7CXk23fgYQ0kKvQCgvwHN/49O7YTVr1qJvP1vwTOuB48vExfxeIic2aBwwXpPrDkJ7jiChNNDg4+6fEK7YuTEWQXkAaUnsB94oBRWutNJzXKNkQEWRAEQWhNKkQ6h1LSdDFHdBEpmPS7vGqTuP3xIEb50g0jzt2UH13xwUd5NHD19otTO9mvk/jFtZM88gnAn8GW/vRSPfFQJ19ezV5ymLTMb935yTF0CZvUeH6yowQ2fwCJ30FAtMlNjuxb1Z+eAi8/DLs2giUWPl1iqlt89hkMHnzSYxXaDycryGdrrdc18QYeQBkwUgRZEARBEE6eAu0gRds5QpERZ20nhSJKcVUeE44XMcqXHiqA3gTQQwXg1QzBbG201hzSR9jm2kEmWfjgzQBLP/qr3nie5IQ+rTV5hVs5VpmfPNidn9xI3nPqdlj9GtgzYcBlMPS6quWqXU74+L+w8C2wRcD3O6CwEN56C2688aTGKbQfTkaQZwNPaq0PNvEGCngHmKm17nDlGUSQBUEQhPaMS2uyKK0hzke0nVSKAZPnHIsfvVQgvVQAvVQgocqrbQfdBLTWpOljbNM7SdGp2LDRT/VmoKUfvurkln92usrIzF1BVt5qAMKDxhAePBaLpQHxdhTDpvdhz/cQ2NXkJkf0rurf8CO8/DcoKYfEUti8Df74R3juObC1zAIpQusjk/QaQQRZEARB6IjYdTkHdAF7dT77yOegLqTMHWkOxbOGMHfDD2s7zmfO1Nn84tpBkk5GoeitEhhsGUCgCjip65n85B/It2/HwxpIdPglBPr1bfiko9tMNLk4GwZcAUOvAatbrNOS4fn7ICkR7BHwzY8wdizMnw/R0Sc1RqFtEUFuBBFkQRAE4XSgXLs4gp19uoB9Op99uoAcygDwwlKZktFLBdJbBbbLtIx8XcAvrp3s1QfQaOJVLIMtAwlXJ1ePuKjkMEczv6a0LJ0A3/5Eh09quCxcWRFsfA/2LTWr8I25C8J7mb7SYnj7H7D8CyAKvloNgYHwv/+ZVfiEDkWzBVkp1RW4HIgBvGt1a631g80eZRsigiwIgiCcjmityabMLctGmJOxowEPFL1UIINUMINUCN3wbVcVM4p0MTtcu9mt9+LAQYyKZphlCJEq/ISvpbWTzLzVZOQsR2ElKvQCQgJHNvy8R7fA6jegOBeGT4H+l4JSpsrF0k/hnafA6QVrjkLKUZNu8cc/IqvvdRyaJchKqRuAuZj66Rng/ihahdZaJ7TEQNsKEWRBEAShs1Csy9mvC9ihc9mhc0nBLFwVhI2BKpiBKoSBKpgA1T5ya8t0Gbv0Xna4dlFCKbEqhuGWoYSpkBO/liObo5lfYy8+gI9XN7pGXI63Z2QDJ9hNObjk9dBtpIkme/mbvv074IX74VgapHjC6g0wZQq88Qb4+Z3k0wqtSXMFeT+wFvid1jr/FIyvzRFBFgRBEDorObqUHTqX7TqXnToXO+UooDv+DFTBDLIEk0AAHs2oV9wSOLSDHTqR7a6dlOGgh4pjmGUIwSe4Op+pdvELaVnf4XSVEB48hojg8fUvMqI17P4WNr0HPiEw7l6I6GP6CnLhvzNg8wooCYevfoKBA00puN69676e0G5oriDnAldrrX84BWNrF4ggC4IgCIKplpFEoRFmVw4HKMAF+GClnwqqTMcIV7WzLVuPUl3KdtcuduhEnDjppXpwhmUwAcr/hK5T7iziWNb35BZuwdMjlOiIS/H3aeAX4pn74OeXwJ4Fw2+E/peZdAqXCz57HT55FZwB8MMucDrhvffM4iJCu6W5gvwOkKK1fvRUDK49IIIsCIIgCMdTpMvZpfPYoXPYrnPJcq8fFosfIyxhjFBhdFW+bTK2Yl3CNtcOdus9uND0Vb0YahmE3wmOp7D4IKkZX1FWnk2Q/1C6hF2Mh7Wea5TZTZWLw+sgZjiMvRu83FU2tqyE/zwIeUWw0w679sBjj8Hjj0tecjuluYLsC7wN2IGlQG7tY7TW3zR/mG2HCLIgCIIgNIzWmmMUs03nsNGVxX4K0EBXfBihwhhhCW+TiX52XcRW13YS9T4sWOinejPEMhCfE4hyu1wOMnJ/JjN3JVaLF13CJhLkP6TuZ9HarL638V3wDobxf6pKucg4Ci/8GfZsgww/+Gk9TJ0Kb74p9ZLbIc0V5DOAT4Ee9Ryitdbtr07MCSCCLAiCIAgnRo4uZbPOZqPOIlHnoYFIvN2yHEY8/q0qywW6kM2uX9ivD2LFykDVl0GW/nidwIIpJWXpHM34kuLSI/j59CA6/DK8bPWUl8vaDz+9aFIuht1gVuFTFnCUwdynYfHHkOkFP26CiRPhk0/A/8TSQIRTS3MFebP75cPAPo6vYoHW+lBzB9mWiCALgiAIwslToB1s0lls0lns0nk40YTiVSnLPQnA0kqynKvz2OzaxkF9GE88GWTpz0DVF1sTq3Jorckp2MCxrB/QOIkIOZfwoDGouiYpltlh9etweO3xKRc/LoA3n4RsBUu2wdCh8PXX0KVLCz6t0ByaK8h24Cqt9XenYnDtARFkQRAEQWgZ7LqcLTqbjTqTHTqXcjRB2Iwsq3D6qMBWkeUsncMm11aSdQreeDPSMpReKgFLE6txOMoLSM1cREHRLny9uxMTMRlPW/DxB2oNexbDhnfBO8hUuYh0r9i3bbVZfS+9GJbtNXL87bfQp0/LPahw0jRXkH8EPtBav3kqBlfrXl2Al4DzAQ9gM3Cf1nqru/9WYCYQDfwC3K213ljt/JHAK8AgIBWYqbWe19h9RZAFQRAEoeUp1uUmZ1ln8YvOoQwXoXgxVkUwxhJJpPI55WNI15msc24knUzCCOUs6wi6qAZqH1fDlITbRmrmN6AU0eGXEuw/uO6Dsw64Uy4yYNiNVSkXSYkw625ITofVR8Biha++grPPbrmHFE6K5gryMGAO8Cz1T9IravYozb0+AwKA64FC4B/AjUAcMBb4Dvg1sBy4F/gz0Ftrna+UCsKkgDwHvAiMBxYAF2mtVzd0XxFkQRAEQTi1lGonW3Q2q3Q6O3QuGuhDIGMtkYxU4XifwmWvtdYc0IfY4NqMnSJ6qO6MsgzDXzVtQY8yRw5H0hdQXJpMkP9gosMuwWqtYxJgWRGseR0OrYGYYTDmbvAOhMw0mHUXJCbC5hzIzoOPP4bLL2/hJxVOhOYKssv9st6DW2qSnlJqG/BfrfUb7q/7AruBCODfgEVrfYu7TwGHgEe11nOVUtOAJ4Du2v1gSqn3gHKt9bSG7iuCLAiCIAitR7YuZbXOYKXrGMcowQsLI1Q4Y1XkKU3BcOhyfnHt5Be9E4AhlgEMVgPwUB6Nnqu1i4zcn8nIWY7NI5CYiF/j59O9rgOrUi58QuD8ByE4Fuz58O/7YOMq2FEMB5Ph1Vdh+vSWfkyhiTRXkG+jATkG0FrPPenR1bzXLcDNwE2YCPI/gTO11uOUUluAOVrrF6sd/zmwX2t9v1LqRSBeaz25Wv99wC1a6+F13Gs6MB0gLi5uxKFDHXqeoSAIgiB0OLTW7KeAla501ulMSnASjhdjLJGMUZFEnKIFSQp1IetdmzmoD+OHL6Msw+ihujep6kZRyRFS0j+jrDyH8OBziAiZgKWu6HfmPlj2LJSXwvj7oOtQKHfAq4/Cj1/A3nLYvldqJbchzRLk1kQplQC8AVwAOIFk4Fda693uJa//obWeXe34uYBDa/1bpdTbgIfWemq1/mnA37TWvRq6r0SQBUEQBKFtKdVONussVup0drnLxvVTQYxVkYxQYXidghSMNH2MNc6NZJNDFBGcZR1JuKqnrFs1nK4y0rK+JbdgM96e0XSLvAovz/DjD7Rnwo9PQ+4ROPN26HORWXnv4//AZ29Ckgs27ILbb4fXXpNaya1MfYLctouq10KZ+ilLgD1AEOCLiSD/rJSKAgrc+6sTDOS7XzfWLwiCIAhCO8VLWTnbEsmfrYN42jqSyZY4snUpb7v2cp9zHXOcezmsC1v0nl1UFFdYJzHGciZ55POFcxErnGso1iUNnme1eBITcQWxUdfhKM9lf8obZOdv4LjAo184THwSoofA2rdM2gXAjffC9McgwQPGDIB33oErr4TCln0+4eSoN4KslNqESU3Y0aQLGbndAEzRWu8+qcEoFQ5kAAO01ruq7c8Cbgeuco/5Vvd+BSQBj1XLQX5ca9292rnvAk7JQRYEQRCEjofWmn0UsMJ1jHU6kzJc9CGQCyzRDFNhWFswLaFUl7HF9Qs7dSIeeHCGZTADVB+sjUSuHeUFpGQsxF58gADfvnSNuBwPa63Jfy4nbJhrVuDrNhLOuQds3rBxGbz0ACTb4afdMHy4qZUc2bQqG0LzOOEUC/fEvAkY6W0KHpjqFiO11ptObpiglEoEFgMzgFLgVuA1oB/QFfgWuBL4Gfgj8FeqqlgEA3uBZzCl4sYBnyNVLARBEAShw1Oky/lZH2OpK5VMSgnFi/MtXRinovBv4iIgTSFX57HWtZEUnUogAZxlGUGsJabBc7TWZOev5VjWEqxWH7pGXEmAbx3ZnbsXGVEOiYfzHgDfUNi3HZ7+PSSlw4okiIkxtZJ7NZgdKrQAJyvIJ5OgPKqZgtwfU07ubMCGKdv2pNb6c3f/rcDjVNVBvqtWHeRRwMvAYEwd5MekDrIgCIIgnD64tGaLzmaJPkqizscTC6NVBBdYuhKjfFvsPsmuFNa6NpJPAT1VPGdZRuLdyLLVJaXHOJL+GaWOdEIDzyQq9EIsllryfmQT/PwSePrCeQ9CaDykJZsycLv3wZpU8PQykeRRo1rseYTjORlBnlpnR+N8obXOOclz2wwRZEEQBEHoeCRrOz+4jrJGZ+LAxQAVxAWqK0NUSIuUinNqJ1tdO9iqt+OFF6MtI4lXcQ1Wu3C5yjmWvYTs/LV4eUYRF3U9nraQmgdlJ8GPz0BZoVl5r9sIyM+BZ++BTethYzaUOGDpUpN2IZwSOkwVi7ZCBFkQBEEQOi4F2sFP+hg/ulLJoYxIvDnfEs05KhKfJtQ4bowsncMK5xqyyKa7imW0ZRS+jawCWFC0lyPpn6FQdIu6Fn+fHjUPKMoxkpxzEEbcCv1+BY5S+L8ZsGwRrMkALLB8OQwc2OxnEI5HBLkRRJAFQRAEoeNTrl1s1tkscR1lHwV4YeEcFcUFlmiimrmstUu72K53sdm1DSsenGUZQS/Vo8Focqkjm+S0jyh1ZNIl7GJCA8+qebyjBFb+F5LXQ9+JMNL9C/x3noLP3oM1x8DLF376CXr3btb4heMRQW4EEWRBEARBOL1I0gUscaWyTmfiQnOmCucySyxdm5mnnKfz+dm5hnQyiFHRjLWc1eCS1U5XKSnpCyko2k2w/1Ciwy/DYqkW1dYu2PQB7PwSup5hUi48vOHNJ2HB+7A6DYJC4OefoXsdK/cJJ40IciOIIAuCIAjC6UmeLmOx6yhLdSoOXIxS4VzeTFHWWrNL72GDazOgGGUZRj/Vu95ostaajNzlZOQsx8crhtio67B5BNY8aM8SWPc2BHUzy1P7hMArj8IX/4NVqRDVxUSSYxquqCE0HRHkRhBBFgRBEITTmwLt4DtXCkt1KmW4GKnCuMwSS7cGor+NX7OQFa61pOo0uhDJOdazCVQB9R6fb99NSvoCLBZPYqOuw9c7tuYBqdtg+QtgtRlJDu4O/5kBXy2AVUehe7zJSZY6yS2CCHIjiCALgiAIQuegQDtY7DrKD/oopbgYocK4ohmirLVmr97POtcmnLgYbhnKQNUXi6p7weKSsnSS0z7CUZ5HdPilhATWqlKRewSWzjIVLi74m6mZ/MKf4ZuvYHUK9OtvqluENr4kttAwLSLISqk+QDfAu3af1vqbZo2wjRFBFgRBEITORWGlKKdSgpPhblGOPUlRtusiVrnWkaxTiCCMc6xnE6KC6zy23FnMkfRPsRfvJyRwJNFhk1DVV+yzZ8L3T0JJvpHk4O7w3L3w3bew6ggMGw7ffw+BgXVeX2gazRJkpdQA4GNgAFBXco3WWje8DmM7RwRZEARBEDonhdrBEtdRluhUinEyTIVyuSWW7sr/hK+lteaAPsRa1wbKcDDMMpghamCduclauziWvYSsvNX4encnNuramktU15bkoFh4+g+weAmsPgJjxsCiReB38ikinZ3mCvLPQCTwALATKKt9jNb6UAuMs80QQRYEQRCEzo1dl7PEdZTv9VGKcXKGCuWKkxTlYl3Catd6kvRhYlQ04y1j8FHH/QIegNyCbRzN/BIPqx+xUdfj4xVdbVC1JDkgBv51FyxZBmuOwAUXwBdfgHfd1xYaprmCXAjcoLX+6lQMrj0ggiwIgiAIAkCRLud7fZQlrqMUuSPK11riT7iOstaaRL2Pta4NeOHFBOs5dFF1T64rLj3K4bSPcbqKiIm4giD/wVWdtSXZvyv8czos+RnWJsPll8Onn4LNVue1hfqpT5Drzh4/nv3UkXcsCIIgCIJwuuGrPLjSEsfT1pFMtsSxU+fyqHMz/3MepEiXN/k6Sin6WXpzmXUiHniwyLmEra7t1BWc9PHqSs+YO/Dx6sqR9M9Iy/oerV2m0y8cLnoMvAPhh39C4VF46FU4bzSMioUvv4SbboLypo9NaJimRpAvBJ4BrtFaHzjlo2oDJIIsCIIgCEJd5OkyPnMdYqVOxx8bv7bEMU5FYWlgBb3alGkHK11rOKgP0011ZbxlNN51pFy4tJO0rG/Jyd9AgG8fukVeW7WoSO1IsncEPPkbWLoONh2BW2+F2bPB0tT4p3DCKRZKqfVA9c7uQAiQBOTWPl5rfWZLDLStEEEWBEEQBKEhknQhHzkPspd8YvHjRmsP+qqgJp+vtWa33sta10Z88GaCdSxR9aRcZOWtIy1rEX4+PYiNugGrxdN01JZkz1B44nZYugG2HoXf/Q5eeQVOQN47MycjyHOoKcgNorWedtKjaweIIAuCIAiC0Bhaa9brTOa7DpFNKSNUGNda4omoZwJeXWTqbH50/kwhdkZYhjJYDaizykVuwRZSMr7Ax6sb3btMwWp138OeCYufhNJ8uPARsAbBzKnw41bYfhTuvx+ee04kuQnIQiGNIIIsCIIgCEJTKdNOvtNH+cZ1BBeai1UMl1hi8FEeTTy/jBWuNSTpZGJVV8ZZxuCtvI47Lq9wJynpn+LlGUX36JvxsLqXx64tyfjD41Nh6TbYlQb//Cc8/HALPvHpiQhyI4ggC4IgCIJwouToUj5xHWKNziAIG1dZujNGRTYpP1lrzS69h3WuTfjgzXnWc4hUEccdV1C0l+Rj/8PTI4Tu0bdg83AvZV1bkl0+JpK8eBsczDQ1kidObOlHPq0QQW4EEWRBEARBEE6W/bqAj5wHOEAh8fhzg7UHvVXTVrnL1Fksda7Ajp3/b+/O4+SqyvyPf55b1dX73p2lsyeEnSQknbAEZHNBUFAERQVEREYZl3EUx+2njDKMC+PoOM7IIijuIszoyI7sBEhIAgl7FrJ3d7rT6fS+1T2/P6o6qa5UdfVevXzfr1e96tY95577nHMv9HNubt2q9E7keDv6sFsuWtreYkf17wgG8pgz/QpCGUXRgrgkuSsDvn45/GU9+Bmwbh3MnTu8nZ1AhvqYNxERERFJYoHl89XAIq72FnKATr4b3sjN4Teodx0pty2zUi4MvJvZNpM1/joe8Z+gI2673Ox5zJl+Bd1+G9v23EFH575oQRm885uQWQCP3AAZXfCV/4BTZkFrC1x8MbS3j0SXJzRdQY4a7SvIr3V08EpH6v9oUrGeV3SmaQlemB1aji+L2T7h+riyw/Y7wP3Htz2Y/Sfqx3D0P2nZMMSQ6MsXIiIyMXW4MPf7u3jA7SGA8SFvLqfb1JR/C5xzvOreYI2/nhyyOStwOuVW2qtOW0c126t+hZkxZ/rlZIWmRgrirySvWwvf/jI8swOuvhpuvXWkujuu6RaLFEY7Qf7v/fv5z/37R21/MjYkSsAPru9rkhBTTrI6fU0SErR/2ISnH/tP1P7BfiRovz/7iO9D0v33J4Y++tlXHxPFkHA/g+xjf/oXu33CvifYP0m2T3WeDGcMwzkGscdwMPtPdg7E7z/p2IxiDH3tPzaJSrQfTbjHj1rXzi/8zbzuDnC8FfEx7whKEnwR7/Dt6ngs/DRttHOWdxqzvZm9yjs669hWdSfOdTNn+mVkZ1ZECmKT5HO+Dnf9En5yM7xeB7fdBp/4xEh0c1xTgpzCaCfIzb5Ps+8PqQ0XfeHcweXYV6TI9f6cYDlZndgykpUPYP9J2xnA/hPVGY7+99X+cMYQu34oY3BYPEn2Ed92qmNw2L4TjE/CfvSx/0RxDDiGJH3sK4YBj2OKfqbs4yD3H9sHkYFIlqCn/BcxEk+k+pugE7P9QCZJPWVeTLkXfWF2cLlnGy+6jp56ZoTMmB4MUhEMMiP6Pj0YpCwQGLMTB985HnfV3OVvI4DxYW8+p1p5ynjbXBsPhx9nH/s51VvBUd4Rvco7u/azrepOwuFWZk/7CLnZcyIFLXXw4PXgd8Hbr4cbPwe/fhjq2+GZZ2DZspHp6DilBDkFfUlPRNItURKeKgGH5Ek48W2lSNIHO9lLVmcg+z/4l2gI++/PZDfp2IxiDH3ufwgTteHaf6Jy4rYfyHmYrNwBfs97dF385551PZ/bnKOqu5vGuAtMmTGJc69XRgYVwSBTAgECaU6g97o2bg9vZhONLLZirvCOoMhCfW7T5bp41H+K3a6KJXYCJ3on9Eqsu7ob2VZ1J13dB5g99VLychZECvbvgAe+AYUzYek18I8fgr9uhJIpsHYtlJYm2ePkowQ5BSXIIiIi40Oz77Onu5s9XV2R97jXvnC4V/0gMDsjg9NycjgjJ4elWVmE0pAw+87xN1fF3f52Qnh81JvPCivr82qy73ye8Z9nk9vKkbaAU70VeHboGQvd4Ra2V/2Kjs46Zk69mILcoyMFO1bDE/8G88+A3MVw3cfgsW3wjnfAX/8KgcAI93Z8UIKcghJkERGRiaHd96nq7mZ3NGGu6u7m5Y4O1rS10QXkmnFKdjZvy8nhbTk5lAf79+Mew6XatfLz8Ga20sQyK+Uybz4FfVxNds6xzn+Jl9wrzLIKzvROJyPmB0nC4Ta2V/+Gto49zJzyfgrzTogUvHQXbPgTVF4Br1fBDd+CdXvgW9+C668f4V6OD0qQU1CCLCIiMrG1+j7PtbXxZGsrT7W2Uh290nxsKHQwWT4+M3NUbsfwneNBt5v/9XeQRYDLvQVUemV9bvOa/ybP+S9QSgnvCJxJdszPW4f9DnZU/47W9u1UlL2X4oKl4Hx44oeway2c/RX4ze3wX7+Cbfvh3nvhvPNGuptjnhLkFJQgi4iITB7OOd7s7OTJ1laeaG3lpY4OfKDY8zg9J4fTc3JYmZ1N4QjfirDHtfLz8Ca20cwKK+Oj3nzyLCNp/e3+Th73nyGXHN4ZOIsCyz9Y5vtd7Kz5I81tm5ledj4lBZXQ1Qb3fwPaGuCsr8N3PgO/exK6g5EfEZk/f0T7N9YpQU5BCbKIiMjk1RAO80z06vLTra00+D4BYElWFufk5HBhfj5FI5Qsh53jfreLv/g7ySXIx7wFLPGSf5Guxu3l4fATeHi8M3AmZTHPSvZdmJ01f6C5dTOzp32E/JwjoKka7vsa5JTACVfBP34Y7n8NjjoWnn0WsrNHpF/jgRLkFJQgi4iICEQS1o0dHTzR2sqTra283tlJyIx35ebywYICTszMHJHHyu10Lfw8vImdtHCKlfNhbz65lvj+6AZ3gAfDj9JBJ2d7pzPTqzgUv9/Jtj130NlVz7wZnyArNAX2bIBHb4SZlZB5Anz57+Dp7XDllXD77TBGH5M30pQgp6AEWURERBJ5s7OTPzQ28n9NTbQ4x5GhEB/Mz+e9+fnkeV7qBgag2/nc6+/iXreLYkJ8NnAMMy03Yd1W18pD4cfZTwOneSez0Dt0u0RXdyNbd9+Gmce8iqvJCObBq/fC2jth8SWwcQfc+F14tRZuvhmuuWZY+zFeKEFOQQmyiIiI9KXF97mvuZk/NDbyWmcn2Wa8Jy+PDxUUcExm6l/IG4itromfhl+njW6u8hYm/QJfp+vkb/6TVLkalnmLWWTHHby63daxh7f2/IKs0BTmTv8YngVh1U9h61Nw+hfgl7fAzXfDvjZ4+hlYvnxY+zAeKEFOQQmyiIiI9IeL3oLxx8ZG7mtpocM5FmVm8sGCAs7NzSV7mK4qN7hOfhp+na008V6byQXebLwEt0KEXZin/GfZ6rZztC3kZK/y4LOSG1teZ2fNHyjIPZaZUy7Gwl3w0PVwYDecdh186+/hT6uhsBTWvwhlfT9JY6JRgpyCEmQREREZqAPhMH9pbuaPjY1s7eqiwPO4MC+PDxYUMD/U9y/l9UeX8/m1v4Wn3V6WWAlXewvJTnBfsnOONf56XnavMcdmcYZ3KsFovbqGVdTUP0xZ0elMLTkbWvbBfV+FjCw46nL4x8vhkc1w5lnwwAOT6kdEkiXIw3vjzDAxs7eb2XNm1mxmdWb2XzFlV5jZFjNrNbPnzWxZ3LaVZrY6Wr7FzC4b/R6IiIjIZFAYCHB5YSF/mTmTX0yfzsrsbH7X2Mh7d+3iyj17eLSlhaFcjMwwjyu9I/iIN58Nrp5/CW+gxrUdVs/MWBFYykneMra7nfzNfxLfRX6Su7TwFIryT6Su4Skaml6C3FI444vQUge7H4Iv3QBLpsIjj0R+RETGXoJsZmcCfwJuAkqBmcBt0bLTgP8GPg0UA3cD95lZQbS8ELg/ur4Y+BTwMzM7ZVQ7ISIiIpOKmbE8O5ubpk7lb7Nn8w8lJezp7uazNTV8eM8enms7PKkdSNvneNP5onccTXRxQ/glXvb3J6x7nHc0K72T2O2qeMp/DuccZkZF2fnkZs1jT+1faGnbDlOOghVXQ9UGKGqHT14D84rhX/4F7r9/0LFOFGMuQQb+FfiZc+5PzrkO51y7c25dtOyTwD3OuYeccx3AD4AO4P3R8ouANuD70W0fBv4HmJxfzRQREZFRVxYM8smiIu6bNYvvlJVR293NJ6qq+ERVFRvb2wfd7tFeEd8ILKaETH7kv8oD/u6EV6eP8o5gqbeILe4tXvDXA2AWYNbUS8jIKGZnzR/o6KqHhWfDUe+CV/8PzjwdLnkHFGbBtZ+Gjo5BxzkRjKkE2cxygRVAu5mti95e8biZ9dwbshhY21PfRc6K9dH1PeXrXO+zZV1MuYiIiMioCJpxUUEB982axT+VlvJGRweX7tnD56ur2dLZOag2yy2LrwUWscxKucvfxm3+Jjpd+LB6i+14jraFbHSv8bL/GgCBQDZzpn0EgB3Vv6U73AaVV8DUY2HN7XD1F+Ck+bBtO/z4x4Pv+AQwphJkIrdFeESuFF8JVAAPEbmNogjIBw7EbdMAFESXU5X3YmbXmNkLZvZCbW3t0KMXERERiZPpeVxRWMiDs2fz98XFPNvWxvt27eLre/eyp6tr4O1ZgE95R3GRN5vnXS3fDW+k3vW+4mtmnOxVMtdmsdpfxxZ/GwChjBJmTbuUrq4GdtX8Ed8M3vYFyC6CDXfC574E0/Phn/8ZamqGoffj01hLkJui73c45zY45zqJ3HKRAZwaLS+M26YIaIzZvq/yXpxztzjnKp1zleXl5UOPXkRERCSJXM/j2uJiHpw9m8sLC7mvpYXzdu7kX+vq2Bc+/CpwX8yM871ZfMY7hhra+U74JTa53umOZx5v81YyjSk85T/Lbr8qEkfWbCrKL6ClfRtVdX/FZebDmV+CjiawLXDuSdDeBl/72rD1fbwZUwmyc+4AsA1I9HVPB7wELO1ZYZEnYS+Jrif6fmLcdifGlIuIiIikVXEgwJdLS7lv1iwuyM/nd42NvGvHDn5SX0+T7w+orSVeCV8PLCKbAD8Iv8wTfnWv8qAFOCdwBoUU8Df/SercPgCK8hdRXvQ2GppeZN+BVVAyF069FvZtgveeBQtK4I474MUXh6nX48uYSpCj/gv4uJkda2ZB4DqgHVgF3ApcZGbnmFkI+CKQReSLeETfc8zsOjMLmdk5RL64d8uo90JERESkD9ODQb5dXs6fZ87kjJwcftbQwLt27OCOhgbaB5AoV1gO3wgs5hgr5E5/C78Kb6HbHdo+00K8K3AWWWTyUPhxGl3kH+zLi8+kIPd4auofobHlNZh7Ciw4A1pfhg+/F0IB+OzfwyT8zYyxmCDfBNwOPArUAe8G3u2cO+Ccexq4lkiifAD4IHCec5F/U3DONQDnAZdEy28FPuWce3a0OyEiIiLSH/NCIf5t6lTumjGD4zMzuam+nvN27uT/mpr6/QzlHAvyee9YzrUZPO6q+U//dbpikuQcy+FdgbNxOB4MP0qba8PMmFF+IdmZM9m19x7aOvbAiR+BQCYsngaLpsPTq+B//3eEej526Zf0ovRLeiIiIjIWrGlr44f19Wzo6ODsnBy+WVZGefDwX89L5gm/mjv9LSy2Yq71jiZoh66H7nV13B9+hCIKOS/wdjIsg+5wC1t334Zz3cyfcTUZm1bBC7+Apjnw1f+AwnLYtBkyM0egt+k1rn5JT0RERGSyWp6dza8rKriupISno0+8eKC5ud/bn+FN46PefF5y+7nZf7PX7RZTrIyzvdOpZz9/858k7MIEA7nMnvZhfL+L7dW/I7zwDCiaBaVNcNpRsHMX/Pu/j0RXxywlyCIiIiJjTMCMK4uKuHvGDGZlZPDFvXv5Yk0N+/v5tIuzvelc6s1jndvHbf4mwjF3DMzyZnCadxJ7XDVP+c/inCMrNIVZUy+ho3Mv1fsfgRVXQWd95F7kinz4zrehurqPPU4sSpBFRERExqj5oRC/rqjg88XFPNLSwoW7dvFoS0u/tn2HV8El3lzWuDpu9zfhxyTJC70FLPOWsNVtZ7Uf+cHivJwFlBaeQkPTeloKcmHuSuAteM+p0N4OX/nKSHRxTFKCLCIiIjKGBc24priYP86YQXkgwGdravjq3r009uNq8rneDC7yZvOcq+UX/uZeSfIiO5Zj7Sheca+z0X8VgCnFZ5ARLKSq7l78pR+BQBDOOgGOKIE774T160esn2OJEmQRERGRceCozEx+N2MGny4q4t7mZt63axdPtbam3O58bxYXeLN4xu3lV/6Wg0mymXGSt4x5Nps1/no2+1vxvBDTy86jo6uWfZ2vwaIPQNcO+Gj0sW/XfnpSPPZNCbKIiIjIOBEy4zMlJfx2xgzyPI9PVVfzrdpaWlI8N/kCm8X5NpMnXQ2/9bcefHycmfE271Sm21Se8p9jl7+H/JwjKcg9htqGJ+lccAoUVMARIVgyA557Hu65ZzS6mlZKkEVERETGmeMzM7lrxgyuKizknqYm3rdrF8+3tSWtb2a835vNuTaDx1w1v/ffOpgkByzAOd4ZFFPEo/5T1Ll6ppWei+FRtf8h3PIrIdwAl18IhZnwuc9G7kmewJQgi4iIiIxDmZ7HF0tL+VVFBRlmXFVVxY11dbQmuZpsZlzszeHtNp1HXBV3+dsOJskhy+CdgbPIJMST4WfwArlMKTmb5rbNNBYEYfZJUFAHZxwPe6rghz8cza6OOiXIIiIiIuPYkqws7p4xg8sKCvhNYyMf2LWLjUmu8JoZl3rzOMum8aDbw//4Ow4myTmWzaneChpo5EV/IyUFy8kKTad63wOEl34IAh68b2XksW83fAeqqkazm6NKCbKIiIjIOJfteXy1rIw7pk+nG7iyqorHkzwOzsz4iDeft9lU7nW7+IvbebBsljeDI2w+G9yr1NNARfl76A63sLfjJTjh/RCoggtPh44O+PJ1o9S70acEWURERGSCWJGdze8qKlgQCvHZmhruamxMWM8z43JvASttCn/xd3KvfyhJPslbShaZPBV+jszQNEoKllPfuIbW+UuhYBqcVAELS+E3v4W1a0era6NKCbKIiIjIBFIWDHLH9OmszM7m+ro6/rO+/uBtFLE8M670juBkK+cefwcP+LsByLRMTvFWUM9+NrhXmVJyNsFAPlX7H8RVfgy8RvjouyHkwaf/bkI+9k0JsoiIiMgEk+t5/GTaNN6fl8d/NzTwzbo6upIkyVd5C1luZdzlb+Nhfw8Ac71ZzLPZvOhvpNHamV56Lu2d1ezL64SZy2B2GJbNgTVr4a67Rrt7I04JsoiIiMgElGHGd8rL+VRREfc0NfG56uqET7gImHG1t5ClVsrv/bdY5e8F4GRvORlk8HT4WXJzjiIvZyG1+x+j88SLIGhw6dlQlAX/8Dno4xFz45ESZBEREZEJysz4bEkJ3yor4+m2Nj5eVcW+BD9RHTSPv/OO5CgK+LW/hb2ujWzL4mSvklr28RpvMr30PByO6vYX4PgLoOAAnH0CVNXATT9IQ+9GjhJkERERkQnugwUF/HjqVDZ3dvLR3bvZ3tV1WJ2geVwdOBIP47bwJsLOMd/mMNtmstZ/ibagx5TiM2lqfYPGuUdD/hR413EwowBuvBH27Bn9jo0QJcgiIiIik8DZubncPn06Tb7PZbt3syHBs5JLLJPLvQVsoYn73C7MjFO95QTweCb8PCUFJ5EZmkJ1w98IL/sIZLfC+06Hzk740hfT0KuRoQRZREREZJJYnJXFbyoqyPE8rqqq4onW1sPqnOSVc7KV8xd/B1tdEzmWwwpvGdXs5Q22UlH2HrrCjdTmNkLFYliSB0eVw+9+D+vXp6FXw08JsoiIiMgkMjcU4tcVFczLyOCz1dXcneBZyR/15lNEJreG36TdhVlo86mwaazx1+NnFlOcv5R9jatpX3Q+ZHtwyZmRDX//69HtzAhRgiwiIiIyyZQHg/yiooKTs7P5Zl0d/7V/f69nJedYkKsDC6mlnT/4b2FmnOadBMAz/mqmFJ9DIJDDns41uGPOjzz2LTMA69ekq0vDSgmyiIiIyCSU63n8dNo0LsjL46f793N9XR3dMUnyUVbIuTaDJ10N6/195Fkey70T2e2q2OpVMa3knbR17Gb/nLlQVA4FWbBte/o6NIyUIIuIiIhMUhlm3FhezjVFRfypqYnP1dT0elby+7zZzCKXX/ibOeA6OdoWMo0prPbXkpG7gNzs+dQ0PkX3kouhJBuq9qaxN8NHCbKIiIjIJGZmfL6khP9XVsZTra18Ze/eg7dbBM3jmsCRdOBzh78ZgJWBkwjj86x7gWml78bRTVXOPijPheZ2aG5OZ3eGRTDdAUxWT1V383hV95DbMcAs7j1umdjlg2WWtH78MvHt9fqcpJ2U+4+vb4fKk8UQsz5pOzFtxcdDovaSxhPtWx/75LD6yZat935j20u57eFjnbIfCdo4bF89G4uIiERdWlBAh+/z/fp6fn7gAFcXFQFQYTlc4s3lt/5WHnPVnO1NZ6m3mDX+OvYE51BWdDq1+x/Hn1WKx3bY8BKcujK9nRkiJchp8nxtN//xSseQ2nA9Lxd5Fxmofk8YIMVkw/qR4Pcx2Tn4OfXEjQRtx7eZrB2StE2SffVqJ2XsfYwT9J4EpupXH8tY8rHudzw9Y93XOCTdf/x66+dYxsTQx1j3HrPE+09Y57CxTDBRTnIck5YnaydRG30dhwTjnXQcUh7LxBPuxLGn6Hf8mKesnyzmnohlIriisJCNHR38uL6e4zMzOTk7G4CzbRobrJ4/+ts42go5zo7iLbbznP8C7ys8lwPNG+mcV0IWwNrnlCDL4HzphCy+dELWsLfrnOuVNMcn0L3KktVx4HAJ6x7WRoJleq1z/dvWxZTHf45rqz8x9O5L7/Je49RHPIe1kXAcko+3i1bqT18PG++EY5moHwnaidv/4fvqXx8PO3ZJxzJF+32MdcLYktVJ1M5A+5dgrImp7/cs+4fGgD7jSvHfm+vdft/j4FLH359x6sdx68//H0QGqs8JV7LlwxLxBBOJuPr9S+L7mHDFt5+k7Z5JRLJ+xE4W++prT5lH3xPuRMv5GUZ5ljE122NKtjE1y6M825iabRSHbEQmJ2bGP5eX82ZnJ1+qqeGPM2dSEQxiZlzlLeSb4fXcFn6TrwUWcXrgZP4cvp/n3YssLzuPtmNWRRLkDS8Oe1yjTQnyBHPwKtWQ/5vRFQGRyS7VBDDVBJqDy4kn3P2ZcPVe5xLWSTyhSrXeJYwnWfu91vcqd4n3k2x8ksbu+m4jurKviVLS8Y5pN9X4JWwrVQzJ+ppkvJOPZZJjFttuqosbKcb7UFyJJ9wJ2+jzPHAJ+srBCffBMh/86OAl3lfiY9bU5djb5ug49J25gzI8KM8ypmRFk+fsQ8tTsj1OmxqgLGtwXzXL9Tx+PHUqH9q9my/U1PCrigpCZhRaiI95R/BT/3X+7O/kA4E5LPGOZ52/gXmhWYTmlFOcHYTXXxvUfscSJcgiIpKQJtwi6eeco7ELatp8atsdNW0+e9scNW2Ove2R5eo2nw31jtp2RziayednwOePy+TTx2SSExz4f4PzQiFunDKFz9fUcGNdHdeXlwOw1CvldDeV+90ujndFLLLj2MZOnnUvcFrxNCjMgh07h3MI0kIJsoiIiMgYZWYUhqAwFODIwr7rhn1HfYdje7PPj17p4IYXO7j9zU6+sSSLD87LIOANLFF+e24uVxcVcVtDA4syM7mooACAS715vBE+wG3hTfxzYAmnB07mL+EH2J+fy6ySHNi6f7DdHTP0mDcRERGRCSDgGeXZHpXlQX59Zi5/fWcu07I9rl3Vxpn3NQ/q6VmfKy7mlOxsvrNvH690RB4ukGUBrg4cSQMd/MbfSqmVsMiOpSoX/Kl50N4F+/YNd/dGlRJkERERkQlo5dQgD787l9tOy6ax0/H+R1q4+G8tvLo/3O82AmZ8f8oUSjyPf6ipoSEc2XaB5fMebxbPuVqe92tZ5B1PW0aI8MyiyIYvrh2BHo0eJcgiIiIiE5RnxgfmhVh9YT7fWZbFC3XdnH5vM597tpWq1gTf/kugJBDgR1OnUtvdzXV79xKOfuPxPTaL+eTxa38LTYTxMoroXFgW2Wjt6pHq0qgYswmymXlmtsrMnJnNjFl/hZltMbNWM3vezJbFbVdpZquj5VvM7LLRj15ERERk7MgMGJ85NpN178vn744K8futXVT+bxPffamd5i6XcvsTsrL4elkZq9ra+On+yD3GATOuDhxJGMfP/U1kBctoO3Z65Hu54/xRb2M2QQa+ALTGrjCz04D/Bj4NFAN3A/eZWUG0vBC4P7q+GPgU8DMzO2UU4xYREREZk0oyPW5cns3zF+TxzpkZfG9DB5V/buKXmzrp9vtOlC/Oz+ei/HxubmjgsZYWAKZaNh/25vO6O0CDF+LArDLIyYBNb4xGd0bMmEyQzexI4FrgS3FFnwTucc495JzrAH4AdADvj5ZfBLQB33fOdTjnHgb+B7hmdCIXERERGfvm5Qe44205PHhuLnPyPP7huTZO/2szD+3uwrnEibKZ8Y3SUo4LhfhqbS3bu7oAOM2mcKKV8BJtNBbmQ1E27Nwzmt0ZdmMuQTYzD7gduA5oiCteDBy869tFjuD66Pqe8nWu95FdF1MuIiIiIlEryoM88K5cfvm2HDp9+NCjrVy3uj1pkpzpefz71Kl4wD/U1NDq+5gZl3sLaCKD/YX5UJqNqzvQ+xdqxpkxlyADnweqnXP3JCjLBw7ErWsACvpZ3ouZXWNmL5jZC7W1tYMOWERERGS8MjMumJPBs+/N4++PCfHzNzv58prkSfKMjAx+MGUKmzo7ub6uDucchRYii1xasrLwpxVgXWGoqhrlngyfMZUgm9kRwBeBzySp0gTEPya7CGjsZ3kvzrlbnHOVzrnK8ugvxIiIiIhMRqGA8Z1lWXzm2BC3vdHJV/pIklfm5PDZ4mLubW7mt42RNKvC8mkO5tA9pyhSae3zoxT58BtTCTJwGlAOvGxmdURujwDYYGbXAi8BS3sqm5kBS6Lrib6fGNfmiTHlIiIiIpKEmfHtpVn8/TEhbnmjk6++kDxJ/mRREWfl5PD9fftY297OTHKozcih/cipkQrj+FFvYy1B/iOwgEjSuwQ4L7r+ncCdwK3ARWZ2jpmFiFxtziLyRTyi7zlmdp2ZhczsHCJf3Ltl1HogIiIiMo6ZRa4kX3tMiJtf7+RrSZJkz4wby8upCAb5x5oaClwWDYFsmo+rAM/g5fF7fTKY7gBiOedaiXm0m5n1xFftnGsGno5eSb4VmA5sBM5zzjVGt28ws/OAnwLfBqqATznnnh3FboiIiIiMa2bGDcuy8B387PVOPIMblmUR+cf7QwoCAX48bRof2b2bX+9rZ0ZpBg0VpVTkhWDL5jRFP3RjKkGO55zbRuRx07Hr7iRyNTnZNmuAFSMbmYiIiMjEZmbcWBlJkv/rtUiS/O2lhyfJR4ZCXF9Wxldq91JSEqKhsACKsmB3TZoiH7qxdouFiIiIiIwRZsZ3l2dx9VEh/vPVTr61LvHtFufl5ZFjHp1+Lg35ubiyXFx9M/j9+znrsUYJsoiIiIgkZWZ8f3kWVx0Z4ievdvLP6zsOS5I9M47KzKS1K5PmUA7+9Hws7MO2t9IU9dCM6VssRERERCT9zIwfrMjCAT9+pQMDvnliZq/bLY4JhXi63aM2L5eueaUEANY8C/MXpCnqwdMVZBERERFJyTPjphVZXLkwxI9e6eCGF3tfST4mM5MDXUHqAzm0HRN91Nu6NWmKdmh0BVlERERE+sUz499OysI5xw9f7sAz+NriyJXko0Mh2jqDtFqIhmNnUhwweGV8PupNCbKIiIiI9Jtnxg9PzsYHbtrYQcDgK4uzWBAKEQ4HaHeZNEwvhfxM2Lot3eEOihJkERERERkQz4wfnZyN7+B7GyL3JP/T4iyOCGXS0p1DfUE+FGfjqut6P693nFCCLCIiIiID5pnxH6dk44DvbujADI6eFuKNzkz25+TjpuTCtgbo7obg+Eo59SU9ERERERkUz4z/ODmbS+dn8K8vdUBbkPqOAPUZ2YRnFGHOwRuvpTvMAVOCLCIiIiKDFvCMn5ySTU4A9tYHaO0KUhvMp3N+aaTCmmfTG+AgKEEWERERkSEJesaS0gBb9xqtXRk0e5k0H1sRKVz3fHqDGwQlyCIiIiIyZMvKArxc55juhWjzs9h33EzI8OC1V9Id2oApQRYRERGRIassC9Lpw3QyqO/MYf+UEijIxG3bme7QBkwJsoiIiIgMWWV5AIBgewa1nSFqCwpwJdmwtz7NkQ2cEmQRERERGbKKHI+KHKOhIfJFvbpQHm5KPjS2Q3t7usMbECXIIiIiIjIslpUF2LLXaO3MYF9GPt2zCiM/FLLxxTRHNjBKkEVERERkWCwrC7Jtv5EdDtHgZdOxcEqkYPX4etSbEmQRERERGRaVZZH7kMtdJs3hLBoWz44UrB9fj3pTgiwiIiIiw2JxSQDPINAepLojh9qjZ0JmAN58I92hDYgSZBEREREZFnkZxjFFHvsbAjR0ZVJbUgyFWbgdu9Md2oAoQRYRERGRYVNZFmRzTeQX9WpyCyOPeqs7kO6wBkQJsoiIiIgMm8qyAA1NHuGuDOpCefjTCrCWTmhsTHdo/aYEWURERESGzbKyAGAUdmWxz8uja05xpGDd6rTGNRBKkEVERERk2BxZ4JGfAdaeQX13Dq1HT4sUrF6V3sAGQAmyiIiIiAybgGecWBpg//4ANZ057Fs6L1Kw4YX0BjYASpBFREREZFgtKwuyvdbjQFcW1QtmQHYQt2lzusPqNyXIIiIiIjKsKssCtLcEaO8Msqe4FFeUBbur0x1WvylBFhEREZFhtawsAM4IdWRRk5UPpbmwryndYfVbMN0BTFYPP+Hz0GP+kNsxA8MwAwwsui522exQ3fh3i9brz7Y9ywm37VmO+5wormiVpO0nbTtuP332Ofoh0bYH+zGAPieKqacfqfo76LEcQH8TxXJ4XHbYtsn62VdcifY9mPHrWT5YPojxs7jtE/W/ZzeHtxdzEEREZFhNzfaYlWv4LZnUBAoJVxQQfLkG6mqhrDzd4aWkBDlNDPCG+Afa4XAOnHM4wDnAcXDZuWi9mOXYz45o/QTL8W0R955s2+ji4dv31Osjtn7HReK2Y/t4aG8ifev/BKPvxJ24z/2Z+CSarCSdFPQjpp5+DHbyGVns34Q7VduxsQxoApXm45Bw0pdqwj2MxyH2eA+l7YPNDWYsE8TSnwl3sslrX30e6WMaO16pzsM+204RtybciVWWBVndEGBvdw6dc0sJsinyJIvzLkx3aCkpQU6Tt5/h8fYz0h3FxOSc6zPpj67u14QiadKepO2RmqjExjzw9lzCOIn7POBJUz/iStXHZMeiv8elZ/tkMSfqY2TZ9aozkD4niiv+fTAT0EGfi3F9Hdzk08XUSz3hTtXfRMcl2aQ7dtuD6wZ5HAb731fsck8omnDLYAwmue7PBDJp2/1sN1U78fWDQcjIgIygEQpBRs/nnlf0cygDMjKMjCCEQpCVBaefZGRlRXa2rCzAvZsC1Hbm0nTcDHIAnn9KCbJIOpgZmszH0mCIDFaqCXf8Mgx+ohJdHHrbSdrv7wSj17rD2ko84U44ARpAHIliGuj4JRuvoU6oYsci5YR7wNumPqajPpZAdzd0dUFnp6O5JbLc1bOuq/fnmLMOgIYGj0suPJQgd70YoKEzi6rKBUw1cC+/OC7+KilBFhERSUIT7ngaDDnEOXcwUe7qgm/cGOaZ1T6XXBh5BsTikgAB38N1ZLJrzgyW5IZgy9Y0R90/Y+opFmb2PTN7xcwazWyPmd1qZiVxda4wsy1m1mpmz5vZsrjySjNbHS3fYmaXjW4vRERERCY+MyOUYeTmGEWFxmkne7yxGWr3Ra4qZweN44o92g9ksbMg+qi36n1pjrp/xlSCDISBy4BSYDEwE7ijp9DMTgP+G/g0UAzcDdxnZgXR8kLg/uj6YuBTwM/M7JRR7IOIiIjIpLNyReRfGJ5dc+i2i2VlQfbtD7InUIgrz4X6lt73foxRYypBds59zTm33jnX5ZyrBf4TODOmyieBe5xzDznnOoAfAB3A+6PlFwFtwPedcx3OuYeB/wGuGbVOiIiIiExCs2YYs2bAqtWHEuDKsgDNjUGq/AK6ZxRiXWHYtSONUfbPmEqQEzgH2BDzeTGwtueDc84B66Pre8rXRdf3WBdT3ouZXWNmL5jZC7W1tcMauIiIiMhks3KFsfFVx4HGSCpWWRagqyXA3s5c2hdMiVR69sk0Rtg/YzZBNrMPELli/PmY1fnAgbiqDUBBP8t7cc7d4pyrdM5VlpeP/YdWi4iIiIxlK1d4+A6eWxtJkBcUeOTj0dieQ/2iWZFKa55OY4T9MyYTZDO7BLgVuMA5ty6mqAkojKteBDT2s1xERERERsiCeTCl/NBtFp4Zy0ozaG3MZMfyI8HAvfpymqNMbcwlyGb2ceBm4L3Oucfiil8ClsbUNWBJdH1P+Ylx25wYUy4iIiIiI8TMWLnCWLfB0dIaSZKXlQWoqc1ka8VMyM+EbboHeUDM7HPATcC7nHPPJKhyK3CRmZ1jZiHgi0AWkS/iEX3PMbPrzCxkZucQ+eLeLaMQvoiIiMikd+oKj+5uWLP+UILc2ZzB9sxSXHE21OxPc4SpjakEGfgxkfuFHzOz5p5XT6Fz7mngWiKJ8gHgg8B5zrnGaHkDcB5wSbT8VuBTzrlnR7UXIiIiIpPUMUdCceGh2yyWlQXoagmywy/Gn5oHDa3g+2mOsm9j6pf0nHMpf6LHOXcncGcf5WuAFcMZl4iIiIj0T8AzTlluPPqUo6PTUZblMTMQZG9HLl0zi8lavQs2vQ5HHZvuUJMaa1eQRURERGScO3WF0d4B6zdEriIvLwuyrzGPlqOmRiqsiv+a2diiBFlEREREhtWi44y8XHhm9aHnIdfty6LmxPkA+C+sSmd4KSlBFhEREZFhlRE0TlpqPPeCo7vbsawsQEtDJm8sOxYCBq+/mu4Q+6QEWURERESG3aknGc0tsPFVxwklAWjN4LWSmbiCTNhZle7w+qQEWURERESG3dJFRmZm5DaLzIBxdGaIrZRCSQ5WO7Z/w00JsoiIiIgMu6xMo3KJsWqNw/cdlaVBdrUUEZ5WAAfaoKsr3SEmpQRZREREREbEyhXG/gZ4fVPkechV9Xl0zCnBHLBhfbrDS0oJsoiIiIiMiBVLjWAAnlntU1kWoGFfDg3HzgDAPfNwmqNLTgmyiIiIiIyI3BxjyQnGqtWOuXlGqC3E9qVHAuCvfy7N0SWnBFlERERERszKFUb1Xnhrh7EoO5O1ixdBhgebNqc7tKSUIIuIiIjIiDm50vAMVq32Oakkg40ZM3GFWdju2nSHlpQSZBEREREZMUWFxnHHwKrVkR8M2dxUhivNwfY1pzu0pJQgi4iIiMiIWrnCY9tOmNYVYE9dAd0VhdDcAW2t6Q4tISXIIiIiIjKiTl1uALzyImS35NI6ryzyqLfnn0xrXMkoQRYRERGREVVeZhy5IHIf8pKsLPYeNweA8DMPpTmyxJQgi4iIiMiIW7nC443NcEIog5dPWgyA27guzVElpgRZREREREbcqSsit1kE9gR4fPZSXGYAtu5Ic1SJKUEWERERkRE3s8KYPRN2rTM2tMyEoiy8qvp0h5WQEmQRERERGRUrVxivvgbuQAl+eR5Wr6dYiIiIiMgkduoKD99BWV0uXTOKsNYu2D/2riIrQRYRERGRUbFgLkwth+DmEI0LpgHgP3VfeoNKQAmyiIiIiIwKM2PlCqN+Y4DtixYC0PXs39Ic1eGUIIuIiIjIqDl1hUegyXh88akAuNdeTnNEh1OCLCIiIiKj5pgjoaTIeDywApcdJLCtKt0hHUYJsoiIiIiMGs8zTlluVK2dCsXZBPY2pjukwyhBFhEREZFRtXKFkflqHuEp+dj+tnSHcxglyCIiIiIyqhYda0xrC9A+uwRr78btfCvdIfWiBFlERERERlUwaKxc4lG7cBYAbY/dleaIelOCLCIiIiKjbuUK4+WFiwHoXPtUmqPpTQmyiIiIiIy6pYuNPxeeB0DGG7rFQkREREQmucyQkdGwFJcXInNnXbrD6UUJsoiIiIikxcXHZ+GKswnUNqc7lF6C6Q5gJJhZAPgucCWQBTwE/J1zbsxMT568AR7/1jA0ZGDW+x0OX9frPcl2/X1P2X6q9yHuP119GOi4TtRxGO3zoT/tT9ZzcihtHKwrIpJGJy/16J5WQMaLuyEchkAg3SEBEzRBBr4CXAicBOwDbgd+Bbw7nUHFmn0anPa1ITbiwLne73D4ul7vSbbr7/twbD/Ybftqw/cH0MZg+pBimwH3a4jj2FcbIgM2nAl7kvbGyoRrxNtIsW2/+xfXzrgbh8Ecy772N8hYhnw+DUcb6T4WI3ROJq0zQLk5xr45Uyhds5OqDY8x/cS3D66hYTZRE+RrgG8757YCmNmXgc1mNtc5ty2tkUXNPTPyEhkJ6U7Sh23C1J/2RzqG0WhjtI9FirqD7s8A+5H2cUjyntbzeSjjPwLjKDJQ5oEFwAtE3s07tOwFEpef9d4jKGUtVff9rxLkkWJmhcBsYG3POufcFjNrBBYB22LqXkMkmWb27NmjG6jICDo4o093ICIy7o3KZKmP/YxWDCPRxqBjH6Y2Um7bz/4NKAYf/DC4cPTdj1kO9y7vKduacy7HLnyErqzcpOfhaJtwCTJQEH0/ELe+IaYMAOfcLcAtAJWVlQ4RERHpRRNuGXlXwr9cyUnpDiPGRHyKRVP0vTBufRHQOLqhiIiIiMh4M+ESZOdcA7ADWNqzzszmE7l6vCFNYYmIiIjIODHhEuSoW4B/MrN5ZlYAfA94cKx8QU9ERERExq6JeA8yRJ6BXAysATKBh4HL0hqRiIiIiIwLEzJBds6FgS9FXyIiIiIi/TZRb7EQERERERkUJcgiIiIiIjGUIIuIiIiIxFCCLCIiIiISQwmyiIiIiEgMJcgiIiIiIjGUIIuIiIiIxFCCLCIiIiISQwmyiIiIiEgMJcgiIiIiIjHMOZfuGMYEM6sFto/ybsuAulHe52SnMU8Pjfvo05iPPo356NOYj76JNuZznHPl8SuVIKeRmb3gnKtMdxyTicY8PTTuo09jPvo05qNPYz76JsuY6xYLEREREZEYSpBFRERERGIoQU6vW9IdwCSkMU8Pjfvo05iPPo356NOYj75JMea6B1lEREREJIauIIuIiIiIxFCCLCIiIiISQwmyiIiIiEgMJchpYGYBM/uBmdWaWZOZ3W1mZemOa6wys1+YWZeZNce8ro2rc4WZbTGzVjN73syWxZVXmtnqaPkWM7ssrnyKmd0TPR61ZvY9M/Niyif0MTOzS83sKTNrNLPuBOVpH99UMYw3fY25mV1pZn7cOf+7uDoa8wGKjsEr0THfY2a3mllJXB2d68Mo1ZjrXB8ZZvYvZvZWdNz3mtmfzGx2TLnO81Scc3qN8gv4OvAmMB8oBO4G7k93XGP1BfwCuK2P8tOAFuCdQCbwZaAGKIiWFwK1wD9Fy98BNAOnxLTxMHBPtO786PH5p8lyzIB3AR8GrgK6x9r4pophPL5SjPmVwOY+ttWYD27MbwROBDKAcuB+4M8619M65jrXR2bcjwYKo8s5wA+BVTrPBzCG6Q5gMr6I/KT1J2I+LwAcMDfdsY3FF6kT5F8Cv4r5bMAO4GPRzx+PfraYOr8C7oguz4uO/4KY8k8Ab022YwacyeHJWtrHN1UM4/mVZMyvpO+kQWM+PGN/PnAg5rPO9dEfc53rIz/mucBNwL7+9FdjHnnpFotRZmaFwGxgbc8659wWoBFYlK64xoEPmFm9mb0Z/WebvJiyxfQeTwesj67vKV8XXd9jXVz5gehxiC2fa2YFOmZjYnxTxTARzTKzajPbaWa/N7N5MWUa8+FxDrAh5rPO9ZEXP+agc31EmNlHzOwAkau/nweujxbpPO8HJcijryD6fiBufUNMmfT2EyL/XFQGvB84A7g1pjyfvsdzsOVE60z2YzYWxjfVPiaaJ4ETgApgOdAOPGxmudFyjfkQmdkHgE8SSRx66FwfQUnGXOf6CHHO/dY5VwhMJ5Icb4wW6TzvByXIo68p+l4Yt76IyOxK4jjn1jrnapxzvnPuFeALwMVmlhmt0kTf4znY8p6yyX7MxsL4ptrHhOKc2+qcezN6zlcTSSoqgJOjVTTmQ2BmlxCZZF/gnFsXU6RzfYQkG3Od6yMvOq63An+NfkFS53k/KEEeZc65BiL32SztWWdm84nMmuL/2UkS86PvFn1/id7jacCS6Pqe8hPj2jgxrrwwehxiy7c55w7omI2J8U0Vw0Tnoq/Yc15jPghm9nHgZuC9zrnH4op1ro+AFGMeT+f6yAgSuRe5Ap3n/ZPum6An44vItzvfIHKjewFwF/BAuuMaqy/gUqAourwQWAXcHVN+GpF7rM4BQsCX6P2N3CIi38i9Llp+Dom/kfun6PGYFz0+X5ksxwwIAFlEvlHcHV3OIvJHKu3jmyqG8fhKMebnAzOjyyVEkovtQJ7GfEhj/jlgH7A8SbnO9dEfc53rwz/mHvAZYEr080zgf4C3iCTKOs/7M47pDmAyvoj8YbwJqCPyzwz3AGXpjmusvoDHgXoij4R5i8jjagri6lwBbAXagNXAsrjy5dH1bdF6l8WVT4keh6bocfk+4E2WY0bkm+QuwWvuWBnfVDGMt1dfYw78ANgTPeeron+IjtSYD3nMHdBF5A/zwddA+qxxH94x17k+ImPuAfcBe6Pjuhv4Db2fOqHzPMXLokGKiIiIiAi6B1lEREREpBclyCIiIiIiMZQgi4iIiIjEUIIsIiIiIhJDCbKIiIiISAwlyCIiIiIiMZQgi4hMEmY218xczKtoGNq8Mqa9F4YhTBGRtFOCLCIy+XwJOIXIA/yH6t5oW/cNQ1siImNCMN0BiIjIqHvDOffccDTknKsFas2sFpg6HG2KiKSbriCLiIwzZlZkZrvM7M649X8xszfNLGeA7fXcenGpmd1hZo3R9i+Lln/ZzPaYWa2Zfc/M9LdDRCY0/U9ORGSccc41AJ8ALjez9wGY2ceB84ErnXOtg2z6e0AV8AHgKeCXZvZvwArgKuBHwJeBDw4hfBGRMU+3WIiIjEPOuQfN7BbgZjPbDvw7cJNzbtUQmn3UOfc1ADN7HrgYuAA42jkXBh4wswuB9wO/H1oPRETGLl1BFhEZv74ItADPAruAbw6xvb/1LDjnGoFa4IloctxjMzBjiPsRERnTlCCLiIxTzrlm4K9AJvBz51zHEJtsiPvcmWRd1hD3IyIypilBFhEZp8ysEvg0sB74hplNS3NIIiITghJkEZFxyMyygDuBB4HTgHrglrQGJSIyQShBFhEZn24ApgGfjD614mPA+WZ2ZVqjEhGZAJQgi4iMM2a2EvgC8BnnXBVA9OkVPwR+ZGYz0xmfiMh4Z865dMcgIiKjwMzmAm8BFwL3Oee6h6FNAwLAz4HjnHOVQ21TRCTddAVZRGTy+TPQZWZFw9DWx4Au4IphaEtEZEzQFWQRkUnCzELAophV6+OecTyYNkuBedGPLc6514bSnojIWKAEWUREREQkhm6xEBERERGJoQRZRERERCSGEmQRERERkRhKkEVEREREYihBFhERERGJ8f8BiRN6hsRvzjcAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 720x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from scipy.sparse import diags\n",
    "\n",
    "##### Set up the domain #####\n",
    "nnx          = 40 #number of grid cells\n",
    "deltax       = 8000.0 #distance between grid cells [m]\n",
    "xgrid        = np.linspace(0,(nnx-1)*deltax,nnx) #grid in x-direction\n",
    "nnt          = 4000 #number of time steps\n",
    "deltatYr     = 0.5*0.5 #time step size in years\n",
    "secpyr       = 365.25*24*3600 #number of seconds in one year\n",
    "deltatSec    = deltatYr*secpyr #conversion to seconds\n",
    "### Initialize the glacier ###\n",
    "hvar         = np.array([np.linspace(10,0,nnx)]).T #from 10 to 0 m thickness (this only enforces small slope for numerical stability)\n",
    "### Accumulation rate ###\n",
    "bdotYr       = 2.0*np.ones((nnx,1)) #constant accumulation rate over the domain [m yr-1]\n",
    "bdotSec      = bdotYr/secpyr #conversion to [m s-1]\n",
    "\n",
    "##### Parameters #####\n",
    "nGlen        = 3.0 #Glen exponent\n",
    "bigA         = 3.5e-24 #ice rheology parameter [s-1 Pa-3]\n",
    "rhoi         = 917.0 #ice density [kg m-3]\n",
    "grv          = 9.81 #gravitational acceleration [m s-2]\n",
    "\n",
    "##### Functions #####\n",
    "def diffParamOnStaggeredGrid(hongrid,deltax,nGlen,bigA,rhoi,grv):\n",
    "    '''\n",
    "    Computes the diffusion coefficient on a staggered grid from the\n",
    "    ice thickness values on the regular grid (hongrid)\n",
    "    '''\n",
    "    h_m   = hongrid[0:-1] #h at i indices\n",
    "    h_p   = hongrid[1:] #h at i+1 indices\n",
    "    dhdx  = (h_p-h_m)/deltax #evaluate gradient in h on the staggered grid (at i+1/2 indices)\n",
    "    hstag = (h_m+h_p)/2 #average h between grid cells to get value on staggered grid (at i+1/2 indices)\n",
    "    bigDstag = (2*bigA*(rhoi*grv)**nGlen)/(nGlen+2) * (abs(dhdx))**(nGlen-1) * hstag**(nGlen+2)\n",
    "    return(bigDstag)\n",
    "def getAlpha(bigDstag,deltax,deltat):\n",
    "    alpha0 = (deltat/deltax**2)*bigDstag\n",
    "    alpha  = np.append(0,alpha0) #first entry of alpha is a dummy value\n",
    "    return(alpha)\n",
    "def getGamma(bigDstag,deltax,deltat):\n",
    "    gamma0 = (deltat/deltax**2)*bigDstag\n",
    "    gamma  = np.append(gamma0,0)  #last entry of gamma is a dummy value\n",
    "    return(gamma)\n",
    "def getBeta(bigDstag,deltax,deltat):\n",
    "    bigD_m = bigDstag[0:-1] #D at i-1/2 indices\n",
    "    bigD_p = bigDstag[1:] #D at i+1/2 indices\n",
    "    beta   = 1+(deltat/deltax**2)*(bigD_m+bigD_p)\n",
    "    return(beta)\n",
    "\n",
    "##### Choose plotting time steps #####\n",
    "nplots    = 10 #number of time steps for plotting\n",
    "StepsPlot = np.linspace(0,nnt-1,nplots).astype(int)\n",
    "listPlot  = [] #list to save the h values at the time steps of StepsPlot\n",
    "\n",
    "##### Loop over time steps #####\n",
    "inith  = np.copy(hvar) #the initial condition\n",
    "for tt in range(nnt):\n",
    "    ### Get parameters ###\n",
    "    bigD  = diffParamOnStaggeredGrid(hvar,deltax,nGlen,bigA,rhoi,grv)\n",
    "    alpha = getAlpha(bigD,deltax,deltatSec)\n",
    "    gamma = getGamma(bigD,deltax,deltatSec)\n",
    "    beta0 = getBeta(bigD,deltax,deltatSec)\n",
    "    ### Right-hand side source vector ###\n",
    "    rhs   = hvar+(bdotSec*deltatSec)\n",
    "    ### Boundary conditions ###\n",
    "    beta      = np.append(np.append(1,beta0),1)\n",
    "    alpha[-1] = 0\n",
    "    gamma[0]  = 1\n",
    "    rhs[0]    = 0\n",
    "    rhs[-1]   = 0\n",
    "    ### Construct matrices ###\n",
    "    coefMat = diags(diagonals=[-1*alpha[1:],beta,-1*gamma[0:-1]],offsets=[-1,0,1]).toarray()\n",
    "    ### Compute new h variable ###\n",
    "    hvarnew = (np.linalg.inv(coefMat))@rhs #right-multiply by the inverse of coefMat\n",
    "    hvar    = np.maximum(hvarnew,0) #remove negative ice thickness values\n",
    "    ### Save for plotting ###\n",
    "    if(tt in StepsPlot):\n",
    "        listPlot.append(hvar)\n",
    "\n",
    "##### Plotting #####\n",
    "ls_cols = plt.cm.rainbow(np.linspace(0,1,len(StepsPlot)))\n",
    "fig = plt.figure(figsize=[10,6])\n",
    "ax  = plt.subplot(111)\n",
    "for ii in range(len(StepsPlot)):\n",
    "    ax.plot(xgrid,listPlot[ii],c=ls_cols[ii],label=f'{int(StepsPlot[ii]*deltatYr)} yr')\n",
    "ax.tick_params(axis='both',which='major',labelsize=13)\n",
    "ax.set_xlabel('x [m]',fontsize=15),ax.set_ylabel('h [m]',fontsize=15)\n",
    "ax.legend(loc='best',ncol=2,fontsize=13)\n",
    "fig.tight_layout()\n",
    "\n",
    "    \n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b1f7f79c",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.9.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}