{
"cells": [
{
"cell_type": "markdown",
"metadata": {
"editable": true,
"slideshow": {
"slide_type": ""
},
"tags": []
},
"source": [
"# 2. Slug Test for Confined Aquifer - Multi-well"
]
},
{
"cell_type": "markdown",
"metadata": {
"editable": true,
"slideshow": {
"slide_type": ""
},
"tags": []
},
"source": [
"### Import packages"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt\n",
"import numpy as np\n",
"import pandas as pd\n",
"\n",
"import timflow.transient as tft\n",
"\n",
"plt.rcParams[\"figure.figsize\"] = [5, 3]"
]
},
{
"cell_type": "markdown",
"metadata": {
"editable": true,
"slideshow": {
"slide_type": ""
},
"tags": []
},
"source": [
"## Introduction and Conceptual Model\n",
"\n",
"A well (Ln-2) fully penetrates a sandy confined aquifer of 6.1 m thickness. Additionally, a fully penetrating observation well (Ln-3) is located 6.45 m away from the test well. The slug displacement is 2.798 m. Head change was recorded at the slug well and the observation well. The well and casing radii of the slug well are 0.102 and 0.051 m, respectively. For the observation well, they are 0.071 and 0.025 m, respectively."
]
},
{
"cell_type": "markdown",
"metadata": {
"editable": true,
"slideshow": {
"slide_type": ""
},
"tags": []
},
"source": [
"
"
]
},
{
"cell_type": "markdown",
"metadata": {
"editable": true,
"slideshow": {
"slide_type": ""
},
"tags": []
},
"source": [
"### Load data"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
"data1 = np.loadtxt(\"data/ln-2.txt\")\n",
"t1 = data1[:, 0] / 60 / 60 / 24 # convert time from seconds to days\n",
"h1 = data1[:, 1]\n",
"\n",
"data2 = np.loadtxt(\"data/ln-3.txt\")\n",
"t2 = data2[:, 0] / 60 / 60 / 24\n",
"h2 = data2[:, 1]"
]
},
{
"cell_type": "markdown",
"metadata": {
"editable": true,
"slideshow": {
"slide_type": ""
},
"tags": []
},
"source": [
"### Parameters and model"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [],
"source": [
"# known parameters\n",
"H0 = 2.798 # initial displacement, m\n",
"b = 6.1 # aquifer thickness, m\n",
"rw1 = 0.102 # well radius of Ln-2 Well, m\n",
"rw2 = 0.071 # well radius of observation Ln-3 Well, m\n",
"rc1 = 0.051 # casing radius of Ln-2 Well, m\n",
"r = 6.45 # distance from observation well to test well, m"
]
},
{
"cell_type": "markdown",
"metadata": {
"editable": true,
"slideshow": {
"slide_type": ""
},
"tags": []
},
"source": [
"Converting slug displacement into volume"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Slug: 0.02286 m^3\n"
]
}
],
"source": [
"Q = np.pi * rc1**2 * H0\n",
"print(f\"Slug: {Q:.5f} m^3\")"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"self.neq 1\n",
"solution complete\n"
]
}
],
"source": [
"ml = tft.ModelMaq(kaq=10, z=[0, -b], Saq=1e-4, tmin=1e-5, tmax=0.01)\n",
"w = tft.Well(ml, xw=0, yw=0, rw=rw1, rc=rc1, tsandQ=[(0, -Q)], layers=0, wbstype=\"slug\")\n",
"ml.solve()"
]
},
{
"cell_type": "markdown",
"metadata": {
"editable": true,
"slideshow": {
"slide_type": ""
},
"tags": []
},
"source": [
"### Estimate aquifer parameters\n",
"The hydraulic conductivity and specific storage are calibrated."
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"...........................................\n",
"Fit succeeded.\n"
]
}
],
"source": [
"# unknown parameters: kaq, Saq\n",
"cal = tft.Calibrate(ml)\n",
"cal.set_parameter(name=\"kaq\", layers=0, initial=10)\n",
"cal.set_parameter(name=\"Saq\", layers=0, initial=1e-4)\n",
"cal.seriesinwell(name=\"Ln-2\", element=w, t=t1, h=h1)\n",
"cal.series(name=\"Ln-3\", x=r, y=0, layer=0, t=t2, h=h2)\n",
"cal.fit()"
]
},
{
"cell_type": "code",
"execution_count": 31,
"metadata": {
"editable": true,
"slideshow": {
"slide_type": ""
},
"tags": []
},
"outputs": [
{
"data": {
"text/html": [
"
\n",
"\n",
"
\n",
" \n",
" \n",
" | \n",
" optimal | \n",
"
\n",
" \n",
" \n",
" \n",
" | kaq_0_0 | \n",
" 1.166109 | \n",
"
\n",
" \n",
" | Saq_0_0 | \n",
" 0.000009 | \n",
"
\n",
" \n",
"
\n",
"
"
],
"text/plain": [
" optimal\n",
"kaq_0_0 1.166109\n",
"Saq_0_0 0.000009"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"RMSE: 0.010 m\n"
]
}
],
"source": [
"display(cal.parameters.loc[:, [\"optimal\"]])\n",
"print(f\"RMSE: {cal.rmse():.3f} m\")"
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAdcAAAFBCAYAAAAlu+WfAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8fJSN1AAAACXBIWXMAAA9hAAAPYQGoP6dpAAB+70lEQVR4nO3dd3hTZRvA4d9Jmu5NJ1Batuwie4Oyt6BWQJYIiiAgshUBlaUCRRwIKiCCoFAXIoKVwoci08GeLWW0hQLdK03O90dtILRAAykdPDdXLpqTM54kp+fp+553KKqqqgghhBDCajRFHYAQQghR2khyFUIIIaxMkqsQQghhZZJchRBCCCuT5CqEEEJYmSRXIYQQwsokuQohhBBWJslVCCGEsDJJrkIIIYSVSXIV90VRFGbOnGnxdlFRUSiKwsqVK60e04MWFBTEkCFDijoMiz2o72DmzJkoimK1/a1cuRJFUYiKirLaPh8W9/OdR0REoCgKERERVo+rNJLkWgrkXmwURWHXrl15XldVlYCAABRFoXv37kUQ4b3L/YXOfWi1Wnx8fHjyySc5duxYUYeXr6NHjzJz5swivfj/+OOPtGnTBh8fHxwdHalUqRJPP/00W7ZsKbKYhHiY2BR1AMJ67O3tWbt2LS1btjRbvmPHDi5cuICdnV0RRXb/xowZQ6NGjdDr9fz7778sXbqUiIgIDh8+jJ+fX1GHZ+bo0aPMmjWLtm3bEhQU9MCP/9577zFx4kTatGnD1KlTcXR05PTp0/z666+sW7eOzp07AxAYGEh6ejo6ne6BxyhEaSfJtRTp2rUr33zzDe+//z42Nje+2rVr19KgQQPi4+OLMLr706pVK5588knT8+rVqzNy5Ei++OILJk2aVISRFS/Z2dm89dZbdOjQga1bt+Z5/fLly6afFUXB3t7+QYb3UEhNTcXJyamowxBFTKqFS5F+/fpx9epVtm3bZlqWlZXFhg0b6N+/f77bpKam8uqrrxIQEICdnR3Vq1fnvffe49bJkjIzM3nllVfw9vbGxcWFnj17cuHChXz3efHiRZ577jl8fX2xs7OjVq1afP7559Z7o+QkW4AzZ87c07GXLFlCrVq1cHR0xMPDg4YNG7J27VrT60OGDMm31Hm3+4crV67kqaeeAqBdu3am6uzc+1T79++nU6dOeHl54eDgQMWKFXnuuecsffu3FR8fT1JSEi1atMj3dR8fH9PP+d1/GzJkCM7Ozly8eJHevXvj7OyMt7c3EyZMwGAwmO3r6tWrDBw4EFdXV9zd3Rk8eDD//PNPge/pffnllzRo0AAHBwc8PT155plnOH/+/D297++//55u3bpRtmxZ7OzsqFy5Mm+99ZZZzDNmzECn03HlypU8248YMQJ3d3cyMjJMy37++WdatWqFk5MTLi4udOvWjSNHjphtl/t5nTlzhq5du+Li4sKAAQNuG2fu+XPy5EmeffZZ3Nzc8Pb2Zvr06aiqyvnz5+nVqxeurq74+fmxYMGCPPu4fPkyw4YNw9fXF3t7e+rVq8eqVavyrJeQkMCQIUNwc3MzfT8JCQn5xnX8+HGefPJJPD09sbe3p2HDhvzwww+3fR/i7iS5liJBQUE0a9aMr776yrTs559/JjExkWeeeSbP+qqq0rNnTxYtWkTnzp1ZuHAh1atXZ+LEiYwfP95s3eeff57Q0FA6duzIvHnz0Ol0dOvWLc8+4+LiaNq0Kb/++iujR49m8eLFVKlShWHDhhEaGmq195p7P9PDw8PiYy9fvpwxY8ZQs2ZNQkNDmTVrFsHBwezZs+e+42rdujVjxowBYNq0aaxevZrVq1dTo0YNLl++TMeOHYmKimLKlCksWbKEAQMG8Oeff973cXP5+Pjg4ODAjz/+yLVr1+5pHwaDgU6dOlGmTBnee+892rRpw4IFC1i2bJlpHaPRSI8ePfjqq68YPHgws2fPJiYmhsGDBxfoGLNnz2bQoEFUrVqVhQsXMm7cOMLDw2nduvVtE8CdrFy5EmdnZ8aPH8/ixYtp0KABb7zxBlOmTDGtM3DgQLKzs1m/fr3Ztrl/gPbt29dUkl+9ejXdunXD2dmZ+fPnM336dI4ePUrLli3z3EvPzs6mU6dO+Pj48N5779G3b9+7xhsSEoLRaGTevHk0adKEt99+m9DQUDp06EC5cuWYP38+VapUYcKECezcudO0XXp6Om3btmX16tUMGDCAd999Fzc3N4YMGcLixYtN66mqSq9evVi9ejXPPvssb7/9NhcuXMj3+zly5AhNmzbl2LFjTJkyhQULFuDk5ETv3r359ttvC/T5i3yoosRbsWKFCqj79u1TP/jgA9XFxUVNS0tTVVVVn3rqKbVdu3aqqqpqYGCg2q1bN9N23333nQqob7/9ttn+nnzySVVRFPX06dOqqqrq33//rQLqSy+9ZLZe//79VUCdMWOGadmwYcNUf39/NT4+3mzdZ555RnVzczPFFRkZqQLqihUr7vjetm/frgLq559/rl65ckW9dOmSumXLFrVKlSqqoijq3r17LT52r1691Fq1at3xuIMHD1YDAwPzLJ8xY4Z6669NYGCgOnjwYNPzb775RgXU7du3m6337bffmr6nwvTGG2+ogOrk5KR26dJFnT17tnrgwIE86+X3HQwePFgF1DfffNNs3fr166sNGjQwPd+4caMKqKGhoaZlBoNBfeyxx/Ls89bPLCoqStVqters2bPNjnHo0CHVxsYmz/Jb5Z7vkZGRpmW53+3NXnjhBdXR0VHNyMgwLWvWrJnapEkTs/XCwsLMvq/k5GTV3d1dHT58uNl6sbGxqpubm9ny3M9rypQpd4w5V+5nMWLECNOy7OxstXz58qqiKOq8efNMy69fv646ODiYnVuhoaEqoH755ZemZVlZWWqzZs1UZ2dnNSkpSVXVG7/b77zzjtlxWrVqlef7efzxx9U6deqYfU5Go1Ft3ry5WrVqVdOy3N/FW89rkT8puZYyTz/9NOnp6WzatInk5GQ2bdp02yrhzZs3o9VqTSWtXK+++iqqqvLzzz+b1gPyrDdu3Diz56qqsnHjRnr06IGqqsTHx5senTp1IjExkYMHD97T+3ruuefw9vambNmydO7cmcTERFavXk2jRo0sPra7uzsXLlxg37599xTLvXJ3dwdg06ZN6PX6QjvOrFmzWLt2LfXr1+eXX37htddeo0GDBjz66KMFbmH94osvmj1v1aoVZ8+eNT3fsmULOp2O4cOHm5ZpNBpGjRp1132HhYVhNBp5+umnzb4nPz8/qlatyvbt2wv4Tm9wcHAw/ZycnEx8fDytWrUiLS2N48ePm14bNGgQe/bsMbudsGbNGgICAmjTpg0A27ZtIyEhgX79+pnFp9VqadKkSb7xjRw50qJ4n3/+edPPWq2Whg0boqoqw4YNMy13d3enevXqZp/75s2b8fPzo1+/fqZlOp2OMWPGkJKSwo4dO0zr2djYmMWl1Wp5+eWXzeK4du0av/32G08//bTpc4uPj+fq1at06tSJU6dOcfHiRYvem8ghDZpKGW9vb9q3b8/atWtJS0vDYDCYNQS62blz5yhbtiwuLi5my2vUqGF6Pfd/jUZD5cqVzdarXr262fMrV66QkJDAsmXLzKoQb3ZzgxpLvPHGG7Rq1YqUlBS+/fZb1q1bh0Zz429DS449efJkfv31Vxo3bkyVKlXo2LEj/fv3v+19Smtp06YNffv2ZdasWSxatIi2bdvSu3dv+vfvf8eW3ImJiaSnp5ue29ra4unpecdj9evXj379+pGUlMSePXtYuXIla9eupUePHhw+fPiODZns7e3x9vY2W+bh4cH169dNz8+dO4e/vz+Ojo5m61WpUuWOcQGcOnUKVVWpWrVqvq/fS+vlI0eO8Prrr/Pbb7+RlJRk9lpiYqLp55CQEMaNG8eaNWt44403SExMZNOmTbzyyiume+mnTp0C4LHHHsv3WK6urmbPbWxsKF++vEXxVqhQwey5m5sb9vb2eHl55Vl+9epV0/Nz585RtWpVs3Mf8v+d9ff3x9nZ2Wy9W39nT58+jaqqTJ8+nenTp+cb6+XLlylXrpwF706AJNdSqX///gwfPpzY2Fi6dOliKjEVNqPRCMCzzz5723tvdevWvad916lTh/bt2wPQu3dv0tLSGD58OC1btiQgIMCiY9eoUYMTJ06wadMmtmzZwsaNG/noo4944403mDVrFsBtGy3d2qjHEoqisGHDBv78809+/PFHfvnlF5577jkWLFjAn3/+medCmGvs2LFmDVbatGlT4I78rq6udOjQgQ4dOqDT6Vi1ahV79uwxldLyo9VqLXpfljIajSiKws8//5zvsW73OdxOQkICbdq0wdXVlTfffJPKlStjb2/PwYMHmTx5suncgJw/Erp3725Krhs2bCAzM5Nnn33WLD7Iue+aXzevm1viA9jZ2eVJdneT3/u+3eeu3tK40Jpy3+uECRPo1KlTvusU5A8mkZck11LoiSee4IUXXuDPP//M03jjZoGBgfz6668kJyeblV5zq9ECAwNN/xuNRs6cOWP2l++JEyfM9pfbkthgMJgSYWGZN28e3377LbNnz2bp0qUWH9vJyYmQkBBCQkLIysqiT58+zJ49m6lTp2Jvb4+Hh0e+DWtySwZ3crfRiJo2bUrTpk2ZPXs2a9euZcCAAaxbt86sqvBmkyZNMrv439yIyxINGzZk1apVxMTE3NP2NwsMDGT79u2kpaWZlV5Pnz59120rV66MqqpUrFiRatWq3XcsERERXL16lbCwMFq3bm1aHhkZme/6gwYNolevXuzbt481a9ZQv359atWqZRYf5DQOK+zz2FKBgYH8+++/GI1Gs4Se3+9seHg4KSkpZn+s3Po7W6lSJSCntqC4vdeSTu65lkLOzs58/PHHzJw5kx49etx2va5du2IwGPjggw/Mli9atAhFUejSpQuA6f/333/fbL1bW/9qtVr69u3Lxo0bOXz4cJ7j5dcF4l5VrlyZvn37snLlSmJjYy069s3VbJBTzVqzZk1UVTXdC61cuTKJiYn8+++/pvViYmIK1Hoyt4/jrcn5+vXreUohwcHBQE5Xp9upWbMm7du3Nz0aNGhw23XT0tLYvXt3vq/l3kO/tWrwXnTq1Am9Xs/y5ctNy4xGIx9++OFdt+3Tpw9arZZZs2bl+TxUVc3z/dxNbonv5n1lZWXx0Ucf5bt+ly5d8PLyYv78+ezYscPsDxfIeW+urq7MmTMn33vj1jyPLdW1a1diY2PN/mjOzs5myZIlODs7m2okunbtSnZ2Nh9//LFpPYPBwJIlS8z25+PjQ9u2bfnkk0/y/aOrKN9rSScl11KqIF0ievToQbt27XjttdeIioqiXr16bN26le+//55x48aZ/oIPDg6mX79+fPTRRyQmJtK8eXPCw8PzLaXMmzeP7du306RJE4YPH07NmjW5du0aBw8e5Ndff73n7iH5mThxIl9//TWhoaHMmzevwMfu2LEjfn5+tGjRAl9fX44dO8YHH3xAt27dTCX4Z555hsmTJ/PEE08wZswY0tLS+Pjjj6lWrdpdG2UFBwej1WqZP38+iYmJ2NnZ8dhjj7F27Vo++ugjnnjiCSpXrkxycjLLly/H1dWVrl27WuUzSUtLo3nz5jRt2pTOnTsTEBBAQkIC3333Hf/73//o3bs39evXv+/j9O7dm8aNG/Pqq69y+vRpHnnkEX744QfTZ3yn0nvlypV5++23mTp1KlFRUfTu3RsXFxciIyP59ttvGTFiBBMmTChwLM2bN8fDw4PBgwczZswYFEVh9erVt61O1el0PPPMM3zwwQdotVqzxkGQU5X+8ccfM3DgQB599FGeeeYZvL29iY6O5qeffqJFixZ5/iB9UEaMGMEnn3zCkCFDOHDgAEFBQWzYsIHff/+d0NBQ0/nbo0cPWrRowZQpU4iKiqJmzZqEhYWZ3X/O9eGHH9KyZUvq1KnD8OHDqVSpEnFxcezevZsLFy7wzz//POi3WTo8+AbKwtpu7opzJ7d2xVHVnG4Hr7zyilq2bFlVp9OpVatWVd99913VaDSarZeenq6OGTNGLVOmjOrk5KT26NFDPX/+fJ6uOKqqqnFxceqoUaPUgIAAVafTqX5+furjjz+uLlu2zLSOpV1xvvnmm3xfb9u2rerq6qomJCQU+NiffPKJ2rp1a7VMmTKqnZ2dWrlyZXXixIlqYmKi2b63bt2q1q5dW7W1tVWrV6+ufvnllwXqiqOqqrp8+XK1UqVKqlarNXVfOHjwoNqvXz+1QoUKqp2dnerj46N2795d3b9//x0/A0vo9Xp1+fLlau/evdXAwEDVzs5OdXR0VOvXr6++++67amZmpmnd23XFcXJyyrPf/N73lStX1P79+6suLi6qm5ubOmTIEPX3339XAXXdunV33FZVc7rztGzZUnVyclKdnJzURx55RB01apR64sSJO77H/Lri/P7772rTpk1VBwcHtWzZsuqkSZPUX3755bZdR/bu3asCaseOHW97nO3bt6udOnVS3dzcVHt7e7Vy5crqkCFDzL6v231et5P7WVy5csVs+e3206ZNmzzdxuLi4tShQ4eqXl5eqq2trVqnTp18f4+uXr2qDhw4UHV1dVXd3NzUgQMHqn/99Ve+v3dnzpxRBw0apPr5+ak6nU4tV66c2r17d3XDhg1mn8ftPk+Rl6KqhXi3XAjxUPnuu+944okn2LVrV6G3vr4f//zzD8HBwXzxxRcMHDiwqMMRpZDccxVC3JObuwfBjXt6rq6uPProo0UUVcEsX74cZ2dn+vTpU9ShiFJK7rkKIe7Jyy+/THp6Os2aNSMzM5OwsDD++OMP5syZYzaoQ3Hy448/cvToUZYtW8bo0aNlgH1RaKRaWAhxT9auXcuCBQs4ffo0GRkZVKlShZEjRzJ69OiiDu22goKCiIuLo1OnTqxevTrPACpCWIskVyGEEMLK5J6rEEIIYWWSXIUQQggrK9IGTTt37uTdd9/lwIEDptFvevfufcdtIiIiGD9+PEeOHCEgIIDXX3+dIUOGFPiYRqORS5cu4eLictdh6oQQQpReqqqSnJxM2bJlLR4f+m6KNLmmpqZSr149nnvuuQI1iY+MjKRbt268+OKLrFmzhvDwcJ5//nn8/f1vO+j0rS5dukRAQMD9hi6EEKKUOH/+vMUzG91NsWnQpCjKXUuukydP5qeffjIbO/aZZ54hISGBLVu2FOg4iYmJuLu7c/78+TxTR4nSQ6/Xs3XrVjp27HhPU5iJ0kPOBZHr1nMhKSnJNESom5ubVY9Vovq57t69O8/MDZ06dcozaffNMjMzzQZFT05OBnImVy6uffHE/bOxscHR0REHBwe5oD7k5FwQuW49F3InZiiMW4QlKrnGxsbi6+trtszX15ekpCTS09PzTZZz5841zdF5s61bt+aZ6FmUPtu2bSvqEEQxIeeCyJV7LqSlpRXaMUpUcr0XU6dOZfz48abnudUAHTt2lGrhUkyv17Nt2zbTJOHi4SXngsh167mQlJRUaMcqUcnVz8+PuLg4s2VxcXG4urretorXzs4OOzu7PMt1Op38oj0E5HsWueRcELlyz4XCPB9KVHJt1qwZmzdvNlu2bds2mjVrVkQRCSGKK6PRSFZWFnq9HhsbGzIyMjAYDEUdlniAdDodWq22SI5dpMk1JSXFbMLtyMhI/v77bzw9PalQoQJTp07l4sWLfPHFFwC8+OKLfPDBB0yaNInnnnuO3377ja+//pqffvqpqN6CEKIYysrKIjIyEqPRiKqq+Pn5cf78eenb/hByd3fHz8/vgX/3RZpc9+/fT7t27UzPc++NDh48mJUrVxITE0N0dLTp9YoVK/LTTz/xyiuvsHjxYsqXL8+nn35a4D6uQojST1VVYmJi0Gq1pj7tKSkpODs7W32gAFF8qapKWloaly9fBsDf3/+BHr9Ik2vbtm25UzfblStX5rvNX3/9VYhRFYA+A3Z/AG4B4B4AbuXBxR+0915/H5OYTmR8KhW9nPB3ky5CQtyr7Oxs0tLSKFu2LI6OjqbqYXt7e0muD5nctjiXL1/Gx8fngR67RN1zLTaSLsJvb5kvUzQ5CdatfE7SdStv/rN7ANjn30l5/b5opoYdwqiCRoG5feoQ0qjCA3gjQpQ+ufdVbW1tizgSURzkdrnU6/UP9P6rJNd7oWhIq/k0+qvROGXEYpNyCQxZOUk36SKc35P/dnauNyXdnMR73daXjd/H4qd6EYcHBlXLtLDDtK7mLSVYIe6D3F8VUHTngSTXe7D+jJapf/W+UdJ8ohYhNRwg8QIkRv/3f+7jfM7/aVchMwkuH815/McD+Pq/P7D1qpZo1YezalmMv2yHKnWgTFUoUwWcvEAuFkIIUSJIcrVQTGK6qQoXwKjCtG+P0rp6O/zLN4DyDfLfMCsVEi/eSLb/PTKvnSM2+hT+XMVWMVBZiaEyMXD0ABy9aXt7t5xE6/Vfsi1TJednz8qgsy/09y2EKFoRERG0a9eO69ev4+7ubpV9BgUFMW7cuDsOISvujSRXC0XGp5oSay6DqhIVn3bnalxbJ/CulvO4iR3w575oXgv7F2/1GlU0MYyrDw2crsLV03D1FCSch4xEuLg/52FGIdu1PMnOldCVrYNzYH3wrZ2TfLXy9QohCldERASLFi1i7969JCUlUbVqVSZOnMiAAQOKOrQiJVdfC1X0ckKjYJZgtYpCkNe9j1Mc0qgCrat5ExWfRpCXY94krU+Ha2ch/lROso0/fSPxZiRik3Qej6TzcGkH5OZerR34PAK+dcCvdk7C9a0Fjp7SMlkIYTV//PEHdevWZfLkyfj6+rJp0yYGDRqEm5sb3bt3L+rwiowkVwv5uzkwt08dpoUdxqCqaBWFOX1q33eS8ndzuP0+dA45idG3ltnimIQ0es7/jiBiqKa5QA3lHDU10QTbXUSjT4OYf3IeN0mz9+VYalkOqUF8rlame5ce9G4ZfF+xC1FaPcg/RDMzM5k4cSLr1q0jKSmJhg0bsmjRIho1amS23u+//87UqVM5efIkwcHBfPrpp9SuXRuAc+fOMXr0aHbt2kVWVhZBQUG8++67dO3atUAxKIrC8uXL+emnn/jll18oV64cCxYsoGfPnrfdZtq0aWbPx44dy9atWwkLC5PkKixz15LmAxJ5NY0rqhtXcGOf4RHT8q8GNaaZRxLEHYbYwxB3BOIOQUI0jhlxPKaN4zH+6yv863tk762ATUBDKN8QyjUE/7qgc5ASrnioPegucpMmTWLjxo2sWrWKwMBA3nnnHTp16sTp06fx9PQ0rTdx4kQWL16Mn58f06ZNo0ePHpw8eRKdTseoUaPIyspi586dODk5cfToUZydnS2KY9asWbzzzju8++67LFmyhAEDBnDu3DmzGO4mMTGRGjVqWHTc0kaS6z26Y0nzAbltFbW3M7h5Q5nKULOX6bU9xyJ574swamrOUVdzlnrKGapoLmGTFA1HouFIWM6KGhuuOVflt+tl2Weozn71EV7u00763oqHRr4NFwuxi1xqaioff/wxK1eupEuXLgAsX76cbdu28dlnnzFx4kTTujNmzKBDhw4ArFq1ivLly/Ptt9/y9NNPEx0dTd++falTpw4AlSpVsjiWIUOG0K9fPwDmzJnD+++/z969e+ncuXOBtv/666/Zt28fn3zyicXHLk0kuZZgllZRVyjrxwEeySnl/jd+ubuSRng/V8ok/AsXD8KF/ZB6Gc+kYwzQHmOANhyAC5u8SDvbFscqrSGwBZSpTExShpRsRal0zw0X79GZM2fQ6/W0aNHCtEyn09G4cWOOHTtmtu7NE5V4enpSvXp10zpjxoxh5MiRbN26lfbt29O3b1/q1q1rUSw3r+/k5ISrq6tpCMFatWpx7tw5AFq1asXPP/9stu327dsZOnQoy5cvp1Yt89tYDxtJriWcJVXU+SXjqX2aUKZuBeC/8ZlVlf2HDvH5um+orzlNI80JaiuRlFfi4diGnAeQbleGg2lV+NNYgz/UOozo3YGQxoEP4B0LUfgKo+Hig/D888/TqVMnfvrpJ7Zu3crcuXNZsGABL7/8coH3ces0bIqiYDQaAdi8eTN6vR4gzzSfO3bsoEePHixatIhBgwbd5zsp+SS5lgKWVFHfNRkrCuUCq7JFbcrm7KYAOJJBQ80pPmyVgUvsXtQL+3HIvEo37VW6aXNGo7r4kxdp0R1xfKQ9VGxLTLajlGpFiVVYDRdvp3Llytja2vL7778TGJjzR6per2ffvn15+qD++eefVKiQc4vm+vXrnDx50uz+ZkBAAC+++CIvvvgiU6dOZfny5RYl1zvJje1WERERdO/enfnz5zNixAirHKukk+T6ELpbMr71wpKpONDtif64/HfP9c+TF1m4Yh1NNMdorjlCA81JyinxcHgtHF6LikK8MYi/jXUINQbzZK8neLqJ5fd+hChKD7LhopOTEyNHjmTixImmKTffeecd0tLSGDZsmNm6b775JmXKlMHX15fXXnsNLy8vevfuDcC4cePo0qUL1apV4/r162zfvr3QGxZt376d7t27M3bsWPr27UtsbCyQM7azJY2gShtJriJfd7qwBPl6mu7dfmB4AgcyaKY9wftNErA7twPd1ePU0URSRxPJS/zA9c0LSY/shEOtrlDlcWKy7KVUK0qEB9lwcd68eRiNRgYOHEhycjINGzbkl19+wcPDI896Y8eO5dSpUwQHB/Pjjz+aJikwGAyMGjWKCxcu4OrqSufOnVm0aFGhxr1q1SrS0tKYO3cuc+fONS1v06YNERERhXrs4kxR7zTnWymUlJSEm5sbiYmJuLq6FnU4Jdb6fdF5qsxCGlXgjzPxjF2+hZaaw7TR/kNbzT+4K6mm7YyKln2GqmwzNGCrsRGj+jxeKK2Q9Xo9mzdvpmvXrnnuIYnSLSMjg8jISCpWrIi9vT1Go5GkpCRcXV1lyrmH0M3ng1arNbsuFGY+kJKruCe3K9lW9HLiquLBt8ZWfGtshRYDDTSn+bz5Vewit6G7eoImmuM00RznddZw+Mcgkq8/g0v9PsTYVpASrRCiVJA/48Q983dzoFnlMmaJMPd+rTZ3Bh/Fhr5PPIlz99ns67qZlpmhzNAP5g9DTQyqQm1NFC5/zIMPG5OyoAH7Vkxk4PzVrN8XDeT0N/zjTDwxielF8RaFEOKeSMlVWN2dSrWX8GGVoROrDJ3wJImO2oNMq3gS+/P/o6rmImM1YYwljCM/fsCfx3sx4WgVLqheMom8EKJEkeQqCkV+DUFubYWcqLhRv/fLHPZ05IXl22mvOUAP7W5aaQ5RS3MOzrzPLjvYa6zORkNr5oSl07qaN4BUHwshijVJruKByq9UG5OYTqriaLpP60ESXbX76K7ZTRPNMRprTtBYc4JZ6krOrGrHO7H1+Z+xDiga5vapQ+tq3pJshRDFiiRX8cDdWqq9tUSbpLhRoeMoBmx5HB/1Kr21v9NXm1NtXOvaVlbZbuWiWoZvDG14f2M8UxWvBza4uhBCFIQkV1Es5FeidXfUMS3sMEsNPVlu7MnrweloDq2jp3Y35ZSrjLMJY4z2W3YY67LO8Bi/Gh9lWthhHvFzISkti4TMon5XQoiHlSRXUWzcWqK9NeECtPjbgTnZA+ik2c8z2t9orj1KO+0/tNP+wyXVk7XZj/P8RwlcUd1Q0KKrcIH+TSsW1VsSQjykJLmKYu3WhJtbffyDsTmbjM2pkB1LiHY7T2l3UFa5xgTdN4xRw9hibMzn2V14/XuFWuXcSc0yyD1ZIcQDI/1cRYkS0qgCu6a046vhTfl96mOM7NOB9wz9aZ65hFeyXuKAsSq2ioGe2t18Z/cGG3RvsHzpAgYt/50W835j/b5o6TsrSoyIiAgURSEhIeG+9pOWlkbfvn1xdXU17S8oKIjQ0FCrxFkQbdu2zTMJQWkmJVdR4txcmr256tjRtg1PfNSSGkQxRLuFnto/eFRzmkc1S7hoU4bPs7vw1sZ00hQHaQAlip22bdsSHBxslvCaN29OTEwMbm5u97XvVatW8b///Y8//vgDLy+v+95fUfnnn3+YN28eu3btIj4+nqCgIF588UXGjh1b1KHlYXFyzc7O5siRI6aZD/z8/KhZs6aM3yqKzM3JNqfaWGFi9ou8kx3CAJtwntX+SjnlKtN1XzLGJozVhg6szO5MvOrG1I2HcLKzoUGgh1QZi2LH1tYWPz+/+97PmTNnqFGjBrVr17ZCVEXnwIED+Pj48OWXXxIQEMAff/zBiBEj0Gq1jB49uqjDM1PgamGj0cjrr7+Ot7c39evXp0uXLnTp0oX69evj4+PD9OnTTRPqClFUcquNv3yuIQNru/C+4UlaZL7PFP3znDH646akMdrme3bZjWGGzSp8ucrotX+ZqoyFKApDhgxhx44dLF68GEVRUBSFqKioPNXCK1euxN3dnU2bNlG9enUcHR158sknSUtLY9WqVQQFBeHh4cGYMWMwGAxATol4wYIF7Ny5E0VRaNu2bb4xREdH06tXL5ydnXF1deXpp58mLi4OgMTERLRaLfv37wdy8oGnpydNmzY1bZ+b8AoqKCiIOXPm8Nxzz+Hi4kKFChVYtmzZHbd57rnnWLx4MW3atKFSpUo8++yzDB06lLCwsAIf90EpcHKdMmUKy5YtY968eZw9e5bU1FRSU1M5e/Ys8+fPZ9myZUydOrUwYxWiQPzdHGhS0ZNAF3i7V02yFTvWGR6jY9a7jMh6hb+MVbBX9Ay1+YUdduOYa7McH/Uq08IOy33Y0khVQZ8GWakP/lHASccWL15Ms2bNGD58ODExMcTExNw2UaWlpfH++++zbt06tmzZQkREBE888QSbN29m8+bNrF69mk8++YQNGzYAEBYWxvDhw2nWrBkxMTH5JiKj0UivXr24du0aO3bsYNu2bZw9e5aQkBAA3NzcCA4ONk0hd+jQIRRF4a+//iIlJQWAHTt20KZNG4u+mgULFtCwYUP++usvXnrpJUaOHMmJEycs2kdiYmKxnDe2wNXCX3zxBatXr6ZTp05my4OCghgxYgSBgYEMGjSI+fPnWz1IIe7VUw3K066Gn6k7z86T9XgyrBFNlMOM1n5Hc+1R+tls5wntLlYaOnHhYnX83SoXddjCmvRpuH9YuBOG39a0S2DrdNfV3NzcsLW1xdHR8a7VwHq9no8//pjKlXPO0yeffJLVq1cTFxeHs7MzNWvWpF27dmzfvp2QkBA8PT1xdHS8YxVzeHg4hw4dIjIy0pTUv/jiC2rVqsW+ffto1KgRbdu2JSIiggkTJhAREUGHDh04fvw4u3btonPnzkRERDBp0iSLPp6uXbvy0ksvATB58mQWLVrE9u3bqV69eoG2/+OPP1i/fj0//fSTRcd9EApcck1OTqZs2bK3fd3f35/U1NTbvi5EUbl59p6cauPH6B8ykGezX6dv5gz2GB/BXtHzos0mGnzfDnYtAr2UYEXx5OjoaEqsAL6+vgQFBeHs7Gy27PLlywXe57FjxwgICDArLdesWRN3d3eOHTsG5Ex+vmvXLgwGAzt27KBt27amhHvp0iVOnz592yrn26lbt67pZ0VR8PPzM8XdpUsXnJ2dcXZ2platWnm2PXz4ML169WLGjBl07NjRouM+CAUuubZt25YJEyawZs0avLy8zF6Lj49n8uTJFn+wQhQFfzcHutdzIDUrm2lhCiFZ03lc+w8LPL/FPfkU/DoT9iyDtlMgeABopVF9iaZzJGHUMVxdXB78ZOk6R+vv8pbGo4qi5LvM2m1gWrduTXJyMgcPHmTnzp3MmTMHPz8/5s2bR7169ShbtixVq1a1aJ93ivvTTz8lPT093/WOHj3K448/zogRI3j99dfv410VngJfNZYuXUrXrl3x9/enTp06+Pr6AhAXF8ehQ4eoWbMmmzZtKrRAhbA28xGgHsfdZQoc+gZ+exsSz8OPY2D3B/D4G/BId8ido1aULIqSk+RsneBBJ1cL2NramhohPWg1atTg/PnznD9/3lR6PXr0KAkJCdSsWRMAd3d36tatywcffIBOp+ORRx7Bx8eHkJAQNm3aZPH91rspV65cvsuPHDnCY489xuDBg5k9e7ZVj2lNBT7TAgIC+Oeff/jhhx/o0aMHFSpUoEKFCvTo0YMff/yRv/76y6KWYrk+/PBDgoKCsLe3p0mTJuzdu/eO64eGhlK9enUcHBwICAjglVdeISMjw+LjCgG3TPiu0UK9Z2D0fug0Bxw8If4krH8WPusAUb+btpOBKIS1BQUFsWfPHqKiooiPj3+gvS/at29PnTp1GDBgAAcPHmTv3r0MGjSINm3a0LBhQ9N6bdu2Zc2aNaZE6unpSY0aNVi/fr3Vk2t+Dh8+TLt27ejYsSPjx48nNjaW2NhYrly5UujHtpRF9V0ajcbUBcca1q9fz/jx41m6dClNmjQhNDSUTp06ceLECXx8fPKsv3btWqZMmcLnn39O8+bNOXnyJEOGDEFRFBYuXGiVmIRAZw/NRkH9Z+GPJbD7Q7iwD1Z2haod2eI3gpd+zZKBKIRVTZgwgcGDB1OzZk3S09OJjIx8YMdWFIXvv/+el19+mdatW6PRaOjcuTNLliwxW69NmzaEhoaa3QJs27Yt//zzzwO5LbhhwwauXLnCl19+yZdffmlaHhgYSFRUVKEf3xKKqhasrfi///5boB3efIP6bpo0aUKjRo344IMPgJzm4AEBAbz88stMmTIlz/qjR4/m2LFjhIeHm5a9+uqr7Nmzh127dhXomElJSbi5uZGYmIirq2uBYxUli16vZ/PmzXTt2vX+BzhJjoUd78DBVWDMxqgqfGtswaLsp7igeqNVFHZNaSeDUBQTGRkZREZGUrFiRezt7TEajSQlJeHq6vrg77mKInfz+aDVas2uC4WZDwpccg0ODkZRFHJzsfLf/aebc7OiKAW+Z5CVlcWBAwfM+sZqNBrat2/P7t27892mefPmfPnll+zdu5fGjRtz9uxZNm/ezMCBA297nMzMTDIzb8w9lpSUBORcfPV6fYFiFSVP7ndrle/Yvgx0mg+NRpDw00y8o3+mr3YX3TV/sszQnQ+ye3MmLgkvR2n4VBzo9XpUVcVoNGI0Gk3XqNxl4uGSew7o9XrT92/V68NtFPhqcHMVhaqq1K5dm82bNxMYGHhPB46Pj8dgMJgaRuXy9fXl+PHj+W7Tv39/4uPjadmyJaqqkp2dzYsvvsi0adNue5y5c+cya9asPMu3bt2Ko6P1W/KJ4mXbtm1W3V+Ccz82ZjZnks06WmqP8LLNd/TQ7ObkrsFsPlayh5YrLWxsbPDz8yMlJYWsrCzT8uTk5CKMShSVrKws0tPT2blzJ9nZ2cCN60JaWlqhHbfAyfXWJKooCuXLl7/n5HovIiIimDNnDh999BFNmjTh9OnTjB07lrfeeovp06fnu83UqVMZP3686XlSUhIBAQF07NhRqoVLMb1ez7Zt2+jQoYPVx73WVbjAoO8r0cGwn5m6VQRp4giKfgej61MY2r8JTt5WPZ6wTEZGBufPn8fZ2Rl7e3tUVSU5ORkXFxdTjZt4eGRkZODg4EDr1q3RarVm14XcmszCUGT1WF5eXmi1WtPYlbni4uJuO4rI9OnTGThwIM8//zwAderUITU1lREjRvDaa6/lez/Fzs4OOzu7PMt1Op1MNvAQKIzvuX/Tiv+N+tQMxeUlOLAA9nyC5vA3aE5vg45vQfCzxbrbR2lmMBhQFAWNRoNGozFVBeYuEw8XjUZj6gus1WqBG9eFwswBRXam2dra0qBBA7PGSUajkfDwcJo1a5bvNmlpaXl+OXI/rAK2yxLCKnK78Pj5+ECX+TA8HPzqQEYC/PAyrOwGl2/c3pCuO0I8XO6r5Hq/VSzjx49n8ODBNGzYkMaNGxMaGkpqaipDhw4FYNCgQZQrV465c+cC0KNHDxYuXEj9+vVN1cLTp0+nR48epiQrRJEo1wCGR8CepbB9NkT/AUtbQstxfOMYwuTvT0rXHSEeIgVOrvXr1zdLpunp6fTo0QNbW1uz9Q4ePFjgg4eEhHDlyhXeeOMNYmNjCQ4OZsuWLaZGTtHR0WYl1ddffx1FUXj99de5ePEi3t7e9OjRo1iP0iEeIlobaD4aavaEzRPh5BbY+S6NjF/QTHmO39U6GFWYFnaY1tW8peuOEKVYgZNrr169zJJrr169rBLA6NGjbzvJbe70RrlsbGyYMWMGM2bMsMqxhSgU7hWg3zo49iOZP04gKD2ONbZz2WhoxSz9IJJUJ6Li0yS5ClGKFTi5Tpo0SbquCFFQigI1e3Ldqym/LBnFQM02+mr/RxPNMSbqRxLk9VhRRyiEKEQFbtDk5eVF9+7dWbZsGbGxsYUZkxClhp+PD/Y9F/C0fibnjD6UV+JZa/s2/nvnQXbW3XcgxH8iIiJQFIWEhASr7TMoKIjQ0FCr7U/cUODkeuzYMTp16sTXX39NUFAQTZo0Yfbs2Rw6dKgw4xOixAtpVIElk18gtv+vpNXuj4IKv4fCp4+btSgWoiQ6ceIE7dq1w9fXF3t7eypVqsTrr7/+0I+AV+DkGhgYyMsvv8yvv/5KXFwc48aN49ChQ7Rq1YpKlSoxbtw4fvvttyKbMkmI4szfzYEmjwTi+OTHEPJlzow7sf/CsjY5c8f+15VMuuyIkkan0zFo0CC2bt3KiRMnCA0NZfny5Q9925h76ufq5uZGv379WLduHVeuXGHp0qUYDAaGDh2Kt7c3a9assXacQpQeNXrAS7uh8uOQnQE/T4Q1T/L9/w7QYt5v9F++hxbzfmP9vuiijvShF5say96YvcSmFv6tsMzMTMaMGYOPjw/29va0bNmSffv25Vnv999/p27dutjb29O0aVMOHz5seu3cuXP06NEDDw8PnJycqFWrFps3by5wDIqi8Omnn/LEE0/g6OhI1apV+eGHH+64TaVKlRg6dCj16tUjMDCQnj17MmDAAP73v/8V/M2XQvc9iIROp6Njx44sWbKEc+fOER4eTrVq1awRmxCll4sfPLsRurwLNvZw+lda/dqTDkrOxTS3y46UYItO2KkwOm3sxLCtw+i0sRNhp8IK9XiTJk1i48aNrFq1ioMHD1KlShU6derEtWvXzNabOHEiCxYsYN++fabuiLlVsKNGjSIzM5OdO3dy6NAh5s+fj7Ozs0VxzJo1i6effpp///2Xrl27MmDAgDwx3Mnp06fZsmXLA5nftTi7p0EkEhIS2Lt3L5cvXzabZUJRlDvOUCOEuImiQJMRULE1KV8NxfP6UT6xXcTq7Pa8lT2QLFUnXXaKSGxqLLN2z8Ko5lzfjKqRWbtn0bxsc/yc8h+e9X6kpqby8ccfs3LlStN82cuXL2fbtm189tlnTJw40bTujBkz6NChAwCrVq2ifPnyfPvttzz99NNER0fTt29f6tSpA+SUKi01ZMgQ+vXrB8CcOXN4//332bt3L507d77jds2bN+fgwYNkZmYyYsQI3nzzTYuPXZpYnFx//PFHBgwYQEpKCq6urmZ9XyW5CnEPfB4heeAWvlw4ihHaTQy0+ZU6mrO8rB9HkJd0fysK0UnRpsSay6gaOZ98vlCS65kzZ9Dr9bRo0cK0TKfT0bhxY44dO2a27s3Dw3p6elK9enXTOmPGjGHkyJFs3bqV9u3b07dvX4vm2AbzObmdnJxwdXXl8uXLANSqVYtz584B0KpVK37++WfTuuvXryc5OZl//vmHiRMn8t577zFp0iSLjl2aWFwt/Oqrr/Lcc8+RkpJCQkIC169fNz0sqToQQtzg7+mGR685DNNP4rrqTLDmLL86T8f/8q6iDu2hVMG1AhrF/PKoUTQEuAQUUUQF8/zzz3P27FkGDhzIoUOHaNiwIUuWLLFoH7cOZq8oiqmGcvPmzfz999/8/ffffPrpp2brBQQEULNmTfr168e8efOYOXPmQ93A1eLkevHiRcaMGSMDSghhZSGNKjBn8iuc7fMzWb71sNMnwpqn4LfZYHx4L1JFwc/JjxnNZpgSrEbRMKPZjEIptQJUrlwZW1tbfv/9d9MyvV7Pvn37qFmzptm6f/75p+nn69evc/LkSWrUqGFaFhAQwIsvvkhYWBivvvoqy5cvt1qcgYGBVKlShSpVqlCuXLnbrmc0Gs0mJ38YWVwt3KlTJ/bv339PdflCiDvzd3PAv15dqLUNtkyF/Z/Bznfgwj7o+xkx2Y5ExqdS0ctJ7sUWsj5V+9C8bHPOJ58nwCWg0BIr5FS/jhw5kokTJ+Lp6UmFChV45513SEtLY9iwYWbrvvnmm5QpUwZfX19ee+01vLy86N27NwDjxo2jS5cuVKtWjevXr7N9+3azxFsY1qxZg06no06dOtjZ2bF//36mTp1KSEjIQz2tZ4GS681Nsbt168bEiRM5evQoderUyfPh9ezZ07oRCvEwsrGD7gshoAn8OBbObidtSTNGJb3EQWNVmV3nAfFz8ivUpHqzefPmYTQaGThwIMnJyTRs2JBffvkFDw+PPOuNHTuWU6dOERwczI8//miaQMVgMDBq1CguXLiAq6srnTt3ZtGiRYUat42NDfPnz+fkyZOoqkpgYCCjR4/mlVdeKdTjFneKWoCJUAs6wbCiKMW+jj0pKQk3NzcSExNxdXUt6nBEIdHr9WzevJmuXbuW/L+e446Sve5ZbK6fIUvV8mb2IL40dECrKOya0k5KsLfIyMggMjKSihUrYm9vj9FoJCkpCVdXV5ks/SF08/mg1WrNrguFmQ8KdKYZjcYCPYp7YhWiRPKtyYGOG9lsaIytYuBt3QresvkcRdUTFZ9W1NEJIfIhf8YJUQJUKOvH6OyxzNc/g1FVGGjzK6ts36Gic2ZRhyaEyIckVyFKAH83B+b2qcsyYy9e0L9CqmpHC81h/NZ3gysnARmXWIji5J5GaBJCPHghjSrQupo3UfFNSFO74/TjYLh2Fj5tT0S9d3jufy4YVaSxkxDFgJRchShB/N0caFa5DN5VGsDw7RDQFDITabXnRQZptgCqjEssRDEgyVWIksrZGwb/wOXKfdEqKjN1XzDb5nO0GDCoqjR2EqIISXIVoiSzscPQ4wNmZw/AqCoMsAlnmW4hzkqmjEssRBG6p+Tarl07hgwZYrZs8ODBPPbYY9aISQhhAX93R6r0msLI7PGkq7Y8rv2LnT7v4a9NLurQhHho3VNyDQoKomzZsmbLypUrR2BgoFWCEkJYJqRRBWZOmsjpLl9htPfEM/EIfNYB4k9LK2IhisA9JdcVK1YwZ84cs2Vz5sxhxYoVVglKCGE5fzcH6jRtj+b5beARBNejyPzkcUbPX0r/5XtoMe831u+LLuowhQUiIiJQFIWEhIT72k9aWhp9+/Y1TROakJBAUFAQoaGhVomzINq2bcu4ceMe2PGKmtxzFaK08aoCw34lyzcYO30Ca3Sz6aDZL62Ii7n8kk/z5s2JiYnBzc3tvva9atUq/ve///HHH39YZX9F5erVq3Tu3JmyZctiZ2dHQEAAo0ePJikpqahDy6NA/Vzff//9Au9wzJgx9xyMEMJKnL052G41KWsG0V77F0t1i3gtexjrDI8RFZ8m4xGXELa2tvj53f/EAWfOnKFGjRrUrl3bClEVHY1GQ69evXj77bfx9vbm9OnTjBo1imvXrrF27dqiDs+cWgBBQUFmDycnJ1VRFNXDw0P18PBQFUVRnZyc1IoVKxZkd0UqMTFRBdTExMSiDkUUoqysLPW7775Ts7KyijqUInMpIU2tMuV7dc1rvVR1hquqznBV5732gnopIa2oQytU6enp6tGjR9X09HRVVVXVYDCo169fVw0GQxFHdnuDBw9WAbNHZGSkun37dhVQr1+/rqqqqq5YsUJ1c3NTf/zxR7VatWqqg4OD2rdvXzU1NVVduXKlGhgYqLq7u6svv/yymp2draqqqrZp08Zsv23atFFVVVUDAwPVRYsWmWI4d+6c2rNnT9XJyUl1cXFRn3rqKTU2NlZVVVVNSEhQNRqNum/fPlVVcz5TDw8PtUmTJqbtV69erZYvX/6277FNmzbq2LFjTc8DAwPV2bNnq0OHDlWdnZ3VgIAA9ZNPPrH4s1u8ePEdj3vz+XDrdaEw80GBqoUjIyNNj9mzZxMcHMyxY8e4du0a165d49ixYzz66KO89dZbhfMXgBDCYv5uDrzdJ5jphuF8lJ0zFeRkm6/w3zsPVPWhaeikqirG9HSMaWkP/KHefdIxABYvXkyzZs0YPnw4MTExxMTEEBAQkO+6aWlpvP/++6xbt44tW7YQERHBE088webNm9m8eTOrV6/mk08+YcOGDQCEhYUxfPhwmjVrRkxMDGFhYXn2aTQa6dWrF9euXWPHjh1s27aNs2fPEhISAoCbmxvBwcFEREQAcOjQIRRF4a+//iIlJQWAHTt20KZNG4u+mwULFtCwYUP++usvXnrpJUaOHMmJEycKvP2lS5cICwuz+LgPgsXDH06fPp0NGzZQvXp107Lq1auzaNEinnzySQYMGGDVAIUQ9+7GkIlNSYqsi+uut+H3UM6cv0CnU73JVjWlfrhENT2duHaPEVcEx65+8ACK4937G7u5uWFra4ujo+Ndq4H1ej0ff/wxlStXBuDJJ59k9erVxMXF4ezsTM2aNWnXrh3bt28nJCQET09PHB0d71jFHB4ezqFDh4iMjDQl9S+++IJatWqxb98+GjVqRNu2bYmIiGDChAlERETQoUMHjh8/zq5du+jcuTMRERFMmjTJos+na9euvPTSSwBMnjyZRYsWsX37drP8kp9+/frx/fffk56eTo8ePfj0008tOu6DYHGDppiYGLKzs/MsNxgMxMUVxekrhLiT3CETXdtPhB6LUVGoHL2BUJsP0JEtDZ1KGEdHR1NiBfD19SUoKAhnZ2ezZZcvXy7wPo8dO0ZAQIBZablmzZq4u7tz7NgxANq0acOuXbswGAzs2LGDtm3bmhLupUuXOH36NG3btrXovdStW9f0s6Io+Pn5meLu0qULzs7OODs7U6tWLbPtFi1axMGDB/n+++85c+YM48ePt+i4D4LFJdfHH3+cF154gU8//ZRHH30UgAMHDjBy5Ejat29v9QCFEFbUYAgnEzVU3DGO7to/cSadF/SvkKnaltqGToqDA77bf8PVxeWBT5auOFj/89TpdObHUJR8lxmNRqset3Xr1iQnJ3Pw4EF27tzJnDlz8PPzY968edSrV4+yZctStWpVi/Z5p7g//fRT0tPT813Pz88PPz8/HnnkETw9PWnVqhXTp0/H39//Pt6hdVmcXD///HMGDx5Mw4YNTW84OzubTp06FcuiuRDCnGuDpxj+6wU+1oXSVvsPK3iHEdkTS+1wiYqioHFwQOPo+MCTqyVsbW0xGAxFcuwaNWpw/vx5zp8/byq9Hj16lISEBGrWrAmAu7s7devW5YMPPkCn0/HII4/g4+NDSEgImzZtsvp9z3LlyhVovdxknJlZvOY2tji5ent7s3nzZk6ePMnx48cBeOSRR6hWrZrVgxNCWJ+/mwNdn3iWId/a85nuHZprj7Ldawnedu2JSYTI+FQqejmVylJscRYUFMSePXuIiorC2dkZT0/PB3bs9u3bU6dOHQYMGEBoaCjZ2dm89NJLtGnThoYNG5rWa9u2LUuWLOHJJ58EwNPTkxo1arB+/Xo+/PDDQo9z8+bNxMXF0ahRI5ydnTly5AgTJ06kRYsWBAUFFfrxLXHPf8ZVq1aNnj170rNnz/tKrB9++CFBQUHY29vTpEkT9u7de8f1ExISGDVqFP7+/tjZ2VGtWjU2b958z8cX4mEU0qgCiyePJKrLGox2rnhfP8jVpd3oPO9HGc2piEyYMAGtVkvNmjXx9vYmOvrBff6KovD999/j4eFB69atad++PZUqVWL9+vVm67Vp0waDwWB2b7Vt27Z5lhUWBwcHli9fTsuWLalRowavvPIKPXv2ZNOmTYV+bEspakHbit/kwoUL/PDDD0RHR5OVlWX22sKFCwu8n/Xr1zNo0CCWLl1KkyZNCA0N5ZtvvuHEiRP4+PjkWT8rK4sWLVrg4+PDtGnTKFeuHOfOncPd3Z169eoV6JhJSUm4ubmRmJiIq6trgWMVJYter2fz5s107do1z/0acYtLf2H84gk0Gdc5ZAxiYNZUEnBBqyjsmtKuxJVgMzIyiIyMpGLFitjb22M0GklKSsLV1bVYVwuLwnHz+aDVas2uC4WZDyyuFg4PD6dnz55UqlSJ48ePU7t2baKiolBV1dTAqaAWLlzI8OHDGTp0KABLly7lp59+4vPPP2fKlCl51v/888+5du0af/zxh+mCWdyqAoQoccrW59/2X1L+x37U0UTxle3bDMh6jWuqa6lt5CREYbM4uU6dOpUJEyYwa9YsXFxc2LhxIz4+PgwYMIDOnTsXeD9ZWVkcOHCAqVOnmpZpNBrat2/P7t27893mhx9+oFmzZowaNYrvv/8eb29v+vfvz+TJk9Fqtfluk5mZaXajO3cMSr1ej16vL3C8omTJ/W7lOy4Yz6B69NNP50vdbGpozrPGdjbP6l+jnJttifsM9Xp9zsARRiNGo9E0kEPuMvFwyT0H9Hq96ft/ENcHi5PrsWPH+Oqrr3I2trEhPT0dZ2dn3nzzTXr16sXIkSMLtJ/4+HgMBgO+vr5my319fU0NpW519uxZfvvtNwYMGMDmzZs5ffo0L730Enq9nhkzZuS7zdy5c5k1a1ae5Vu3bsWxAJ27Rcm2bdu2og6hxGgQ5Ef/yNdYa5uTYDc5zebvHQb+snEp6tAsYmNjg5+fHykpKWa3rZKTZX7bh1FWVhbp6ens3LnTNEZD7nUhLS2t0I5rcXJ1cnIynbD+/v6cOXPG1ME3Pj7eutHdwmg04uPjw7Jly9BqtTRo0ICLFy/y7rvv3ja5Tp061ayDcVJSEgEBAXTs2FHuuZZier2ebdu20aFDB7nnWkBdgZjEDC5E1cHrt0H4p53HL+5jLvb4isg0OwLLOOLvZl/UYd5VRkYG58+fx9nZGXt7e1RVJTk5GRcXFxRFKerwxAOWkZGBg4MDrVu3RqvVml0XCnM2HYuTa9OmTdm1axc1atSga9euvPrqqxw6dIiwsDCaNm1a4P14eXmh1WrzjOoUFxd32yG6/P390el0ZlXANWrUIDY2lqysLGxtbfNsY2dnh52dXZ7lOp1OLroPAfmeLVPBS0cFr2YQ+BOs7IZy+TCJy7vzctY0khSXEjFMosFgQFGUnP6tGo2pKjD3uXi45J4LN+eO3OtCYV4bLD7TFi5cSJMmTQCYNWsWjz/+OOvXrycoKIjPPvuswPuxtbWlQYMGhIeHm5YZjUbCw8Np1qxZvtu0aNGC06dPm903OXnyJP7+/vkmViHEPfKuxpW+G7iiulFLc44vbefirKaUiGEScy+gt/ZkEA+n3KrfB/1HtsUl10qVKpl+dnJyYunSpfd88PHjx5tGe2rcuDGhoaGkpqaaWg8PGjSIcuXKMXfuXABGjhzJBx98wNixY3n55Zc5deoUc+bMkTlkhSgEp9RyvJH1Gl/Zvk1tTRSrbN/h2aypxb4FsY2NDY6Ojly5csV0Qc3KyiIjI0NKrg8RVVVJS0vj8uXLuLu7o9VqH2iDNouTK+QM5LBhwwbOnDnDxIkT8fT05ODBg/j6+hZ4yCqAkJAQrly5whtvvEFsbCzBwcFs2bLF1MgpOjra7JchICCAX375hVdeeYW6detSrlw5xo4dy+TJk+/lbQgh7qCilxNnKc+zWdNYZ/s29TWnWWH7Lg5KM/44E19sR3FSFAV/f38iIyM5d+4cqqqSnp6Og4OD3HN9CLm7u1tlwnlLWTyIxL///kv79u1xc3MjKiqKEydOUKlSJV5//XWio6P54osvCitWq5BBJB4OMoiEdazfF820sMPU4AxrbefgqqSxy1ibYVkT0Cu2xfoerNFoJCsrC71ez86dO2ndurWcCw+ZW9vo3HpdKFaDSIwfP54hQ4bwzjvv4OJyo4l+165d6d+/v1WDE0IUrRvzwTYh7noNtD/2p6XmMB/pFvOi/hWmhR2mdTXvYlmC1Wg02Nvbo9Vqyc7Oxt7eXpKreGAsvgGxb98+XnjhhTzLy5UrR2xsrFWCEkIUH7nzwV7xCOa5rElkqDoe1/7FIt2HqKqBqPjC6ysoREllcXK1s7PLt2/QyZMn8fb2tkpQQojip6KXE/uowQj9eLJULd21e5ir+wxHncIfZ+KLfStiIR4ki5Nrz549efPNN03DRimKQnR0NJMnT6Zv375WD1AIUTz4uzkwt08dfleDGaN/GYOqEKLdzv5lL9F/+Z8yk44QN7E4uS5YsICUlBR8fHxIT0+nTZs2VKlSBRcXF2bPnl0YMQohiomQRhXYNaUdg4eN4WLrdwAYZvMzY7TfYlQpEf1ghXgQLG7Q5ObmxrZt29i1axf//vsvKSkpPProo7Rv374w4hNCFDP+bg74uznwB334/LdDzNR9wXjdBpJwZKWhc7HvByvEg3BP/VwBWrZsScuWLa0ZixCiBKno5cQXxs646NN4VbeBmbovSMSFIK/Hijo0IYpcgZPr+++/X6D1ZLQkIR4Oufdgp4WBZ3YyQ21+YYHtJyRENuMPl2bFdpAJIR6EAifXRYsWmT0/f/48/v7+2Njc2IWiKJJchXiImPrBXmlM+v7JOBzfiMO3z7Egayp/Ub1YDzIhRGEqcHKNjIw0e+7i4sKOHTvMxhoWQjx8cu/BxniG8ueRM7TT/s3ntu/ydNYbTAtTiu0gE0IUJhnFWghhFZHXsxipH8t+YzXclDRW2c7HT70sg0yIh5IkVyGEVVT0ciJLsWNY1gROGsvhp1znC9v5uBgTZJAJ8dCR5CqEsIrcBk4piguDsqZwSS1DZc0l9KufZtjyHTLIhHioFPie661DHiqKQkpKSp7lMtOMEA+vGwP9p5GUVgOHDb2prznNR7rFDNe/WqwH+hfCmgpccnV3d8fDw8P0SElJoX79+qbnua8LIR5uuQP9X3OsyLCsCaSrtrTT/sNsm88xqEa5ByseCgUuuW7fvr0w4xBClDIVvZz4m2q8rH+ZT3QLCbGJIA5PrqY+SkxiupReRalW4OTapk2bwoxDCFHK3BhkQuH17OeYq/uMMTZhTF3vwRjj49IHVpRq9zz8oRBC3E3uPdgDUcG8/811xtiE8bbN58TpPaQPrCjVpLWwEKJQ+bs54Olsy8LsvqzPbotWUflAt4RanOanf2Oki44olSS5CiEKXUUvJzSKwmvZzxFhqIejksnntu/yxeYI6aIjSiVJrkKIQpd7/1VVdIzSj+GwMQgvJYmVuvm4qUkyD6wodSS5CiEeiNyJ1l/p9ihDsyZyQfWikiaW5bYLsVEzpYuOKFUK1KCpT58+Bd5hWFjYPQcjhCjd/N0c6FbXnzmbPRicNZkw2xk01JwkVPcRV1MaSRcdUWoUqOTq5uZmeri6uhIeHs7+/ftNrx84cIDw8HDc3NwKLVAhROmQW0UcRXlGZL1Klqqli3Yvl76ZJPdfRalRoJLrihUrTD9PnjyZp59+mqVLl6LVagEwGAy89NJLMvShEKJAbu6iM+nra4TqPmSEzU9cUL2ki44oFSy+5/r5558zYcIEU2IF0Gq1jB8/ns8//9yqwQkhSq/cLjrfGVrwjv5pAGbYfEFbZb900RElnsXJNTs7m+PHj+dZfvz4cYxGo1WCEkI8HHK66MBHhl6s+68P7BLdB3y3+SepIhYlmsXJdejQoQwbNoyFCxeya9cudu3axYIFC3j++ecZOnRoYcQohCilcu+/ahUNr2c/x05Dnf/6wL6Hv3pFuuiIEsvi4Q/fe+89/Pz8WLBgATExMQD4+/szceJEXn31VasHKIQo3XLvv/70bwwv/TSWb5RZ1NCc5zPb93gqawY//RtDt7r+cg9WlCgWl1w1Gg2TJk3i4sWLJCQkkJCQwMWLF5k0aZLZfVghhCio3C46aYojz2VNIk515xHNeT7ULWbeT4ekiliUOPc0iER2dja//vorX331FYqiAHDp0iVSUlKsGpwQ4uGRW0V8WfFiWNYE0lQ7WmsP8ZbNCoyqKlXEokSxuFr43LlzdO7cmejoaDIzM+nQoQMuLi7Mnz+fzMxMli5dWhhxCiEeAjeqiB/h5Z+vs0y3kH4224lS/fjE0IMDUdfpXk+qh0XxZ3HJdezYsTRs2JDr16/j4HDjJH/iiScIDw+/pyA+/PBDgoKCsLe3p0mTJuzdu7dA261btw5FUejdu/c9HVcIUfzkVhFvVxvwZvYgAKbqvqKLZg9j1v0l1cOiRLA4uf7vf//j9ddfx9bW1mx5UFAQFy9etDiA9evXM378eGbMmMHBgwepV68enTp14vLly3fcLioqigkTJtCqVSuLjymEKN5yq4hXGzqxIrsTAIt0H1GPU1I9LEoEi5Or0WjEYDDkWX7hwgVcXFwsDmDhwoUMHz6coUOHUrNmTZYuXYqjo+MdB6QwGAwMGDCAWbNmUalSJYuPKYQo/kIaVeD9/vV5K3sgvxrqY6/oWWa7AH/iZJAJUexZnFw7duxIaGio6bmiKKSkpDBjxgy6du1q0b6ysrI4cOAA7du3vxGQRkP79u3ZvXv3bbd788038fHxYdiwYZaGL4QoQRoEeoCiYYz+ZQ4bg/BWklihe5f3f9ovLYhFsWZxg6YFCxbQqVMnatasSUZGBv379+fUqVN4eXnx1VdfWbSv+Ph4DAYDvr6+Zst9fX3zHQUKYNeuXXz22Wf8/fffBTpGZmYmmZmZpudJSUkA6PV69Hq9RfGKkiP3u5XvuGTzcrTh7V41ef37owzLmsB3dm9QVXORj3ShDNFPZmrYIZpV9MDfzf62+5BzQeS69VwozHPC4uRavnx5/vnnH9avX88///xDSkoKw4YNY8CAAWYNnApDcnIyAwcOZPny5Xh5eRVom7lz5zJr1qw8y7du3Yqjo6O1QxTFzLZt24o6BHGfnIAZ9eGvq24Mi57AN7azaKk9wtvq50zJHs7n32+nvpd61/3IuSBy5Z4LaWmFN4ewoqrq3c/KQpKVlYWjoyMbNmwwa/E7ePBgEhIS+P77783W//vvv6lfv77ZYBW54xlrNBpOnDhB5cqVzbbJr+QaEBBAfHy8zOJTiun1erZt20aHDh3Q6XRFHY6wgpjEDNou2Elb5SDLdQvQKirz9M+wzNiTt3vV5KkG5fPdTs4FkevWcyEpKQkvLy8SExOtng8sLrlqtVpat27Nxo0b8fT0NC2Pi4ujbNmy+TZ2uh1bW1saNGhAeHi4KbkajUbCw8MZPXp0nvUfeeQRDh06ZLbs9ddfJzk5mcWLFxMQEJBnGzs7O+zs7PIs1+l08ov2EJDvufSo4KVjbp86TN0Ib2YPYpZuFVN064jO8mH69wrtavjdcYhEORdErtxzoTDPB4uTq6qqZGZm0rBhQ3788Udq1apl9pqlxo8fz+DBg2nYsCGNGzcmNDSU1NRU0yQAgwYNoly5csydOxd7e3tq165ttr27uztAnuVCiNInpFEFnOxsGL0WgpRYhtr8wkLdx/TLKkNUfBMZf1gUGxa3FlYUhY0bN9KjRw+aNWtmVnWbOxSiJUJCQnjvvfd44403CA4O5u+//2bLli2mRk7R0dGmCQKEEKJBoAcaBbMuOp/aLqCy7kpRhyaEyT2VXLVaLYsXL6ZWrVqEhITw+uuv8/zzz99zEKNHj863GhggIiLijtuuXLnyno8rhCh5cgeYmBZ2mDH6l/laeZPamij4YSAM2woOHkUdohCWJ9ebjRgxgqpVq/LUU0+xc+dOa8UkhBB3lDsGcVR8Gt729WFdN4g/CV8PggEbiUk1EBmfSkUvJ7wc7+syJ8Q9sfisCwwMNGut265dO/7880969Ohh1cCEEOJO/N0c/rvHWgYGfA2fd4bInZxdOZz2Z57CqCpoFHi7V02cijpY8dCx+J5rZGQkZcqUMVtWpUoV/vrrL86ePWu1wIQQosD86sCTK1AVDZUufMeLmpy2IEYVXv/+KAmZd9leCCu7p/lc82Nvb09gYKC1dieEEJap1pHIRjMAmKT7mh6aP4CcBHslw/LGlkLcjwIlV09PT+Lj4wHw8PDA09Pztg8hhCgqDi1e4LPsLgC8p/uEhspxFMBWU2Rj5YiHVIHuuS5atMg0483Ng/YLIURx4u/mgEvPuWz98QodtftZbruQJ7JmseiwH15VLtC/acWiDlE8JAqUXAcPHpzvz0IIUdw83bgi/3qv5p8V3ainOctK3Tv0yZrFa98dpVY5d+oFSFcdUfgKVC2clJRU4IcQQhS1FKMtw7Imct7oTZAmjuW2C7Ali94f/SHT1IkHokDJ1d3dHQ8Pjzs+ctcRQoiiVtHLiWuKG0P0k0hUHWmgOcUi3UegGpkWdlgmWheFrkDVwtu3by/sOIQQwmpyR3GauhFGZL3KF7Zz6ardyzR1LbOzn+VA1HW615NxiEXhKVBybdOmTWHHIYQQVhXSqAKP+LnQ+yOYqH+B920/ZLjNZi6qXoxZB6lZ2YQ0qlDUYYpS6p7HBUtLSyM6OpqsrCyz5XXr1r3voIQQwhrqBXgwu1dNpn2nUk5/lcm6dbxhs5oYfRmmboRH/FykgZMoFBYn1ytXrjB06FB+/vnnfF+3ZD5XIYQobE81KM/JI4f4+FQPyilXeNYmnMW6D+if9Rq9P4J5fepICVZYncUjNI0bN46EhAT27NmDg4MDW7ZsYdWqVVStWpUffvihMGIUQoj7UtFFRaMozMgectM0de8RRAxTNx6SBk7C6ixOrr/99hsLFy6kYcOGaDQaAgMDefbZZ3nnnXeYO3duYcQohBD3xd0uZwB/FS0v61/mb2MlPJUUVunm4UkiK3ZFFXWIopSxOLmmpqbi4+MD5AyFeOVKzgTFderU4eDBg9aNTgghrOSpBuX5dlRz0rFnWNZEooy+VNBcYYXtfNb+7wj/nL9e1CGKUsTi5Fq9enVOnDgBQL169fjkk0+4ePEiS5cuxd/f3+oBCiGEtdQL8GBEq4pcxY3B+snEq67U0UTxkS6Upz7aKQNMCKuxOLmOHTuWmJgYAGbMmMHPP/9MhQoVeP/995kzZ47VAxRCCGsa2rIiGgXOqX4My5pAmmpHa+0h3rVZymth/8r9V2EVFrcWfvbZZ00/N2jQgHPnznH8+HEqVKiAl5eXVYMTQghruzHAxCH+UaswUj+OT3Xv0Uv7B9dUF6KuNP1vEnYh7t19z+fq6OjIo48+KolVCFFihDSqwLejmqMosMNYj1f1LwIw1OYXap/5pIijE6WBxSVXVVXZsGED27dv5/LlyxiNRrPXw8LCrBacEEIUlnoBHszrU4dpYYf5wdiCMvoUZuhW4bL7HfD0hUbPF3WIogSzOLmOGzeOTz75hHbt2uHr64uiKIURlxBCFLqQRhVoXc2bqPg0grwegwO+sPMd+GkCOHhC7T5FHaIooSxOrqtXryYsLIyuXbsWRjxCCPFA+bs53LjH2m4apMXD/s8hbAQ4uEPlx4o0PlEyWXzP1c3NjUqVKhVGLEIIUbQUBbq+B7WeAKMe1j1L/Inf+eNMvLQiFhaxOLnOnDmTWbNmkZ4uJ5oQohTSaOGJT6BSW9Cnoln7NNM/DaPFvN+kH6woMIuT69NPP83169fx8fGhTp06PProo2YPIYQo8WzsiO3yGf/8N0ziatu5lOWKTLQuCszie66DBw/mwIEDPPvss9KgSQhRap1NglFZk/jG9k2qaC6xRjebkKzpRMWnST9YcVcWJ9effvqJX375hZYtWxZGPEIIUSxU9HIiUXFlQNY01tu+RZAmjrW2c3ByaAOUKerwRDFncbVwQEAArq6uhRGLEEIUG7kjOcUrZeif9RoXVC8qaWLw/fZpSI0v6vBEMWdxcl2wYAGTJk0iKiqqEMIRQojiI6RRBXZNaceC4d2xfW4TuJSFK8dhdW9Iu1bU4Yli7J7GFk5LS6Ny5co4Ojqi0+nMXr92TU44IUTpcaMfbBkY/AOs6Aqxh+DLvsT2XsfZZC0VvZzkPqwwY3FyDQ0NLYQwhBCiBPCqmpNgV3aDSwe58EF3ns+aQoZiz9w+dQhpVKGoIxTFhEXJVa/Xs2PHDqZPn07FihULKyYhhCi+fGpw5Yn12H7Zk4aak3xu+y7PZU1kWthhWlfzlhKsACy856rT6di4caPVg/jwww8JCgrC3t6eJk2asHfv3tuuu3z5clq1aoWHhwceHh60b9/+jusLIYS1ndJUZGDWVJJUB5pqjrHKdh4OaipR8WlFHZooJixu0NS7d2++++47qwWwfv16xo8fz4wZMzh48CD16tWjU6dOXL58Od/1IyIi6NevH9u3b2f37t0EBATQsWNHLl68aLWYhBDiTip6OXGYygzMmkqi6kgjzUm+tJ1LReesog5NFBMW33OtWrUqb775Jr///jsNGjTAycnJ7PUxY8ZYtL+FCxcyfPhwhg4dCsDSpUv56aef+Pzzz5kyZUqe9desWWP2/NNPP2Xjxo2Eh4czaNAgC9+NEEJYLrebzrQwhf5Zr7Padg7BmjPw7VMw8Htwkn6wDzuLk+tnn32Gu7s7Bw4c4MCBA2avKYpiUXLNysriwIEDTJ061bRMo9HQvn17du/eXaB9pKWlodfr8fT0zPf1zMxMMjMzTc+TkpKAnPvHer2+wLGKkiX3u5XvWBTWudAn2J9mFT2IvtaATLUx6g/9UGIPoa7sysUeXxGZ7kRgGUf83eytelxx7249Fwrz+mBxco2MjLTawePj4zEYDPj6+pot9/X15fjx4wXax+TJkylbtizt27fP9/W5c+cya9asPMu3bt2Ko6Oj5UGLEmXbtm1FHYIoJgrzXLgKOFd4lean5+Nw5TiZn3ZlQtZrXMaDkEpGmvmqhXZsYbnccyEtrfDukVucXG+mqjknTFGNLzxv3jzWrVtHREQE9vb5/3U4depUxo8fb3qelJRkuk8rI02VXnq9nm3bttGhQ4c8fbHFw+VBngtXolty7YveVNbE8LXtmzyrn8rXkb681Ke1lGCLgVvPhdyazMJwT8n1iy++4N133+XUqVMAVKtWjYkTJzJw4ECL9uPl5YVWqyUuLs5seVxcHH5+fnfc9r333mPevHn8+uuv1K1b97br2dnZYWdnl2e5TqeTi+5DQL5nketBnAtR+DMxczprbOcQpIkjzHYmQ7ImczGxKRW8XAr12KLgcs+FwjwfLG4tvHDhQkaOHEnXrl35+uuv+frrr+ncuTMvvvgiixYtsmhftra2NGjQgPDwcNMyo9FIeHg4zZo1u+1277zzDm+99RZbtmyhYcOGlr4FIYQoFBW9nIhRvHkq6w2OGgPxVhJZZ/sW1dIOFnVo4gGzuOS6ZMkSPv74Y7OWuT179qRWrVrMnDmTV155xaL9jR8/nsGDB9OwYUMaN25MaGgoqampptbDgwYNoly5csydOxeA+fPn88Ybb7B27VqCgoKIjY0FwNnZGWdnZ0vfjhBCWM2NVsSHCcmaznLbhTTVHIXv+oPyCTEBXYiMT5XhEh8CFifXmJgYmjdvnmd58+bNiYmJsTiAkJAQrly5whtvvEFsbCzBwcFs2bLF1MgpOjoajeZGAfvjjz8mKyuLJ5980mw/M2bMYObMmRYfXwghrCmkUQVaV/MmKj6NQPcOED4Gjn6PuuE5PtEPYqWhExoFGS6xlLM4uVapUoWvv/6aadOmmS1fv349VatWvacgRo8ezejRo/N9LSIiwuy5zMYjhCjubgz2Dzy5gtTvX8XpnxXM1K3CR7nOu9lPy3CJpZzFyXXWrFmEhISwc+dOWrRoAcDvv/9OeHg4X3/9tdUDFEKIEk2j5Z86r/H7/lQm6r7mJZsfqKjEMF4/kqj4NEmupZTFybVv377s2bOHRYsWmYZBrFGjBnv37qV+/frWjk8IIUq8it7OPGvszaWsMszTLaeLdh8BSjze9vUBGc2pNLqnrjgNGjTgyy+/tHYsQghRKt08XGJ0lg/LbBdRWxMJX3UhvsdKTtpUlUZOpcx9DSIhhBCiYG40dGqCQdcNfhwMV47htLYna/Qv8rPaVBo5lSIF7ueq0WjQarV3fNjYSK4WQojb8XdzoFnlMvhUqE7sUz+w3RCMg5LFh7bvM0W7hulh/xCTmF7UYQorKHA2/Pbbb2/72u7du3n//fcxGo1WCUoIIUq7s8kahuknMEldx4s2mxhh8xP1NGe4eL4qECj9YUu4AifXXr165Vl24sQJpkyZwo8//siAAQN48803rRqcEEKUVhW9nEDRMC+7P38Zq/Ce7hOaaI6T+kM3nkt5iT3GR6Q/bAlm8fCHAJcuXWL48OHUqVOH7Oxs/v77b1atWkVgYKC14xNCiFIpt5GTVlH4xdiYJ/Rvc82pMk5Z8azRvc0L2h9RVSPTwg5LVXEJZNFN0sTERObMmcOSJUsIDg4mPDycVq1aFVZsQghRqt08mlOQlyMnYzoQt+ZFemn/YKruK1pp/uVV/UgORF3H01mqiUuSAifXd955h/nz5+Pn58dXX32VbzWxEEIIy5iN5oQP/bNH8YexFjNsvqCl9gg/a6Yw5esRbDU0lGriEqTAyXXKlCk4ODhQpUoVVq1axapVq/JdLywszGrBCSHEwySnqrgu08I07Muqzvu6D6itiWKZbiFrlXbMyR4gwyaWEAVOroMGDSqySdGFEMVHbGos0UnRVHCtgJ+Tn0XPgXv62cHGgfTs9Nsu83O68/zPJcnNVcXRSR34feM0XrD5if4222mtPcRU/fMciAqWauJirsDJdeXKlYUYhhCiOLg5Efo6+mJMTSP20kliLp3E1+jM8egD/PTv1zhkqDjoFWo5VuLi5TPY6lUOZysE2PlwLSkOm2yVS6pCGZ07SekJgMqVm4+jKBg1YNConNUqZNtAlg0ctIFMnUKmTiXdTiHdFlLtIdUu5/8UBw0pDpDoqGLQ5vyxr1E0zGg2gz5V++R5Hw42DiRnJnMh+wL74vZRyaNSiUjEuVXFMYmOjDYMYLuxPu/YfEIFzRVW285j/YY/map/llTFUaqJiykZ9UGIUu7WkqRqNGK4epXYs4e5fO4YnqkaHBMzOHNmP1FnDuKWqpKVBlczbNDoswFwAlKA8sALZns/RU3TzyoQS3mz59fxyjcq9TY/3/z81uUAN/rSp9hDoiNcd1GI/2E6Zxscwr1cJQ4YI/nk0jdccVG57gyq5r8at3BQUBhcazAdAztyMeUiAME+wcU24d48bGLnrPlMslnPEJtfCNFG0FbzN2/pBzItDB7xcyE1yyAl2WJEUVU1vzO41EpKSsLNzY3ExERcXV2LOhxRSPR6PZs3b6Zr167odLqiDueByE2iAS4BeGfoyDp/nj/2fcsf+8LwTjDilQTVsjywi09G1esLvN8sG0h2uLkEqZBmB+m2kGkLGTqFTNuc9bJsQK+FbC0YNGDUgFHJeeTSqDceWgPoDGBjAFs92GaDnR4cslTss8AxExwywSlTxSkDXNJzHpoCXrWyNRDvCpfdFS67Q5y7QqwHxHrk/J9pmxPYkFpDGFBjQLFNsjGJ6UTFp3E1NZNVX33FfN0yKmliAfifoTYzDUM4YywrDZ7u4tbrQmHmAym5ClGC3FwK9TY4khUZSWZkJIcObOHMvzvxva5ikwAJWTnrlwOeMtvDNVRykt01Z7jmAtedFRKcFa47w3VnSHCCJEeFJEdIcryRgIoLRVVxSge3NHBPVfFIhjIpCs/59SblUjRnTu/DMwU8k8HGCH4J4JeQf2k43gUulVG4WOYz3vFeQYe2z9Gx3TC0bm6mdW4t+ReFG9XE6YzhEbpkzeMF7SZG2XxPK+1hftZM5jNDVz7K7sXUjYdwsrOhQaCHlGKLkCRXIYqRPFW4qorh6lUyT59m7+8bOLzvZ8peNWK4BtdTb2zn9d8jlxEwertz3CHRVGKLd4N4V4U+LUcwL3KZ6Z5lLgUF9abko/z37+aqWAUFRVEwqkY0iobulbqz6eymAj3P2ReoqBb/fPPxURRSHFVSHOGil2K651qpah9iU2MZuLFTzvGMKp4p4JMAPgkqPgkqftfB77qK/3VwzgCvZPBKVqkbBWCALcs5yXJsfHywq1KFaG+FDfo/ifKGiz4apraeme+93QeVeG9UEx/mfUMffjC2YIbNStpp/2GkzY88pd1BaHZfxq41oCpaKcUWIakWFqVSSawW/u7vtXy1aS7lrhipEA8tMwNxvZiA4fr1226jeHuhL+/NbxznYpmcqs44D4UrbjD38QVM2jkJo3ojOWoUDV92+ZJnf342z/Jxj44j9GCoKTHOaDYDgFm7Z5kta162OeeTzxPgEmBqDVzQ58A9/WyvtSfDkHHbZTcntrBTYaaY78Q5TcX/GpS7plIuXiXgCgTEq3gn5b++UYFLnlCx4WN41KnPfpervH1tDQmOar6NqgpTbjWxo62GJz76nXbKQabZrKWyJgaAU8ZyLMh+im3Ghizu30BKsf95kNXCklxFqVTckuutrXD1Fy+ScewYmcdPkHHiOGlHj2C8FJv/xoqCsawP+50uc94LLpZRuFRGIcYTPuy5ggCXADr9V1rLpVE0/NL3F/649Eee5Ninah+zBHTz8lsTY27sty4r7nJjttfak5KZwsE9B0kLSOPL41/eMelqFA1bOoXhHpPM8f1bCd+xkoArUOGyinta/ttcc4azfgqR/hqG9HkL/wYtuOpkfGAl2vX7opkWdhhF1dNfG844m414KikAHDYG8V72U+xUg5nbpy6tq3k/1BMCSHItRJJcHw7FJbmqWVls3r6MrVuXEhhrJOgyVLtqizYtM9/1410g2lvhvDec91Z4rudM6jfqzmVj4m0TqJ+T322TJdw+OZbEpHkvbj4XrmZdNSXdS6mX2HNpDxtPb7zt52b6zFUV91SodFlhptcQko/8y5V/9+J3Lf8B2nMT7tmyGpo+9iyPd3oBGw+PQqtGjklM50DUdcas+wtnNZVhNpsZpv0ZZyUDgH+NFfkouxfb1IYYVM1D2/BJkmshkuT6cCjs5JrfRdKYkUHmyZNkHD1KxpGjOf+fPAn5tcy1scGualXsq1fHvsYjpAb6EHJiMsn2N34db06ewB0TaG5MD0OytNTdzoU7fW53KuF32tgJXaaBoDioFKtSORZapZRFjbqQb8LN8vXggGcip/zhbFkNz/R6jSfqPGPVhJtbijWoKh4k8YLNJgZrt+Kg5LRwO2UsxzJDN34wNEePLe/3r/9QVRlLci1EklwfDoWZXMNOhTFvx0wqxBmpFAtPGOvhez6FzDNnwGDIs36KPUT6KkT5QpSPwjlfhZnPLKdRQLM8+71T8gRJoPfifs+F233m+X1f5Z3L89Km5wiKg8oxqulRNp/b5kYF0gO82eNxldP/JdzBPd6gT42n7ivh3lyKNargSRJDbLYwRLsVVyWnbjtedWWN4XHWZLfnCh4Mb1WRoS0rlvokK8m1EElyfTgU9IJakIuYMTWVjGPHyDhyhOv/HODsnm2UvZp/daDW0xP7mjWxr1UL+5o1Sa7oTdfdgzFy+xLprfFI8rSuwvxDK7/GW/lV378bPIOVYdOpHJOTdKvEqJRJzru/LBvIrFSWXa6x/yVcLSO6vkHz8i0tTrY3l2I1gBNpPKP9jcE2WymvxAOQrWr41diAtYbH+EOtw8TONalT3q3U3pOVfq5CFKLchHrk6pE8rWN7+Xe8Ua175AgZR4+SFRkJN/0NmjsCkakhix907zKWOi16YePrazYGtyswQ5mZp4Rzuwukn5OfJNUS5Nbvy8/JjxnNZuT5vuuWbc6RSjYcqngj6XomQ6UYI1ViVCpfgsqxKs4ZYHvyEl1MaxlJXTmd//kpnPGHs/5aenQdS/eWw+461vut09ntPHmFaWFOfJ7ZhY6a/Qy12UJjzQk6a/fRWbuPS6onYdta8ZqhNefwZ3LnR0p1oi1sUnIVpU5saixnr5/l9L7T9Ovej6tZV01/9d/cetYxQyUoTqVSbM49s0qxKmWv5b9PGz8/7GvVIrtqBWZeW80ZX5VE5xtj296uJHpzTFIiLRpF0bgtv+/71mrkm7s+Qc7gGP7XlBsJN0alYlzOqFW3Ul2dMVQLxLl2PcoENyK5og8XXfRUcAu84/l1a5VxVeUC/bXh9NH+DzflRnPoQ8YgfjA0Z7OhCZfwZnirinSr61/ih1iUauFCJMm15LvTrCs3J08Fhe4Vu/NT1E+4pBioGAsV43ISasVYFb+E/Pdv4++Pfa2aONSqlVO9W6sWNmXKmF4vyL1RUXwUl5bjkDfp3i3hag0q5eOhSkzOH3+VYlQCL+eMPHWrVDs456vgUedR3Os8ik9wE7RBFTifHpOnOnn9vmimbjxkGh7Ejiwe1xzkSe1OWmv+xUa5cYBDxiB+MTQi3Pgox9QKaBSlxJZqJbkWIkmuJcutifTWi9GtIwApqorP9ZwEGnRZJSgWgi7njNSTn8tuEOmnmPopzn/+G8qWf6RAcUlJtGQoTsk1P3dKuPmNUmWTrVLhSk5tS8X/alwqXMkZo/lW2Rq44AXRPhqqN+pI5UfbEetnR/mKdbmSkslnu/ex6aAeg/7GcI+eJNFVu4fu2j9ppBxHq9w49kW1DDsNddlprMvvxlok42xWqnWy1Rbr0q0k10IkybXo3a4R0d0S6a1/1Ttk5FxUKlxRCbyck1ArXAH7fHq+GIEYT4jyVf5LpjkteFMdbj9tmSgdintyzc/NCffW2hg1n9mCcku4QXGqqXYm8DI45d+dmhT7nH7UF7zggpeCb/V22JXtTujBq6i2VzFmeaFmu1GGRNprD9Jec4CWmsOmLj0ARlXhiBrIn8aa7DU+wn5jNa6Tc01VoFhWJUtyLUSSXIvW7apU75RIddkqZa9ChSsQcMWYM2LOldsPU5dlA9FeEOWnEOWjEOmrEO2TMwD9rWPjjnt0HLW9aksJtBQricn1VjePOHXr0JW3pap4J+bU3FS4nPM7U+FyzrjKt5tVKM02Z4jHmDIK0bbVuOgQwEVnLy65OaFXvWhmvEgrzSFaa/6lquZinu3PGv34S63K38bKHDJW4phagUxsUYB+jQNoXsWLAA+HIku4klwLkSRX67lTN5b8XrtdN4XcsW5tMwyUuwrlrqqUvwrl4nPGfPVLuP3FIN41Z0Sjcz5wzkch2kdDbBmFbEU1u+d6p7FxRelWGpLrze50a6Qgcv9YDbiiUv5qTom3fLyK73XQ3iEbJDpCjLMrl+wDiXFyJ8UtGxfHdDzdrtDU9iLNlEt5tslWNZxRy3JUDeSEMYDjagCn1fJcVMsAmnyrlAuzalmSayGS5Fpwd0qed2rUc7vX9sbsZcTPz+GTkDMzSdlrUPaaSpOsALLPRd/2vijcqMbyqFmPzeohznnDBR8N7Wr1MJuFJTd5Rl6P5NS+U6bWwpJMH16lLblC/v1rzyef53D8YVONz91mGLqV1pAza1C5qzkTGpS9quJ/XcXvGrcdVzmXQYFrTnakOTiS7WTAwTETB+cUXB3TyHACO4ds0h2hgiEbP4OBNNWOKNWPs6o/kaof51Rfzhl9Oa96cxkPjGjylHajr6WhKIrZz5aOLiXJtRBJci3YwAl3G6v2duPcqqpKyBed8Ew04p2o4puQM59mJ5t6GC9eIjs27o4TXSc4wYUyCpe8FKrWf4w1aRFEe6kkOWuY0XxmvoPL59e4qDReUMW9edjOhdvNQlSQe7e345CZU7L1u57zO+2boOKdAD6JKl5J+TemupURSHICjb0RR/tsNPZG9A4qNvYGshxUbO0MZDiAVqdywcGNZFy5pHVDn+VNqt6XWNWTy6o7l1V34nEjGxsUYF7fgo+R/NANIvHhhx/y7rvvEhsbS7169ViyZAmNGze+7frffPMN06dPJyoqiqpVqzJ//ny6du36ACMuXiwZKq2gQ+zdPGWXUTUya/csmvk3w1tvz4Wjv1PnTDbuKeCVBF5JKt6JBq6u7ocSd5VPsvK2KMrmLyBnVKMMHcR6QIynhqr1WlP30c78ro3i7YufkWx/Y/qurlX78Gg+iTO/jvtSIhUiR36/HwB9qvYxuyUCsObYGr44+sVtWybnSrdTiPLLacdwK0VVza4FZf773zMZPJNz/vdIyalydk8FUjWALQB2/+0jt+zp9N//VYAMXQrlHVJIcbiIq0029rZGjLYqTjoDifaQbmvDp2WeYWoYtK7mXSwaTN2syJPr+vXrGT9+PEuXLqVJkyaEhobSqVMnTpw4gY+PT571//jjD/r168fcuXPp3r07a9eupXfv3hw8eJDatWs/sLitPbvFve7Pkj6X+SXNt3fNpKntI3hm2GC4fh3D9WvERP3NEwezcU1T8UgB95ScrizX3+1Igj4bJ+C1/I+Qs18gwRmuuOVM0n3ZQ2FQx8l4V6mNbUAAVxyyyU65QP2bEmZX4NHUkLsmUiHEvbv19+nVhq8yoMaAPCXdm6uX70ZVFK67wHUXOFUu/1GjFKOKazp4JIN7as51xS0V3FJzZhtyTQO3NBWXNHBJz+nHa6/PeeQ0XDRPVR7/PTIGfY9W0REV37TYJdcirxZu0qQJjRo14oMPPgDAaDQSEBDAyy+/zJQpU/KsHxISQmpqKps2bTIta9q0KcHBwSxduvSux7NGNYC1BxGwdH+qqqJmZRF7/Tz9wvqg1Rux14NDJjhlKbwdPBUXvQ3GlFSMKSkYUpIxpqQSHx/NsegDOGTmDLPmmnb7pvp3onV3J9XdjuOay1xzhng3DS0b9KZpcA905cqyKXk3s/bPLtJBFh62qkBxe3Iu3JubWyhvPbe1QCVcq1BVHDMxJVqX9JzrlXM6OGWAc4aKU0bOz5931HDFVcPazj9Q1y/orrt+aKqFs7KyOHDgAFOnTjUt02g0tG/fnt27d+e7ze7duxk/frzZsk6dOvHdd98VZqgmsamxvPfbTDocMaAxgkY1cuDPN3i0ZhTONo6oBgMYjKjGm/7PNqAajWAwmC83GEnPTOHaud8Ya1TRGEGXbcD45eucdFuDTbaKmpmJmpmJ8b//1cxM1Kwbfc3y+3Midf1bpOazXAfUze9NKQpad3e0Hh5oPT2w8fAgWpPAjuSDXHOCBBcNfZoNo/2jT6H19kZjm1Ol43+bgRT6EEDzwNbSiEiIEuzmUm4d7zr5lnBvTbxWoSik2UOaPcTlLLjz6qhkcgUIss7xraRIk2t8fDwGgwFfX1+z5b6+vhw/fjzfbWJjY/NdPzY2Nt/1MzMzycy8UTxLSsrpHKnX69HnN8/mXZy9fhaXVCPDtpqfSOm/LSfd4r3laJHPMgNHKUAbAYxAli6nb2eaXc69kSrl6mDn6oHW2QWNkxOKszNaZ2c0zk4cTD3BVxe/J8VOJcVRw/CWr9K9bgiKVmu2X1+gQlqc6ZfJ19H3v7jA8N/nVsa2DGX+Gxbw1s/yTq89CLnHLIpji+JFzgXruPl3GjD9/Ij7I4RUzbmlY29jT0Z2Rp7/98TuYcWRFRixUgK+iUbR4O/gX6Dv99ZzoTDPiSK/51rY5s6dy6xZs/Is37p1K46OjhbvL9GYSLqdwu5HFIxKzpyMqgZq2NXBVmOHqtGARkFVbvpfqzF/rlFy1lM0pCtZ/Ja1HYMCRk1OktTbKHRz7Y2Dzg1VZ4Oq02G0sUG1sUG10eUss7HBqNOxX3+Q7zN+MDW37+XQC71dw9u/AddGtPOtxlXDVcpoy6CNcePnmF/u+J4vc9niz6m42LZtW1GHIIoJOReKTkUq8qrrq1w1XMVWsSVLzbrr/2eyz7Azc+cdq58VFHra9+RAxAGL4sk9F9LS7tLH6D4UaXL18vJCq9USFxdntjwuLg4/v9tMyeXnZ9H6U6dONatGTkpKIiAggI4dO95zHbvDGQfednnbdE/x9cav06xy73vaF0DSme94e6/5/roUcH9d6cELaS/mKWE+7PR6Pdu2baNDhw5yn+0hJ+dCyRX3X+2ZvY09MSkxAPg7+5t+rutd16Jr3q3nQm5NZmEo0uRqa2tLgwYNCA8Pp3fv3kBOg6bw8HBGjx6d7zbNmjUjPDyccePGmZZt27aNZs2a5bu+nZ0ddnZ2eZbrdLp7/kV76pGnaBXQymr3FO93f+XdylPerfzdV3wI3c/3LEoXORdKnpuvbfWpb1p+88/3IvdcKMzzocirhcePH8/gwYNp2LAhjRs3JjQ0lNTUVIYOHQrAoEGDKFeuHHPnzgVg7NixtGnThgULFtCtWzfWrVvH/v37WbZs2QON29pdRKTLiRBClB5FnlxDQkK4cuUKb7zxBrGxsQQHB7NlyxZTo6Xo6Gg0Go1p/ebNm7N27Vpef/11pk2bRtWqVfnuu+8eaB9XIYQQ4k6KPLkCjB49+rbVwBEREXmWPfXUUzz11FOFHJUQQghxbzR3X0UIIYQQlpDkKoQQQliZJFchhBDCyorFPdcHKXco5cLs3ySKnl6vJy0tjaSkJOl+8ZCTc0HkuvVcyM0DhTHE/kOXXJOTkwEICAgo4kiEEEIUB8nJybi5uVl1n0U+K86DZjQauXTpEi4uLijKnQeEtrZGjRqxb9++Yrnve9m+oNsUZL27rXOn1/N7LXckrvPnz1t9tgtreFjPhYKsK+dC8dq3pfsoSeeCqqokJydTtmxZsy6f1vDQlVw1Gg3lyxfNaEZarbbQfrnvd9/3sn1BtynIendb506v3+k1V1fXYnlBfVjPhYKsK+dC8dq3pfsoaeeCtUusuaRB0wM0atSoYrvve9m+oNsUZL27rXOn1wvzcy0sD+u5UJB15VwoXvu2dB9yLuR46KqFxcOhMCdBFiWLnAsi14M8F6TkKkolOzs7ZsyYke+kDeLhIueCyPUgzwUpuQohhBBWJiVXIYQQwsokuQohhBBWJslVCCGEsDJJrkIIIYSVSXIVQgghrEySqxBAWloagYGBTJgwoahDEUUkISGBhg0bEhwcTO3atVm+fHlRhySKyPnz52nbti01a9akbt26fPPNNxbvQ7riCAG89tprnD59moCAAN57772iDkcUAYPBQGZmJo6OjqSmplK7dm32799PmTJlijo08YDFxMQQFxdHcHAwsbGxNGjQgJMnT+Lk5FTgfUjJVTz0Tp06xfHjx+nSpUtRhyKKkFarxdHREYDMzExUVS2UqchE8efv709wcDAAfn5+eHl5ce3aNYv2IclVFGs7d+6kR48elC1bFkVR+O677/Ks8+GHHxIUFIS9vT1NmjRh7969Fh1jwoQJzJ0710oRi8LyIM6FhIQE6tWrR/ny5Zk4cSJeXl5Wil5Y04M4F3IdOHAAg8Fg8TSlklxFsZaamkq9evX48MMP8319/fr1jB8/nhkzZnDw4EHq1atHp06duHz5smmd3Htotz4uXbrE999/T7Vq1ahWrdqDekviHhX2uQDg7u7OP//8Q2RkJGvXriUuLu6BvDdhmQdxLgBcu3aNQYMGsWzZMsuDVIUoIQD122+/NVvWuHFjddSoUabnBoNBLVu2rDp37twC7XPKlClq+fLl1cDAQLVMmTKqq6urOmvWLGuGLQpBYZwLtxo5cqT6zTff3E+Y4gEorHMhIyNDbdWqlfrFF1/cU1xSchUlVlZWFgcOHKB9+/amZRqNhvbt27N79+4C7WPu3LmcP3+eqKgo3nvvPYYPH84bb7xRWCGLQmKNcyEuLo7k5GQAEhMT2blzJ9WrVy+UeEXhsca5oKoqQ4YM4bHHHmPgwIH3FIckV1FixcfHYzAY8PX1NVvu6+tLbGxsEUUlioI1zoVz587RqlUr6tWrR6tWrXj55ZepU6dOYYQrCpE1zoXff/+d9evX89133xEcHExwcDCHDh2yKA4bi9YWohQbMmRIUYcgilDjxo35+++/izoMUQy0bNkSo9F4X/uQkqsosby8vNBqtXkancTFxeHn51dEUYmiIOeCyFVczgVJrqLEsrW1pUGDBoSHh5uWGY1GwsPDadasWRFGJh40ORdEruJyLki1sCjWUlJSOH36tOl5ZGQkf//9N56enlSoUIHx48czePBgGjZsSOPGjQkNDSU1NZWhQ4cWYdSiMMi5IHKViHPhntoYC/GAbN++XQXyPAYPHmxaZ8mSJWqFChVUW1tbtXHjxuqff/5ZdAGLQiPngshVEs4FGVtYCCGEsDK55yqEEEJYmSRXIYQQwsokuQohhBBWJslVCCGEsDJJrkIIIYSVSXIVQgghrEySqxBCCGFlklyFEEIIK5PkKkQxFxERgaIoJCQkPPBjK4qCoii4u7vfcb2ZM2cSHBxs9jx329DQ0EKNUYjiSJKrEMVI27ZtGTdunNmy5s2bExMTg5ubW5HEtGLFCk6ePGnRNhMmTCAmJoby5csXUlRCFG8ycL8QxZytrW2RTpvm7u6Oj4+PRds4Ozvj7OyMVqstpKiEKN6k5CpEMTFkyBB27NjB4sWLTVWqUVFReaqFV65cibu7O5s2baJ69eo4Ojry5JNPkpaWxqpVqwgKCsLDw4MxY8ZgMBhM+8/MzGTChAmUK1cOJycnmjRpQkRExD3FOm/ePHx9fXFxcWHYsGFkZGRY4RMQovSQ5CpEMbF48WKaNWvG8OHDiYmJISYmhoCAgHzXTUtL4/3332fdunVs2bKFiIgInnjiCTZv3szmzZtZvXo1n3zyCRs2bDBtM3r0aHbv3s26dev4999/eeqpp+jcuTOnTp2yKM6vv/6amTNnMmfOHPbv34+/vz8fffTRfb13IUobqRYWophwc3PD1tYWR0fHu1YD6/V6Pv74YypXrgzAk08+yerVq4mLi8PZ2ZmaNWvSrl07tm/fTkhICNHR0axYsYLo6GjKli0L5NwX3bJlCytWrGDOnDkFjjM0NJRhw4YxbNgwAN5++21+/fVXKb0KcRMpuQpRAjk6OpoSK4Cvry9BQUE4OzubLbt8+TIAhw4dwmAwUK1aNdP9UGdnZ3bs2MGZM2csOvaxY8do0qSJ2bJmzZrdx7sRovSRkqsQJZBOpzN7rihKvsuMRiMAKSkpaLVaDhw4kKeR0c0JWQhhHZJchShGbG1tzRohWUv9+vUxGAxcvnyZVq1a3de+atSowZ49exg0aJBp2Z9//nm/IQpRqkhyFaIYCQoKYs+ePURFReHs7Iynp6dV9lutWjUGDBjAoEGDWLBgAfXr1+fKlSuEh4dTt25dunXrVuB9jR07liFDhtCwYUNatGjBmjVrOHLkCJUqVbJKrEKUBnLPVYhiZMKECWi1WmrWrIm3tzfR0dFW2/eKFSsYNGgQr776KtWrV6d3797s27ePChUqWLSfkJAQpk+fzqRJk2jQoAHnzp1j5MiRVotTiNJAUVVVLeoghBDFk6IofPvtt/Tu3fuetg8KCmLcuHF5Rp0SorSTkqsQ4o769etn8TCGc+bMwdnZ2aolbyFKEim5CiFu6/Tp0wBotVoqVqxY4O2uXbvGtWvXAPD29i6ycZGFKCqSXIUQQggrk2phIYQQwsokuQohhBBWJslVCCGEsDJJrkIIIYSVSXIVQgghrEySqxBCCGFlklyFEEIIK5PkKoQQQliZJFchhBDCyv4PkTKayO90IlIAAAAASUVORK5CYII=",
"text/plain": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"hm1 = w.headinside(t1)\n",
"hm2 = ml.head(r, 0, t2, layers=0)\n",
"plt.semilogx(t1, h1 / H0, \".\", label=\"obs ln-2\")\n",
"plt.semilogx(t1, hm1[0] / H0, label=\"timflow ln-2\")\n",
"plt.semilogx(t2, h2 / H0, \".\", label=\"obs ln-3\")\n",
"plt.semilogx(t2, hm2[0] / H0, label=\"timflow ln-3\")\n",
"plt.xlabel(\"time [d]\")\n",
"plt.ylabel(\"Normalized Head: h/H0\")\n",
"plt.title(\"Model Results - Single layer model\")\n",
"plt.legend()\n",
"plt.grid()"
]
},
{
"cell_type": "markdown",
"metadata": {
"editable": true,
"slideshow": {
"slide_type": ""
},
"tags": []
},
"source": [
"### Comparison of results\n",
"\n",
"The values in `timflow` are compared with the KGS analytical model (Hyder et al. 1994) implemented in AQTESOLV (Duffield, 2007), and to the MLU model (Hemker & Post, 2014). The parameters of `timflow`and AQTESOLV are very similar. MLU differs a bit, this is likely due to the added skin resistance (skin factor of 2.42). "
]
},
{
"cell_type": "code",
"execution_count": 18,
"metadata": {
"editable": true,
"slideshow": {
"slide_type": ""
},
"tags": [
"hide-input"
]
},
"outputs": [
{
"data": {
"text/html": [
"\n",
"\n",
" \n",
" \n",
" | | \n",
" k [m/d] | \n",
" Ss [1/m] | \n",
" RMSE [m] | \n",
"
\n",
" \n",
" \n",
" \n",
" | timflow | \n",
" 1.17 | \n",
" 9.38e-06 | \n",
" 0.010 | \n",
"
\n",
" \n",
" | AQTESOLV | \n",
" 1.17 | \n",
" 9.37e-06 | \n",
" - | \n",
"
\n",
" \n",
" | MLU | \n",
" 1.31 | \n",
" 8.20e-06 | \n",
" - | \n",
"
\n",
" \n",
"
\n"
],
"text/plain": [
""
]
},
"execution_count": 18,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"t = pd.DataFrame(\n",
" columns=[\"k [m/d]\", \"Ss [1/m]\", \"RMSE [m]\"],\n",
" index=[\"timflow\", \"AQTESOLV\", \"MLU\"],\n",
")\n",
"\n",
"t.loc[\"timflow\"] = np.append(cal.parameters[\"optimal\"].values, cal.rmse())\n",
"t.loc[\"AQTESOLV\"] = [1.166, 9.369e-06, \"-\"]\n",
"t.loc[\"MLU\"] = [1.311, 8.197e-06, \"-\"]\n",
"\n",
"t_formatted = t.style.format(\n",
" {\n",
" \"k [m/d]\": \"{:.2f}\",\n",
" \"Ss [1/m]\": \"{:.2e}\",\n",
" \"RMSE [m]\": lambda x: \"-\" if x == \"-\" else f\"{float(x):.3f}\",\n",
" }\n",
")\n",
"t_formatted"
]
},
{
"cell_type": "markdown",
"metadata": {
"editable": true,
"slideshow": {
"slide_type": ""
},
"tags": []
},
"source": [
"### References\n",
"\n",
"* Butler Jr., J.J. (1998), The Design, Performance, and Analysis of Slug Tests, Lewis Publishers, Boca Raton, Florida, 252p.\n",
"* Hemker, K. and Post V. (2014) MLU for Windows: well flow modeling in multilayer aquifer systems; MLU User's guide. https://eur03.safelinks.protection.outlook.com/?url=https%3A%2F%2Fmicrofem.com%2Fdownload%2Fmlu-user.pdf&data=05%7C02%7CMark.Bakker%40tudelft.nl%7Cad7f16364d2d4fd55dbf08de73832eaa%7C096e524d692940308cd38ab42de0887b%7C0%7C0%7C639075204580287861%7CUnknown%7CTWFpbGZsb3d8eyJFbXB0eU1hcGkiOnRydWUsIlYiOiIwLjAuMDAwMCIsIlAiOiJXaW4zMiIsIkFOIjoiTWFpbCIsIldUIjoyfQ%3D%3D%7C0%7C%7C%7C&sdata=OBoe8seXZUfoat89Dfr4g6lF%2Bn1FdtXqtp%2F18BMXCn0%3D&reserved=0\n",
"* Duffield, G.M. (2007), AQTESOLV for Windows Version 4.5 User's Guide, HydroSOLVE, Inc., Reston, VA.\n",
"* Hyder, Z., Butler Jr, J.J., McElwee, C.D. and Liu, W. (1994), Slug tests in partially penetrating wells, Water Resources Research 30, 2945–2957."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.13.3"
}
},
"nbformat": 4,
"nbformat_minor": 4
}