{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"\n",
"# Removing `if` Statements from Expressions \n",
"\n",
"## Author: Patrick Nelson\n",
"\n",
"## This notebook introduces a strategy for avoiding `if` statements in numerical programming, utilizing the absolute value function to boost the efficiency of NRPy+'s automated code generation. This approach serves crucial applications such as improving execution speed, allowing automated management of complex expressions, and implementing limit constraints in GRMHD simulations.\n",
"\n",
"### NRPy+ Source Code for this module: \n",
"* [Min_Max_and_Piecewise_Expressions.py](../edit/Min_Max_and_Piecewise_Expressions.py) contains functions that can be used to compute the minimum or maximum of two values and to implement piecewise-defined expressions.\n",
"\n",
"## Introduction:\n",
"\n",
"Conditional statements are a critical tool in programming, allowing us to control the flow through a program to avoid pitfalls, code piecewise-defined functions, and so forth. However, there are times when it is useful to work around them. It takes a processor time to evaluate whether or not to execute the code block, so for some expressions, performance can be improved by rewriting the expression to use an absolute value function in a manner upon which we will expand in this tutorial. Even more relevant to NRPy+ are piecewise-defined functions. These inherently involve `if` statements, but NRPy+'s automatic code generation cannot handle these by itself, requiring hand-coding to be done. However, if it is possible to rewrite the expression in terms of absolute values, then NRPy+ can handle the entire thing itself. \n",
"\n",
"The absolute value is a function that simply returns the magnitude of its argument, a positive value. That is, \n",
"\\begin{align}\n",
"|x|&= \\left \\{ \\begin{array}{lll}x & \\mbox{if} & x \\geq 0 \\\\\n",
"-x & \\mbox{if} & x \\leq 0 \\end{array}. \\right. \\\\\n",
"\\end{align}\n",
"\n",
"In C, this is implemented as `fabs()`, which merely has to take the first bit of a double-precision floating point number 0, and is thus quite fast. \n",
"\n",
"There are myriad uses for these tricks in practice. One example comes from GRMHD (and, by extension, the special cases of GRFFE and GRHD), in which it is necessary to limit the velocity of the plasma in order to keep the simulations stable. This is done by calculating the Lorentz factor $\\Gamma$ of the plasma and comparing it to some predefined maximum $\\Gamma_\\max$. Then, if\n",
"$$\n",
"R = 1-\\frac{1}{\\Gamma^2} > 1-\\frac{1}{\\Gamma_{\\max}^2} = R_\\max,\n",
"$$\n",
"we rescale the velocities by $\\sqrt{R_\\max/R}$. In NRPy+, we instead always rescale by\n",
"$$\n",
"\\sqrt{\\frac{\\min(R,R_\\max)}{R+\\epsilon}},\n",
"$$\n",
"which has the same effect while allowing the entire process to be handled by NRPy+'s automatic code generation. ($\\epsilon$ is some small number chosen to avoid division by zero without affecting the results otherwise.) See [here](Tutorial-GRHD_Equations-Cartesian.ipynb#convertvtou) for more information on this. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Table of Contents\n",
"$$\\label{toc}$$\n",
"\n",
"This notebook is organized as follows\n",
"\n",
"1. [Step 1](#min_max): Minimum and Maximum\n",
" 1. [Step 1.a](#confirm): Confirm that these work for real numbers\n",
"1. [Step 2](#piecewise): Piecewise-defined functions\n",
"1. [Step 3](#sympy): Rewrite functions to work with symbolic expressions\n",
"1. [Step 4](#validation): Validation against `Min_Max_and_Piecewise_Expressions` NRPy+ module\n",
"1. [Step 5](#latex_pdf_output): Output this notebook to $\\LaTeX$-formatted PDF file"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 1: Minimum and Maximum \\[Back to [top](#toc)\\]\n",
"$$\\label{min_max}$$\n",
"\n",
"Our first job will be to rewrite minimum and maximum functions without if statements. For example, the typical implementation of `min(a,b)` will be something like this:\n",
"```python\n",
"def min(a,b):\n",
" if a\n",
"\n",
"## Step 1.a: Confirm that these work for real numbers \\[Back to [top](#toc)\\]\n",
"$$\\label{confirm}$$\n",
"\n",
"For real numbers, these operate exactly as expected. In the case $a>b$,\n",
"\\begin{align}\n",
"\\min(a,b) &= \\tfrac{1}{2} \\left( a+b - (a-b) \\right) = b \\\\\n",
"\\max(a,b) &= \\tfrac{1}{2} \\left( a+b + (a-b) \\right) = a, \\\\\n",
"\\end{align}\n",
"and in the case $a\n",
"\n",
"# Step 2: Piecewise-defined functions \\[Back to [top](#toc)\\]\n",
"$$\\label{piecewise}$$\n",
"\n",
"Next, we'll define functions to represent branches of a piecewise-defined function. For example, consider the function \n",
"\\begin{align}\n",
"f(x) &= \\left \\{ \\begin{array}{lll} \\frac{1}{10}x^2+1 & \\mbox{if} & x \\leq 0 \\\\\n",
"\\exp(\\frac{x}{5}) & \\mbox{if} & x > 0 \\end{array} \\right. , \\\\\n",
"\\end{align}\n",
"which is continuous, but not differentiable at $x=0$. \n",
"\n",
"To solve this problem, let's add the two parts together, multiplying each part by a function that is either one or zero depending on $x$. To define $x \\leq 0$, this can be done by multiplying by the minimum of $x$ and $0$. We also will need to normalize this. To avoid putting a zero in the denominator, however, we will add some small $\\epsilon$ to the denominator, i.e.,\n",
"$$\n",
"\\frac{\\min(x,0)}{x-\\epsilon}\n",
"$$\n",
"This $\\epsilon$ corresponds `TINYDOUBLE` in NRPy+; so, we will define the variable here with its default value, `1e-100`. Additionally, to get the correct behavior on the boundary, we shift the boundary by $\\epsilon$, giving us\n",
"$$\n",
"\\frac{\\min(x-\\epsilon,0)}{x-\\epsilon}\n",
"$$\n",
"\n",
"The corresponding expression for $x > 0$ can be written as \n",
"$$\n",
"\\frac{\\max(x,0)}{x+\\epsilon},\n",
"$$\n",
"using a positive small number to once again avoid division by zero. \n",
"When using these for numerical relativity codes, it is important to consider the relationship between $\\epsilon$, or `TINYDOUBLE`, and the gridpoints in the simulation. As long as $\\epsilon$ is positive and large enough to avoid catastrophic cancellation, these functional forms avoid division by zero, as proven [below](#proof).\n",
"\n",
"So, we'll code NumPy versions of these expressions below. Naturally, there are many circumstances in which one will want the boundary between two pieces of a function to be something other than 0; if we let that boundary be $x^*$, this can easily be done by passing $x-x^*$ to the maximum/minimum functions. For the sake of code readability, we will write the functions to pass $x$ and $x^*$ as separate arguments. Additionally, we code separate functions for $\\leq$ and $<$, and likewise for $\\geq$ and $>$. The \"or equal to\" versions add a small offset to the boundary to give the proper behavior on the desired boundary. "
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:15:07.589772Z",
"iopub.status.busy": "2021-03-07T17:15:07.583863Z",
"iopub.status.idle": "2021-03-07T17:15:07.781320Z",
"shell.execute_reply": "2021-03-07T17:15:07.781822Z"
}
},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEGCAYAAABiq/5QAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAAsTAAALEwEAmpwYAAA5PklEQVR4nO3deVxU1f/H8deZQRaXEv2akriUFq6ISyqmBiLhli1qZpqUu0luWZqa+5ZpaZa5pJmpaWKLX5dSUdRkUjGX1DL9llvhkpFLiixzfn/MwE8NFGRmLjCf5+MxD2Dmztw3d4bPHM49c47SWiOEEMJ9mIwOIIQQwrWk8AshhJuRwi+EEG5GCr8QQrgZKfxCCOFmPIwOkB3/+c9/dMWKFY2OIYQQ+cqePXv+1FqXuvX6fFH4K1asSHx8vNExhBAiX1FKncjseunqEUIINyOFXwgh3IwUfiGEcDP5oo9f5F0pKSmcPn2apKQko6MIA3h7e+Pv70+hQoWMjiJyQAq/yJXTp09TrFgxKlasiFLK6DjChbTWXLhwgdOnT/PAAw8YHUfkgHT1iFxJSkqiZMmSUvTdkFKKkiVLyn97+VCBLvwWC0yebPsqnEeKvvuS5955nFm/CmxXj8UCISGppKSY8PY2ERMDwcFGpxJCiDuzWCAsDJKTwdMTh9evAtvij42F5GSF1iaSkzWxsUYnEs5iNpsJCgqiRo0adOjQgatXrxIfH0///v1dmmPOnDksXrw414/z2muvUb16dV577TUHpLLZt28f69aty/h59erVTJkyxWGPLxzLVr8gLc321dH1q8C2+ENCwMsLrl9PwWSCkBAZdVBQ+fj4sG/fPgA6d+7MnDlzGDx4MPXq1XNpjj59+jjkcebNm8dff/2F2Wx2yOOBrfDHx8fTqlUrANq2bUvbtm0d9vjCsUJCbC399BZ/SIhjH7/AtviDg2HzZhP+/vMpXbozDRpYjY4kXKBJkyYcO3aM2NhY2rRpA8A///xDt27dqF+/PrVr1+brr78GIC0tjSFDhlCjRg0CAwOZNWsWAHv27OGxxx6jbt26REREkJCQwLlz56hbty4A+/fvRynFyZMnAahUqRJXr15lzJgxTJs2DYD33nuPatWqERgYyHPPPXfbHDdq27YtV65coW7duqxYsYIXX3yR6OjojNuLFi0KQGxsLCEhIbRv354qVarQuXNn0lfT2717N40aNaJWrVrUr1+fixcvMmrUKFasWEFQUBArVqxg0aJFREVFAXD8+HGaNWtGYGAgYWFhGb/Xiy++SP/+/WnUqBEPPvjgTTmEcwUH27p3xo93fDcPFOAWP0CjRoqpU315/vmVrF8fSevWrY2OVKANHDgwo+XtKEFBQcyYMSNb26amprJ+/XpatGhx0/UTJ06kWbNmLFy4kL///pv69evTvHlzFi9ezPHjx9m3bx8eHh789ddfpKSk8Morr/D1119TqlQpVqxYwYgRI1i4cCFJSUlcunSJ7du3U69ePbZv307jxo257777KFy48E37nDJlCr/99hteXl78/ffft81RpEiRjPutXr2aokWLZhzH9evXZ/n77t27l0OHDnH//ffz6KOPsmPHDurXr0/Hjh1ZsWIFjzzyCJcuXaJw4cKMGzeO+Ph43n//fQAWLVqU8TivvPIKkZGRREZGsnDhQvr3789XX30FQEJCAt999x0///wzbdu2pX379tl6LkTuBQc777xkgW3xp2vfvj3+/v68++67RkcRTnLt2jWCgoKoV68e5cuXp3v37jfdvmHDBqZMmUJQUBAhISEkJSVx8uRJNm3aRO/evfHwsLV/SpQowZEjRzh48CDh4eEEBQUxYcIETp8+DUCjRo3YsWMH27ZtY/jw4Wzbto3t27fTpEmTf2UKDAykc+fOLFmyJOPxs8pxt+rXr4+/vz8mk4mgoCCOHz/OkSNH8PPz45FHHgHgnnvuydh/ViwWC88//zwAL7zwAt99913GbU899RQmk4lq1apx9uzZu84q8pYC3eIHKFSoEFFRUQwbNowDBw4QGBhodKQCK7stc0e7sY8/M1prVq1aRUBAwB0fS2tN9erVsWQyhq5p06Zs376dEydO8OSTT/LWW2+hlMr0P8m1a9eybds2/vvf/zJx4kR+/PHHHOVI5+HhgdVq66a0Wq0kJydn3Obl5ZXxvdlsJjU1NduPm1037iO9K0nkfwW+xQ/Qs2dPChcubFhhEsaKiIhg1qxZGYVr7969AISHhzN37tyMgvnXX38REBDA+fPnMwp/SkoKhw4dAmznD5YsWcJDDz2EyWSiRIkSrFu3jsaNG9+0P6vVyqlTpwgNDeWtt97i4sWLXLlyJcsct1OxYkX27NkD2LqBUlJSbrt9QEAACQkJ7N69G4DLly+TmppKsWLFuHz5cqb3adSoEcuXLwdg6dKlmf4HIwoWtyj8JUqUIDIykqVLl3Lu3Dmj4wgXe/PNN0lJSSEwMJDq1avz5ptvAtCjRw/Kly9PYGAgtWrVYtmyZXh6ehIdHc3QoUOpVasWQUFBxMXFAbYirLWmadOmADRu3JjixYvj6+t70/7S0tLo0qULNWvWpHbt2vTv35/ixYtnmeN2evbsydatW6lVqxYWi+Wm8wGZ8fT0ZMWKFbzyyivUqlWL8PBwkpKSCA0N5fDhwxknd280a9YsPv74YwIDA/n000+ZOXNmto+tyJ9Ufvj3rV69ejq3C7EcOXKEKlWqMGbMGEaPHu2gZOKnn36iatWqRscQBpLXQN6llNqjtf7XuGa3aPGD7V/g1q1bM3v2bJlbRAjh1tym8INtuOG5c+cy+jOFEMIduVXhDwsLo2bNmrz77rsyQkEI4bbcqvArpRg4cCAHDhxgy5YtRscRQghDuFXhB3j++ecpVaqUDO0UQrgttyv83t7e9O3blzVr1nD06FGj4wghhMs5rfArpRYqpc4ppQ7ecF0JpdRGpdRR+1ff2z2Gs7z88ssUKlRIxiu7sd9++40GDRpQuXJlOnbseNMnYtNduHCB0NBQihYtmjGhGfz/J1jHjBlz08+3k5qaSuvWrfnPf/7DwYMHb7rttddeo0qVKgQGBvL0009nzO2THSNGjKBcuXIZk7cJkR3ObPEvAlrcct0wIEZr/RAQY//Z5UqXLs3zzz/Pxx9/TGJiohERhMGGDh3KoEGDOHbsGL6+vixYsOBf23h7ezN+/PiMGTfTbdiwgREjRnD16lU++uijbHUb9u3blypVqvDVV1/RsWPHjPl/wPYJ4oMHD3LgwAEefvhhJk+e/K/7jxkz5qaJ1dI98cQT7Nq1686/sBA3cFrh11pvA/665eongU/s338CPOWs/d/JwIEDM/5whWs5ckm5UaNG3VR4R4wYccf/5LTWbN68OWOmycjIyIzZKG9UpEgRGjdujLe3903XR0REEBERwcyZM7lw4QKDBg3ixIkTPPTQQ/z5559YrVaaNGnChg0bABg7diz33nsv06dPp3Hjxnz00Ud06tSJixcvAvD4449nTKTWsGHDm94U7qRhw4b4+flle3shwPWTtJXWWifYvz8DlM5qQ6VUL6AXQPny5R0epFatWoSGhjJr1iwGDRp0xxkMhWM4ekm5bt268cwzzzBw4ECsVivLly9n8+bNBAUFZbr9smXLuO+++yhevHjGc+7v78/vv/+e7X1u3LiR2NhY+vfvT8mSJZk5cyYDBgxg6NCh9O3bl/r161OtWjUef/xxgH99Ujw4OJjt27dn+tgLFy6kY8eO2c4ixN0wrNpprbVSKsvOUa31PGAe2KZscEaGQYMG0bZtW1atWiV/bC6S2ZJyuSn8FStWpGTJkuzdu5ezZ89Su3ZtKlSocNvZOv/888+73yHQvHlzwsPDGTNmDD169Mjo4+/RowcrV65kzpw5d7UuwcSJE/Hw8KBz584A/Pjjj7zwwgsAnDlzBk9Pz4z/bmJiYihZsmSufg/hvlxd+M8qpfy01glKKT/A0BnTWrduTeXKlZkxY4YUfhdxxpJyPXr0YNGiRZw5c4Zu3bpx+fLlLGeYXLZsGVWrVuXvv/8mNTUVDw8PTp8+TdmyZbO9P6UU8P8nd9N/vnr1akY3zZUrVyhWrFi2H3PRokWsWbOGmJiYjMerWbNmxhvImDFjqFixIi+++GK2H1OIrLi68K8GIoEp9q//XnvOhUwmEwMGDOCVV17h+++/p2HDhkbGcQvpS8rFxtqKviNWGHr66acZNWoUKSkpLFu2DLPZfMcWd2hoKNHR0Tz33HN88sknPPnkk7nOMXToUDp37kyFChXo2bMna9asydb9vvnmG6ZOncrWrVv/tZKXEE6htXbKBfgMSABSgNNAd6AkttE8R4FNQInsPFbdunW1s1y+fFkXL15cP/vss07bR0F2+PBhoyNorbXu3bu3Hjp0aLa3/9///qcfeeQRXalSJd2+fXudlJSktdb666+/1m+++WbGdhUqVNC+vr66SJEiumzZsvrQoUOZPl5sbKxu0KCBTk1N1Vpr/fTTT+uFCxdmK0ulSpW0v7+/rlWrlq5Vq5bu3bv3v7YZPXq0/vjjj/91/WuvvabLli2rlVK6bNmyevTo0dnapyPlldeA+DcgXmdSU91mWubbef3113nnnXf49ddfnXIiuSDLC1PyWq1W6tSpw8qVK3nooYcMzeKO8sJrQGTO7adlvp30D+ekL0Qt8o/Dhw9TuXJlwsLCpOgLkU1S+LENF23Xrh3z5s3jypUrRscROVCtWjV+/fVXpk+fbnQUIfINKfx2gwYN4uLFi5l+OlIIIQoSKfx2DRs2pEGDBsycOROr1Wp0HCGEcBop/DdIn7tl7dq1RkcRQginkcJ/g3bt2lGuXDneffddo6MIIYTTSOG/gYeHB1FRUWzZsoX9+/cbHUfkEwMHDmTbtm133G7KlCksXbqURYsWUapUKYKCgggKCsqYKPD8+fO0aHHrhLZCOJ4U/lv07NmTwoULywpdIlsuXLjA999/T9OmTe+47bfffpsxcVvHjh3Zt28f+/bto0ePHgCUKlUKPz8/duzY4dTMQkjhv4Wvry8vvfQSy5Yt48yZM0bHKZAspyxM3j4Zy6ncz8t8N9MyA+zZs4fHHnuMunXrEhERQUJCAhcvXiQgIIAjR44A0KlTJ+bPnw9A0aJFGTRoENWrVycsLIzz588DsGrVqoxW+u3uf+nSJZKTkylVqtRtcz311FMsXbo0ZwdBiJzK7OO8ee3izCkbMvPLL79oaKibNduo4+Jcuut8J6cf1487Gad9Jvho81iz9pngo+NO5u4A//bbb7p27dpaa63T0tL0gw8+qI8fP54x/cGtl0OHDunk5GQdHBysz507p7XWevny5fqll17SWmu9YcMG3bBhQ/3ZZ5/piIiIjP0AesmSJVprrceOHav79euntda6a9euevXq1RnbZXX/VatWZUwF8fHHH+syZcromjVr6nbt2umTJ09mbHf69Gldo0aNXB0TV5MpG7IWF6f1pEnasDpCFlM2GF7Us3NxdeGPi9PaZErSkKJ9fKxS/G8jp3/0k7ZN0uaxZs0YtHmsWU/aNinXGZo3b65/+OEHvX79et2uXbs7bv/jjz/qYsWKZbwZ1KhRQ4eHh2fc3rNnT12iRAl96tSpjOtMJpNOSUnRWtvm+alVq5bWWuvw8HBtsVhuevzM7t+zZ08dZ38h/fnnnxlzA82ZM0eHhoZmbJecnKxLlCiRwyNgLCn8mYuL09rHR2uz2fbViDqSVeGX1UcyERsL4Akorl+3EhurHDKLpICQiiF4mj1JTkvG0+xJSMWQXD9mTqdl1lpTvXp1LJksAWa1Wvnpp58oXLgwiYmJ+Pv7Z/o46VMn+/j4kJSUdMf779q1iw8//BDgpnn0e/Toweuvv57xc1JSEj4+Pjk8AiIvcvTaE44kffyZCAkBLy8FpGK1Xic4+LrRkQqM4HLBxHSNYXzoeGK6xhBcLvd/CU8//TTffPMNu3fvJiIigmLFimWcOL31Uq1aNQICAjh//nxG4U9JSeHQoUMAvPvuu1StWpVly5bx0ksvkZKSAtgKenR0NGB782jcuDEAVatW5dixYxlZMrv/oUOHqFKlCmazGYCEhISM7VevXn3TBGe//PILNWrUyPUxEcZLX3vCbHbc2hOOIi3+TKTPGb9w4Qk++qgLhw51ISSkn9GxCozgcsEOKfjpPD09CQ0NpXjx4hnF9U7bR0dH079/fy5evEhqaioDBw7Ew8ODjz76iF27dlGsWDGaNm3KhAkTGDt2LEWKFGHXrl1MmDCB++67jxUrVgC2xXzmzp1Ljx49OHLkSKb3L1as2E3DNN977z1Wr16Nh4cHJUqUuGmakC1bttC6dWuHHRthHGesPeEwmfX/5LWLq/v401mtVt24cWPt7++f0ScrbpYX+nfT0tJ0rVq19C+//OK0fRQpUiTL2x599FGdmJiY5e3NmzfXf/zxR7b206RJE/3XX3/lNJ6h8sJrQGSOLPr4pavnNpRSjBo1itOnT8vkbXlUXpiWefr06Zw8eTLL2zdu3Iifn98dH+f8+fMMHjwYX19fR8YT4l9kIZY70FrTqFEj/vjjD44ePYqnp6chOfIqWYRDyGsg75KFWO5Seqv/5MmTfPrpp0bHyZPyQ+NBOIc89/mTFP5saNGiBfXq1WPixIkZozyEjbe3NxcuXJAC4Ia01ly4cAFvb2+jo4gcklE92ZDe6m/bti3Lli0jMjLS6Eh5hr+/P6dPn86YwkC4F29v7yw/6yDyLunjzyatNXXq1OGff/7h8OHDeHjIe6YQIm+TPv5cSm/1Hz16NGMMtxBC5EfS4s8Bq9VKUFAQKSkpHDx4MFsfFhJCCKNIi98BTCYTb775Jj///HPGx/eFECK/kRZ/DlmtVmrWrInJZGL//v2YTPLeKYTIm6TF7yAmk4mRI0dy8OBBvvrqK6PjCCFEjkmL/y6kpaVRrVo1fHx82Lt3b8YUvUIIkZdIi9+BzGYzI0eOZP/+/fz3v/81Oo4QQuSIIYVfKTVIKXVIKXVQKfWZUirfffSvU6dOVKpUiXHjxsmnVoUQ+YrLC79SqizQH6inta4BmIHnXJ0jtzw8PBgxYgR79uxh/fr1RscRQohsM6qrxwPwUUp5AIWBPwzKkStdunShYsWK0uoXQuQrLi/8WuvfgWnASSABuKi13nDrdkqpXkqpeKVUfF6dB6ZQoUIMHz6cnTt3snHjRqPjCCFEtrh8VI9SyhdYBXQE/gZWAtFa6yVZ3Sevjeq5UXJyMpUrV6ZcuXJ89913MsJHCJFn5KVRPc2B37TW57XWKcAXQCMDcjiEp6cnw4YNIy4uji1bthgdRwgh7siIwn8SaKiUKqxszeMw4CcDcjhMt27duP/++xk3bpzRUYQQ4o6M6OPfCUQDPwA/2jPMc3UOR/L29mbo0KFs3bqVrVu3Gh1HCCFuSz656yDXrl3jgQceoEaNGmzatMnoOEIIkaf6+AskHx8fXn/9dWJiYtixY4fRcYQQIktS+B2od+/elCpVivHjxxsdRQghsiSF34GKFCnCkCFD+Pbbb9m5c6fRcYQQIlNS+B3s5Zdf5p57IoiM/BmLxeg0Qojcslhg8mRc/vdsOWVh4raJWE45fseyYriD/fhjUa5d+y9HjihCQ9PYssVMcLDRqYQQd8NigbAwSE4GT0+IicElf8+WUxZCF4VyPfU63oW82Ry5meByjtuxtPgdLDYWrFYPwIPr1zVbtuT9UVNCiMzFxtqKflqa7WtsrGv2G/NrDNdTr4MJUqwpxB537I6l8DtYSAh4eipMJiuQjMm0zehIQoi7ZPt7BrPZ9jUkxDX7PfP9GUgDEyY8zZ6EVHTsjmUcvxNYLLB5cxofffQCHh67OXToEJ6enkbHEkLcBYvF1tIPCXFNN8/x48epXr06dZ6sQ8u+LQmtGHrX3TxZjeOXwu9E69evp1WrVsyYMYMBAwYYHUcIkcdprWnTpg1bt27l8OHDlC9fPlePJx/gMkCLFi0IDw9n3LhxJCYmGh1HCJHHrVq1inXr1jFu3LhcF/3bkcLvREop3n77bRITE5k4caLRcYQQedjFixfp378/QUFB9O/f36n7ksLvZLVq1eKll15i1qxZ/Prrr0bHEULkUSNGjODMmTPMmzcPDw/njrSXwu8C48ePx8PDg2HDhhkdRQiRB+3atYvZs2cTFRXFI4884vT9SeF3gfvvv5/XX3+dlStXEhcXZ3QcIUQekpqaSq9evfDz82PChAku2acUfhcZMmQIfn5+vPrqq7IwuxAiw8yZM9m/fz+zZs3innvucck+pfC7SJEiRZgwYQLff/89K1euNDqOECIPOHHiBKNGjeKJJ57g6aefdtl+ZRy/C6WlpVGnTh0uX77MTz/9hJeXl9GRhBAG0VrzxBNPEBsb65Ax+5mRcfx5gNlsZvr06fz222+8//77RscRQhjoiy++YO3atU4fs58ZafEboFWrVlgsFo4dO0bJkiWNjiOEcLGLFy9StWpVSpcuze7du502fFNa/HnI22+/zaVLlxg3bpzRUYQQBhg5ciRnzpxh7ty5Th+znxkp/AaoXr06PXv2ZPbs2fzyyy9GxxFCuNDu3bv54IMP6NevH/Xr1zckg3T1GOTs2bNUrlyZ8PBwvvjiC6PjCCFcIDU1lUceeYRz585x+PBh7r33XqfuT7p68pjSpUszbNgwvvzyS7Ztkzn7hXAH7733Hvv27eO9995zetG/HWnxG+jq1asEBARQpkwZdu7cickk78NCFFQnTpygWrVqNGvWjNWrV6OUcvo+pcWfBxUuXJhJkyYRHx/PZ599ZnQcIYSTaK2JiooC4P3333dJ0b8dKfwG69y5M3Xq1GH48OFcu3bN6DhCCCf48ssvWbNmDePGjaNChQpGx5HCbzSTycT06dM5efIkM2fONDqOEMLBLl26xCuvvEKtWrXyzEp8hhR+pVRxpVS0UupnpdRPSikXrGSZd4WEhNC2bVsmTZrEuXPnjI4jhHAQyykLLSa24A/zHy6ZZz+7jGrxzwS+0VpXAWoBPxmUI8+YOnUq165dY+zYsUZHEUI4gOWUhdBFoVi8LXh08yDNL83oSBlcXviVUvcCTYEFAFrrZK31367OkdcEBATQp08f5s6dy08/uf37oBD53jdHvuF66nUwgTZpYo/HGh0pgxEt/geA88DHSqm9SqmPlFJFbt1IKdVLKRWvlIo/f/6861MaYNSoURQpUoSePRcweTJYLEYnEqJgsFhw6d+U1prti7dDGpgw4Wn2JKRiiGt2ng0uH8evlKoHfA88qrXeqZSaCVzSWr+Z1X0K6jj+zPTrt4TZs5/BZPLGy8tETAwEu/UZECFyx2KBsDBITgZPT1zyN/Xee+8xYMAAXp70Mv6N/QmpGEJwOdf/IWc1jt+IMw2ngdNa6532n6MBWYzWrkyZjoDCajWRnKyJjVVS+IXIhdhYW9FPS7N9jY11buHfuXMnQ4YM4YknnmDW0Fl58oOZLk+ktT4DnFJKBdivCgMOuzpHXtW8eSG8vBSQglIphIQYnUiI/C0kxNbSN5ttX535N3XhwgU6dOhA2bJl+eSTT/Jk0QfjRvW8AixVSh0AgoBJBuXIc4KDYcsWM/XqrcZqbYaX1w9GRxIiXwsOtnXvjB/v3G4eq9XKCy+8wNmzZ1m5ciW+vr7O2ZEDyFw9eVRiYiLVqlXDz8+PXbt25Znxv0KIzE2cOJGRI0cye/Zs+vbta3QcQObqyXd8fX15//332bt3L++8847RcYQQt7FlyxZGjRpFp06d6NOnj9Fx7uiOLX6l1CvAEq11omsi/Zs7tvjTPfPMM6xfv54DBw7w0EMPGR1HCHGLhIQEateuja+vL7t376Zo0aJGR8qQmxZ/aWC3UupzpVQLZfS0cm7m/fffx8vLi169epEfuuWEcCepqak899xzXL58mejo6DxV9G/njoVfaz0SeAjbJ21fBI4qpSYppSo5OZsA7r//ft5++21iY2NZsGCB0XGEEDcYOXIk27ZtY+7cuVSvXt3oONmWrT5+bWtqnrFfUgFfIFopNdWJ2YRdjx49CAkJYciQIfzxxx9GxxFCAGvWrOGtt96iV69edOnSxeg4OXLHwq+UGqCU2gNMBXYANbXWfYG6QDsn5xOAUor58+dz/fr1jMUchBDGOX78OF27dqV27dr5cjr17LT4SwDPaK0jtNYrtdYpAFprK9DGqelEhsqVKzN27Fi+/PJLVq1aZXQcIdzW9evX6dChA1arlZUrV+Lt7W10pBzLTh//aK31iSxuk2kkXWjw4MHUrl2bqKgoEhMNG2QlhFsbPHgw8fHxLFq0iEqV8uepThnHn494eHiwYMECzp8/z2uvvWZ0HCHczvLly5k9ezavvvoqTz31lNFx7poU/nymdu3aDBkyhAULFrB582aj4wjhNn7++Wd69OjBo48+yuTJk42OkysyZUM+dO3aNQIDA7Farfz4448ULlzY6EhCFGj//PMPDRo04Ny5c+zdu5eyZcsaHSlbZMqGAsTHx4f58+fz66+/MmbMGKPjCFGgxZ2Mo/EbjTl06RDLli3LN0X/dqTw51MhISH07NmT6dOns2fPHqPjCFEgWU5ZCPk4hH2++/Do7kGRgH8tFpgvSeHPx6ZOnUrp0qXp3r07KSkpRscRosD58JsPSbGm2NbNVXlr3dzckMKfjxUvXpwPPviA/fv3M23aNKPjCFGgbN26lRVTVqC0wqzMeW7d3NyQk7sFQPv27VmzZg0HDhzg4YcfNjqOEPnenj17CA0NpVy5ckxbMY19ifsMWzc3N7I6uSuFvwBISEigWrVqBAYGsmXLljy73JsQ+cHPP/9MkyZNKFKkCDt27MjXJ3NlVE8B5ufnx7Rp09i2bRvz5883Oo4Q+daJEycIDw/HbDazadOmfF30b0cKfwHRrVs3mjVrxuDBK3njjYtYLEYnEsKxLBaYPBmnvbbPnj1LeHg4V65cYcOGDVSuXNk5O8oDpKunAFm58jTPPlsC8MLHx0RMjHLawtJCuJLFAmFhkJwMnp6OXzT977//JjQ0lCNHjrBp0yYaNWrkuAc3kHT1uIFjx/xRyhswk5RkJTbW6ERCOEZsrK3op6XZvjrytX316lWeeOIJDh06xJdffllgiv7tSOEvQEJCwNtboVQaWl+nWDH5YJcoGEJCbC19s9n2NSTEMY+bnJxM+/bt2bFjB0uXLiUiIsIxD5zHSeEvQIKDISZGMWpUKuXLd2PChNYkJCQYHUuIXLO9tmH8eMd186SlpdG1a1fWr1/P3Llz6dChQ+4fNJ+QPv4C6uDBg9SvX5/69euzadMmPDw8jI4kRJ6htaZv377MnTuXqVOnFthpzqWP383UqFGDDz/8kK1btzJ69Gij4wiRpwwfPpy5c+cybNiwAlv0b0cKfwEWGRlJ9+7dmTRpEuvXrzc6jhB5wtSpU5kyZQq9e/dm0qRJRscxhBT+Am7WrFkEBgbSpUsXTp48aXQcIQw1f/58hg4dSseOHfnggw9QShkdyRCGFX6llFkptVcptcaoDO7Ax8eH6OhoUlJSePbZZ0lOTjY6khCG+Pzzz+nduzctW7Zk8eLFmM1moyMZxsgW/wBAFmt3gYceeoiFCxeyc+dOhg4danQcIVzu3ZXv0ml2J2q2qkl0dDSenp5GRzKUIYVfKeUPtAY+MmL/7qh9+/b079+fGTNmsGrVKqPjCOEyU5ZMYfD+wVgfs3I0+Cj7L+w3OpLhjGrxzwBeB6wG7d8tvf3229SvX59u3bpx7Ngxo+MI4XSzZs3ijXlvgBkwQXJacoFZTCU3XF74lVJtgHNa69t+rFQp1UspFa+Uij9//ryL0hVsnp6efP7555jNZjp06MC1a9eMjiSEU1itVoYMGUL//v1p4t8EH0+fAreYSm4Y0eJ/FGirlDoOLAeaKaWW3LqR1nqe1rqe1rpeqVKlXJ2xwKpQoQKffvop+/btY8CAAUbHEcLhkpKS6NSpE9OnTycqKootn24hpmsM40PHE9M1Jt8tpuIUWmvDLkAIsOZO29WtW1cLxxo2bJgG9OLFi42OIoTDXLhwQTdu3FgDetq0adpqtRodyVBAvM6kpso4fjc1fvx4HnvsMfr06cOhQ4eMjiNErv322280atSIXbt2sXz5cl599VW3Had/J4YWfq11rNa6jZEZ3JWHhwefffYZRYsWpUOHDly5csXoSELctfj4eIKDgzl37hybNm2iY8eORkfK06TF78b8/Pz47LPPOHLkCL17907vfhMiX1m3bh2PPfYY3t7e7NixgyZNmhgdKc+Twu/mmjVrxtixY1m2bBnz5s0zOo4QOTJv3jzatm1LlSpV+P7776latarRkfIFKfyC4cOHExERQVTUUqKifpf1eoVTOHLNXK01I0aMoHfv3jz++ONs3bqVMmXK5P6B3YTMxy8AWL/+b1q39kLrQnh7m9i82STr9QqHceSaucnJyXTv3p0lS5bQo0cPPvzwQ1lvIgsyH7+4rX37itvX6/UgKcnKunX/GB1JFCCOWjP34sWLtGzZkiVLljBhwgTmzZsnRf8uSOEXgG0NUy8vhclkBZJZtaq/jPQRDuOINXO/2vMVAd0D2Pq/rSxevJgRI0bIcM27JG+VAvj/NU1jY02YTPGMGPEJTz11grVr1+Ll5WV0PJHP/f/ry1b0c9rNM2/dPPpY+qCra7xqeVE5pLIzYroNKfwiQ3Bw+h9kU/z8FhIZGcnzzz/PihUr5N9pkWv///rKvtTUVCZNmsSYTWPQoRpMkKpTiT0eK1Mv5IJ09YhMde3alRkzZvDFF1/Qq1cvGeMvXO7EiROEhoYyevRowh8Ol4nWHEiacSJLAwYMIDExkbFjx+Lr68u0adOkT1W4xOeff06vXr2wWq18+umndOnSBcspC7HHYwmpGCKt/VySwi9ua/To0fz111+88847lCxZkuHDhxsdSRRgV65cYcCAASxcuJAGDRqwdOlSKlWqBEBwuWAp+A4ihV/cllKKGTNmkJiYyIgRI/D19aVv375GxxIF0J49e+jUqRPHjh1jxIgRjB49mkKFChkdq0CSwi/uyGQysXDhQi5evEi/fv0oXrw4nTp1MjqWKCCsVivTp09nxIgRlC5dmi1btvDYY48ZHatAk5O7IlsKFSrEihUraNq0KV27dmXdunVGRxIFQEJCAhEREbz++us88cQT7N+/X4q+C0jhF9nm4+PD6tWrCQwMpF27dmzfvt3oSCIf++9//0tgYCA7duxg3rx5REdHU6JECaNjuQUp/CJH7rnnHr755hsqVKhAmzZt2Lt3r9GRRD5z7do1oqKiaNu2Lf7+/vzwww/07NlTRoy5kBR+kWOlSpViw4YN3HvvvURERPDLL78YHUnkA5ZTFgZGD6RGixp88MEHDB48mO+//54qVaoYHc3tyOyc4q4dOXKExo0bU7hwYXbs2IG/v7/RkUQetfXXrYR/Gk6KNQXS4J1a7zCowyCjYxV4MjuncLiAgAC+/fZbEhMTadz4NUaO/Efm8i/g7mZO/bVr1/LM4GdsRd8EZk8zSWWSnBdS3JEUfpErderUYcqUrZw4sYCJE71o1kxL8S+g0ufUf/NN29c7Pc9Hjx6lTZs2tGnThmJ/FsPLw0umXMgjpPCLXLt4sTYmU/pc/mmsWHHW6EjCCbI7p/6VK1d44403qFGjBtu2bWPatGn8svkXtry4hfGh44npGiOfwDWYfIBL5JptLn8TycmatLQU5s/vTETEq7Rs2dLoaMKB0ufUT19F69Y59bXWLF++nNdee43ff/+dyMhIJk+ejJ+fHyBTLuQl0uIXuZY+1/r48Yovv7zEww9foE2bNsyYMUNm9SxA/v95/vfSifv37yckJITnn3+eMmXKEBcXx6JFizKKvshjtNZ5/lK3bl0t8o8rV67op59+WgO6Z8+e+vr160ZHEk5y4cIF/fLLL2uTyaRLliyp582bp1NTU42OJeyAeJ1JTZUWv3C4IkWKEB0dzRtvvMH8+fOJiIjgwoULRscSDpSWlsbcuXN5+OGHmTNnDi+//DJHjx6lZ8+emM1mo+OJO5DCL5zCZDIxadIkFi9eTFxcHA0bNuTnn382OpbIJcspC32W9KHa49Xo06cP1atXZ+/evcyaNQtfX1+j44lsksIvnOqFF15g8+bNXLx4kYYNG7Jx40ajI4m7tHjLYpp81IS5R+dytNFRxi4cS2xsLIGBgUZHEzkkhV843aOPPsquXbsoV64cLVu2ZPbs2UZHEjmwc+dO2rZtS+SoSNJIAxOYCpkoVLmQzK+TT7m88CulyimltiilDiulDimlBrg6g3C9ihUrEhcXR8uWLenXrx9RUVGkpqYaHUtkQWtNbGws4eHhNGzYkB07dtCjeQ9Z97agyOyMrzMvgB9Qx/59MeAXoNrt7iOjegqO1NRUPWTIEA3o8PBwnZiYaHQkcQOr1arXrVunGzVqpAFdunRp/fbbb+tLly5prbWOOxmnJ22bpONOxhmcVGQHWYzqMXyoJvA1EH67baTwFzwLFizQhQoV0gEBAfro0aNGx3F7aWlpetWqVbpOnToa0OXKldPvv/++vnr1qtHRRC5kVfgN7eNXSlUEagM7M7mtl1IqXikVf/78eZdnE87VrVs3Nm7cyPnz56lfvz7vv78nx5N/ibtz40RrqampLF26lJo1a9KuXTsuXbrEggULOHbsGP369cPHx8fouMIJDJuWWSlVFNgKTNRaf3G7bWVa5oLrf//7H2FhIzlxYgFKeeHtbSImRt30qVDhOBYLhLxgIaXsFkwnm1Am5T1+/z2aGjVqMHz4cDp06ICHh8zkUlDkqWmZlVKFgFXA0jsVfVGwVapUicjIjwEvtDZz7VoaH3xwSKZ6cJJ567aQ/FwYOmQUaZ0juH5feb788kv2799Pp06dpOi7CSNG9ShgAfCT1vodV+9f5D0tWnjj42PGZLKiVApLl/agZcuWsrKXg2j7CJ3OnTuzeHs3MCeDKQ1MybR7tSRPPfUUJpOM7HYnRjzbjwIvAM2UUvvsl1YG5BB5RPrkXxMmmNi6tRAzZjyHxWKhZs2aDB8+nH/++cfoiPnSuXPnmDp1KgEBAYSGhrJ27VqeDqqLp7kQCjNehTyJbBpqdExhhMzO+Oa1i4zqcT8JCQm6a9euGSNMoqOjtdVqNTpWnpeWlqa//fZb3a5dO+3h4aEB3bhxY/3JJ5/of/75R2stQzLdCVmM6pE1d0We9t1339GvXz8OHDhAeHg4s2bNIiAgwOhYec7vv//Oxx9/zIIFCzh+/DglS5YkMjKSHj16ULVqVaPjCYNkdXJXCr/I81JTU/nwww958803uXr1KoMHD2bkyJEULVrU6GiGsZyysPnXzXj+4cn2z7azdu1arFYrYWFh9OzZk6eeegovLy+jYwqDSeEX+d7Zs2cZNmwYixYtwt/fn+nTp9OhQwe3mi9Ga82S2CV029qNVJ0KaVBiTQl6t+pN9+7dqVSpktERRR4ihV8UGHFxcURFRbF3717CwsLo3v0jjh+vSEgIBXL8v9aa3bt3s2rVKr744guOlTkGzQATKEyMCxnLyMdGGh1T5EFZFX4ZtCvynUaNGrF7927mzp3L0KFfERNzH5CGlxds3Khp0iT/v6zT0tKIi4vLKPanTp3Cw8ODsLAwGtTuztK0caCT0VZP7rsaZnRckc/I4F2RL5nNZl5++WUGDPgCpbwAM9evW2nV6i0GDRpEfHx8vvsQWEpKChs3bqRPnz6ULVuWpk2bMmfOHGrXrs0nn3zCuXPn+Oabb6h+zzBMn8bAlvGYPo3hwr4C+G+OcCrp6hH5msUCYWGQnKwxm9MIDn4Ti+UdkpOTCQgIoEuXLnTu3JkHHnjA6Kg3sZyyEHs8luD7g7l8+DKrVq1i9erVJCYmUqRIEVq3bs0zzzxDq1atKFas2M33zfidwdPz3wufC5FO+vhFgWWxQGwsGX38iYmJrFq1iiVLlrB161bAthhM586defbZZylZsqSheb/e8zXPrn2WZGsypAKfQPErxWnbti3PPPMMjz/++B0nR7v1dxYiM1L4hVs6ceIEn332GZ9++imHDx+mUKFCtGzZki5dutCmTRv27fNxegH9+++/iY2NZdOmTWzatIkjpY5knJxFK14sH8ncLnPx9PR0TgDhtqTwC7emtWb//v0sXbqUZcuW8ccff1C4cBjXr69F60IUKqRZuTKR1q1L5HremqSkJCwWS0ahj4+Px2q1UqRIER577DEKB9QnuvBbYEoGqydzH42hV0tptgvHk8IvhF1aWhqxsbG8/noiP/zwFLbBbSnAKDw936F8+fKUL1+eChUqZHxN/75cuXIZH4yat97Cqj2xPFO7KY/4+WQU+u3bt5OUlITZbKZBgwY0b96c5s2b06BBAzw9PZk8GUbOsWAtH4vpZAgT+gTzxhsGHhBRYEnhF+IWtpOkmuRkMJutREV9idm8i5MnT3LixAlOnjxJQkLCv0YHlSlTBu+H6nG8aYxtpss0T/jkQTh9iOrVq2cU+qZNm3LPPfdksV85OSucT8bxC3EL26ygyt7HbyY4uD3Q/qZtrl+/zunTp296Mzhx4gRfXfj7/6c31snUfCKUb9/ciJ+fXzb3KydnhXGk8Av35m+BxrHgHwL8uwJ7eXlRqVKlf02F0GC9hd471oO29dNHPfF8top+uuBgKfjCOFL4hduynLIQtjiM5LRkPM2exHSNIbhc9qqx7WRsDKv2xNKuboicnBX5ihR+4bZij8eSnJZMmk4jOS3Z9oGqbBZ+sBV/KfgiP5IpG0S+ZzllYfL2yVhOWXJ0v5CKIXiaPTErM55mT0IqhjgnoBB5jLT4Rb6Wm+6a4HLBxHSNIfZ4LCEVQ3LU2hciP5PCL/KE9LlrclqAc9tdE1wuWAq+cDtS+IXhctNqT++uSb+vdNcIcWdS+IVD3G2LHXLXapfuGiFyTgq/yLXctNgh96126a4RImek8IsMRvazS6tdCNeRwl/A3G3xNrqfXVrtQriOFP48yIjiLf3sQrgPKfxZyM3Jytze14jiLf3sQriPAl34jWg55/ZEp1HFW1rtQrgPQwq/UqoFMBMwAx9prac4eh9GtZxze6LTyOItrXYh3IPLC79Sygx8AIQDp4HdSqnVWuvDjtyPUS1nR3SZSPEWQjiTES3++sAxrfWvAEqp5cCTgEMLv1EtZ0d0mUjxFkI4k8uXXlRKtQdaaK172H9+AWigtY66ZbteQC+A8uXL1z1x4kSO95Wbk6xCCJHf5bulF7XW84B5YFtz924eQ1rOQgjxb0bMx/87UO6Gn/3t1wkhhHABIwr/buAhpdQDSilP4DlgtQE5hBDCLbm8q0drnaqUigK+xTacc6HW+pCrcwghhLsypI9fa70OWGfEvoUQwt3JmrtCCOFmpPALIYSbcfk4/ruhlDoP5Hwgv81/gD8dGMdRJFfOSK6ckVw5U1BzVdBal7r1ynxR+HNDKRWf2QcYjCa5ckZy5Yzkyhl3yyVdPUII4Wak8AshhJtxh8I/z+gAWZBcOSO5ckZy5Yxb5SrwffxCCCFu5g4tfiGEEDeQwi+EEG6mQBR+pVQHpdQhpZRVKVXvltveUEodU0odUUpFZHH/B5RSO+3brbBPHufojCuUUvvsl+NKqX1ZbHdcKfWjfbt4R+fIZH9jlFK/35CtVRbbtbAfw2NKqWEuyPW2UupnpdQBpdSXSqniWWznkuN1p99fKeVlf46P2V9LFZ2V5YZ9llNKbVFKHba//gdksk2IUuriDc/vKGfnsu/3ts+LsnnPfrwOKKXquCBTwA3HYZ9S6pJSauAt27jkeCmlFiqlzimlDt5wXQml1Eal1FH7V98s7htp3+aoUiryrgJorfP9BagKBACxQL0brq8G7Ae8gAeA/wHmTO7/OfCc/fs5QF8n550OjMrituPAf1x47MYAQ+6wjdl+7B4EPO3HtJqTcz0OeNi/fwt4y6jjlZ3fH3gZmGP//jlghQueOz+gjv37YsAvmeQKAda46vWU3ecFaAWsBxTQENjp4nxm4Ay2Dzi5/HgBTYE6wMEbrpsKDLN/Pyyz1zxQAvjV/tXX/r1vTvdfIFr8WuuftNZHMrnpSWC51vq61vo34Bi2pR8zKKUU0AyItl/1CfCUs7La9/cs8Jmz9uEEGctlaq2TgfTlMp1Ga71Ba51q//F7bOs2GCU7v/+T2F47YHsthdmfa6fRWidorX+wf38Z+Ako68x9OtCTwGJt8z1QXCnl58L9hwH/01rf7YwAuaK13gb8dcvVN76GsqpDEcBGrfVfWutEYCPQIqf7LxCF/zbKAqdu+Pk0//7DKAn8fUORyWwbR2oCnNVaH83idg1sUErtsS8/6QpR9n+3F2bx72V2jqMzdcPWOsyMK45Xdn7/jG3sr6WL2F5bLmHvWqoN7Mzk5mCl1H6l1HqlVHUXRbrT82L0a+o5sm58GXG8AEprrRPs358BSmeyjUOOW55devFWSqlNQJlMbhqhtf7a1Xkyk82Mnbh9a7+x1vp3pdR9wEal1M/21oFTcgEfAuOx/aGOx9YN1S03+3NErvTjpZQaAaQCS7N4GIcfr/xGKVUUWAUM1FpfuuXmH7B1Z1yxn7/5CnjIBbHy7PNiP4fXFngjk5uNOl430VprpZTTxtrnm8KvtW5+F3fLzjKPF7D9m+lhb6nd9VKQd8qolPIAngHq3uYxfrd/PaeU+hJbN0Ou/mCye+yUUvOBNZnc5JTlMrNxvF4E2gBh2t7BmcljOPx4ZSI7v3/6Nqftz/O92F5bTqWUKoSt6C/VWn9x6+03vhFordcppWYrpf6jtXbqhGTZeF6MXIK1JfCD1vrsrTcYdbzsziql/LTWCfZur3OZbPM7tvMQ6fyxndvMkYLe1bMaeM4+4uIBbO/cu27cwF5QtgDt7VdFAs76D6I58LPW+nRmNyqliiiliqV/j+0E58HMtnWUW/pVn85ify5fLlMp1QJ4HWirtb6axTauOl7Z+f1XY3vtgO21tDmrNytHsZ9DWAD8pLV+J4ttyqSfa1BK1cf2N+/UN6RsPi+rga720T0NgYs3dHM4W5b/dRtxvG5w42soqzr0LfC4UsrX3i37uP26nHH22WtXXLAVrNPAdeAs8O0Nt43ANiLjCNDyhuvXAffbv38Q2xvCMWAl4OWknIuAPrdcdz+w7oYc++2XQ9i6PJx97D4FfgQO2F94frfmsv/cCtuokf+5KNcxbH2Z++yXObfmcuXxyuz3B8Zhe2MC8La/do7ZX0sPuuAYNcbWRXfghuPUCuiT/joDouzHZj+2k+SNXJAr0+flllwK+MB+PH/khtF4Ts5WBFshv/eG61x+vLC98SQAKfba1R3bOaEY4CiwCShh37Ye8NEN9+1mf50dA166m/3LlA1CCOFmCnpXjxBCiFtI4RdCCDcjhV8IIdyMFH4hhHAzUviFEMLNSOEXQgg3I4VfCCHcjBR+Ie6CUuoR+8R23vZPqh5SStUwOpcQ2SEf4BLiLimlJmD7xK4PcFprPdngSEJkixR+Ie6Sfd6e3UASto/2pxkcSYhska4eIe5eSaAottWvvA3OIkS2SYtfiLuklFqNbTWuB7BNbhdlcCQhsiXfzMcvRF6ilOoKpGitlymlzECcUqqZ1nqz0dmEuBNp8QshhJuRPn4hhHAzUviFEMLNSOEXQgg3I4VfCCHcjBR+IYRwM1L4hRDCzUjhF0IIN/N/XuNmwldFI74AAAAASUVORK5CYII=\n",
"text/plain": [
"