{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "24c87f06-8222-4ff4-ad7c-c8e93259a7f3",
   "metadata": {
    "tags": [
     "remove-cell"
    ]
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Exception reporting mode: Minimal\n"
     ]
    }
   ],
   "source": [
    "%load_ext jbmagics"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "6e619510-18bc-407e-a5ea-227cddc0b4aa",
   "metadata": {},
   "source": [
    "# Example: Head Controlled Well\n",
    "\n",
    "This example demonstrates how `pymf6` can be used to control the rate of an\n",
    "extraction well so that the well water level does fall below under a given water\n",
    "level."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b29f3fac-76bd-4242-a2d6-96588cbc625a",
   "metadata": {},
   "source": [
    "## Setting up a model with `flopy`\n",
    "\n",
    "First, we create a simple model with `flopy`.\n",
    "You may use any GUI tool such as ModelMuse to create such a model.\n",
    "This is not specific to `pymf6` but helps to better understand the problem.\n",
    "You may skip ahead to the section with use of `pymf6`.\n",
    "\n",
    "### Creating the input data\n",
    "\n",
    "First, we import a few need libraries:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "758f26a5-c0ae-484e-a2d9-c2f4bc5c0fda",
   "metadata": {
    "tags": [
     "remove-input"
    ]
   },
   "outputs": [
    {
     "data": {
      "text/markdown": [
       "```python\n",
       "\"\"\"Create, run, and postprocess a MODFLOW 6 model with flopy.\n",
       "\"\"\"\n",
       "\n",
       "import flopy\n",
       "from flopy.utils.postprocessing import get_specific_discharge\n",
       "from matplotlib import pyplot as plt\n",
       "import numpy as np\n",
       "```\n"
      ],
      "text/plain": [
       "<IPython.core.display.Markdown object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "%%include make_model.py\n",
    "end_at = 'import numpy'\n",
    "import_module = True"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d36f10db-1db7-42ac-8a63-f34b759825a5",
   "metadata": {},
   "source": [
    "The function `make_input` creates our model with help pf `flopy`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "ba22b0f1-6d9e-448c-92b9-ea8e6bc02285",
   "metadata": {
    "tags": [
     "remove-input"
    ]
   },
   "outputs": [
    {
     "data": {
      "text/markdown": [
       "```python\n",
       "def _make_wel_stress_period(gwf, wel_q, wel_coords):\n",
       "    \"\"\"Create stress period data for the wel package.\"\"\"\n",
       "    period = flopy.mf6.ModflowGwfwel.stress_period_data.empty(\n",
       "        gwf,\n",
       "        maxbound=1,\n",
       "    )\n",
       "    period[0][0] = (wel_coords, wel_q)\n",
       "    return period\n",
       "\n",
       "\n",
       "def make_input(\n",
       "    wel_q,\n",
       "    wel_coords,\n",
       "    model_path,\n",
       "    name,\n",
       "    exe_name='mf6',\n",
       "    verbosity_level=0):\n",
       "    \"\"\"Create MODFLOW 6 input\"\"\"\n",
       "    sim = flopy.mf6.MFSimulation(\n",
       "        sim_name=name,\n",
       "        sim_ws=model_path,\n",
       "        exe_name=exe_name,\n",
       "        verbosity_level=verbosity_level,\n",
       "    )\n",
       "    times = (10.0, 120, 1.0)\n",
       "    tdis_rc = [(1.0, 1, 1.0), times, times, times]\n",
       "    flopy.mf6.ModflowTdis(\n",
       "        sim, pname=\"tdis\",\n",
       "        time_units=\"DAYS\",\n",
       "        nper=4,\n",
       "        perioddata=tdis_rc,\n",
       "    )\n",
       "    flopy.mf6.ModflowIms(sim)\n",
       "    gwf = flopy.mf6.ModflowGwf(sim, modelname=name, save_flows=True)\n",
       "    flopy.mf6.ModflowGwfdis(gwf, nrow=10, ncol=10)\n",
       "    flopy.mf6.ModflowGwfic(gwf)\n",
       "    flopy.mf6.ModflowGwfnpf(\n",
       "        gwf,\n",
       "        save_flows=True,\n",
       "        save_specific_discharge=True,\n",
       "        icelltype=[0],\n",
       "        k=[0.5],\n",
       "        k33=[0.1],\n",
       "    )\n",
       "    sy = flopy.mf6.ModflowGwfsto.sy.empty(\n",
       "        gwf,\n",
       "        default_value=0.2\n",
       "    )\n",
       "    ss = flopy.mf6.ModflowGwfsto.ss.empty(\n",
       "        gwf, default_value=0.000001\n",
       "    )\n",
       "    flopy.mf6.ModflowGwfsto(\n",
       "        gwf,\n",
       "        pname=\"sto\",\n",
       "        save_flows=True,\n",
       "        save_specific_discharge=True,\n",
       "        iconvert=1,\n",
       "        ss=ss,\n",
       "        sy=sy,\n",
       "        transient={0: True},\n",
       "        )\n",
       "\n",
       "    stress_period_data = {\n",
       "        1: _make_wel_stress_period(gwf, wel_q / 10, wel_coords)[0],\n",
       "        2: _make_wel_stress_period(gwf, wel_q, wel_coords)[0],\n",
       "        3: _make_wel_stress_period(gwf, wel_q / 10, wel_coords)[0]\n",
       "    }\n",
       "    flopy.mf6.ModflowGwfwel(\n",
       "        gwf,\n",
       "        stress_period_data=stress_period_data,\n",
       "    )\n",
       "    flopy.mf6.ModflowGwfchd(\n",
       "        gwf,\n",
       "        stress_period_data=[\n",
       "            [(0, 0, 0), 1.],\n",
       "            [(0, 9, 9), 1.]],\n",
       "    )\n",
       "    budget_file = name + '.bud'\n",
       "    head_file = name + '.hds'\n",
       "    flopy.mf6.ModflowGwfoc(\n",
       "        gwf,\n",
       "        budget_filerecord=budget_file,\n",
       "        head_filerecord=head_file,\n",
       "        saverecord=[('HEAD', 'ALL'), ('BUDGET', 'ALL')])\n",
       "    sim.write_simulation()\n",
       "```\n"
      ],
      "text/plain": [
       "<IPython.core.display.Markdown object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "%%include make_model.py\n",
    "start_at = 'def _make_wel_stress_period'\n",
    "end_at = 'sim.write_simulation()'\n",
    "import_module = True"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4b81de4e-1305-44c1-84b7-5c79064dc564",
   "metadata": {},
   "source": [
    "We specify some data for the model:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "6872ac8b-6ebd-4bc4-b593-33dc57cfebd6",
   "metadata": {},
   "outputs": [],
   "source": [
    "model_path=r'models\\mf6'\n",
    "name='headconwell'\n",
    "wel_coords=(0, 4, 4)\n",
    "wel_q = -0.5"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4276b017-42ed-4b5c-b5ff-5d761aac42f9",
   "metadata": {},
   "source": [
    "and create the input files."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "09ad303d-cec3-4230-96bb-c5d502e06c92",
   "metadata": {},
   "outputs": [],
   "source": [
    "make_input(\n",
    "    wel_q=wel_q,\n",
    "    wel_coords=wel_coords,\n",
    "    model_path=model_path,\n",
    "    name=name\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9b356fd7-1fd3-45bd-a095-e669f7363992",
   "metadata": {},
   "source": [
    "### Running the model\n",
    "\n",
    "The function `get_simulation` retrieves a `flopy` simulation object:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "eeff0f50-8067-4be0-be44-437e1ef0c173",
   "metadata": {
    "tags": [
     "remove-input"
    ]
   },
   "outputs": [
    {
     "data": {
      "text/markdown": [
       "```python\n",
       "def get_simulation(model_path, exe_name='mf6', verbosity_level=0):\n",
       "    \"\"\"Get simulation for a model.\"\"\"\n",
       "    sim = flopy.mf6.MFSimulation.load(\n",
       "        sim_ws=model_path,\n",
       "        exe_name=exe_name,\n",
       "        verbosity_level=verbosity_level,\n",
       "    )\n",
       "    return sim\n",
       "```\n"
      ],
      "text/plain": [
       "<IPython.core.display.Markdown object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "%%include make_model.py\n",
    "start_at = 'def get_simulation'\n",
    "end_at = 'return'\n",
    "import_module = True"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3963453e-e8de-4673-9787-0f5ac3048814",
   "metadata": {},
   "source": [
    "This function runs the model:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "fc2741eb-7925-4b81-8f4c-28652e126d7c",
   "metadata": {
    "tags": [
     "remove-input"
    ]
   },
   "outputs": [
    {
     "data": {
      "text/markdown": [
       "```python\n",
       "def run_simulation(model_path, verbosity_level=0):\n",
       "    \"\"\"Run a MODFLOW 6 model\"\"\"\n",
       "    sim = get_simulation(\n",
       "        model_path,\n",
       "        verbosity_level=verbosity_level)\n",
       "    sim.run_simulation()\n",
       "```\n"
      ],
      "text/plain": [
       "<IPython.core.display.Markdown object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "%%include make_model.py\n",
    "start_at = 'def run_simulation'\n",
    "end_at = 'sim.run_simulation()'\n",
    "import_module = True"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c0391751-0c29-49bf-9c10-e7686f55f886",
   "metadata": {},
   "source": [
    "We run the simulation:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "e6e6aa07-019a-4f53-a8eb-e4f861eef6cb",
   "metadata": {},
   "outputs": [],
   "source": [
    "run_simulation(model_path=model_path)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ba8555ae-c57a-47f2-aab1-58aa249a0df9",
   "metadata": {},
   "source": [
    "There is visible output from MODFLOW 6.\n",
    "Setting the `verbosity_level` to `1` or higher in `run_simulation` will show\n",
    "the simulations steps."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cd19fa5d-5097-4b09-8f98-7cfb0f172f07",
   "metadata": {},
   "source": [
    "### Looking at the model results\n",
    "\n",
    "Let's look at the model results.\n",
    "This function visualizes the head at the end of the last time step: "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "4a4c13b2-3e7f-4396-8a43-8b5eb663424d",
   "metadata": {
    "tags": [
     "remove-input"
    ]
   },
   "outputs": [
    {
     "data": {
      "text/markdown": [
       "```python\n",
       "def show_heads(model_path, name):\n",
       "    \"\"\"Plot calculated heads along with flow vector.\"\"\"\n",
       "    sim = get_simulation(model_path, name)\n",
       "    gwf = sim.get_model(name)\n",
       "\n",
       "    head = gwf.output.head().get_data()\n",
       "    bud = gwf.output.budget()\n",
       "\n",
       "    spdis = bud.get_data(text='DATA-SPDIS')[-1]\n",
       "    qx, qy, _ = get_specific_discharge(spdis, gwf)\n",
       "    pmv = flopy.plot.PlotMapView(gwf)\n",
       "    pmv.plot_array(head)\n",
       "    pmv.plot_grid(colors='white')\n",
       "    pmv.contour_array(\n",
       "        head,\n",
       "        levels=np.arange(0.2, 1.0, 0.02),\n",
       "    )\n",
       "    pmv.plot_vector(qx, qy, normalize=True, color=\"white\")\n",
       "```\n"
      ],
      "text/plain": [
       "<IPython.core.display.Markdown object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "%%include make_model.py\n",
    "start_at = 'def show_head'\n",
    "end_at = 'pmv.plot_vector'\n",
    "import_module = True"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5e6c8ef1-d22b-4141-9454-3984a0342e53",
   "metadata": {},
   "source": [
    "We look at the heads:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "e4117cb5-06be-4239-9428-f6d4adfb07fc",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAQEAAAD8CAYAAAB3lxGOAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAAsTAAALEwEAmpwYAACIYUlEQVR4nOydZ3hUVdeG7zMz6b333iCBAEF6L9IRUKQoRcXesL+W194LVmyoqIiigCIivffeAiSk997LTKaf70cIBkhIMnMi4XtzXxcXkMxZszOZec7ea6/9LEEURTrppJP/XWTXegCddNLJtaVTBDrp5H+cThHopJP/cTpFoJNO/sfpFIFOOvkfp1MEOunkf5wWRUAQhKWCIBQLgnC20ddcBUHYKghCyoW/Xdp3mJ100kl70ZqZwA/AuMu+9iywXRTFCGD7hf930kkn1yFCa4qFBEEIBv4WRbHbhf8nAcNFUSwQBMEH2CWKYlS7jrSTTjppFxQmXuclimIBwAUh8GzugYIg3AvcC+Dm5to7OKjZh5qGYI1RFFAbdJKFtJZbICKg1kkXE8DawgJBBLVW4riWFgCoNRLHtboQVy3ha2tdH1MjYUwAq3aNK14X47WytuDkyZOloih6tOU6U0Wg1YiiuARYAnBD767ikfUGSeMLrl9xotaDB458K1nML/vejb7Wgtt/XS1ZTIAfZk1HroZ7PpE27pKF05Fr4ZGXf5M07mevzkTQizzx5C+Sxfxw0W0IRpGnHv5JspgAHyyeC6LIM/f9IGnc976+oz7uPRLHXXIHAvDM3d9LGPNOevQJzWrrdabuDhRdWAZw4e9iE+N00kkn1xhTReAvYP6Ff88H1koznE466eTfpjVbhCuAg0CUIAi5giAsAN4BbhQEIQW48cL/O+mkk+uQFnMCoijObuZboyQeSyeddHIN6BgVgzJfkAdc61F00kmHxsJSQf9h0u/EdwgREByeRnD9qVMIOvl/hyAIksRxcLLh7a/mY2dvLUm8xnQIERCrXwRDkWRCICDwdd97uC98NI4WNhKMsB4veztGhIUgl+gX24CHkx0KmfS/CmdHW8ljSvyj/7/E2dWOMTf15IX3ZtB/WKTZ8Xz8Xfho2T10iwvi+MFUCUZ4Ke1eJ9AqxFrEigUILt8huP6EWD4XxDowVgBtryuQCQJ7ihOZEzKUmcEDWZV1kF8y91OlUzHQPZIDpckmDXNAUCDvTRhLqVLFmnMJ/H7mHOnlFfTx9yOvupr86hqT4r4+bxxR/h7sOpPG9pMpHE7KQac3MKFPF3acSkWt07c5pqWlglVf3kNeUSV7Dqew53AKqZklyGUC40bEsH772ZaDNEH37gE8/9xkjh5N5/DhdE6czESl0hIR4YWFhZyEhHyT4vbuG4qNjQWpKUUUFVTSUMg6blJPTh7LoKiwyqS4TdGlmz/RPQI4cyKLtORCjAaj2TEDgt0ZNrYbfYdEEhXjB8CKb/dwcFeSWXG7xgbw8se34exqR2piPpXlSrPHejkdQwTgCiFAsxtRdw7qVl79MlGHTLz0w2cQjfyevZO/cvYzyX8gc0OGMDNoIKuyD3Fb8CBePbOSHYWn0Yvyq8aWCfVvDqNYf5f+81wih7NzmBoTzfTuMdzXrw8n8vLJqaxiSEgAT63fyp6MzDb/6C/9tIlRPSMY3SuCT+6filKjZfeZdII8nVkwri8v/riZhOyiNsXU6/U8+fpqhvaPYNywGO6aMZD8okr2HE5l4qhuDOsXybtfbqasom1vqoKCStb+dYJ+fcMYNy4WURQ5czaXrMxSbropjrVrT/Dtd7upq9M2ukoErj6FGDKiC+Mm9kSukKFUashILSYtpZCAIHcefXo8G9ed5Jcf9lNW2jah9fByxNbOCjt76wt/W+HobMtdD4/GwkKOSqkh8UwuZ05mEX8sk3Ons9sUv4Hy0lr6D4sioqsvAHu2nOXHL3aYFKsBR2dbptzWH2dXOwCOHZB+FgAdSQTgHyFwW41gOxsMRYh1fwHqph8uitTVLMJXV4VC6HPxQ20t0/Fc0BbURgVvZxn5I+cwUwL6sCBsJBYyBY+EaBnrtJ1XUoZQY2h+uTCuy0l6+mWyaNdkNHpLAApqavny0BG+PHSEPv5+3NI9hsldItlXdA9PjhpIr/MD+Gz/IYxtMHCNiD1Iqe1OFnzUBzdHO0b2CGd0rwi6B/sA8MNTM/lmw2GWbjmCwdhyXAfPWvrMjufE6hgW/5DH4h92ER7swdB+EYwaFIWDnTUDe4fy00d3sOibbWzf3/q7VUlJDStWHGLFikPY2VnRu3cwffuGMnx4V2QygWnTejNoUAQff7KZQ4fSsPasI/LuVNJ+DkGZZd9s3I/f3cDnH20mKNiDsAgvQsM9CQ33okuMHxYWcm66+QbGTezJujXH+fWnA1S2QrxkMoHl65+45GtarR5Vrebi/23trHBxs0en1VOYX9Hq16ExXj7OPPrfyYRGelNTXUdeVhnvv7QGc018qytVJJzKZuDIrhTklHP8f0IEAKxGgKz+fIEg90K0nQeqJVc8TBRF6kQteu0xqhWPoRfLgH8EIMy2hE9zhgMCGqOe1JpCFIIMjT6Lgqp38La9gZuDxvJj+p4mh2EQ05gQfYITuaEXBeByjubm4e1gT7/gcso1Zwh1nMYjA/vTy9eHJ9ZtpLyursUfV6HQE9k1m6wMb0CgrFrFqr3x+Lg60ivMDwuFHAu5nAcnD2RwtxBe/HET2SWVV43pF1uErYuausp/kkipmSVk51cQ28UPX09nLCzkODrY8OoTkxnWP5JFS7ZRVdPyeBujVGrYsyeJ8+cL6NkjECsrBQaDEStrC55+agLbtp9jY8FyrFw1aMqsWoyn0xpITS4kNbkQgKAQd97/bA7FRVWUFtdQWlyNRqPjhn6h7Nx6DkML03ijUeTB27+iTqlFWatGpdSg0xlwdLLlkx/uZv/ORLZviCcjtW2zrAYEQWDyzL7c9cgoCnIrWDjvGwaN6MraXw+j1bR9CXc57l6O3PHIaFYs2c2RvUlkpZWYHbMpOp4IqNchaveD3QOINjPBeGVF8kUBwIi13d1UG6YB314hAEeqgy9ec6I8g7E73uDp4J34Whk4UTWJCAdfJvj2YkP+ycvi66kzLkKptea3k4OuOtx1iefx9dyCi50dC34twtH6F5ytrQh1c6U8N6/FHzcsMg8rKz0JZ4Mv+fqna/fx6dp9CAIo5HIs5DIsFPJW3F1E/GILKU13QaO89IOn1ep57NVVF/8vCCCXy5DLZBhbMcNojuLiaubM/fqKrwsyIz1eLqYywRl9rUWb4+ZmlzNj0scmjwsgLanwiq9pNDoW3PKZWT9zQLA7j710E5ExfvyyZDcrf9yHQW8k7XyBWXEb8/BzkygprGLVD/vQ66U9c9OYjicCAMZyxJo3ofZjEJVgORC0B4DLBECwwMpmHNReOQNoLAANDHROIsgml5/y+7G/Mh6Ib/Lp1ZrVGEllxYkbUWqvviUT4l5EsHsJf53qQ3alEmjbGju6Wybl5Q4U5rs1+X1RBJ3egE5vgFacEnQJqMLORU3S9tAWHyuKoNcb0WN+YqwpnLpWY+moJ/NI0z9bS7R0pzcVKU7u+Qe7IZMJPDjrS3IySi9+XSoBGDw6mv7Du/D4/G/aVQCgo4pAA+KFD5S+EBAQReMlAmApKC4+riUB8LOq4GbPEyTU+rC/MrzZp/SyrKJOsw6FMIxTeSEtDnFEl7PUqK05mhHW5h/Pzb0SL+8K9u/uTkuJs9bi36MQvVZGQWKbTpO2C+59S9HVKKhKcLrWQ5Gcg7uSOLQ72ex1f1PYOVjz4H8msm7lERJP50ge/3I6RJ1AixjTmxUAo7EOB82TLQrA40HbUBkt+Sm/P8194GQYucXnBIJgh7XwcIvDCnAtIcKrkL3J0eiNbdfT6G5Z6PUyks9LUyQlyIz4diumKMkdg/ba6rvCTodzTBWlx9wQjdfH26yttFf3LnWdllU/7uP7T7e2S/zL6dgzgQtcvgS4KACiSFHVq8iNKXzUggBoRTkfZd5Ihd6u2ecZ4pqCn3UldjbPI9O0fPca0eUsKo0lh9Mj2vwzKSz0RHTJIS3FD42m6cRjW/EIK8fKTkduvLck8czBrXc5MrlI6RH3az2U6w6D3sia5Qf/tefr8BJ9VQEwqNDoUlBavtoqASjROTT7PJ6W1Yx0TyS+2g9Li6snAwF8nMrp6pPH/tQuaA1tT3qFR+Riaakn8bKEoDn49yhEq1JQnGLaGlw6RDz6lVKbZUtdoXQVm520Dx1aBFoUAIx4Oj6JTjHiimvbIgAyjNzifRyNwYJ1xT1aNbbhXc6i1llwIK2LST9bdLdMysscKCxwNen6yzGgwbtrCfnnPBEN1/bXauuvwta3rnMWcJ3QYUWgNQLgIbPGznrgFde2RQAABrum4G9TyV9FPVEZWt7P9nCooptfNgfTolDr2j6Vd/OoxNO78sK2oDQJwWJOobA0knv62i8FPPqWYtQJlJ2URuA6aV9a5TYsFTfc0Es8euCrFh8nijpUNYvQa49hbXd3/TYg9UnAwqpX0ehS8HR8sl4AFF2pNViQXFNQf7EhDeoeB8EKbD4Cmf9Vn8tgyKaq9hEsFP1wsHsegEgHH0SDQGJx08UZpfqPqRMP4qtYglxofea7q6cHghG25LxHGXuIZTEKofkcRWuJ8vdgQ9EbFKnSGSy+jSCRtoeHeCCIkJraevc4AzrWCq/gTRcGinOvjBnuCSKkpZhWoNMcYRFeAKQlX1kXYFbcSG9AbLLewKy4UfVinS5h3NAobxwcbY+LonhDW677VxODRlHgRG3LW1fO+l9w0x+jRPFYfSFQLSAqcdC8iNxYnwOo0I+AWoh0sEA0CujUCgQxFYXhCcAanexT0F5dAETRgEr3MQK2WMoeRaeufzlEOwFEEJoo+tKJhajYgwOTUBjauPUlgs5YR7m4H1f6YYFdfVm9mdQZqsmpO0WgOApBFJAkKI3DtD5eHmfRCnWEGvu0cN2/d/PpyHSEV+FfFQG1QdcqV2CFYKCX/UiO1pTRuBLQwbbkwi5AGpAG1DsD69QK/rPrc17tuwatQcFLR8dSqNrf4vOMDYpnemQSX8cP51jR5otf/3XSTGRqgTu+u9IVePLAQ/SMgJdXKqipa5tr8A8LplOu3YlRUPP9coHivFUtX9QK/vu0B0YM/PhpJTWF0jkOf/LmTASDyJOP/dzqayLuT8XG24L3Xz0K4rErvr/o49sRjCJPPySt2/D7n88FkfZzG26PuCCpi/F739xh0nUdMiegF+Ucrak/jtmaSkBBTG0kANMoVDm3+BzedpVMCTvJ8aIgjhW1XBQE4GirpGdEOieSw6mpM+2sfgk7KC9xoDhPus5theJh7ERfagrNX1qYg4WzFqcu1ZQdcQOx03jgeqFDikADrREADGkoDAvbJAACRu6M3ovGoODn8wNpbXJuUPcEBEFk/5noNv0cDajEDFRCBudPBbf6OVvC1rGOKtLwEftKFtNU3PuUI8ig1MQy4U6uDR1WBFojAH5WFfVJQKxaLQAANwadI9S5hBVJ/anRtm4f2866jt5RqZxODaWytvkjsVej1LgTQbQg9Zx0NmqBXeoTS95iH8limoaIe78yqlPs0ZS2vMPSScehw4iAp5UjD0SMQUBotQA8HrQNBCt08k9bLQDetpVMDTvByeJAjhS2fMimgQHdziOXGdkXH9PqaxpjodBRYTyAK/3RatpeXNQcgV0LcSIMG67tnrx9qBJrTw2lhztnAQDWNtJUgTbQe0Dbz6a0lg4jApZyBXNDhjDSK4aRLkmtEgCtKK/fBhSa3gVY0L03gQ7OF/8vYGR+zD60RjnLE1u/DLCx1NC3axLnMgMpq3Y04aeDyJA8jKjxYKRJ1zeFk3sNTh61eAv9JItpKu79yjCoZVSccr7WQ7nmxPYOZszknpLFc3Kx457Hx0oW73I6jAjkqsrZXHCaBeEj2VTWjRfTJrWqFLi5OgAruZxJoV1YMXnmRSEYHZhAuHMxK84PoFrb+sRev5jzWFno2XO6G/cO60OEV9vvdt0jMrEmADuuPMFobanAz6PtJ+0CowsxGgQ86d3k9/19nPF0v3qhlCl4eV86VpmlAde4CspPumDUXt2y7f87sb2Def3T28nKkM4A5PZ7h+HjL10i+XI6jAgAfJ+2kxB7T4Z6xZChvnJ625ZKQI3BwLwNqyhW1rJi8kz6+9QxNfw4p4oDOdyGZYCVhZb+0UkkZgVQVePO4IhgVj90O0+MHYyNRet2WD1dK/Fyr8RdNgKhidnHhIHRrH7zDp6bNxovl9bmG0QCuxRSlOWGpdD063D79H78/NXdPHjncBwdpLGqtrCU88VXd7L4i/n07Vv/Ojp1qUZuZaTsqOkVgn4BrkR08ZFkjNeKBgGwtrYkNbFAkpi+Aa5MvOUGrG0sJV9iNNChThFmq8rYUnCau8NGsLsoARu5BSpDvWFlW0uBAWp0WuZtWMV3Y7txT/e/KK2z51z5dKCy1WPq0yUFGyste051Q6M3MP+7VUyLi+GpcUMY3z2SN9bt5HxBCT0Dfdh8NqXJGN0iM9Hp5bhYD2rSPPnP3WdQqbXcO2UgkwZF8/uueH5Yf4T+3YI4k1pATvGV43X3r8TWUU38nnCamFwA8MHnW0hILuDOWYOYeGN3fvn9MKvXnWDaxF6sWX8SjbbtFlg6rYGnn/iZufOH8PZ7s0hIyGNZ4scY1DJuHjKOH1Nbrs9oihvHx3LbnUMoKariwN5kDuxJIv5EFuNv6kVpSTWH9jX92prC1Fn96D80CpVSg0qpQVlb//epI+mcPp5pUszGApCTWYqytmlfzLZyx0OjUFjUz66cXe0ozNO2cEXb6VAiALA0bRe/DV7IE10nUqVV8W3aDrra5bPAb3+bBKCBUOcUCpVLsbNwp5/vV/TxcUMUt7M7N7PFa+1t6hjQLYGUXB/yy+qXAKIIfxw/x87EdJ4aP4Qv503ldHYBXXw8KK1RcTzrUksxGysNUSG5JGf60bdr0/v4RlFk06HzbD2azKSB0dx9U3+mDu1OfkkVjvbWPPTBajLyyy+5JqRbHnqdjPy05iswDUaRdZvj2borkek3xXH79P7cPCkOgH5xITz7+h8mCUF6egmvvvwHYeGe3DavP07R1TjqY5g3fxh2djZ8+fk22lqN/sOSXezceo4BQyMZODSKKdP7UFujpqS4mpAwT7asP82XH29G2cgktK04u9rRrVcQgaEe9Or7z2ywvLSW7z/fRvyJNnf1BkCukBEZ7YtOo8fa2pKURNNs1y8nvIsPASH/zIjrRcA0M9Sr0aGWAzIExvjEUqlVMTNoIGH2DszxOcjCoB3UGKzbJAB2Fmru7raLh3tup0Rlzy/nZ+JjF0Kosyv/6TusxZSgQazgjvFbsVTo2X685xXfr1DV8cLvW3jsl7/pEeiDlYWCxXNuItTjnymxIBgZP/QoMsHIyYTm3YwuPqfByNq9Z7n5ue/5Y9dpwvzd8XC256tnZhDu/8+bwTesmOBuBaSf9sega1nH1Rody1cdZva935BXUImHmwO9ewTxzos3Y2Vp+n0gLa2IFVkfUakto5/rGABumd6XJ5+agEzW9pqFrIwSfv1xP48uWMptN33M91/vxMfXGYAxE3uw5Of76d2v9Us5gBHjuvPYC5P59veH+W3L0zz35i2EhHlRUliFXm9g1bL9LLjlM7asO2WySYhBbyQlMR8HJ1tW/riP5HMte0u2htTzBXy9aBMALz6yHFEi67LL6VAzASMiu4sSmB00iKq6nbganmGgczGbSmP4uyS2xT4BDfTyyOT2rgexs1CzNq0XmzJ7MjXcn1K1Cj97R6LdPbkpvCtrUxObvN4gVlDMSzjZKVm+ZSQFZU0nAgUBQjxcSC4sJdLbHSdba76eP5XZX/9GaY2SQXEJBPiUsmVfHGWVrd9VkAkClhYKUnNLCfZxxdXRli+fuZVHFv1OTmUmfcafo7zQgTN722Zm4mBvRXWNmuLSGjzd64Xg7f/ezHNvmDIjEAmcnoNLjyqMh3sg9PKm3KYWV1d7xk/siY2tFW+/uRa93jSfwNKSGkqKqog/mYVMJkMmF5DJZMyePxgXV3u2bWzaH/JybprRF02djp2bznD2RBbnz+ah0eh46JkJrP3tMLlZZSaN73IWPHIjB3ef57tPt2JjK93aPTLGj4Lcco5IuBy6nA4lAgDJNQV8kfAMo13+xEoRwZtpvUmta13Cyd5CzewuB+nrnUFWtRsfnxhDbm39B3h18jnWpZ1nbnQvHu7Vn6duGMzG9GS0xksX6Q6WdRRrXsRAKcu3jCSryKvZ5xNF+HrXEb7edQR3BzsGhAUyKDyIN6bdyOL9X9E7JpXT50NITA9s02ug0el5b3l94wpLhZxQPzciAz0ZFhdMhsPfABxc1wNjG30Dcgsq+e/bfwLg5mpHl3BvukR4M21CL35be7RNU3jv0UV4DSmlYJsXuX/JOf7zCgCcnGwJCfUgNNST/gMi2LfX9A48h/almJ0LePyu75r8+ufvbTArbmMGj4omvKsvD8z+EoA6lXTr9qgYP5Ikmlk0R4cTAYC/Cq1xsbiR+VGfQfq3QMvHLeM8M7m9ywFsLbT8mRrHpsxYDOKlHxKNwcC3Z46xMukM9/foyy2RMaw4/88dxcGyjvvidmIQNXjyEllFTc8UmqK0Rsm6U4msO5WIm3MVMyccJb/YlT3Hurc6RlNo9QbOZxVzPquYXqPOEx5Wxf41PVBVmefYU1auZP+RNPYfSWvztW59ygi4KZ+yoy7krvO95HtVVSpOnczi1EnT1tfXG3KFjDsfGsW2v0+Rldb6I9etJSrGjzW/HJI8bmM6pAhoRQXfZHuRqf2LYHsPUmubFwGjsYo63Vc80GMPWdVufHhiHHm1V585VGs1vHd0LzaKfyr3GgTAxVqJu+WrWOtigNaLQANWllomDT+CRmvB+t19MUpksukfVUh4rxySjgaRn+YpSUxTcIyqJvi2LKqTHMj4Jeh//qDQuClxeHo78ezXOyWP7epuj4e3E0lncyWP3ZgOKQINbC2Mx17R/P52N/s8qmofQBSVrEnpzeas7lfc/a9Gnb7ef76xAHx7ahhvD44Bk6zpRcYOPo6DnYrftwxGVSfN3ry9i5IbxiZQmufEmb0tJxjbCxs/FeEL0lEX2pD6Xeg1tzHrCBzZl0xVpYqSomrJY3v5OlNZoSRFopqD5ujQIgBQq79yv9VOrmGy1ym6O+Qjk0VgJXuKDZmHTYp/uQBkVJp+l+3f4zwh/kXsONSDghJpauhlCgMDJsdj1Ms49HfsNbPvtnTREHl/KnqVnOSvwjCo/7crAxsoKaqmpCihXWInxucyc+R77RK7MdeZlIt0d8hlYfA2utoVsrkkGke7RchlwSZFk1IAQvwL6NcjiXOpgZxJNm08TdFrZBLOnrUc3tCNuhppZhZtRW6rJ/KBVGQWIslfhaOrap/KtU6uDWbNBARBeBy4m3qXpDPAnaIoSlMqdRl2cjU3eZ2mm0M+uXUurC6Mo0TryOxw0+5IUgqAs0MtYwcfp6jUmZ2HeiDVuf6g6HxCY/NIOBhCUea1OSVoQEfEPWlYuWtJ/iIcdaeF+P87TBYBQRD8gEeBaFEU6wRBWAnMAn6QaGwXEIl1yGWSZzxWMj2bSmLYXx6O0YxJjJQCYKHQM2nEYQxGGet398VglGaa7OhWS9yNiRRnu3DuQPsdI70aRtHIQeEXHMKUpH4fQk2q9IeROrn2mOw2fEEEDgE9gGrgT+BTURS3NHdNXO848ZutrfflMxrLUdZ9gU5/ELk8Cnubx5DLL91zj3TwQTQKJJS1bnvGIFbU1wGIpbhbvoi1/Ep/gGg3TzDC+YKrnwQzihrSDIuoFRMJlz+Lg+zqXgNdfOrde5Ozrx7XIGo4Kr6FDiV9hRexasHRODLQA4yQmiHdFpUoiqQ7buZY9RZ6GacQxVBJ4oaHe11wG5bYvTei3r33enMbljJu2L/tNiyKYp4gCB8A2UAdsKUpARAE4V7gXoCevXq1qkW1KIroxR2oxS8QUWMl3I2leAtinZzL69pE23pnYJm65Sl4QyWggVI8eal+G7CpXQAjCCLIrrJDYERDOouoJZFg7sPJENPk4aDGCCL1Y22hOC9JWI1SKCDO+Bg2tOKIsVi/ABEM0pWVnmcnp6u30M9xPKHlQ5Dcwbgj2Oy2GqG+PFTqmCBxXNNimTMTcAF+B2ZSfyxvFbBaFMXlzV3TNTZW1M25s8XYfQNTuLPfTtLLPFl2dDhFNc7NPvbnWdMR9DTpDPwPIsHeRUwaeORiKfDVKgF/WDAdmQ7u/bDpmAqFnnFTD+HjX8quzb1JSWydXdiSJ6Yj18EjLzXvCOzdtZi+t58hdW8gCZtbVxb82WszEfRGnnz8l1Y9viVce5cTNi+DQLEnAw238cxDrXcbbokGV+CnH1gmWUyA97+cB8Az9/8oadz3vprfrnGlfB3e/3IePW8I+Vf7DowGMkRRLAEQBOEPYCDQrAi0lhO5oVjIDRzIiEQ0awNDJMI/n6E9zhDoVUqNyrpFAWgJUwWgNVg7qek5LZGKXAcSt12bPIBDRDUht2VSnWJPv7DZkjUy6aTjYo4IZAP9BUGwpX45MAq40mjeBPRGOfszTOvxByAg0iUoh6E9zuDrXkFlrR1/H+jDyZRw9AbTE3ftKQAIInHTzyGTi5xY2e2aFOLY+KoIX5CGusSK1O/CkL+l4Dqbt3diAubkBA4LgrAaOAHogZPAEnMHZKWQo9G3sLhuBplgJCYki6E9zuLpUkVZlQNr9g7gTFqw2Vn7dhUAIGJYJu4hlZxYHY2y3LSeBuZg6awl8r5UjBo5KV9FYKjr8HVkZmFhKUenNe199m/GbM+4DZh1uxFF8WVRFLuIothNFMW5oiia7vgAfH/rNB4dNMCEceiIi0jl4Vv+Yvrw/YjAql2D+GrtFPr5zcbX2bzGmO0tAC4BlUSNyCD3tBe5p6RvKBoYePXqRbmNnoj7U5BZGUj+KgJtZcvFQAqFjIlT4yTPl0mef2uG6XOubGTbEWPa2VsxcWrTHpJS0aEWfEdy8pgWE428le8EhUxPjWE9eTzAlCGHUGstWbFtGF+umcTZ9BAMRoExMRH8d7LpDr+XC0BtcQzfPTUDP/e2G4M2Gd9aR+8Z5xDV9tjmjkfqBiI9ewby7dK76dkrqMnvCwoj4XenYe2hIfXbcOoKWlcM5OPnwgOPj2XaTOmcjr18nPhy+X2ERpies2kNoRFe3DSjr6Qxbe2smH3nEBwcpa3qjO4ewJjJPSSNeTkdSgT+PJeAu50tQ0ODr/o4S7mOIREJPDP+TyqM36DAg582j2TJX+M5nx2AeOGDZDCKvLJ2G4PCgxjXLbLN42lqBlCjVCOXCyx5croEQiDS46bzWDtq8C6dyZtP3MqUMdL+wk+dymb79nO89PJUvLwuG68gEnp7Jo7htWT8HNymYqCcrDJ+XLKLu+4fSWCwNNWM5aW1lJXU8M7iOQSGtF+F5Lgpcbi62ePpLY2QA3Tt7o+VtQXDxnSTLCZATI8AwqN8CAlvv5OjHUoECmpq2Z+ZxS3dmi66sbbQMqLLGf4zfg0TYk9QWOWMp/wNvHiL1DxfmrqLnskt4tcjp3lu0nDsrVpf825E0+QSoFat5eFP1lBSWWu2EATEFeAXW8z57aF8900GXy3fy9P338i9tw82OWZTfLRoE4WFVbz2xi1YWf2z1g+YkotrXAU5a/0oP9H2JdPvKw6RlJjPMy9NQS43/62k0xl49T8ryUgt5t3P5+IbYN4yriksrRSMHF/v8RAV4ydZ3OjY+vfHmEk9JYsJ9SIAMHpC+80GOpQIAKw6c46R4aG42tjgZlufILOzVDMm5hT/Gb+GMTGnySr34IudY1m6bzTWsm4ILSwfPt5S74C7cMwgIr3csWtBDIyihlQWNZsDaEoIbh/Vq00/p527ku6TkihJdyF1b/1UfcXao7y8aB2zbrqB/z46nogQT0YOjGpT3KbQavW8/OLvuLnZ89TTEwkKdqfrzUq8RxRTtNuDwh2mTb+NRpH3X19LQJA7s+cPpnuvQCwszUvAajV6Xn7yV/JzK3j/i7l4+zozZYZ0LdYGj+yKg2P9kkdaEfC/GDMotHnz17agUMjocmGMo8Z3RyZvn4RJh0oBC8C2lDTqdDo+mjyBxJJkzpR9T9+QFBRyA2fzAtl1vhsFVW27Q9RqtLy9fhcfzJjAwLBAPti0l53n05t8rJ1NHamGd1CSetUkYIMQLF44je+fmYmrgy2HE7NJzW/Zs87SVkuf2Wcw6OScXBV9iTHH9v1JlFeqeOvZKfTrFYLBaOTIqUxqVWblXCkpqeHVV9bw4Gux9PcrIl17nvLTzmSvCcCcPERhfiVLPtvKw0+OZ8zEWBYv2sSRA6lmjVWt1vHi47/wzmdz+ejbO3F2tuPowTTyc8pbvrgFuvcKIiezFFt7K8mm2DKZgKe3M3UqLdmZJfQfGkVWuvnNR8KjfCgtrsYv0I3sjFJ6xAVz8miGBCO+lA41E5AJAl/ffBMI5VjbLcfJ5QUGhCVxNi+Qj7dOZsXhoW0WAAAHayvGXsgJhHi4MjC8ac8/X88yZk/cRZ2YTSiPtLgLoJDLSM0rw93JDplM4OFpLU/jtdQw8K6T2LnWcfy3bqibOB5cVFrNmfP5uDjZ4u5iL9nywGtoMUuzXqPYmMFg5ynkr4w02xnI1s6K7j0DkStkePu6MGiY6fUdjbGytiAnqxQ3dwfkChlz75bm7MInb/1NblYpxw+msujVtZLElCtkPDJ/CUUFlRzclcTKH/dJEjcvt5wn76uvVFz6xXZOHZNeAKCDiYBBFHlmwxZKVXlkVv+FJcNZtOUmVh0bREmN6WvvGrWGxdsOUFhVA8CA8Msz5SKxUencPGYfWr2CSMWruNBy9rhSqWbz0STOZdYfAhkaG0pcRPNTTK1Yw3HhQ+zcVBxe3oPS9KYFrbishp0HksjILgVg6tieREeYs3Uo4js+n6qIA7jVRfGfLl8wKXAuI4aZn8RSKTX88sM+8nLqZ0ADhkRKMm2tKFOyY9MZ0lOKABgxtrtk02wnFzsqK1SUl9VKEk+nNaBSaqmtrsPByabNPReao6aqjpqaOgDs7a0li3s5HUoEAEqUSl7alMn4wA34Wj1KhVKa46tpJeXM/vpXzhcUE+rhirdTfbsvudzAjYNOMKJfPFl5nvy6fhg2zTQ4bYqjSTnMfXsFz3z9N9lFFTx6c9N3bSsbLSeNH6KiuF4A0pqf0ej1RjbuPMe8x3/g6Tf/4HRiLk/fPwa5CV7+CPXW4H7jCig55Mam5+147aV1VFQouemmuLbHa4LszFIeWbCUowdTcXaxo1ts29yVm+P44XQenLuED99cR2WFknn3DpckrpOzLVUVKkliNaamRo29RO3eGtBpDWjUOsnjNqbDiQDA4ZxcPj+QgLONtD94SY2Sud+sYn9KFgPDg3CwU3HruD1Eh+Vw6FQX1u3sj1ZnmmvOthMpTH9lGX8fTCQ29NKeelY2WkbcfIw6iukpPnRVAWiMKMLB4+k88uJvfPD1VoID2mZZJsiNhM7LwGtICQXbvchcEQRGgQP7U7j7rm8pLKwiKkqa/n+1NWpefOpXfl22n0HDzE9mNmA0imxae5I7b1lMVnoJARJsRzq72lFZoZRgdJdSU13XLh/W2lo1du0oAh0qMdiYJYePEuAs3T5uA0qNlgd/+pMxcTB74i5kMiN/be9PRp75lXp6o5HVey5titEgAPbOKmJlj+JmjKbehKltnEtum9mkzNJA+F3pOHWtJmetH4U7Lv35KitVvPLyH7i6Nt0azRSMRpGlX+64mCmXkjqVlmVLdpkdx8JCjp29NVXtIAK11Wq8L3RMkjquvf3/oAiIwNs7d7dL5NguSYRGnaO8yoH1u/pRWdPaTsBto7EA7P2rF2NujW6X57kcua2eyPtSsQtUkvFLEKWHm797lpdL/2FIiG9fi2xzsLRScPJoOsVFVZLHzkwrvqQOQyrOnsqmolya/EVTdFgRAFBqTfL9bhYLhZ7RA08QGZxPcqYf2w70Qqdvn5fgcgEoypHGfbglLJy0RD2QgpW7htSlYVSecf5Xnvd6QVmr4dkHf2qX2Jv/Osnmv05KHveTd9ZLHrMxHVoEpEIuMxAelE+f7sm4ONaw91gMJxLCkbpOv4FrJQBWHmqiHkxBYasn+auITk/ATlrF/2sRcHGsoVtEJl3DcrCx1lJZbcef2waRUyjNVlNTXCsBKCeHrgvr+/6d/ywKVe6/fxy5k+uT/3ciIJcZCAssoHtkBv7eZRiMAmnZPpxNCSanwIP2uvvDtROAIjGFfcJ3GLUykr6MQFNybfoTdHJ9YrLHoCn06t1b/GrNn5LG7OrpASKcyj9NmXEXZcY9GKjFEk/cZMNxkw3DogW33stpcAVOyml96adWrOGk8UPqKCZW9hCuwpVJwKiA+rGmZphfUtpArhjPIfEnXCy8GKhZgG1rjElbSXs4A7evK/D/dtywSG8cHG3+VY/BNiOIIJewNYlR1JFVu4vU2vUU6+MRRDnOxOHOSByJQTDIWnQAbm6ciCBvZV5SSw2nhA8v1gHUbwM2gVgfW9Ab2z6oJkjjMMeEVfhahzHL80kKUmsw6QduDmO9jbFglPJGIdKes7F2C21Koda1jNsG/tWZQHS3WNFqWMtuwy3h5FRL126ZREVnY2OjxVL0YP8ed5LPBFKnNH8q/NUztyLTX90VuAErOw0D7jz1TynwVQqBLroCP2a+e6/3qEICpuRTlejIgqhnsdQpeOau78yO25j3li4AQeCZu7+XLua3dyIKQvu49wrw9EPSZv7f/3wuCAJPPWK2f+4lfPDZHACeXChd3EWfzKFnXHDHngmYg0xuICS0gOhumfgFlGI0CmSmezM2fB6OxPDFod//9TF5hJfR65YELKz1LQqAdIj4T8nDZ1QxZcddyFgehGKRFZLOADr5n6LDi4CjUy3Rje761VW2HD7QlaSEIFQqa25d2L1++v4vIsiNdL0xjfDB2VQX2XHwh17UFLVPwdElyESCZ2Xj0b+Moj0eZP/ub/YpwE466VAi4OPqiL+7E8dTswgJK6Brtwz8G931E84Gk5vtSVsXfrbWFqjU0hQe2bmp6D3jLM5+NWQc9uPcxgiM+vZv0y0ojITdkYFLbBV5G73J3+hDW14HVw8H+g3rwsbVR9tvkEBIhBcZF07+/S9haalAq22htZQJODraUF1dJ3ncxnSoA0ST+nbl1blj8PUv4cbxR3FyUnLkQFeWfz+GLRv6kZvtRVsF4NFbh/LeQzdJMDqRgF75DHvwCA7uGk781pMz67pIKgCOjk2bfMqtDUQ+kIpztyqyVvmTv7FpK7Wr4e3nwsKXp+Ll5yLBSJvG3dORxSvup99Q6Q4QSYWrm327OfMAzJfI7+By5t05pF1KkRvToURg28kUfFwdcRK68/eagfzy4xhOHItCpTS9HfauEyn0iwkiLsr0Qy0KKz1xM87R65ZEqgscme33CfePuxOZhJndSTf14qtv7sLnsgMoCgcdUY8mYx9aS/qyYIr3muaGk3Q2F2WNmrgB4RKMtmnG33IDFeVKjh1IabfnMBW5QsZbH96Go1P7tFafPLU3PZpxdDaHqK6+TP5fshzPKConvaCM0b2iyM3xRJRgvRufVsC+0+k8cPMgk653Cahi2MOH8Y0pJnFrKPu+68XbH+5nQFwIj9wx3OzxNbBzewJlZbUs+uj2i0Jg6aqh68JkrD3VpCwJM8kMtAGD3sjpo+ntJgIKhZzxt/Rmw6qjGCTaApWSkqJqbGwt+fz7u4mIkr63AwI88exErK1bbrjbFjw8HZl52wCsbaSN25gOJQJQPxsY3at1jThby1drDtAzwo8B3YJbfY2IkYhhGQy6+zgA+7/tTcruEBAFziUX8PKHf3Pz+F7MnFyv0ua67SqVGp59+teLQhAQq6DrY8ko7PUkLY6gOtH8IqATB1Lp2S9U0hlMA4NHR+PgZMPGP45LFtPbz5kevYMli7dnRyJePs589PWd3DghVrK4AIIg4Ovnwh33DJMspkIhw9XVHhcXO6ZMa9OuX5vocCKw/VQK/u5OdAnwJNJfmhr/pOxith9L5v6bBxId4o2fx9U/UGqxguPCh3S9MZ38c57sXtyPipxLr9l3NI1PvtvBI3eOYOTAKJ68Z5TZ42wQgtTSJELuScDSUk5U1myUmdLsPJw4mIqDky0RMX74BUlT0uzkYkdIpBeTZ/Zl/7YEKsyw7PL1d2HclF48/epUlv/9GD/+uZDo2AC8/ZxxcrHF0sy18b6diUD9ceJHnhrPuMk9zYrXmAbD62m39iW6mzQuxm7uDhcFe+bs/tjammZ40xIdandAIZfh7mRHWbWSD++dzPZTqSzKNd9TwEIhZ+W2U3z97Aw+f/IW3vxxK3klTZ8n9wou5Yj4Gkb0nPy9Kzknm8/C/7HpFJ7uDrz0+EQUchmbdiUQfz7PrLHqbatYU/kRjg62PNnjLex7uLP5tywqzDz3b2mlICLal+pKFc+8fSv7tp3j+4+3mBUTwMvPmdc+nYOzqx0fvfInLm72JglBUKgHN982gMEju17iznPHgyO548F/OkgZ9Ebq6rS89/IaDu9NbtNzFBVWcT4hj4Agd4oLq9i+ue3mLs1RXqbExsaCPTsT0UrUN9DD05GTJzLpFRfMr78cxNvHmfS0YkliN6ZDiYDeYGRITAhujvVuNxYSNLQAsLSQs3BmffbW3taKyAAPth297A0kiMQMTKNr/wys8CNWvI81J68uQG4udvSODUJxYZy3TetD/Numi4CFo5aoB1ORKQSmuzyNh40vAHPnDubTTzabHBfq/fx79gvD0dkWR2dbyRpcuns44nzBnWjmgqEcuHC3bStZ6SV89MZfLH5vPX0HRTBibHf6DY7ks3fXk5yQj42tZaM/ViZvQ677/RhJCfl89PUd3HnfCJZ8ts2kOJfz9KPLmT1vECFhnqRKdB4gLaWI/zzxC2s3PkV5WW27CAB0wOXAJ3/uI6Ow3l/eQiHN9puyTssjH/5Bcnb9ixgZeGmG3cpWw7BbjxM9IIPMs77cIDyLHS0nj8oqlDz9xu9s2HEWgMF9wgn2N22aLbfVE/lgKgo7PYmLQ3nvhd3s2X0egImTe+Inwdbedx9tpvLCXVqvk0YE3DzrPQtUSg2vPvYzNVXm7WnrtAb27zzPG8+uYta4D0hOyCczrZjEM7mcOJzO/p3n2bb+NMWFpjkDbd0QT3ZmKZ+9v5Hptw2gp0Q5h6LCKg7tTyG6m//F5ibmUlenxWgUyUwvITS8/fozdjgRUOv0/PfHTegMBiwlEgGAaqWahxb9TnpeKVGB/+QaPALKGTP/EK4+VRzZGMOxzTHIBatWx62sruOtxZt45KXfyM4rZ9ZNbU/gyCwNRN6XhrWnhpRvwlDl2JKTXcarL//Bww/+wLmzedy5wPyEU211Hd98sBEAnUSFLe6ejhiNRt59bjVZadKdjoR6F6DMdrr77dp2jh2bz/D0izdJZg568lgGer2BPv1CJYnXQHpaMaFh/yO9CBtIyC7imw2HsZBLW4lXWVPHgx/8jlKtw83Rhi790hl263G0agXbf+5H1jlfk2OfPJvD/Md/JLegAsc2mEIKciNhd2VgF6Qk/YdgalIudQNKTMjniceWs2VTPA5mv1lFjij3svnwbnQ6aUTAzdORHxdv5/CeJEni/Zt89sFGRBEefmq8JPE0Gj2njmfRf5B0u1sKex0pOZmEhrafCHSonEBjlm45wpyR0hdJlFUpWfjpcnpPOoGTbyHZid4c29IVg878l0KnN7B8zZHWXyCIhNyehXN0NRkrAqmIb37Kf+RI023TWo+I65wKHEfX8of+O7yc29Y7sTnij2WwZa30vnr/BspaDe+9tpaAYOnMXw4dSOGu+4Yjl8swGNpWLyGzMmAXqMI2UIVdUP0fKzctuSc1uLpNw9mlfazSO6wIGIwiy3dIt+fcgJtvJT0nxWNlq+X4lq6kx/vRrufbm0Uk8OZc3G6oIOcvX0oPtl8r7sYCUL3VHgt/HcV9T+Cc70jln06Y8/NfrwLQQPzJLOJPZkkWb9+u82RnldLSEX0DeirJx2NoCXYXPvTWXmqEC3NzTaklykxbive4o8204rOKzZLlcS6nw4oA1AuBVMgUBsJ75dB9cCqqGmt2/NKXymJHyeK3Fd9xhXgNK6FwuyeF29ov6dNYAKo2OlDxmzPIwe2OcpynVqPw0lO21A1R13kaUQoqK1VUnsy+9IuCiLW3GrsgFbbB9R/41bLnMQoGgmaArlqBMtuW8uPOKLNtUWXboq+9tEJwbdqxdhuzWSIgCIIz8C3QjXrLmLtEUTwowbikQRDxDCwnsGsh/hFFWFgZyE325OimaPTa9ivDbAmPwSX4TSig9LArOWvbcybShAAggAHKvnNFX6jA5dYqFO56ij/1wFjT/qch//8jYummveQDbxugQm5dvzQw1MlQZtsSJQ7FzRjIF68cQ1thwbWZjdZj7kzgE2CTKIrTBUGwBDqAxa2Ii3c1gV0LCYgqxMZei06jIDfZi+xEb4qzXbmWL3gWJwmankPFGScyVgS141iaEYCLCFStd0JXZIH7vWX4vFhI8cee6PKvnThez1i6a/Cbko9DZC0WDvVJV6NOQJVrQ9khV5RZtigz7dAUW4EocN9nkwDQVsRfLey/gsn2YoIgOAKngVCxlUHi4nqLS37906Tna46oC6XFp7LOUigeppAj1FGMgAJ3uuMt9MWNWORC69/ckYHSG4ICCF4FrCxYhJsYxDDxXhRI84ELD/cCo0h6Un2rMhGRvcIGzgpH6CEOZKA4BuEqYlNELhuEFRjQM1acQQBhAIRG+YAA6UnSmWGGRnkDQvsZd0poigqNjFFTmy5OMqAnUdhJgrANATkBYixuBOImBuCED/Jm7rNhF/b9m4tr0ljDvf51o9FQoAT4XhCEHsBxYKEoipekLwVBuBe4F6Bnr17ItWY842VoqOR01SFSlHsoEVNBFHAlkhBxHJ5iLywumZi0QeyM9fdMqQxBAUrJYnfh13hY+DGkdj6WyJDMEswogtGIqKxDRGSf1VbOWp6gh7YvAzRDAPVVf3pP3LhFmMsGm9X8LfuJoZqxROt6gtEIMhkYJLRuak+fUQFEqQ9HXQjXVNwiUjjG79QIxQSIPejFFGwF50uvaymulB4HJoYyZyZwA3AIGCSK4mFBED4BqkVRfLG5a6JjYkXXLvNMG+kFFFZ6fKKL8Y8twj20HEEGDmIgh9ZbUhDvhbq69YU+zfHJmzMRDCJPPv6L2bEArL3r6PpoEs62ztyofZCX71opSdwG3lu6AFFZx1MjX8H9aQVOMxRULtdT9knbagEEO/B60wK7QXIql+v5z5TnkdvY8Z95SyQb67vL7kVUyHj6gWWSxQR4/8t5iHJBEhPXxiz6+HYQBB5/ZsXFrynsdQRNysG9dxnqMisy1wRRldS2U54fvTcbBHjsP79KNtaP351Frx5B/+pMIBfIFUXx8IX/rwaeNSNes8gUBjwjyvDvUYhXZBlyCyPKMhuSdwdz57A7sRO9WbGvZWfga4Glq4aoB1Iw6mUMF+/DhvZpDSYimiUAAKISCp/U4f64iPMcBZt1a7iRGe0w2usUQcSzXwkB43ORWRrJ3epL/g4fRH2HrLlrNSaLgCiKhYIg5AiCECWKYhIwCkiQbGSCiHtwBX6xRfjGFGNho0dTa0HWMV9y472pzHUEBOyG+dCmqf6/iMJeR9QDKcgsjZz/NAr7Z9xoD1fghiWAOQJwEQOUfqBHmy0iPJXKGpYidxExVPxv7xzY+qoIuTkT+yAlVakOZP4RhLqkfVyK/m3M3R14BPj5ws5AOmB+UwHq3XxumHUGG0cNeo2cggQPck97U5rhgmi8PlRXZmUg8v4ULJy1JH8RSV1Be71hLiQBLU+YLwCNqF5pYP5D09lqsxbv1/UUv++GLqt9zrN3ZHSoOcNmui08h16pIPWXUMpOXtsdJqkxSwREUTwFSGZ5ckNsIDGRvqzYsIfKXEfOnfGiKMkdg870u5CDvTVODtbkFlRKNcwWERRGIu5Jxca3jtRvw6nNaC878vptwLNCDj20ffnjkz2SRg8yhDHNsIBfjEvwfrmU0sUu1J34/3H3axkRl56VbBDeo45qig97kLPRH0PdlR8ZRwdrqmskbK3VCH8/F3LzKtoldgMd6rbq5mLPvFv6IzfacXRFLPlnvcwSAICnHhrDswulOSDSgHC1m4BMJGx+Bo4RtWT8HExVgnS9AS/lnzqAHuJABmhGNPkov3BvvENMP3zijjeFL3mgy1Pg8UQ5DhNq6KjLr6vhH9B6f0YrNw0R96URflcGVthxI4+S+UdwkwIAEBXpw/gx3aUa6iU8/vAY5O3cqqxDicCBY2nIZQL9eoZIFvOnlYeIjfZnUN8wSeLNmNmPV1+/pRmrJ5HgmVm4xFaS9XsA5cfbqyvxpYVAV6sDWPD27Sza+apZQmColFP0ujuqo9a4zqnG9a5KkF9fQuDv78rrb07H2bn5ejZBYcRnTAHdnkvAIayW7N/9GSM+hhtXdxGuqVWz8OEbCQ2RvuV9j27+TBzXQ/K4jelQIlCj1HDibA7D+kt3FDM1o5gtuxK4d95QSRT19OlsIiO9WfzFfPz9L727+E/Ow6N/GXmbfCje015HP6+sBLxaIdAHC76kNK/cbCEQtTJKP3Wlaq09DqNVeD5dhmArXR2FIMCg4V2Ijg3A1s78bd7LOXIkjchIH75deg/9+l95Q3CIqCHmP4n4Tyqg8pwTZ96Mpmi3JzJanonW1KixsrLglf9OldQH0MbaAoVCzl3zBmPfDq9JAx1KBAD2HElhQFwoCoV0Q/vu5334+7gwdmQMgFluu0nnC3jgvh+orq7j8y/nX3xDRU+vxWd0EcX7PC50B2o9r3w2l9c+n8eU2wcQGPrP3WTSzH7YXPKmaqkU+EpU1SqeG/8mpXnlvLn9ebxmuuM+ywWZKdZtokDlb06Ufu2MdbQG71dK8Opmx7AJscxdeCMKC9OXbqJYbzH2wps38+fO/7Dsz0d59YOZ3PHACAIkMEU1GkQ2bjyNi6sdb70zk0cfG4uVlQKFvY6QOZl0eSQFQS6S/GUYad+Hoqtq/Ye5VlmfDwjwd+XJhePMHmsD9hd8KZydbJkza4BkcS+nw4nA3iOp2NtZMbRvBEP7STMjKCiq4s+NJ7nrtsFMnxxHWLB507aKCiVPPfELO7Yn8MabtzL9GU/shiRRddqNrNUBtDVzvHtTPMpaNbPvHcGStY/x07ZneOL1m5k4oy9L/nqMQaPqW527zKpskwA0oKpW8dzkN/m8cAniY3pkC0WivzSt/4C9kw0TIicwtm42ToGWeLxUwJw3BxPbN5Thk3rSb0RXuvQMxDfIDXsnG4SrJlAuJTe7jCfu/YG8nHK8fZ0ZMDSKPv3DTLYSu5yNG05jvHAyNbq7Nz1n6en+3wRc4yrI3+TN2bejqTLB2r22VnPx3z17BDJmdIwk43VoZE4zZVJP/HycJYl7OR3qKLFCIWPq2J5UVKl4aeEEft94kj2Hze9m4+hgTWZOGbe6O/DoPaN4//PNpKSbZ1ul1xv55OPNZGYXUTVyP8FWXQhwv5kPxbYbgu5cf5qd608jCAKhUd7EDQwnbkA4IRfq4V/8+HYO7UxkqWoRRo1A9VYH2io0aqGOHF0enjYeVBgrMfTUY2ltiVbdtjru2qo69m48w2jbOKJH9OC46iDpylR6hPbl3ucm4XBZh5+/lh/gy9f/anX8ooIqnrzvB95ZPAdvXxcCQjz4btVDrPhhH5v/OonejFLuosIqjh/LoNI5iSO6NRh7l6NKtidrZSDqYtNdmwwGI2cT8rC3s+Lk6Sy2bDtncqzG2NlZ8cuqw9x2az9efmMttUpNyxeZQIeaCej1RkrKanBxskWhkJvd0KOBmlo1kaH/nNmPCpemA42VlQL7G3JQGqoZ4XULE8fHMWCg6R1+RFEk7XwBq5bu5cTBVGqqVJSX1FCUX4F/iAfz/O5Dhgyvx0oQrNv2YTCWiQTtCaLEUIq9wo5nYx5n6iOmTV0Lc8v4u/h3jqsOYn3WD/GYGzv/Ps2Mvq8yMfp5Zg14nfsmfMjTt3/N3z+3/WR5eWktT933I/t2JDJ/2qfs3ZHI/Y+N4fvfH2bC1Dizloq/xf/IdvV3WBrsmWTzKFnfRZslAA08+99VrN94mpHDo7EwY1nUmHOJeSxZupu8/AqiIrypaqfGpB1qJgDw19Z4hvaNoH9ciGQiIIrw4VdbAZgyviddIqQQARG3kdkk2xWgSXfnTLIO37F1PPnkBBac+4YqM113Vy3dy6qle6/4unU3V7weL8Hj/lKKP/FoU2vyYy/G47zbgX5T+uIe48asZ6ex8dsd1FS0pU+AiPPsapwm11KzxY6sH0TeYQVWF9pkGQ1GqsqVVJnZJ6G6qo5Fr/+FwWDk64+3sOqnA8yYN5AHnxzH8DHdeObBtp49EPGfnI+ibxHlJ5xJWx/AXV/F8dTT7rzx2p9mjRXqG8ds35nA/feMYEC/MPbsa1tPhKZoWLqcScgjJtp0/8uW6FAzgQbe+XIzNbVq5DLphtcgBH9uPEVokAeWlqbrnyATCZ6dhd+4AkoOuXHms0C+/Hw7M29dzDdLdjJuvLQtrhqjPmtD+XIXbHuqcZlV2ebrK7fV8NMjq7iry0IOrTvOjGemtOHqSwWg/Id/rMk0ddK0fm9MY4++8rJavvpoC/OmfcoPX+5oWyBBJHhWNj43FlG8z520H0MoL6njtVf+YMjQKKbf2leS8VZUqjh0JI0xo7tJEq+Bcwl5xHTxu3p9ihl0uJkAQGl5LR9+u50+PaTt8iqK8NGFGUFEiCfnkvLbHENmaSDsznSco6vJ2+hD/qZ/OhRpNHo2S9jVpjlqdjqg8NbjNLYGfZFpv8KirBLeu2Mx/pGtvcM0LwD/JuWltZSXtn7mIiiMhM7LxLVnJXmbvMnf8M/v6+yZXL76YjsL7hnOli1nqDZz9gawZdtZXnp+Ci7OtlRUqsyOB3A2IQ8HB2uCA93JyCqVJGZjOqQIAGzdm4haI/3dpUEI3FzaXsqrcNAReV8qtr4qMlYEUXqoPc1Br07Fr85YeOlwvb2CbFIJwLT+d7nJrRHCjiEAbUVmZSDi7nQco2rI/t2fot1X1kms+eMYhw6mSiIAAAcPp1FaVktEuDdHjpnrEF1PZnYptbVqorv6/m+JANRvF7YHolg/22gL1p5qIu9PQWGvJ+Xb8HYsB24lokDJV+74PF/EloCVTJXNaZ+nuU4FQGGvI/K+NGz9VaQvC6LsWPO1BgUSnivR6QzcNv+ri+t5KRBFuPuhHygurZEsZmM6ZE6go2EXXEuXheeRWRpJWhx57QXgAqJaRtHHHiiwZKPNKuStL49vXXxEDsq2XncCYOmipcvCZGx86kj5NuyqAtAeSCkADRQWV7dLXOgUgRZx7l5Jl4eSMdQpSPyoC8psu2s9pEswlCuYIN5GnVCH9weWtKGDWoscstrFKdn+60oArL3q6PpYEhYOepK+iKDqXMcQ7I5MpwhcBY9BJYTflYYq34bEj6LQlLVf/bY5eOLLKPUkrGIEPF+2kOSz6vqIglOWh+lm7HPdCEAZ2XRZmIwgFzn/aQS16e11hPv/FyZ7DJpCXFxv8btlf0oaMzykwRlYusaVYcEe7KpYyYHKdfiK0QwU56DAfAG43BVYKkKjfMBoZPW53zhkvZPemoH01Q41KZaIWD8DsDzMQMdhTPO8jYzzEroNd/UBQXq3YSGglFVFH2FptGO4eB8OSJO0Db/gCpwqcWPU8AsNRlPNrFy9JGaoJw4O/67bsEkIemlFR7gQTpDIEdeAnnUlX3Omdh9hxv7cYJx24SSZBPEbnHal3vAVAJmMXophVBmrOG51ACeFF13Enm0cnsgh2VZOyepnANM8b0OQyRAlPMyFINS7AkvospvDaQ4V/oKbpQ/9tQuwkTlK5nggAghgVEj7O2uo8ZIybhvqxi7hX50JREfHit4+t0ka88NFtwGiJC6zMmsD4Xel49Slhu6Gcfz4YD5SToPf/3wuglHkmbu/lywmwHvf3gkGsd4VWC7i+Z8yrLtoKHrLHc351s5grtwFeHfZfYgKGc/c94N0Y/36jnpXYImcnD0GlBA0IxsPIZghwt0888I6SeI28PFbMxFlAo+8Iq1D9GevzAABHnxTurhfvDCD3jGBbZ4JdOYELmDhqKXLo8k4RNTQ1ziTbuJorod18BUYBEo/cUVfrMDj8TIU3q3xHLwetwFFfEYXEDwrm6rzjgwX78NS6AANsK5DOkWACxnlJ5Kx9tCQ8nU4oUhTRnqtMCplFL/nBgh4PlWGzO5qh42uQwEQRAKm5uI/OZ+yY66kfhOOgv89E1Sp+J8XAfvQWro+loxMYeT8J5FUn792nYqlRF+soORDVxQeejweK2vGDuw6FACZSMjsLLxHFFO024P05cGIxg4+5g5OhxKBPn1CWP7T/WY5/1xO99gAbpoS1+T3XGIriHooBV2tgoQPo1Dl/v+aTmqSrChb4oJ1jBa3BZWAyLAJDYebrj8BECyMhN+Vhnu/MnLX+5L9R4Dp2TAJCfR1bbfDPX1iAtsncCM6lAjk5JTj6+tM167SHZsMCnLn/gdH4ep6aZGP59Biwu7KQJVrS+JHUWjLO2YNgLko99tS+bsD9sNVuExTcfezE1FYyMwWAAfH9rEet7Zuukmr3Lq+j4NzTBWZKwMp2PLPQaDWEBnmxQ09pT2QdjF2iCdzpki/hLS0kPP0naOwULRv45cOJQKFhVVkZ5fRr2+oZDE3bTxNZYWKWbdd8GgTRPxvyiVoei6VZ5w4vzgCg6r1O6VPvDCZW2b3l3S20t5U/e6Acr8NjrdWkmWZRNRTcrNnAFHd/Hj+ren4+ktbqzz/jiHMnNUPK6t/ficKBx1RjyRhH1JL+rIQSva33R4uLbOE/z45kdumS/9hrahSce/swfSNlVZkrK0sCPRxZfaE3pLGvZwOJQIAR46k01cie3Codyv6efl+Jt8Uh5unLaFzM/EZXUzRHg9Sl4Yi6tr2EmSll3Dn/SN4b/FcvNvJ8+1ybO2sCAhxp1e/ULz9XEyIIFC6xAW7Sje+TV9MbfdMs5cAxw6k4uxmzzerH+KhZybg7CpNOfXvq48yb/4Qfvr5AaZMjcPWQ0fXhUn1Sdsl4ZSfNE10DAYjB46kcd/8Ybz23BRsbKRpCw9QWa1CJhN45bGJ+HhKl1Oysaof411T++Pl1j49LKFDikAakZHeuLhIV6O/dcdJjhXso8cTWbjdUEHBhgCyV/ubtJ78fcUhHrrjW2xsLPnqp3sZf1MvbGwtGS2hkYhMJvDAMxP4Zs0j/LH/Bf7Y/wLfrHmU3gPDKTb1xJtOYKxuFgG2QYz0GIdsY9sNUS9n6WfbUCjk3DSjLz/8uZC59w7HwgyzFoDS0hpWrTyCm5s9s+7rRt+XsrF2Ekn6IoLq8+adA9i1PwmAYQMj+XrRXAIlmsVUVNf7Bjg52PDWU1OwMvM1aMDGygKDqEdQ6Fg4Z7gkMZuiwx0lPh2fQ12dlmnTelNdVcfq34+aFEdmYcSlRwVufctxCK/lr6qTWNpb07VsKtpCyCXb5DFmZZTw6N1Lue3OITz69AQmTo0jONSTnKxSkhLablRyOe5ejsgUMgIuNLPQavV88tpatv992qy4PbpHMcqrD7XVddRN28+vX+00K975s7ns25HA4JHRWFlbkJyQj05rfi/E3349xIjxIfxa9i4iRvJ/6IEy0/weB8dPZ1NdU4ejgw11ai1OEuU1Gnv/qeq09O8ZzG6Tj8GLOLgr8QiqIKRrMouTNhDnOhl359FEh3qTkC5tuTV0MBGwsJDzyivTEEWRObcPZOWqI22MIGIXrMS9XxmucRUobIyoSy2pPOjL8KjhzBg5CYXMgs9ObeH0KdNFAOqnlz99u5vM9GJefHM6AC+9fSsP3/ktFSb46wmCQK/+oUye0Ze+Q6MoLaomI7kQJxc7XntiBefP5Jo1XgtLBe88/guPvHYz+7ecZePKtr62TfP959sJ7+JLeWkNDz87kaz7Sig0s3eeTqbky+SXwaGWWZ7PoL/Lhf++sBq12jyTGYPByP7DqWh1BiaNiZXsaK7BKLJ2Wzze7o5odfo2CsA/H3r3wArcAyuxtq93gNYpbanK9sfPtzvfbjzeLgIAHUwEdDoDq1YdYUD/esfexh5zV6OOarxHF+Letwwbbw0GjYyKU86UHnajJs0eRIHQMXIUo+vXWCGh0rSLkstlBAS6kZleQnCoBx6ejrz41nSeefinNlljDx/fnTn3jcA/2J3jB1N5/clfObIniXE39+bInmRKi6vNHqtOq+f86RyqK1U4ONuabQTaQG5WGW/8ZyV52WW8+dkc3v/6Dp6+7weThUBhpyPq4WR0VlpUf/bho/SjvL/oNt55dybPP7cSlaptFumX881PeykrV+LqbMcLT0xkwcIfqJPAH/H9JVsZ2iec156YjLuLHaUVzb2+IvbuKnLZRQXJTFgYf/FDr6q2oijdldJsF0qyXFBV2gACU4I86R0tsOOI+ealTdGhRADg1Kls1m84zcQJPTAam/8gCXIjzjFV7BG+pYDzBNxkpCbNjoxfvCg/6YJRc+m2ytYtZ7G1teLRx8YSKpEIGAxGfvlhH7/8sI+AIDcGDevCwGFRPPjEOD59b0Or41hbW3DsQAqvLPyZ3Kyyi1/fsPqYJONsTHWFEicXaeshUhLrl0AvPLLcLCFoEAArdw3JS8KpSdEDZTzx2M988OFs3n1/Fs8+8xtKM/z3yy6I3/uLN/P9Z3fwyN0jee+ztveKuBxRhH3H06mqqWPCiG4s++Nww3ewd1fhHlKBe3AFbsEVWNtrOc9BrHBp8kN/OccTshnZN9LsMTZHhxMBgK+/3sGA/mEYmjgZaOOrwqNvKW69y7Gw11MhOtKFEfz2egGakqv7x6/98zi2tpbcNmcgglD/i5OKnKwyfl22n1+X7cfDyxErKwUaTevWx5vWnJBuIC1QXaHEJ6B9nHZUSs0lQvDg7V9R00rvvsYCkLIknJqUf7LsubnlPL6wXgjefmcGCx/9yezfXVV1He9+uon3XpnOwaPp7D1kfpMbg8HI+l1nGTLUhb05ubg1+tAD1FVZUZLmSmmmCw/ddDs2ggcr1q1qMe7xhBzumjYAN2c7yiqlmcE1pkOKQG2thk8/20pwUP2ZcIWtHrfeZbj3LcPOX4VRL1Bx1pnSw+68eM99yBBYVtK6U4QrfjmIrZ0VXl5OFErU3upySorMn763F9WVKqJiA9otfoMQjBwfK4kANFBQUMnjj/1McJC7ZOJ9+HgGf/x9glk39zFbBKwd1XQZkUaO/14KqzTETrr0Q1+a4YKq4p87ve1NrW8OezalAK1OT1zXALYePG/WOJuiQ4oAwJ695ylRyAibn4ZLt0pkChFlri1ZfwRQdtwN/YUCHxky2nrW/7tvdmEt4T7x9UR1hQorm/Y9bKNSavh7det2dVojAA0UF1VTLLHAfvn9bmRm1PwKciNhA7KJGpaOIEB+gic9fPtz9rCMo0erkKIUW6PTc+BUBg627VPV2iFFwD64lvD5aVg669DVKije70HJEXfq8qVby6rboVnG9cDaZftZu2z/tR4G0DYBaC+0ZmxpeoSW0X1CEg4eKgoSPTizKZK6ShtO0OBkLV1V6bMft76fY1vpkCKgLrZGmWNH1ho3Ks85IRo6XE1TJ2bSEQTAVKwd1XQbm4xft2Jqy2w4uLwnxSnXrgeFuZgtAoIgyIFjQJ4oipPMHxLoVQpSlpre2LOTjo2a2utSAC6f+iduDyX1QBBGffse8GlvpJgJLAQSgevjN9nJNUVNLTuFr647AWhu6v//AbM8BgVB8Ad+BN4EnmhpJhAX11tcunSNyc/XFOHhF1xbU4skjOkFIqSlSFuhFRbhDYikJ0kbNzTKG0RIT5TYxVhiZ2A1teyz/ZYKfTFDDHfhRYQkcaH+dyYires0gHeAjAOV35OmOoCN6EEUs3Cnu9lxI4Lra1VSskrMjnUxZpAHjvb/vtvwx8AzQLNHnARBuBe4F6Bnz14IUndRaXAbljKueNnfkiIgSu5AccHBV0pXYJDUGbhhBqDUlzPT+0n0hR6SvrxSuwIbRT3ZbGNnwQZEjITIphAojEEuWGD+KQbqX1cBDBJu1Igm/vpNFgFBECYBxaIoHhcEYXhzjxNFcQmwBCAmOlZ86uGfTH3KJvlg8VxA5OmHpIv7/udzQYSnH1gmWUyA97+cB8Az9/8oadz3vpoPtM94pXAGbpwEHK24D32hBwv/+5tEo6znkzdmYlQIPPhWy8U3LeEZXEbPMedxcFPhLvYgQj6ThR/vAtaaHbuBJY9PR5TBgs9XSxbzu4emm3SdObeOQcBNgiBkAr8CIwVBWG5GvE7+H3L5LoCUSwCpsXFQ029qPENmn0CQiexf2ZMewkPYCNdv5r81mDwTEEXxOeA5gAszgadEUWyf1ridXJdcL9uAgsxIRN8sug7KQBBEzu0JI/lQEEbD9Z31by0drk4gPNKb2lo1hfmVksW0tFJgZ2dl0hHfa4WllQJtK88edESuFwFoPPXPT/bg9LYoVFX/P7L+rUWSTJIoirukqhF44rmJTJjcS4pQF3n57Vu5+6FRksZsb8ZPjcPRSbo3Y2iEl2SxWqIjCsDludimpv4Hf+9pkgDEhftJNMpLifBxIzbIp11iN6bDleKlpxYTGtH6wxWt4ejBVPoNikAmQZbb19+Ftz69XTJPvebo1TeUoaNjJIs35+6hdO8lnRFmRKQ3dnZX1rKbKwDRke3zph8aF87gXvUGtqFxOYy59wA+4SWc2xPG1m8GUJhm+vHyobGhzBjWQ6qhXkQuk/HqrBuxkP8PuQ0DpKUUESbxXevA3mQcnWzpFmu+h7u6ToeHpyOLvpqPm0f7mD/KZALdewUxYmw3yWLm5ZTz+oeziIqWxs69qlLF0u/vYeq03igubE1KMQOICPXii3dvJzZa2rvrsYRsXn9wHHc8XEuvsecpzXFmyzcDOb8/1Oy1f25pFU/PGM7QWOlcsqHebSrM2417bmzfjlgdUgTcPRxxcpbusFBxYRUpSQUMHBZldqzyslqefuBHtDo9i76ej6e3eeaXTRES7oW9gzXdegVJFj8ztRhbOyve/OR2QsLNn2kVF1eze/d5Hnl0DN8tvYeBIwMlWQL8vTUeZ0cbFr99G2+/MI2QQGky82pdHcuSXqTG4SBhinHk7h0p2do/r6QKuUzG2wsmEB0k3Q2s4XTjgtF9iPRtvx2KDicC6Rcq/0LDvST19j+wJ4mBQ6MQBLCzN+9IZmWFimceWEZtjZpFX8/Hx8+F2XcOlmikENs7+OK/R4wzvzoN/nldHZ1seOezOfgFmu+0u+KXA9TVabH3MuI68zj2PgYyl0aZlQMwGIws+WkPAIP6hvP9J3cwc0qbCuCuwMpWw9Dbj1GlOMso7/u5JXIhP74+lzH9zb8pAOSWVgJgY2nBJw9OwddNmhxIQx7DQi7n1VljkLdTr4sOJQIWlnICg9wpKa5m0rQ45tw5RJK4tnZWnDiSjrePM6++P5NuPcxfFtRUq3nmwZ8oLanh0+8XcOcDI4noIs16Niral4K8CgrzKwmP8pYkZk5mKXq9AYPeyJEDqchk5v/qKypUrFi3hcWpz1KpLYGdgyhPND9XsutAMonJ9SXQldUq9h421bm3PgE4fO5RHD1qOfxHTwIYCUBKdgnbJfLsKyivwWA0ojMY2Hc2AwcbKc79i9g7VHCo8Huyaw9xLDWXUK/2cYTqUCKg0xqYc9cQPDwdGTqiKzqdQZK4CoWMl96+FYD+gyIlm2HI5QJZacUXly73LrxRkrifvrOeQ3uSKMqv4K3npako0+uNLPt6F6t+PsgNA8IoMrV/QSMsnLSc9VqBVqejS/YcXnnwfno3msWYw9fL9vD31ngKiqr45I2Z+Ho7tzmGjYOaobcfw8pWx94VvclL9uDvved4+cuNdA3xYvZ4aTr76A1GPli1mzX7ztAjzJfkPFPPA4h4elbQv+85bpu1nd5DVnK66gfSq45RWFlDSkGpJOO9nA4lAgA/frP74r9LJHDZBaiuquOz9zde/L8Ud0Gonw1s3RDP6eOZAPToHcyAoeYbQqqUGupUWmxsrST1Qfz1x/2sXn4AGxtLJt1s3gdAbqsn8oEU5DY6kr+K4LsPTrBt61lee2M6sT3Mty87cSabr5ft4alXVlFaXttmIbgoADY69v3Wi/K8+mu/X3uYTQcS+eSX3Txw62B6RkqTgPxt1yl+3n6CQA8XhnRvfYJQEIz4+pYwZFA88+Zs4dabd9OrZyq1tTbs2tOD75eNpaZgFEOiQyQZZ1N0OBFISsznwN76TjFlJTWSxT2wJ4ltG+MBJNkqbODc6RyefmAZzzz0E4lncrnn0RsvZsvNoV4EpLcBq66q48+VR5gxb9Al/f7agszSQOS9qfWtwb4Jp+S8HFGEjz7cxJ4953nr7RnEdDP/w1VVXYdSpW2zEFwhAPn/XKO9MLtctfUUu46l8PpDE3GRqAlJTkkVO0+nMm/01QVWLjdQKZ4ky/gNd87bxLSb9hPdNZOSEie27Yhj6Y/j+evvQZxLCEGlsmFvQjo3hPlhY9k+lngdTgTgn9lASYm0fnJffLSZspIaszzlmuPU0QwWLljK1x9vIVoCI8+6Oi3W7eQF+PvPB7GysmDSLW1PuAlyI+F3pWMXpCTtxxBqUv/ZJjUaRd5/dz2HDqby+hvTsZFo/I2F4IXHJlz1sVcTgMt569utqDU6Xn1ggmTviZ+2Hicuwp9uwZfmciwsdIRF5jJq/FHm3beRNBZRIR4hJ9eTjZv78N0PE9i4uT9JyYFoNJe+bvvPZyETZPSPbJ825R2ubBjqC4Z2b0+gVMKZAEBtjZoP3/7b7N2Bq3F4n/nW1QB1Sg22du0jAjXVav787TA3z+7Pmt8OY2zC2r1JBJGQOZk4da0m45cgKuOvbI5qNIq8/dY6IiK9qaszr1FIYxqEwNG++bu2WqxotQAAKNVanvtsHbPHxWGhkKPRmV+mHZ9RwKm0fOaM7s0ry/8gKLSQkPB8/AJLUCiM1KksSUvyY3L3GdjLYvhs+58txqyp03A6M5+h0SHsPJtm9hgvp0OKAMC3X+xoFzPQowdTcXC8en+CjkBRYRXnTue0W/zVPx9k49qTrRcARAJvycEtroKctX6UHm5+39pgMHI+0fyejJejVGlRNtOBSE0FJ1jUagFoIDWnlNe/2SLZGO3s61h25FNcfdKZe28RMhnUVNuQEB9CRqoPRQVuiKLAnNgebeqHu3z3Saws2qdysMOKQKEE2evmqKlWt1tsqThzIoszJ7LaLb6yVoOytvWdfHzHF+A1pISCbV4U7pBm21IqrB3VHJMtQktNmwRAShwclQy78SR+AfUZ/IoyB04djSQj1ZfSEtNbwDew/Yzp26Qt0WFFoJOOg+fQYvzGFVBy0I3cde1zWMZUrB3VDLj3JDpEevEYK/OP/+tjCAotYMSYEyCIHNnflYxUXyor2qekvD3oFIFOrkomxwm6JYeKeGcyVwYhpZe+uTQIgKWdll7Gp3BUhAL/ngjIZEb6Dkygxw2plBQ5s3V9H2qq2/dgWXvQKQKdNEu+kMhh4VeqU+xJ+zEEjB1TAA4v7cmE+0Kl8f5rJXb2dYwafxQfv3LOnQ7h4J5uGK5TExKz3IbbSlxcb/H7b/+QNGbDiUMpnYHrXYGRzGX3YtzI6yduiZDBLotvcLfwZZD6HiyQJpkqhSuwmgqOyRaho4ZexoU4E0p4iCcIkJwtnXsvQGSgBwiQlPtP3GrxDBl8gREtQSzAVRjY5rhR/h6IAiSZXF3YREw/D5xs/3234bbTXqLz72mZ+bTHDVUiV2CASvLZLSzFUeHKLJ+nKcjWSPbymusKrBbrdwG01NCLx3BUXJgBXHDvNUr8jhYFQFbvCiyKRgoNf1BkWIu14Eew4lGsZb6YUtwuyurHbJByo8rEqp9/dSYQEx0r+tveLGnM976+A4Bn7vtB+pjt4QosIKkzMtS7I4sygScfa11n5qth5aahy2NJIArc4vgMtnI3HnteOmfgj98y3RX4aoVAXzx/K0YF3L1YOvdegG8fno7RAu7/bTkTBhwj0KuUsxmB7Dgei95guuIsmzcdUQG3rTbfHbmBX6bfyoCAwDbPBDpkxWAn1waFg47Ih1KQKUSSvgjHDvOPG0tFWyoBpabWeJ45Y3fh61bO5iO92HIkziwB6Gj8//lJOjELuY2eqAdTsXDQk7Q4AnVhxzHbvHYCIFIo/kW+bhVanR1/7B5AaZX0JjItjUEmiBhN7SzSCjrkTGDg8C7MmD+oXWLfOKknUTEda6/7WiOzMBJxbxrWXmpSvw1FmdVxtrmulQBYWWmZOPEQ+eJKnGV9+WXrsH9VAAREYn1yeXDQFoaFJbTrc3VIEejWK4g+A9unScXdj95IlAQn3P6/IMhEwu5Mxz5ESfqyYKqTrr0zcAPXSgA8PcuZMWMn/v4lBAjzCVI8jFbf9Ak+f2dHPO2lE02FzECfwFSeHLGO22/YjZO1hiq1dFZ7TdEhRcDO3gqlsvUlra0lIMgdZ1c7zrZjOW5Hw8npKm8gQST49iycu1WTtTKAilNXHgj6N7Gx+ueDprDSMWjGSaxstWYLgLWFglsGdG+Fa69IbGwa06btRRQF/vhjCB7CjQhXOWFYrqxj2fxbCXY1fXwA1gotw8PP8Z9RfzK9x2G0BgVFFfOIsP+Ok7nSGpheTocUAVs7K1RtqGtvLd3igqipriMzTdrOtVJw44RYHnhsDDGxAVd45JvDwEERvPv+LIYO63KZz4FIwLRc3PuUk/u3DyUHWm+57eRow7OPjWfy2FicJDqLD2Brbcnnz03n9gm9GHLrORzclBz8vYfZMwC1To+fmyNrX5jPhN5dmnx9LS11jB17lMGDz5Cd7cWqVcMpKWlZFFU6HRWqOlbcOZPuvm03GXWyVjGh6wmeG72G8V1PUVjtwpKDo/h0z3gKq2IYHhLGA33a1224QyYG7eytycsukzxu915BJJzOwShhB2NPbydcXO1ISjDv1Ny2jfEMGtqFj76+g7KSGvbuSmTvjkTOxZs33o0bTjNgYAQvv3ozlZUqtmyOZ8P60+i7nMV7eAmFOz0p2NK2A0FV1XWs23iaT96exWMP3MipMzns2p/E3oMpVFapTB5rWZWS5euPMua2LNRVJUTo51Ht5cDe/LSLZiCm8v32Y0wfGMvb88Yzf2RvPlm3jwPn62eE7u6VjB17FAcHFQcOxHDqVDhtKeY4nJlDXIAvP86bzsMr13EgPbvFa3RiDlW6NTwzahcyQeR0XiB70qLJr/5nR6ZaU38jfHzAQI7m53E0L69NP3Nr6VAiIJMJiOKFmYCEywEnFzt0Wj2xcUGsXXnE7HgKhYz+Q6MYPzWOuH5hfPXhJtw8HXFyscXZ2Q4nF1ucXGz589cjJJ1r3S9OFOG919fyWfACAoLcmHprX0LDvVj05joK8irMGu+i9zfQNdoXV1d7Zszsj8rvPOcdCyg94krOn36YUr10Limfxd/u4PEHbuSGnkHExvghAH9tOm3WWCscdnC2Kp2B7rcz2HM2g/xrqNPoOHA6w6y4NXUavt16mCenDqOLvyfPTh/Ba79uQ2mxjyFD4lGrLfnzz8EUFrbdzPNQZg4PDOmHnaUlH948gYdXruNYdtO/9yCPYoZEn6PAuBwBSw5nRrA3vQsVdfZXPLZBBOQyGZ+Mn8jkn3+irK6uzeNriQ4lAgDvfDEPb19n6qJ9mT53IKt/OmB2TGdXO15+fyYe3k7Y2lkR3sWH1PMFbY7j4+fCxFtu4MaJPS7pQPTgU+MxGIzUVNVRVamkskJFVYWyzdN6lVLDq8+u5LPvFmBlbUFImOdF52FzqKpS8f6763n73ZkUq/JJctxAXYobmb8E0qZD7Zfx54ZTdIn0YfyobshlMvrGhbD7QDJV1aa9Ud38K4gZmk72WV8m9R6Fwd2Iq5Mdof5uHIzPMLvY9Ne9p7ltWC9c7W2xVMhx8N5P39hTZGd7sm1bb9Rq08xmTubkU6fTodbpWXfmfBMCIBLpm8/QmLMEe5agVFvhJMzCzmICf53b3Gzcao2GtPJywlxdeWP3rnYriu1QImA0iui0elzc7HFxs+f3nw9KEldVq8YvsF7hx02JY9Wy/SbFKS6sYtfmMxTklhMV40eXbn4EBHuw/JvdrFi6R5JlRnZmKR+88Rd+ga44Otrw/Gs3039wJIs/2Ehtjek+CEcOp7F2zXH2GX7GMlLOf4a/zCeHdnPsqHl32A+/2IqnuwPfLd/HC09M5PvFd/LuJxs5fLxtcQXBSM8x51FVWXNiUxdePLyR2ePjyMgr57Hbh9Ezyo/Xvt5MtdL010CrN/D5+oMUV9UydWwhXWNPkpMZyvqN3RHNEEON3sBza7cglwm8N3Ucv586x/miEgTBSLfAbIZGn8PHtYJKpS1/H7uB46nhfHf7rBb1N6eqkhkrf2Xp1JsZFBjIhhRpLNIvp0OJAMCZk1n0GRRBeWktJw6nSxKzsXnG0sXbTF5qGAxGUpMKSU0qZP0f9UdWbe0s8Qt0kzTPsHdnIja2ltSptBw5kMrTL97E1z/dy/uv/8WpC87GpvDj2r+IeCybgj1enAwu4u13Z/H90t2s+PmAyXdZrVbPC2+uoa5Ox4KFP/Dw3SN575XprFl/ki+/34WmlZ2VQ+Nycfaq5eAfsRh0cnKLK/l0xR40Wj3n0gp465FJ/Pj67byweD0J6aYflFp/LJHYHsnkGM7hLhvM5D5PceLwGvLKqkyOCbAxof4DOiOuOy+NH8KHe99lcHQC7o41FFc58vvBAcRnBmMwtv6kYcPU/4/EBJ4cOJDXdu1CY5C+U3WH2x04e7I+qbJryxmMBmkOh6oveN2dP5vLtvXmrVkvR6XUkpLY9qVFS9RdsNE6dTyT++Yu4czpHN75dA7+gSY2oBBEfG/OQF+rIHe9Fx9+sJGPFm1k7rzBPP0f8xpK112wgaur0/H+Z5t54c01jBgSxeMPtK4Pg5WthpihaRRluJKf9E+LNI22/g2fmFHE/Bd/Ji23lCUvzqRHpOn9FH39ihgw4BypqX589J03JVV1fPfIdPzczC8EslToWJ2wiHzjfUztfxiNzoJf9gzls78nczI9rE0C0Ji/k85jrbBgdFj7bBV2uJlAckIeWo2e7evjJYtpNIooazV88f5G/s0DU1JRW6PmnZfX8Mevh8g1cdfEY2Ap9sEq0pYFY6ir/7VvWH+KtNQi7B2k9VzcdyiVhPMFKCxad4/pPjIFuYWBU1u60FySslqp5umP1jJuYFfOppk2E7C1rWP06GNUVDiwY0cv9Ho9D361hi/un8ao2DCW7TxhUlwAd8cq7hq1DUfbOmrrorjB6zneXX2WGo35ZqsVajU7MzK4uWs065OlXxJ0OBHQ6Qxs/uskqUnS3l3/Wnm41Zn6jkqyiTMOhYMO/8n5VCc5UH7s0r3vJIlf5wbKK5WtepybfwVB3Qs4fyCY2vKrV96JImzcn2jSeATByJgxx7CwMLB2bV/0+vq3fq1ay92LV6PVm74F6e5YxYLRWwH4ZssYSqv8+GxGMJ4OmdRoyk2O25g1iQksnjgJZ2trKtXSemR2OBEA+OGL7ZLH/PnbPZLHvF4InJqLzNJI1qoAOpI9WONk4PkD7VsV169fIr6+ZWzd2puKy/z/pBKApdtGU1LtDOi4a7kU5jkivX3TqdbYsjNDxthlP0guANBBRcCcLHhz6LTSJ1SuBxwianDrU0HeJm/UxR3Lav3yZGB7ERRUSFxcCufOBZOSYn5jmAaaFgBpcLGu5bYe+4n1zuZYXgjflPmQUVkpWfzGmCwCgiAEAMsAb8AILBFF8ROpBtaJ+QgKI0EzslGXWlKwtWPZhF9MBqa7XZIMlBoHByWjRh2npMSJffukafMO7ScAAkaGhyQwNfoYAiKrzvZjR3o3SWI3hzkzAT3wpCiKJwRBcACOC4KwVRTF9j332Emr8R5ZhI2XhqQvwxB1HWsj6GIycGsU7bVEkcmMjBlzFEEQ2by5r2RGoO0lAL4O5cztuZdQ12LOFfnzc/wgylTtf6pTMnsxQRDWAotFUdza3GPi4nqL33/zuyTP10B7mGy2uyGohKaocMEYVYDU1KKLX6uhlI3C+/gRwyBxnklxpTAFvSJmiCcFmgT+LHqBIMYTLkyTJG5k4JXGnTnGZZSwhRBhIS5CH5PiRvl5gAwSC+vj6sRciowvAiJestexEExbXnT1rjcwTSgpQRS1aI0/ozWuRMAeK/n9KISRVz292BTRHh44WV8jo1FBEIKBXsDhJr53L3AvQM+evdrBaPRCPEnjirRrAq0dGqI2jisicoI1yJDTiynmPZ8Aoky68RoxsLf8G6xwIVg+wZyq5Uu5YAhqvHAaudJwmBLjFjzk43BS9DHdjvyCIahoJaIz5lKsfgkQ8bR+DYXMH9HUYl6h3stBZnWcWtViDMY8rCxGYmezAJnMCWh7Cz5BZtpYzJ4JCIJgD+wG3hRF8aop0Zjo7qK/lcRGo9/cAbSD0aggtI/RqEzgqUeWSxr3g8/mIMoEHn9mBQCuseVEzE0j889Aiva3/XhrAx+9NxujQuCRV1ZKNVSefsWNZOE39m3oQW6a6WO7nCWPT8dgCfOWrcbZvpbbx+yirMqBlTuHYDSavhRaNm86opXII7uXcH/sTgQBvoofTrHKvOKiFWMnohOXoNFuoUJny98lsaTXtf44d1N81vs+4lzD/92ZgCAIFsDvwM8tCUAn/w5yKwNBN2WjzLWl6GD7JdxMwcpOQxprcRWiyU1rn7Ep5AYmDTyC0Siw/mAfswSgAZ0xV0IBEOntkYlScw8iVeyvCGNXRSR68dpt1Jn8Cgn1C5bvgERRFD+UbkidmIPf2DwsHHRk/BHUoToGAUSPScWIjkjZbNpruTWiVzyeLtVsOtybGpX5tlw6sX4JIIUAuFjV8lC37dwfswtBcMPZ/iO2lUdfUwEA884ODALmAiMFQTh14c8EicbViQnY+inxHlRE8SFPlDlXnk+/lrgGVhLYs4AgxmArSLcMaEy5YS/dw7I4nBBJRoH5W6LujlUXk4DmCICAkZF+CbzWZw1dXQpYmdoHW8tPUSjCzB6jFJgsQaIo7uNfKD+zs7dGWdvxW4lfa4wYCbk5C51SQc7GjmWkKsiMxE48j6rSihCn9rlP1Im55Op/IKfYjQNnu5gd759tQGs8rV+jWGXa8XNfuwrmRe4nzKmEc+W+LE8eSKnagQXd5HCVdKVckGEQjdjItMgQURpN8zpoDR1r8/gyrG0s+eK3+xl/c+9rPZQOTxoHsQ9Ukv1XIAZ1xyoEDemTi5N3LWc3RSJH+jezwkJPOp8iw4oNB/sgmunR37gOwEv2OhayprcBr3YHVMj0TA05zou91+JpU823iUP5OH4MperWtSxXCHKeDPPk+eAN3Op1rK0/Qpvo0CJwx8OjsLa2YP8O0w6NtAZH5/a1c/43qKOGeNZTleJI2amO0zUI6pOBXUamUZzqSkFieyQDRYaOOoWaAoIsHkKpNq80+vJCoKvVAYQ6ufJ8n+E4W136nJFOhbx8w1omBsVztDiUl47ezOGiMFo7cbaVaZjpuYcg+acoZO5sK482+edpDR1WBKJ7BDBlVl++eHcD1ZWmm1c2R0CwOy++P5Nvf38YaxtLyeP/m5wU/sKAjsw1gfwbB4TaUnYQPSYVucJI/IarVwYqZDKig7xwtm+be3HXbllEdMnFl1twkMVc8X2ZIPDwsP4sHD6AfsH+WF7FdrytlYBpVeV42tizZ/p9PBjbH1crI3Mj9/N0r43IBSMfnh7D0vNDqdW1Xph62OfwfPAGejlko1XMZVDQLnq5jW319abQseaNgJW1BUajyOMvTeHwnmR2bzknSdzY3sGkJhVgbW3J3PuGM/amXqSnFPHO86svmo60BZlMIDjck+jYAGQyGX+ZYWAaEubJ6LHdEASh0QdMoKy0htW/XlF/dQkOUTVkyVKJEcdwsORSL8JePQK5d8Hwev/DWjW1tWpqatUkJOSzbafp1d03jYpl3i39qKlVU1ldR0W1isoqFTsOJnMm6R/X5YZkYPKeYJRlVz8mrDca8XZ14NsnZqA3GMgtrSKvtIrckip+3xtPbumVzj/uHpUMHB5PTqYnccE3NbnCNooi3x88zrJ5t/Lg0P5o9HpO5RZwKCOHlSfOUKqsv8GYWgq86MReJgRH8VBsIH3dX0RjKGVLTnfWZvRE24YWyfZyNbd6HqeXQw45ahe+yB3B5IBR9BcsebrrzSj1anYUS+ex0ZgOJwJPvjKV0uJq3DwceFaiYp3Y3sG8/untbP/7NKMm9aC8tJZ3X/yDPVvPtclkJCDIjeFjuxMd60+Xbv7Y2tWvb7etP824Kb2orVZTU1NHTXUdtdVqKspq0bXCKjsjrZgTRzN4/D8T8fL+JwO94qf9ODrZUF3VtHGnIDcSOCsHe9GNaEYBl3bkPXk6m7fe/5tXXphKTHR9stBoFDkdn4MgmF5kuXZbPOk5pbzx5GQiguun+OWVStZua/QmFcSLycDkPSGtirvjZCr3V67mowduokuAJ10CPIlPz+ebDYeueKyFhY7RE4+grrNkx+beTL2/+UmtUqvjnl/WsOLOmQS7udAvOID8qpqLAuBqX2PyWYCc2ip+Td5CsN07gI4q/Wvsyi9Ea2ydnwJAN7s8bvM+jLWgZ11JLNsrumBERoayvgxcJsiY6j+AxOocCtTmmc42RYcSAWdXO4aNrT8xte3vU1hYmj+8BgGwtrZk1MQefPfJVjb8cRy9CWfIc7LKOH4oDTcPB6Ji/C9+vVuvIG4YGI6Dgw3yRg0+3np+Nbu3tm4mc/xoBvfMXcLdD4zgpptvQKXScPOMvsy8fSCJZ3M5dCCVQwdSyEz/pzbeLlSJtZeGHsYZyIWm22Tl5JTz4MJlLHx4DOPHdKeqWsUr/51KQWElW7adY8u2s+QXVLb5tTiTlM9dzyzn9Scn06OLH4Ig8PNHd3LiXA7rtsVzKvcoTt61nFgT3aZjwvEZBdzx/q989vA0vFwcCPFxY+Nbd7N67xl+3XmSkqr6D1doRD5OzirWrR6Euq7lZGO5qo4FP//BijtnATCtRzRO1la8uXkXg7rtxUJh4OtNY004DCSi0ixCb1OHj8MXRLhFsjj+11Zf3dM+m/k+B8nTOLO8sD+F2n9uAhnKIk5WpBHtGMi2olPtIgDQwUSgz6B/+g+Gd/FBNNO8s7EANKDR6E0SgAYS4nNIiM/hyw82MnBEV8ZM6sH2jWcuehfa2lni4GiDvYMNxUVtM6+sq9Py2Yeb2b0jkemz+vHmy2voERdM/4HhTJ7WmwX3jyA9tYj77vgWALug+juZpxh+1VSARqPnvUUbiD+Tg6OjDfsPpDBmdDfGjenO/DmDiD+bwy+/HuLw0bYZu5ZVKnn01ZU8Mm84W/YmYmNtwU2jY3n+oXEcKKrkcNVhDCVBtLUOPqekijve+5VnZ43kjZ+3MXVQN24bGcecUXFsPHqe5duOExJ+kJpqG/Jz3VsdN7eymrt/+YM+gf4kFhbz0viR/Hb3CHbkf8n+xBiTTgMO9k4hzCmPk6XTWbhvH8vG+LJi3Cxmb/qV7JrKq14r6nYx3+cgmXVufJU3DI14qZBnKYv5b/xPzA8ZxdzgEazPP4pBlMZ3szEdSgT6DokE4Nj+FN58dpVZrci6dvdnzr3D2fLXKVIS8khKyCcns1Qy81KNRs/OTWfYuekMllb/vIwqpRaVUktRgenutfGnskk8l4dOZ+DIwVSOHEyFRZsIC/fC1/8fezDbIBWaMkusne1bdYxl05YzKBQy9Hoj3y/byw8/7aVXjyDG3tgNe3vTtu70eiMfLd1xMe6xM9k4O9ow4u6zOLr4IepsMOUwTKVSzQvfb8RgFFm+/QS/7jzFjb0jmXtjb965dzSbSr/nXHwIbU2EJhWVklxUigjc/M3PPD4xCVcnCx4Z8DJ5RYc4k1/UYowGHC1VTA87SlKlF4vP2COiZd6WlSwbM6NFIRCMO0D3erMCAKA16tEZdfySuYspfv0Z59Ob9flH2/TztoYOIwJyhYze/cP469fDfPnBJrM/rOfP5kl6qOhqaFtpq90WmsolpKUWkdbouLBdkAplpi30bH1cvf6f11UU4cSpLE6cMr9Ba+O4VbVKtNbZnD7kQ0FxtckxDY1mgnqjkY1Hz7Px6Hl69iih3wgjGammuQ43RLWzrsXe/hgn07uAu4pyZduapswKP4yl3MBPSYMQL4hRja5lIRjgmYrC8CXIY/gqL7JJAYD69uQLvA+iE+X8nT+C+cEj2VRwXPLZQIfZIuwaG8APX+zg83c3SHK3vh5dhduCwl6HlbsWVVbHq3Nw9q1BYWWgJKN9uhy7+2WiUlpRVGBeTcSgrokIgsjW+EieX7eVvKrWC1asWw59PDNZnxVLUd2l5cQNQnCyJB/jZR/YAZ6pPBGzGVGIAZt3WhSAoc5pVBus+TlzNx7WTtzo3avtP2gLdJiZQNLZvP+pluHmYhtYnw9QdkARcA+uT2CVZUkvAnK5gYDgIlKT/M3qGmRrpaZPRArxmcFUKtt2zsJKruP2iAPk1jqzKbtpy7IanZaHd/11ydcaBCCp2ptw9w+waiaZ21gA1pTGsqa0B1DFL1m7oR2akXWYmcD/qhGoqdgFqRCNoMrueCLgFlxBdZEdWqX0RVh+gSVYWurJSPUxK07/yCQsFQb2JlxZYNQS00KO42yl4qfkQRjE1u18NBaAN05NBqHp31vTAlDPN2mb2VRgem+E5ugwItBJ27ANVqEutMaoaT+XXlMQZEbcAispzWyfpUBIeD4ajYL8HNMNOCwVOvpHJZGQ409xlXObrg11LGaEXyK78rqSXt26MujLBUBtaFocryYA7UmnCFyXiPVJwQ64FGjIB7SHCAiCkaDQQrLTvc0yC+kTkYKtlZY959rm4isXDMyN3E+lxpY/Mlp3qK2jCwB0isB1iYWLDgtHfYdMCrZnPsDHrwwbG63JuwIAcpmBQV0SSS/0Ires9TUGAGMDzuJvX8nPKQPQGJpezzdGCgGY5pbIZNekNo2zrUjmNtwa4uLixO+XSOw2HHXBwTdJQrfhKG9AaD+34dTW70U3RQ7x7JP/yBjDQtwIJCzcq95tOE06V2CA8DBPRAFSMktafvAFTvIJaioYwCtNfj8iuN5lNym39TEbyBZ/pJRd9OBL5MKlh3Ki/D0QG7kCN0etcTvl4md4yF7GRmg5097gCny29ARK7f0oZP2xsfxvi9cJxh0oDK8iCjHo5R9ckQOIdvVEkImk1uRio30HS/0G1BZ3orFccPExjrpfcdN/Ra18JCUWL7Z4civc3hdHS9tr4zbcFjqW4dX1SZmQjUyU48xld8T2enFbGdeInkpS8WHAVR8nCtDWI/+iaKRSPIYjscjk1lfkyEWhfpxXc+oSRSPVxjVYCMFYKXq2zulYAAQDGsNHCIIljvb3IpNdvfhJ1O0C3esgj0GwefvCLsCl1wgyETuZBj/jq9Tqt+NsOxMX+ylAKQAa1VrU+p+wsByIr8N9+AmtaEQrNy1P8q/OBGKiu4sBltL4zDfw3rd3IgLP3PODdDG/uQMEgacfWCZZTID3v5wHMoEnF5rnNhz5SApyayOJ70cBsOiTOYhygcf+0/qa9dbw8buzMCoEHnyzdW7DLr5VjLjzKId+707++aYtxL54YQYGS1jw+eomv98cnp4V3HrzbrZujyM5JfCK73/30HQM1nDb6lXNxujqlcsdfXez4sRATuW17lDTL9NvRVSso7buM/4qjuVkTdBVH99wFuBqlYAAn8XdUy8A6u04CZa4yP+p2NSIOjSiHgVybASLVvcfEFx/QmbVv80zgc6cwPWGIGIb2DGTgh5BF/IB2dLnA0JD8jEYBLKyTPcOHBaeQLnKjvj8q3+QG2MUy1Gql5JZ58bJmivFpzGtFQABERvtO5IKgDl0isB1hrWnBoWNEWUHrA/wCKqgqtgOjUrq+gCRsJB88vLd0WhNix3kUkyIawl707tgbMNaRGP4AlHUsq4klquti3wsK5nnc6hVArDA+yCW+g31S4BrLADQKQLXHXZB9UdplZkdSwQEmRG3gEpK22EW4OpSg7OzkvQM03cFhoUnotRacjQ7vNXXxHploRf3YGs9i3Jd81WFMozc5n2EOoMF3xUMblUpsNriTlzsZ1/83rUSAOgUgesO2yAVBo0MdVHHajPu7F2DwtJASTtsDYaG5COKkJFpWpWgh30VMd65HMyIQmdoXS7cWqFldo/9yAjGpoWuWSNdkgiyLmd1cW9qDU3/Xi7fBmy8C3AtBQA6ReC6wy5QVV8qLFkTP2loyAe0x0wgNLSAwiJXVCrThG9oaCI6g5wDmZGtvmZq16M4Wyuxlj+O0EyNP4CnRTXj3c5wqsafk7XNuRI3XwfQpAAI7VNt2RzXhQiERXnz2EtTkEnYGPN6RJAbsfWv69D5AK3E+QAHByUe7lWkp5s2C3CwUhHnn8HR7FCU2taJSKhLEcNCEtiVEYNc1rXZxwkYud37MFpRwariG2gqZ3BVAVCtbXoGIPcD5MC/Y4DboUXAytqCux8fw6c/34eXnwt2DuZNgf/taZbU2PiqkVmIHTcf0C5LgQIAk/MBoyJSkclEDmW2rkRYLhiY03MvlXV2/JnY/E5bkK0nz4crCLEp41DNjRi4sjvR1QTAUfcratVPTS8B9GdBHgIeexEcXwGZeY1KW6LDHCVuoGffUFQqDQ4ONjz8wiTs7K356NW1bFt3yuSYVtYWTL2tP3VKDX/9dgTfQFeiYvzq/3TzY/v6eP5e1XrHFplcICzSm249AonpEUB+TjlLv9hh8vgCAt3wD3DF08sJD08HPD0c8fB0ZP++ZFb/9o/bcIOdWGtnAiOHdSE81BN7e2scGv1Zunxfm63EGrC0kHPbhBsI8nXFxsoCaysLNFaZHDbswFDp33KAZnC1t+XJm4YQ7OlywXVZQCbA0doXyC1zorrm6m7FTWGt0NLL/zxWQn8O3v1fNHo91RoNNRoNh/NyeWnHdgyX1cmMjTiNn2MFnx8ag0bf/J1YqUvDTdiCg/Uonur+BU8LAkXqCjbkH+e79C0tlgK76U9hYTkQa92Jpm9Ogi0yuQvY3gbWU0H1A2Lt55ji0tQSHUoEgsI8eXHRLIoLKgmJ9Gb736dYsmgTVRWm9R2QyWXcOLknc+8fjoeXE1npxcx7cCQOjjZo1DrSkgpIOptHeivLgwOC3Fjw8Gh69QnBxvafN8j5s3k89/rN1NaqUdaoL/59/HA6hfmVLcZVKjX0jAtm8k1xl1iVVVXVYWEhJz2tmIy0YhSOBWirFGjLWjdN3HcwFT8fF6ZNjsPmQm8FjUZH/z6hONhbc/ZcLoVtdP7R6gz8tukEd07tz+wJvbFQyDlZfg5FoRX9goZhWVPE2dQCtK1wWW5Mea2K11Zt44GxA5g/ojdymQydUc2hlDIsdT3xcrKnqKq2TTGDXEtQyA0s2uPEudxdPDN4CB52dnjY2XGupJge3t6cKCi4+HiZYKR/QArH8kKIL7p6LUFfx0y0RvgxP5b/uumwllviZe2CxqDFSmZBuHV2kwLgYaFkrkc8tfIR+DrcDxV3NP0E+nhE9VYE6xtBsELUHqE9BAA6kAi4utvz2mdzsHOwJsTBm6/e38ifPx80OV6/oZEsePRGgsL+Oe6pUMj57pOtJJ/LIzOtGIO+bQ5GOVllvPPSHwwZ2ZXRE3rQo3cwMplAaUk1ggDevs7Y21tj72CNnYM1FeXKVolAeVktX362lZUrDjF7zkAmTOqJRqPDYDAyZlwsfv6uyGQCtTVqUtMLeMFqJZpWWJpptXp++vUgG7acYcH8wYwb3Z2cvApCgz2YMKY7VlYWlJTWcC4xj/2HUtm6o3W9COo0Or74bS/r95zjyfkj6dt9InV54fj29eCeqa5otHrOpuZzPCGH7YeSySoob1Vcjc7Ax3/vY8vpFF6bdSPh3u74q99hWh8fnhrmyPm8YvYmZLAnIYMzWYUYW6h2TSr2482t01Bprfnu5AmOFxTw2YSJ6AwGgp2dWT1zNpmVlaxJTGBNYgK51dW8uXsaCqHl98WGsm4cqgqhXF9A2YlveK/nXSTV5DI3ZCQzAoewInsPb2fbknhZE9MSnR1PZ97IY9FP4ydUXtUiRKz5EBQRYKxAcHoLsXwuGHJa81K2iQ4hAtY2lrz66Ry8fJ0ByEguxMnFFksrhUn+fVbWFtRWq/nxix04u9rV/3Gxw8HZlmMHUikpNN0EVF2nY+v6eLauj8fDy5FR42MpKqhk5+azJsdsoKy0hsUfb2blLweZefsAvv1qJ3V1WqytLQgK8SA01BM/f5dWCcAlcctree+jTfy57iT9bgjlp18PolDIiAz3IqarH92i/YgI82q1CDSQVVDOo++sZkTfCESjyK5jf+PubEdcdABxXQMYPziaorKaVotAAwk5Rcz68BfuHHkDfxw6S5VKTVyIH0OiQxgdG8E9N/Zj7GvfUlBR02IsVaNk4KnCAib9vPz/2jvX2CiqMAw/L61EEVQqaKAgl0iwBCUQL9yiBEzAoKKJGAgqMST8QURiYoA/+E9jUDTGaLhZooSKgNIY5SISEWpAbuHWFkqBtrRSLhYoWqD088dOoZa2292d7uy65/mzM5M5Z9+dnX3nnNlvvo+XBw5k8e5d9Mu4lxezspj88CPMHjacdQX5zF7/E61LbyvO14ZiBw5cOMnM3V/Q6bYOFF0qZ+IDI3mt92im9BpFTslW1pbm8ff1m70eq8kAteKnd/0YVjULrpeizktRxldtYgSBm0C7tHbMmDee8tJz5H6zgz2/F3GuMvyX2xJXaq5xaF+JTwqb58zpi+Rkb/O938rKi3y6cMON9ZqaaxTml1OYX95Cq/AcKTrNEe8JxtraOg4XVHC4oIJvv4ut4OWWnUdvLJ+tuszGvAI25hXE1Gft9ToWb7pZ1WlnUSk7i0r5MHcr3Tp3apUBNMWFKzUs3h36vEfPn+OD7dtYkLed4T170ql99MVSi6pvTiuWFW9iVclvvNRzJFN6PcX+quPsrzoRXce1oTqc9te0NjOC4E1AYuH876mLscaAI3WI1gCao86MbSX+XjSqa2vIPv4zOSd/pabOh7m8Vf/XCM4+B+bPcQjcBGIpBOJwJDq+GEA9nhFw+1jfDAASPE7A4XA0wqrhH38T8zgTcDhSHGcCDkeKE5MJSBonqVBSkaQ5folyOBzxI2oTkJQGfAY8AwwAJksa4Jcwh8MRH6LOMShpGPCumY311ucCmNl7zbUZMmSIZfucbbhvCmYbvqXf+mzDxT5nG+7rZRs+GXlm4Obo1yuUFbjwlH99AvTP7Art4PAZf/sd0LVrKCtwdWwxGo15sGN3OqZduxEH4AvpWbRLuzuu2YYzgYYRC2XAE413kjQdmO6tXhn0WN/YQ+viQxfqU78mB8mkN5m0QnLp7R9pg1hMoKnncm8ZVpjZImARgKRdkbpUUCSTVkguvcmkFZJLr6SIwz9juTFYBjRMpdID8HfM5HA42pxYTOAPoJ+kPpLaA5OA3DBtHA5HghH1dMDMaiW9AWwglAtpmZkdCtNsUbTvFwDJpBWSS28yaYXk0hux1rhWIHI4HImHixh0OFIcZwIOR4oTFxNIpvBiST0lbZGUL+mQpFlBawqHpDRJeyX9ELSWcEi6R9JqSQXeMW65hHGASJrtnQMHJa2UlFAVXyQtk1Qp6WCDbRmSNkk66r2GTQHd5iaQhOHFtcDbZpYFDAVmJLhegFmAj6FnbconwHozewgYRILqlpQJvAk8amYDCd38nhSsqlvIBsY12jYH2Gxm/YDN3nqLxGMk8DhQZGbFZnYVyAEmxOF9o8LMKsxsj7d8idBJmhmsquaR1AMYDywJWks4JN0FPAksBTCzq2ZWFaiolkkH7pCUDnQgweJgzGwr0DiB4wRgube8HHghXD/xMIGmwosT9kfVEEm9gcHAjjC7BsnHwDtAZKmTg6EvcAb40pu+LJEUeUGBOGBmp4AFQAlQAVwws43BqmoV95tZBYQuaMB9YfaPiwm0Krw40ZDUEVgDvGVmkSXnjxOSngUqzWx30FpaSTowBPjczAYDl2nFcDUIvLn0BKAP0B24U9IrwapqG+JhAkkXXqxQBco1wAozWxu0nhYYATwv6QShadZoSV8HK6lFyoAyM6sfWa0mZAqJyNPAcTM7Y2bXgLXA8IA1tYbTkroBeK9hHy2NhwkkVXixQjWhlgL5ZvZR0HpawszmmlkPM+tN6Lj+YmYJe7Uysz+BUkn1T7qNASIrdhA/SoChkjp458QYEvQmZiNygane8lRgXbgGbZ5tOMrw4iAZAbwKHJC0z9s2z8x+DE7S/4qZwArvglAMvB6wniYxsx2SVgN7CP1jtJcECx+WtBIYBXSRVAbMB94HVkmaRsjIJobtx4UNOxypjYsYdDhSHGcCDkeK40zA4UhxnAk4HCmOMwGHI8VxJuBwpDjOBByOFOdfGEJkKQA1KeQAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "show_heads(model_path=model_path, name=name)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fe8d0305-1d97-4475-ac9b-6c9e98077e68",
   "metadata": {},
   "source": [
    "The well is roughly in the middle of the model.\n",
    "The yellow cells have a constant water level and \"feed\" the model.\n",
    "The well causes a cone of depression."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8588027f-8d91-4849-892a-68925fb89259",
   "metadata": {},
   "source": [
    "This function shows the head the well location over time:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "668fe6e6-fcd5-4e20-ab18-4de21c20dad9",
   "metadata": {
    "tags": [
     "remove-input"
    ]
   },
   "outputs": [
    {
     "data": {
      "text/markdown": [
       "```python\n",
       "def show_well_head(model_path, name, wel_coords):\n",
       "    \"\"\"Plot head at well over time.\"\"\"\n",
       "    sim = get_simulation(model_path, name)\n",
       "    gwf = sim.get_model(name)\n",
       "    heads = gwf.output.head().get_ts(wel_coords)\n",
       "    _, ax = plt.subplots()\n",
       "    ax.plot(heads[:, 0], heads[:, 1])\n",
       "```\n"
      ],
      "text/plain": [
       "<IPython.core.display.Markdown object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "%%include make_model.py\n",
    "start_at = 'def show_well_head'\n",
    "end_at = 'ax.plot'\n",
    "import_module = True"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f606676f-d460-46ac-baa9-185bb0c46f16",
   "metadata": {},
   "source": [
    "This is the head at the well over time:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "7b3b60f6-10b0-4e10-b1fc-dd8c57d6873d",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAAsTAAALEwEAmpwYAAAk+UlEQVR4nO3deXxcdb3/8ddntuxb2yRt0xW6t1BaQqsi60VsVUREuUXcAC0ocPEuXrjen9e7+bgq4oKivSBcQEBAAUEpFhdWQehCoftCaZt0SdJmX2cy+f7+yDTEkDbTdtIzy/v5ePSROWdOZj6n38w733zP95xjzjlERCT1+bwuQEREEkOBLiKSJhToIiJpQoEuIpImFOgiImki4NUbjxo1yk2aNMmrtxcRSUmrV68+4JwrHew5zwJ90qRJrFq1yqu3FxFJSWa263DPachFRCRNKNBFRNKEAl1EJE0o0EVE0oQCXUQkTQwZ6GZ2t5nVmtn6wzxvZnabmW03szfNbH7iyxQRkaHE00O/B1h0hOcXA1Nj/5YCPz3+skRE5GgNGejOuReA+iNscjFwn+v1F6DYzMYkqsCBtuxv4ZtPbaQjHB2utxARSUmJGEOvAKr6LVfH1r2LmS01s1Vmtqquru6Y3qy6oZ07X3ybN6obj+n7RUTSVSIC3QZZN+hdM5xzdzjnKp1zlaWlg565OqTKiSMAWLXzSH80iIhknkQEejUwvt/yOGBvAl53UEW5QaaXF7ByZ8NwvYWISEpKRKA/CXw2NtvlPUCTc25fAl73sConlbBmVwORaM9wvo2ISEqJZ9riL4BXgOlmVm1mV5vZtWZ2bWyT5cAOYDtwJ/DlYas25vwZZbR0dfPs5trhfisRkZQx5NUWnXOXD/G8A65LWEVxOGdaKaUFWTy8sooLZ48+kW8tIpK0UvJM0YDfx6cWTOCPm2tZvUtj6SIikKKBDrD07JMoL8ziq796g6aOiNfliIh4LmUDPS8rwG1L5rH7YDufuetV9jV1eF2SiIinUjbQARaeNJJlnz6dbTWtXHDr8/zsxR20h7u9LktExBPWe0zzxKusrHSJugXd7oPtfP2J9Ty/tY6inCCXzKvgwtnlLJg0goA/pX9niYj8FTNb7ZyrHPS5dAh0AOccK3c2cO/LO/nDphq6unsoyglyxqQS5k8sYf6EEk6pKCIvy7PbqIqIHLcjBXrapJuZsWDyCBZMHkF7uJsXth7gT5trWLWrgT9seme+ekVxDlPL85lWXsDJpXmML8ll/IhcRhdlE1RvXkRSWNoEen+5oQCL5oxm0ZzeOeoNbWFer2pg495mttW2srWmlZffOki4+50zTX0GY4pyGFeSQ0VxDmWF2ZQVZFFemE1ZYRblBb1fs4N+r3ZLROSI0jLQByrJC3H+jHLOn1Het6472sO+pk6q6tupbuigqqG97/FfdhykrrWLSPTdw1EF2YHekC/IYkReqO/fyLwQI/KyKMkLMjKv97mS3KDG8EXkhMmIQB9MwO9j/Ije4ZbB9PQ4Gjsi1LZ0UtPcRW1zJ7Ut/b62dLFxbzMH28JHnAdflBOMhX2IkrwQI3JDFOUGKcoJUpwbpDgn1Pe4KCdIUW6QgqwAZoNdxFJE5PAyNtCH4vNZX+97xhBXF4hEe2hoD9PQFuFgWxf1beG/+newLUx9a5jdB9t5s7qRxvYIXd2Hv7CY32e94Z4T/KuwL84JUpQb6v3ab31BdpDCnAAF2UHyQn79MhDJUAr0BAj6fZQVZFNWkA0UxPU9nZEoTR0RGtsjsa9hGjsiNMfWNXaE+56rbwuzo66Npo4IzZ0RjjQxye8z8rMCvQGf9U7QF2YHKcgOUJgdoDDn0OO//mVQmN37NRTQMJFIKlKgeyQ76Cc76Ke8MPuovi/a42jpjPT9MmjsiNDSGaGls5uWzgjNHbGv/Zar6ttp6eymObbd0LX5KMjuH/qxXwJZAQqyA+RnBcnPDryzHPtFMHNMAVkBHTQW8YoCPcX4fUZxboji3BATRx799/f0OFrD3TR39IZ739fOfstdf72+ubObPQ0dtHR109rZTUdk8Pu5XnnmJL5x0ezj3EMROVYK9Azj8xmFsSGYY9Ud7aGtK0pLV4TWrm5aOru5/sE11LeFE1ipiBwtBboctYDfR1Guj6Lcd34pFGYH/2pev4iceDr6JQkRCvgU6CIeU6BLQoQCPsK6x6uIpxTokhBBv3roIl5ToEtCZKmHLuI5BbokRMjvI6JAF/FUXIFuZovMbIuZbTezmwd5vsTMHjezN83sNTObk/hSJZlpyEXEe0MGupn5gduBxcAs4HIzmzVgs68Ba51zpwKfBX6Y6EIluWmWi4j34umhLwC2O+d2OOfCwEPAxQO2mQX8EcA5txmYZGblSMYIBXyDXm5YRE6ceAK9Aqjqt1wdW9ffG8DHAcxsATARGDfwhcxsqZmtMrNVdXV1x1axJKWg33fEK0iKyPCLJ9AHuxbrwK7Yt4ASM1sL3AC8DrzrKlDOuTucc5XOucrS0tKjrVWSWFZAB0VFvBbPqf/VwPh+y+OAvf03cM41A1cCWO/FuN+O/ZMMEfSbxtBFPBZPD30lMNXMJptZCFgCPNl/AzMrjj0H8AXghVjIS4bQmaIi3huyh+6c6zaz64EVgB+42zm3wcyujT2/DJgJ3GdmUWAjcPUw1ixJKOT3E+1xRHscfp/umCTihbiutuicWw4sH7BuWb/HrwBTE1uapJJgoDfEI9Ee/D7d5ELECzpTVBIi5O/9UdJMFxHvKNAlIbJi9yHVTBcR7yjQJSGCsR66ZrqIeEeBLgkRCijQRbymQJeECGnIRcRzCnRJiKAOiop4ToEuCdE35KIeuohnFOiSEFmxHnpEPXQRzyjQJSGC6qGLeE6BLgkR0rRFEc8p0CUhDh0U1SwXEe8o0CUhDh0U1SwXEe8o0CUhsnRikYjnFOiSEO8Muei+oiJeUaBLQrxz6n/U40pEMpcCXRJCJxaJeE+BLglxaAy9K6JAF/GKAl0SIuj3EfAZ7RENuYh4RYEuCZMT8tMRVqCLeEWBLgmTq0AX8ZQCXRImNxTQkIuIh+IKdDNbZGZbzGy7md08yPNFZvYbM3vDzDaY2ZWJL1WSXU7QT0e42+syRDLWkIFuZn7gdmAxMAu43MxmDdjsOmCjc24ucC5wq5mFElyrJLnckJ92DbmIeCaeHvoCYLtzbodzLgw8BFw8YBsHFJiZAflAPaCuWobJCfnp0JCLiGfiCfQKoKrfcnVsXX8/BmYCe4F1wI3OuXdNSDazpWa2ysxW1dXVHWPJkqx6h1wU6CJeiSfQbZB1Ay/Y8UFgLTAWOA34sZkVvuubnLvDOVfpnKssLS09ylIl2WnIRcRb8QR6NTC+3/I4envi/V0JPOZ6bQfeBmYkpkRJFTmhgAJdxEPxBPpKYKqZTY4d6FwCPDlgm93A3wCYWTkwHdiRyEIl+fXOQ9ehExGvBIbawDnXbWbXAysAP3C3c26DmV0be34Z8F/APWa2jt4hmpuccweGsW5JQjnB3oOizjl6j4+LyIk0ZKADOOeWA8sHrFvW7/Fe4MLEliapJifkp8f13rUoO+j3uhyRpDVcnZ64Al0kHrmh3hDvCEcV6CJAfVuYbTUtbK9r5a3aNt6qa+WtulaWnDGe68+fmvD3U6BLwhwK9PZIlBKPaxE5UZxzHGgNs622he21rWyraWVbbQvbalo52Bbu2y4n6OfksjxOn1jClLKCYalFgS4JkxPq/XHSgVFJR845apq7+sJ6W20r22tb2FbbSmN7pG+7guwAU8vyuWBmOVPL85lSls/U8gLGFGbj8w3vsSUFuiRMbvDQkItuciGpraUzwtaaFjbta2Hz/mY272thS00LLZ3vdFaKcoJMK89n8ZwxTCvPZ2pZAVPL8ykryPJsUoACXRIm59CQi3rokiKiPY5dB9vYvL+Fzfua2bS/N8Cr6jv6tinICjBjTAEXnzaWaeUFvT3usgJG5YeSbjaXAl0SJqffGLpIsmnr6mbTvmbW72nq63lvqWmhM3bbRJ/BSaX5zB1XzJIzJjBjdAEzxhQytig76YL7cBTokjD5Wb0/Tm1d6qGLt1q7utmwp4l1e5rYsLeZdXuaeKuuFRe7aMmIvBAzxxRwxcKJzBhdwMwxhUwpy0/52VkKdEmYwuwgAM0dCnQ5cdq6ulm3p4l11b0Bvn5vE28faOsL7/LCLE6pKOIjp47hlIoi5lQUeTrOPZwU6JIwRTmxQO+MDLGlyLFxzrHjQBuv727k9d0NrNndyJb9zfTEwntsUTazK4r42GkVnFJRxOyKQsoKsr0t+gRSoEvCZAd9BP1GU4cCXRKjMxJlza4GVu5s4PWqBl7f3dj381WQFeC0CcV84LwpzJtQwinjihiVn+Vxxd5SoEvCmBmF2UGaFehyjNrD3aze1cCrO+p59e2DrK1qJBJ1mMG0sgIWzxnNvAnFzJtQwpTS/GGf151qFOiSUEU5QfXQJW7RHscb1Y08v6WOF7fV8WZ1E909Dr/PmFNRxFVnTuY9J43k9Eklfcdo5PAU6JJQBTlBmjt1UFQOr66lixe21vHc1t4Qb2yP4DM4dVwxS88+iYUnjeT0iSV9s6Ykfvofk4RSD10Gs7exg6fX7+fpdftYvbsB52BUfhZ/M6Occ6aXctaUUZTk6b7yx0uBLglVmB2gur7d6zIkCTS0hXli7R5+vXYva6saAZg5ppC/v2Aa588oY9aYQo2BJ5gCXRKqUD30jBbtcby0/QCPrKri9xtqCEd7mDWmkH9eNJ3Fc8YweVSe1yWmNQW6JFRRTpDmzojuWpRhOsJRfrm6ip+9+Da769spzg3yqYUTuKxyPLPGvut+8TJMFOiSUIXZQSJRR0ckSm5IP17prrkzwl0vvs19r+ykoT3CaeOL+eoHp3Ph7HKyAql9Gn0q0idOEqrvbNGObgV6GuuMRLn/L7v48bPbaWyPcMHMcq455yQqJ5boLzMP6RMnCXUo0Bvaw4wuypxTrjPJy9sP8LXH17HzYDtnTR3FTYtmMKeiyOuyBAW6JNio/N6pZwdauzyuRBKtuTPCfzy5kUfXVDNxZC73XbWAs6eVel2W9BNXoJvZIuCHgB/4mXPuWwOe/ypwRb/XnAmUOufqE1irpIDSgt5raSjQ08uGvU1c98Aaqho6uO68k7nh/Kkpf6nZdDRkoJuZH7gd+ABQDaw0syedcxsPbeOcuwW4Jbb9RcDfK8wz06FAr2tRoKeLX62u5muPr2NEboiHl76HykkjvC5JDiOeHvoCYLtzbgeAmT0EXAxsPMz2lwO/SEx5kmryswJkBXwcaA0PvbEkNeccP/7Tdm79/VbOnDKS25bMY2SGX80w2fni2KYCqOq3XB1b9y5mlgssAh49zPNLzWyVma2qq6s72lolBZgZpQVZ6qGnOOcc3/7dFm79/VY+Pq+Ce65coDBPAfEE+mBzkNxhtr0I+PPhhlucc3c45yqdc5WlpTqYkq5G5SvQU93tz25n2fNv8amFE7j1srkE/fFEhXgtnlaqBsb3Wx4H7D3MtkvQcEvGKy3I0kHRFPbE2j1895mtXDKvgv++eI7mlaeQeAJ9JTDVzCabWYje0H5y4EZmVgScAzyR2BIl1WjIJXVt2NvETY++yYJJI/j2pafq4lkpZsiDos65bjO7HlhB77TFu51zG8zs2tjzy2KbXgI845xrG7ZqJSWU5mdR3x4m3N1DKKA/1VNFRzjKdQ+soTgnxO1XzFfbpaC45qE755YDywesWzZg+R7gnkQVJqmrojgH52B/UycTRuZ6XY7E6bvPbGHnwXYe/MLCvumnklr0K1gSblxJDgDVDboueqpYs7uBu//8Np9+zwTeN2WU1+XIMVKgS8KNK+ntlVc3dHhcicTDOcc3n9pEaX4WNy+e6XU5chwU6JJwo4uy8Zl66Knij5tqWb2rgRsvmKr7eKY4BbokXCjgY3RhtnroKSDa47hlxRYmj8rjssrxQ3+DJDUFugyLcSW5CvQU8PzWWrbUtPCVC6bq5KE0oBaUYTFuRA676jWDNdnd/5fdlBZk8aFTxnhdiiSAAl2GxZSyfGqau2ju1A2jk1VVfTvPbqllyRnj1TtPE2pFGRZTSvMBeKu21eNK5HB+8dpuDLh8wQSvS5EEUaDLsJhS1hvo2xXoSamnx/HommrOm17G2OIcr8uRBFGgy7CYMCKXkN/H9joFejJ6vaqBmuYuPnraWK9LkQRSoMuwCPh9TB6Vx7YaBXoyWr5uPyG/j/NnlHldiiSQAl2GzayxhWzY2+R1GTKI57bU8p6TR1KQHfS6FEkgBboMmzkVRdQ0d1Hb0ul1KdLPnsYO3qpr4+ypumZLulGgy7A5paIIgPV71EtPJi9t673949nTdNewdKNAl2Eze2whZvBmtQI9mby6o55R+SGmxmYiSfpQoMuwycsKML28gNW7GrwuRfpZs7uBeRNKdGu5NKRAl2FVOamENbsa6I72eF2KAAdbu9h5sJ35E0q8LkWGgQJdhtUZk0bQFo6yaV+L16UIsLaqEYD5E4o9rUOGhwJdhtWCySMAeGXHAY8rEYCNe5sBmB07YC3pRYEuw2pMUQ4zRhfwp821XpciwOaaFsaPyNGNLNKUAl2G3Xkzyli1s4GmDl150Wub9zUzvbzQ6zJkmMQV6Ga2yMy2mNl2M7v5MNuca2ZrzWyDmT2f2DIllf3NjDK6exwvxuY/izc6I1F2Hmxn5pgCr0uRYTJkoJuZH7gdWAzMAi43s1kDtikGfgJ81Dk3G/hk4kuVVDVvQgnFuUENu3jsrbpWoj2OaeUK9HQVTw99AbDdObfDORcGHgIuHrDNp4DHnHO7AZxz+uRKH7/POHdaKc9tqSPa47wuJ2PtPth70+7Jo/I8rkSGSzyBXgFU9Vuujq3rbxpQYmbPmdlqM/vsYC9kZkvNbJWZraqr05/fmeTC2aOpbwvzylsHvS4lY+2u7w30CSNzPa5Ehks8gT7Y6WQDu1kB4HTgw8AHga+b2bR3fZNzdzjnKp1zlaWluo5EJjl/RhmF2QEeXVPtdSkZa1d9OyW5QQp1hcW0FU+gVwPj+y2PA/YOss3vnHNtzrkDwAvA3MSUKOkgO+jnorljeXr9Plp0n1FPVNW3M2GEeufpLJ5AXwlMNbPJZhYClgBPDtjmCeAsMwuYWS6wENiU2FIl1V16+jg6Iz08vW6/16VkpN317UwYqfHzdDZkoDvnuoHrgRX0hvQjzrkNZnatmV0b22YT8DvgTeA14GfOufXDV7akonnjizmpNI+HV1UNvbEkVHe0hz0NHUwYofuHprO4Thdzzi0Hlg9Yt2zA8i3ALYkrTdKNmfHphRP5z99uZG1VI6eNL/a6pIxxoDVMd49jTJECPZ3pTFE5oS47YzwF2QHufHGH16VklEN3jSoryPK4EhlOCnQ5ofKzAnxq4QSeXrePqtg0Ohl+dS1dAJQVZntciQwnBbqccJ9/3yT8PmPZ8295XUrGqD0U6OqhpzUFupxwY4pyWHLGBB5eWdV39qIMr9rm3kAfla9AT2cKdPHE9edPwe8zfvDHrV6XkhFqWzoZkRciFNBHPp2pdcUT5YXZfP59k3j89T19N12Q4VPb0qXhlgygQBfPfPncKZTkhvjGk+txThftGk61LV2UKtDTngJdPFOUG+SfPzidlTsbeGLtwKtJSCIdUKBnBAW6eOqyyvHMHVfEfz+1kfq2sNflpK3G9jAluSGvy5BhpkAXT/l8xrcuPZWmjghff0JXixgOkWgPbeEoRTm6ymK6U6CL52aOKeQrF0zjqTf38Zs3NPSSaM2xe7kq0NOfAl2SwjVnn8Rp44v5+hPrqWnu9LqctNIYC/TiXAV6ulOgS1II+H3cetlcuiI93PCL1+mO9nhdUtpoigV6oXroaU+BLknj5NJ8/ufjp/Da2/V8Z8UWr8tJG00acskYCnRJKh+bV8Fn3zuRO17YwdPr9nldTlrQGHrmUKBL0vl/H57FvAnF/OMv32D9niavy0l56qFnDgW6JJ1QwMeyT59OcU6Qq+9dyd7GDq9LSmlN7Qr0TKFAl6RUXpjN3VeeQVtXlKvuWakbSx+Hpo4IuSE/Qb8+7ulOLSxJa8boQn766flsr23l6ntX0RGOel1SSmrqiKh3niEU6JLUzppayvf+9jRW7qxn6c9X0dWtUD9ajQr0jKFAl6T30blj+falp/LitgNc98DrRDRH/ai0dEYozFagZ4K4At3MFpnZFjPbbmY3D/L8uWbWZGZrY//+LfGlSia7rHI8/3XxbP6wqYYbH1KoH42OcJSckN/rMuQECAy1gZn5gduBDwDVwEoze9I5t3HApi865z4yDDWKAPCZ906iq7uH/35qE12R1dx+xXyygwqqobSHo4wt1v9TJoinh74A2O6c2+GcCwMPARcPb1kig/vCWSfxzUvm8KcttVz5fytp6+r2uqSk164eesaIJ9ArgKp+y9WxdQO918zeMLOnzWx2QqoTGcQVCyfyvcvm8trOej5916s06DrqR9QRiZKrQM8I8QS6DbJu4P3C1gATnXNzgR8Bvx70hcyWmtkqM1tVV1d3VIWK9HfJvHHc/qn5bNjTzKU/fZndB9u9LilptYe7yQ0NOboqaSCeQK8GxvdbHgf81UWrnXPNzrnW2OPlQNDMRg18IefcHc65SudcZWlp6XGULQKL5ozm/i8s5GBbmEt+8mfWVjV6XVLSifY4OiM96qFniHgCfSUw1cwmm1kIWAI82X8DMxttZhZ7vCD2ugcTXazIQAsmj+CxL7+P3Cw/S+54hWc27Pe6pKTSEemdt69AzwxDBrpzrhu4HlgBbAIecc5tMLNrzeza2GafANab2RvAbcASp9u4ywlycmk+j3/5TKaPLuSa+1dz5ws70I9fr/Zw70HjHA25ZIS4Wjk2jLJ8wLpl/R7/GPhxYksTid+o/Cwe+uJ7+IdH1vLN5ZvYuK+Z//n4KRk/rfHQ5RJyM/z/IVPoTFFJGzkhPz+5Yj7/dOE0fr12D59Y9jJ7MvxKje1hDblkEgW6pBUz4/rzp3LnZyrZdaCdj/7oJV7dkbmHcw4FuuahZwYFuqSlC2aV8/h1Z1KUG+SKn73K//357YwcV+8bctEYekZQoEvamlKWz6+vO5Nzp5fxH7/ZyJfuX9N3955MceigqIZcMoMCXdJaYXaQOz97Ov/6oZn8YVMNF/3oJdZVZ85t7Q5NW9SQS2ZQoEvaMzO+ePZJPHzNe+mO9nDpT1/mvld2ZsQQjA6KZhYFumSM0yeW8NTfncWZU0byb09s4LoH19Cc5re26wv0oMbQM4ECXTJKSV6Iuz53BjcvnsGKDTUs/sGLvPZ2vddlDZv2rkMnFqmHngkU6JJxfD7j2nNO5pfXvhe/z1hyxyvcsmJzWt40oz0SJeAzQgF91DOBWlky1vwJJSy/8SwunT+O2599i0/89GV21LV6XVZC6W5FmUWBLhktPyvALZ+cy0+vmM/Og+18+LaX+MVru9PmgGlXdzTjL3+QSRToIsDiU8aw4itnM39iMf/y2Dqu+flqDrR2eV3WcQt3O0J+fcwzhVpaJGZ0UTY/v2oh/+/DM3luSx0Xfv8Fnnpzn9dlHZdwtEfj5xlELS3Sj89nfOGsk/jt372fcSU5XPfgGq57YA0HU7S3Hu6OqoeeQdTSIoOYVl7AY196H1/94HSe2bifC7//Ak+vS73eeri7h2BgsLtISjpSoIscRsDv47rzpvDbG85ibHEOX3pgDdc/uIb6FLopdSSqMfRMopYWGcL00QU89uX38U8XTmPFhv1c+P3nU6a3Hu7WGHomUUuLxCHo93H9+VP5zQ3vZ3RRNl96YA1L71vF/qZOr0s7oq5oD6GApi1mCgW6yFGYMbqQx798Jv+yeAYvbKvjgu89z89f2UlPT3LOW4909xDyaww9UyjQRY5S0O/jmnNOZsVXzua08cV8/YkNfPJ/X2FrTYvXpb2Lpi1mFrW0yDGaODKPn1+9gFs/OZcdda18+LYX+d4zW+jqjnpdWp9wd48OimaQuFrazBaZ2RYz225mNx9huzPMLGpmn0hciSLJy8y49PRx/OEfzuEjp47ltj9tZ/EPk+cKjjoomlmGbGkz8wO3A4uBWcDlZjbrMNt9G1iR6CJFkt3I/Cy+/7ence9VCwh393DZ/77CV3/5hucnJEWiPQTVQ88Y8bT0AmC7c26Hcy4MPARcPMh2NwCPArUJrE8kpZwzrZRn/v5srj3nZB5/fQ/nffc5fv6XXUQ9OmiqHnpmiaelK4CqfsvVsXV9zKwCuARYlrjSRFJTbijAzYtn8LuvnMXssUV8/dfrueQnf+aNqsYTXkuXDopmlHhaerA5TwO7Gz8AbnLOHfFokJktNbNVZraqrq4uzhJFUtOUsgIe/OJCbrt8HvubOvnYT/7M1x5fR8MJOtPUOUe4u4csDblkjHhauhoY3295HLB3wDaVwENmthP4BPATM/vYwBdyzt3hnKt0zlWWlpYeW8UiKcTM+OjcsfzxH8/hqjMn8/DKKs6/9TkeXrl72Oeud8deX2PomSOell4JTDWzyWYWApYAT/bfwDk32Tk3yTk3CfgV8GXn3K8TXaxIqirIDvL1j8zitze8nyll+dz06DouXfYy6/c0Ddt7hrt7b6mnIZfMMWRLO+e6gevpnb2yCXjEObfBzK41s2uHu0CRdDJzTCGPXPNebv3kXKrq27noxy9x06/epK4l8bNhFOiZJxDPRs655cDyAesGPQDqnPv88Zclkr4OzV3/wOxyfvTHbdzz8k6eWrePG86fwufPnERWgq69Eo7d9FpDLplDLS3ikcLsIP/64Vms+MrZLJw8gv95ejMXfv8Ffr+xJiH3NFUPPfOopUU8dlJpPnd9/gzuvWoBQb+PL963is/c9Rpb9h/ftWEO9dCzFOgZQy0tkiTOmVbK0zeexTcumsWb1Y186LYX+bcn1h/zNMe+HrqGXDKGWlokiQT9Pq48czLPffU8PrVgAvf/ZRfn3PIsd76wg87I0V30K6Ix9IyjlhZJQiPyQvzXx+bw9I1nM29CCd9cvokLvvc8T6zdE/f8dY2hZx61tEgSmz66gHuvWsD9Vy+kMDvIjQ+t5eLb/8wrbx0c8nsV6JlHLS2SAt4/dRS/veH9fO+yuRxs7eLyO//C1fesZNsRbqrRFVWgZxq1tEiK8PmMj88fx5/+6VxuWjSD196u54M/eIF/eexNapvffW/TiA6KZhy1tEiKyQ76+dK5J/P8P5/H5943iV+trubc7z7H93+/ldau7r7twuqhZxy1tEiKGpEX4hsXzeYP/3AO500v44d/3MbZ33mWu196m67uqKYtZqC4Tv0XkeQ1cWQet18xn6VVjXxnxWb+87cbueultzmlogiAoHroGUMtLZIm5o4v5oEvvIf7r17IyPwQv9uwH1APPZOohy6SZt4/dRRnTjmTp9fvZ2tNC6PyQ16XJCeIAl0kDZkZHzplDB86ZYzXpcgJpL/FRETShAJdRCRNKNBFRNKEAl1EJE0o0EVE0oQCXUQkTSjQRUTShAJdRCRNWCLuLn5Mb2xWB+wasHoUcMCDcoaD9iX5pMt+gPYlWZ2IfZnonCsd7AnPAn0wZrbKOVfpdR2JoH1JPumyH6B9SVZe74uGXERE0oQCXUQkTSRboN/hdQEJpH1JPumyH6B9SVae7ktSjaGLiMixS7YeuoiIHCMFuohImkiaQDezRWa2xcy2m9nNXtdzPMxsp5mtM7O1ZrbK63riZWZ3m1mtma3vt26Emf3ezLbFvpZ4WWO8DrMv/25me2LtstbMPuRljfEys/Fm9qyZbTKzDWZ2Y2x9SrXNEfYj5drFzLLN7DUzeyO2L/8RW+9pmyTFGLqZ+YGtwAeAamAlcLlzbqOnhR0jM9sJVDrnUupkCTM7G2gF7nPOzYmt+w5Q75z7VuwXbYlz7iYv64zHYfbl34FW59x3vaztaJnZGGCMc26NmRUAq4GPAZ8nhdrmCPtxGSnWLmZmQJ5zrtXMgsBLwI3Ax/GwTZKlh74A2O6c2+GcCwMPARd7XFPGcc69ANQPWH0xcG/s8b30fgCT3mH2JSU55/Y559bEHrcAm4AKUqxtjrAfKcf1ao0tBmP/HB63SbIEegVQ1W+5mhRt6BgHPGNmq81sqdfFHKdy59w+6P1AAmUe13O8rjezN2NDMkk9RDEYM5sEzANeJYXbZsB+QAq2i5n5zWwtUAv83jnneZskS6DbIOu8Hws6dmc65+YDi4HrYn/+i/d+CpwMnAbsA271tJqjZGb5wKPAV5xzzV7Xc6wG2Y+UbBfnXNQ5dxowDlhgZnM8LilpAr0aGN9veRyw16Najptzbm/say3wOL1DSqmqJjb2eWgMtNbjeo6Zc64m9iHsAe4khdolNk77KPCAc+6x2OqUa5vB9iOV2wXAOdcIPAcswuM2SZZAXwlMNbPJZhYClgBPelzTMTGzvNgBH8wsD7gQWH/k70pqTwKfiz3+HPCEh7Ucl0MftJhLSJF2iR2AuwvY5Jz7Xr+nUqptDrcfqdguZlZqZsWxxznABcBmPG6TpJjlAhCbqvQDwA/c7Zz7prcVHRszO4neXjlAAHgwVfbFzH4BnEvvJUBrgG8AvwYeASYAu4FPOueS/mDjYfblXHr/rHfATuCaQ+OdyczM3g+8CKwDemKrv0bv+HPKtM0R9uNyUqxdzOxUeg96+untGD/inPtPMxuJh22SNIEuIiLHJ1mGXERE5Dgp0EVE0oQCXUQkTSjQRUTShAJdRCRNKNBFRNKEAl1EJE38f+oj8vj4JLH3AAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "show_well_head(model_path=model_path, name=name, wel_coords=wel_coords)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d571fc75-8053-4bbd-98b1-297f63c9bc96",
   "metadata": {},
   "source": [
    "There are four stress periods.\n",
    "The first period is steady state.\n",
    "Periods two to four a transient with different extraction rate.\n",
    "Period three has the full rate and periods two and four have only a tenth of the\n",
    "rate. See generation of the input files for details."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cdb69577-5713-4c44-af88-f1abbf29d8e4",
   "metadata": {},
   "source": [
    "## Working with `pymf6`\n",
    "\n",
    "Finally, we will use `pymf6`.\n",
    "The objective is the avoid a drop of the water level at the well below a given \n",
    "limit."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "18003894-6b06-4e43-8239-4bbda08f5dc1",
   "metadata": {},
   "source": [
    "### Creating new input files\n",
    "\n",
    "We generate the input files again.\n",
    "There are no differences in the files to the model above.\n",
    "The only difference is the model directory:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "6df65903-126c-47e8-80cc-22b191a09ba6",
   "metadata": {},
   "outputs": [],
   "source": [
    "model_path=r'models\\pymf6'"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7748d444-caba-45a8-a583-c26be632ba67",
   "metadata": {},
   "source": [
    "Creating the input files:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "96c2ccda-369e-4e07-b944-e02bdddb536e",
   "metadata": {},
   "outputs": [],
   "source": [
    "make_input(\n",
    "    wel_q=wel_q,\n",
    "    wel_coords=wel_coords,\n",
    "    model_path=model_path,\n",
    "    name=name\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a518fc88-a3b3-4edd-bc70-08ce9af05a67",
   "metadata": {},
   "source": [
    "### Run the model with `pymf6`\n",
    "\n",
    "Now run the model with `pymf6`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "6128591d-110f-44bd-b8c5-110dd1e7899b",
   "metadata": {
    "tags": [
     "remove-input"
    ]
   },
   "outputs": [
    {
     "data": {
      "text/markdown": [
       "```python\n",
       "from pymf6.mf6 import MF6\n",
       "\n",
       "def run_model(nam_file):\n",
       "    mf6 = MF6(nam_file=nam_file)\n",
       "    head = mf6.vars['SLN_1/X']\n",
       "    wel_index = 44\n",
       "    tolerance = 0.01\n",
       "    head_limit = 0.5\n",
       "    upper_limit = head_limit + tolerance\n",
       "    lower_limit = head_limit - tolerance\n",
       "    wel = mf6.vars['HEADCONWELL/WEL_0/BOUND']\n",
       "    been_below = False\n",
       "    for step in mf6.steps():\n",
       "        if step < 21:\n",
       "            if head[wel_index] <= upper_limit:\n",
       "                wel[:] = wel[:] * 0.9\n",
       "                been_below = True\n",
       "            elif been_below and head[wel_index] >= upper_limit:\n",
       "                wel[:] = wel[:] * 1.1\n",
       "```\n"
      ],
      "text/plain": [
       "<IPython.core.display.Markdown object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "%%include pymf6_head_controlled_well.py\n",
    "end_before = 'if __name__ =='\n",
    "import_module = True"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cd46d845-8ee1-4bd9-a4f2-c763d8bdd754",
   "metadata": {},
   "source": [
    "The class `MF6` provides the functionality of `pymf6`.\n",
    "The dictionary `mf6.vars` provides all MODFLOW 6 variables.\n",
    "The head has the name `SLN_1/X` and the values of the boundary condition `WEL`\n",
    "is `HEADCONWELL/WEL_0/BOUND`.\n",
    "The name will change depending on model name and the chosen name for the\n",
    "boundary condition.\n",
    "\n",
    "`mf6.steps()` provides a generator that allows to iterate over all MODFLOW 6 \n",
    "time steps.\n",
    "`step` is the current time step in the unit specified in the input files, here\n",
    "`DAYS`.\n",
    "We don't want influence the last, which starts at day 21.\n",
    "We adjust the extraction rate by reducing o increasing the rate by 10%,\n",
    "if the water level at the well, found at index 44, is not within the upper\n",
    "and lower limit."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "aa034ee6-c813-4243-9341-f5f8f60e09e7",
   "metadata": {},
   "source": [
    "Running the model:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "99658a58-9276-478d-a5f5-ae2abed7b6ee",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "run_model(r'models/pymf6/mfsim.nam')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9c11cadb-f236-4e2d-885b-ce1f14358156",
   "metadata": {},
   "source": [
    "### Looking at the results"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d7b7ad46-d58a-41dd-9dc2-8000936524cb",
   "metadata": {},
   "source": [
    "Now the water level stays just around 0.5:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "e44a531e-1d7b-430e-a9b2-3bd96aae50b9",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAAsTAAALEwEAmpwYAAAnuklEQVR4nO3deXxc5X3v8c9vFo00I1mSJXmTvK8YGww4ZksIJCEBUgpJ0xS4NySQlNBCmubVJCW5t9l6k+aWNrdJQ8OlDc0e0oQSlpACN+xJAC94X0A2trUv1r6NNDPP/WNGshCyLZuRzyzf9+ullzTnnJn5HR/rq2ee85zzmHMOERHJfj6vCxARkfRQoIuI5AgFuohIjlCgi4jkCAW6iEiOCHj1xpWVlW7RokVevb2ISFbavHlzu3OuarJ1ngX6okWL2LRpk1dvLyKSlczs0LHWqctFRCRHKNBFRHKEAl1EJEco0EVEcoQCXUQkR5ww0M3sXjNrNbOdx1hvZvYtM6s1s+1mdm76yxQRkROZSgv9e8AVx1l/JbA89XUL8J03X5aIiJysEwa6c+5ZoOM4m1wD/MAlvQCUmdncdBU40b7mXv7u13voi8am6y1ERLJSOvrQq4G6cY/rU8vewMxuMbNNZrapra3tlN6srmOA//vMAfY1957S80VEclU6At0mWTbprBnOuXucc+udc+urqia9cvWEVs0tAWBvc88pPV9EJFelI9DrgfnjHtcAjWl43UlVlxVREgqwt0ktdBGR8dIR6A8BN6ZGu1wAdDvnmtLwupMyM1bNLVELXURkghPenMvMfgpcClSaWT3wRSAI4Jy7G3gUuAqoBQaAm6ar2FGr5szgly83kEg4fL7JenxERPLPCQPdOXf9CdY74La0VTQF5y4s44cvHGJnYzdn1ZSdzrcWEclYWXml6CXLqzCDJ/e2el2KiEjGyMpArygOsW5+GU/sbvG6FBGRjJGVgQ5w7bpqdjX2sPlQp9eliIhkhKwN9A+cV8OMwgB3PVVLshtfRCS/ZW2gR0IBbrtsGU/ubeXh7dM2SlJEJGtkbaADfPStizlnQRmf/cU2fre/3etyREQ8ldWBHvD7+Ncb11NTHuZD332Jf//tayQS6n4RkfyU1YEOUFkc4oE/v4jLVlbx5Yd3c+2//JbnX21Xv7qI5B3zKvjWr1/vNm3alLbXSyQcv9zawJ2P7aOpe4iVs0v4wHk1XLl2DjXl4bS9j4iIl8xss3Nu/aTrciXQRw2NxHl4WyM/fOEQ2+u7AVhTPYOLl1ZywZIK3rJ4JsWhE14gKyKSkfIq0Mc7dKSfR3c08+TeFrbWdTESd/h9xvJZxZw5r5Q11TNYW13KijklzCgMTmstIiLpkLeBPt7gcJwthzt58cARdjR0s6Ohh/a+6Nj6qpIQS6siLK0qZklVMUurIiysiDCvrJBQwH/a6hQROZ7jBXre9D0UFfi5eFklFy+rHFvW0jPEzoZuXm3tY39rHwfa+3lkexPdgyNj25jBrJIQNeVhasqLUl/Jn6vLiphTWki4IG/+GUUkg+V1Es2eUcjsGYW884zZY8ucc3T0D7O/rZ+6jgHqOwep70x+33K4k0e2NxGfMDSypDDA7BmFzEm93pzS0NjPyceFVBaH8OtWvyIyjfI60CdjZlQUh6goDrFh8cw3rI/FE7T0RqnvGKCxe5Dm7igtPUM0dw/R3DPE/v3ttPZG3xD6fp9RVRxi1owQlcUhKosLqCoZ/Tn5VVVSQFVxITOKApgp/EXk5CjQT1LA76O6LNndcizxhONIf5SW7ijNPcmgb0kFfltvdKyr50j/8BuCH6DA76OiuCAV8snwHwv+khBVxSEqiguYGSmgrChIwJ/1lxOISBoo0KeB32fMKilkVkkhayk95naJhKNrcIT2vihtvdFx34fHHrf0DLGrsZv2vsnD3wxKi4LMDCcD/lhfFZEQ5ZEgFZEQRQU6ySuSixToHvL5bCxwV8wuOe6248O/vTfKkf5hOsZ/DQzT0TfMoSMDvFzXRWf/MLFj3AahKOgfe9/ySAEVE8K/PBykLFxAWThIeeq7RvqIZD4FepY4mfCH5MndnqFYKvCjdPSPvO77kf5hOlN/DF5r76Ojb5j+4fgxXy9c4Kc8XEBpUZDySDLwy1OBX1qU/H7JiiqqSkLp3G0ROQkK9BxlZpQWBSktCrK4MjKl5wyNxOkcGKazf4SugWE6B0boGhyma2CEzv7U44FhOgeGaeruoSv1ePSDwA3nL+Br71s7jXslIsejQJcxhUE/c0uLmFt67BO+EyUSjt6hGO/95+foj8amsToRORENj5A3xeczSsNBQgEfsbjucCniJQW6pEXQ72MknvC6DJG8pkCXtAj47ZijakTk9FCgS1oEfGqhi3hNgS5pUeBXH7qI1xTokhbJLhe10EW8pECXtAj4fYyohS7iKQW6pEXQpxa6iNcU6JIWAb8xElMLXcRLUwp0M7vCzPaZWa2Z3THJ+nIze8DMtpvZS2a2Jv2lSiYL+H2MqIUu4qkTBrqZ+YG7gCuB1cD1ZrZ6wmafB7Y6584CbgS+me5CJbMFfaZRLiIem0oLfQNQ65w74JwbBu4DrpmwzWrgNwDOub3AIjObjeSNgN9HTOPQRTw1lUCvBurGPa5PLRtvG/B+ADPbACwEaia+kJndYmabzGxTW1vbqVUsGSnoN0Z0paiIp6YS6JNNbjnxN/frQLmZbQU+AbwMvOHWe865e5xz651z66uqqk62VslgAZ9a6CJem8rtc+uB+eMe1wCN4zdwzvUANwFYcnbj11JfkicCfvWhi3htKi30jcByM1tsZgXAdcBD4zcws7LUOoCPAc+mQl7yRIFGuYh47oQtdOdczMxuBx4D/MC9zrldZnZrav3dwBnAD8wsDuwGPjqNNUsGUgtdxHtTmrHIOfco8OiEZXeP+/n3wPL0libZJODzEUs4nHMke91E5HTTlaKSFkF/MsR1PxcR7yjQJS0C/uR/Jd3PRcQ7CnRJi4BPLXQRrynQJS2Coy10jUUX8YwCXdIikOpD17yiIt5RoEtaBH3J/0qaV1TEOwp0SYuxFrr60EU8o0CXtAhqlIuI5xTokhaj49CHNWuRiGcU6JIWAZ9a6CJeU6BLWgR0paiI5xTokhYahy7iPQW6pMXolaIahy7iHQW6pMXovVw0Dl3EOwp0SYugxqGLeE6BLmmhUS4i3lOgS1rofugi3lOgS1oE1Ycu4jkFuqSF7uUi4j0FuqTFWAtdfeginlGgS1qMjUNXC13EMwp0SQuNQxfxngJd0iKoGYtEPKdAl7QoSLXQh2NqoYt4RYEuaRHw+wj4jGgs7nUpInlLgS5pEwr4iI6ohS7iFQW6pE0o6GdILXQRzyjQJW3UQhfxlgJd0qYw6Ceqk6IinlGgS9qEAj6GRtTlIuKVKQW6mV1hZvvMrNbM7phkfamZPWxm28xsl5ndlP5SJdOFAj610EU8dMJANzM/cBdwJbAauN7MVk/Y7DZgt3PubOBS4B/NrCDNtUqGCwX8GrYo4qGptNA3ALXOuQPOuWHgPuCaCds4oMTMDCgGOoBYWiuVjBcKqoUu4qWpBHo1UDfucX1q2XjfBs4AGoEdwCedc/rNzjOhgJ8hjXIROaH4NN0iIzCFbWySZROreQ+wFXgHsBR4wsyec871vO6FzG4BbgFYsGDBSRcrmS3ZQleXiwgkQ7u+c4D9bX3sb+3nQHvy+/62Pm68cBGffNfytL/nVAK9Hpg/7nENyZb4eDcBX3fOOaDWzF4DVgEvjd/IOXcPcA/A+vXrdRenHFMY8GscuuSdoZE4r7X382prH6+29FLb2sf+tj4Otg8wPO7uoxWRApZWFXP56tmsrZkxLbVMJdA3AsvNbDHQAFwH3DBhm8PAO4HnzGw2sBI4kM5CJfOphS65bGgkzv62Pl5t6ePV1t7U9z4OHelntAfF7zMWzgyzpKqYy1bNYmlVMUurIiypLKY8Mv3jRE4Y6M65mJndDjwG+IF7nXO7zOzW1Pq7gb8FvmdmO0h20fy1c659GuuWDKQrRSUXDA7HqW1NhXaq1f1qax+HOwZwqeAO+IxFlRFWzSnh6rPnsXxWMctnF7O4MkIo4Pes9qm00HHOPQo8OmHZ3eN+bgTend7SJNskhy0q0CU7OOdo6Bpkb1Mve5p62Nvcy57mHg62H21xB/3G4soIa6pLed851SyfVcLy2cUsqohQEMi86zKnFOgiU1EY9DEcT5BIOHy+yc6li3ijPxpjX0sve5t62dvckwzwpl56o0dHVy+sCCdb3GfNY9WcEpbPLmFhRXhsvtxsoECXtBn9qBmNJSgq8O5jp+S3roFhdjb0sKOhm50N3exq7ObQuO6SklCAVXNLuPacalbNLWHVnBmsnFNCcSj74zD790AyRij1ETQaiyvQ5bTo7B9mR0M3O1LBvaOhm7qOwbH182cWsWZeKX90bg2r5s5g1ZwSasqLSF4DmXsU6JI2oeBooKsfXdIvGouzq7GHlw93seVwJ1sPd9HQdTS8F8wMc1Z1GTdsWMja6lLWVM+gLJxfdyBRoEvaFI52uWiki6RBU/cgWw518fLhTrYc7mRnY8/YnLXVZUWsm1/Ghy5Mhfe8UkrDQY8r9p4CXdJmtIWuWYvkVNR1DPD7A0d44cARXjzQMdb6DgV8nFVTykcuWsS5C8o4Z0E5s2cUelxtZlKgS9qE1EKXk9DUPchva4/w+/3JEB8N8IpIARcsqeBjb1vMuQvKOWPujIwcIpiJFOiSNoXBoydFRSYaiSfYcqiTp/a18fS+VvY29wJQHg5ywZIKPv72JVywpILls4pz9qTldFOgS9qMttB1x0UZNTQS5+l9bfxqRxNP72uldyhGwGesX1TO565cxSUrqlg5u0TXLaSJAl3SZnTYoqahy28j8QRP7W3l4e1NPLmnhf7hOOXhIFeumcM7Vs3i4mWVlBTqBOZ0UKBL2oyOPR9UoOel2tZefraxjgdebqC9b5jycJA/XDeP966dxwVLZhLIoisus5UCXdImkrrSbmBYk1XlC+ccz73azj3PHuD52nYCPuOdZ8ziT94yn0uWVynETzMFuqRNcUHyv1NfVC30XJdIOB7e3sh3nt7P3uZeZpWE+OwVK/nj8+ZTVRLyury8pUCXtAmHkl0uA1G10HPZM6+08b9/vZfdTT2smF3MnR84i2vWVWtoYQZQoEvaBP0+CgI++tTlkpOau4f4woM7eXx3CzXlRXzzunVcfdY8jVDJIAp0SaviUIB+tdBzinOOn75Ux9ce3cNIPMFnr1jJR9+62NOJHGRyCnRJq3CBnwH1oeeM/miMzz+wgwe3NnLR0gr+7v1rWVgR8bosOQYFuqRVcShAn1roOeFgez9/+oNN1Lb18VeXr+C2y5apeyXDKdAlrcIFfvrVh571djf2cOO9LxFLJPjhzefz1uWVXpckU6BAl7SKhAL0DinQs9nuxh7+5J7fUxwKcN8tF7JsVonXJckUaZyRpJVOima3Q0f6ufHelygOBfj5rQrzbKNAl7QKFwQYGNZJ0WzUF41x0/c2Ek8k+OFHN1BTHva6JDlJ6nKRtCoO+XVSNAs55/jr+7dzsL2fH3/sArXMs5Ra6JJWkVSXixudYl2yws821vGr7U18+j0ruXBphdflyClSoEtaRUIBYgnHcFz3RM8W7X1RvvboHi5YMpNbL1nqdTnyJijQJa0iqVvo9uvioqzxtV/tYXAkzv+6dq3GmWc5Bbqk1egtdDXSJTtsr+/iP19u4OOXLGXZrGKvy5E3SYEuaTU6E03P0IjHlchUfPvJWmYUBvj425d4XYqkgQJd0qosnAz07kEFeqbb19zL47tb+MjFizUlXI5QoEtajQZ614ACPdN95+lawgV+brpokdelSJoo0CWtyooKAAV6puvoH+ZXO5r44Pr5lEcKvC5H0mRKgW5mV5jZPjOrNbM7Jln/GTPbmvraaWZxM5uZ/nIl04210AeHPa5EjueXLzcwEndct2G+16VIGp0w0M3MD9wFXAmsBq43s9Xjt3HO3emcW+ecWwd8DnjGOdcxDfVKhisM+gkFfHSrhZ7RHtzawJrqGayaM8PrUiSNptJC3wDUOucOOOeGgfuAa46z/fXAT9NRnGSnsnBQXS4ZrKFrkG313Vy1dq7XpUiaTSXQq4G6cY/rU8vewMzCwBXA/cdYf4uZbTKzTW1tbSdbq2SJsqICOgfU5ZKpHt/VDMCVaxTouWYqgT7ZpWPHulHH1cBvj9Xd4py7xzm33jm3vqqqaqo1SpYpDQfp0rDFjPXsK20sqYywuFJTyeWaqQR6PTD+zEkN0HiMba9D3S15rzwcVB96hhqOJXjxtQ7NQJSjphLoG4HlZrbYzApIhvZDEzcys1Lg7cCD6S1Rsk1ZUYFGuWSorXVdDAzHuXiZAj0XnfB+6M65mJndDjwG+IF7nXO7zOzW1Pq7U5u+D3jcOdc/bdVKVtBJ0cy16VCyN/T8xRpVnIumNMGFc+5R4NEJy+6e8Ph7wPfSVZhkr/JIAdFYgr5ojOKQ5lDJJNvqulhcGaEsrIuJcpGuFJW0m1USAqCtN+pxJTLRtrpuzqop9boMmSYKdEm7KgV6RmrpGaK5Z4iza8q8LkWmiQJd0m5WSSEArb1DHlci4+1u6gHgzHm6OjRXKdAl7dRCz0z7W/sAWD5bE0DnKgW6pF1ZUZCg3xToGaa2tY+ZkQJm6u6KOUuBLmnn8xmVxSFaFegZpba1j2VVmmYulynQZVrMKgmphZ5BnHO82trHUs0bmtMU6DItqkrUQs8kXQMjdA+OsLRK92/JZQp0mRZzSgtp6h70ugxJaehKHoua8iKPK5HppECXaVFTHqZrYIS+aMzrUgRoTAX6vDIFei5ToMu0qE4FR0OnWumZQIGeHxToMi1GP9o3dA14XIkANHYPURDwUaEhizlNgS7TojoV6PVqoWeEhq5BqsuKMJtsvhrJFQp0mRZVxSFCAZ+6XDJEY9cg88oKvS5DppkCXaaFmVFdVkRdp7pcMkFz9xBzS9V/nusU6DJtFlaEea1dge415xxH+oapKFb/ea5ToMu0WVJVzGvtfSQSx5pTXE6HvmiM4XhCJ0TzgAJdps2SqghDIwmaenQbXS919Cfnd50ZCXlciUw3BbpMmyWVyfuGHGjr87iS/DYa6Gqh5z4Fukyb0fuGHGjTvOFeGg30cgV6zlOgy7SpKglREgrwamuv16XktSNqoecNBbpMGzNj1dwS9jYp0L3UOdaHrkDPdQp0mVar5sxgb3Mvzmmki1c6+ocpCPgIF/i9LkWmmQJdptUZc2fQF43pFgAeOtI/TEWkQJf95wEFukyrM+YmJyQenXFeTr+O/mF1t+QJBbpMq1VzZhDwGdvru7wuJW/1DI5QWhT0ugw5DRToMq2KCvycOW8Gmw91el1K3uqLxoiEAl6XIaeBAl2m3TkLytlW181IPOF1KXmpLxqjRIGeFxToMu3OW1jO4Ehcwxc90q8Wet5QoMu0O29hOQCbD3V4XEl+6ovGKC5UoOeDKQW6mV1hZvvMrNbM7jjGNpea2VYz22Vmz6S3TMlm88qKmFtayObDXV6XkneisTgjcUexWuh54YRH2cz8wF3A5UA9sNHMHnLO7R63TRnwL8AVzrnDZjZrmuqVLHXuwnK26MToadc3FANQoOeJqbTQNwC1zrkDzrlh4D7gmgnb3AD8p3PuMIBzrjW9ZUq2W7+wnIauQeo1g9Fp1R+NA6gPPU9MJdCrgbpxj+tTy8ZbAZSb2dNmttnMbpzshczsFjPbZGab2traTq1iyUpvW14JwDOv6LifTr3REUAt9HwxlUCf7HrhiTfmCADnAe8F3gP8jZmteMOTnLvHObfeObe+qqrqpIuV7LW0qpjqsiKe2qtAP51GW+gK9PwwlUCvB+aPe1wDNE6yzX855/qdc+3As8DZ6SlRcoGZcenKKn63v51oLO51OXmjL9VCj4R0Y658MJVA3wgsN7PFZlYAXAc8NGGbB4G3mVnAzMLA+cCe9JYq2e7SlbMYGI6z6aBOjp4ufakWeomGLeaFEwa6cy4G3A48RjKk/8M5t8vMbjWzW1Pb7AH+C9gOvAT8m3Nu5/SVLdnooqUVFPh9PL1P58xPl9FRLjopmh+mdJSdc48Cj05YdveEx3cCd6avNMk1kVCA85fM5IndLXz+qjN0O9fToD+qYYv5RFeKyml11dq5HDwywK5G3U73dOhLBXqkQIGeDxToclpduWYOAZ/x8LaJ59VlOvRFY0QK/Ph8+jSUDxToclqVhQu4ZEUVD29rJJHQtHTTbWA4TpGmnssbCnQ57a4+ey6N3UNsOazRLtMtGosTCijQ84UCXU67y1fPoTDo4/4t9V6XkvOisQShoH7N84WOtJx2xaEAV581jwe3No6dtJPpER1JUKgWet5QoIsnbjh/AQPDcR7c2uB1KTktGourhZ5HdKTFE+vml7FqTgk/femw16XktOhIglBAv+b5QkdaPGFm/LfzF7CzoUcTSE8jnRTNLwp08cz7z62htCjIvz57wOtSclY0phZ6PtGRFs9EQgH++wULeGx3M6+193tdTk6KxhIUBtVCzxcKdPHUhy9aRNDn49+eUyt9OkRH4mqh5xEdafHUrJJC/ui8Gn6+qZ6GrkGvy8k5QxqHnld0pMVzt79jGQD//JtXPa4k9yRb6OpyyRcKdPFcdVkRN5y/gJ9vrldfeprppGh+0ZGWjHDbZcso8Pv4P0+84nUpOSMWTxBLOJ0UzSMKdMkIVSUhbrp4EQ9vb2S37pWeFsPxBIBa6HlER1oyxscvWUppUZAvP7wL53Rr3TdraESBnm90pCVjlIaDfPrdK3nxtQ4e3t7kdTlZLxpLThAdUpdL3lCgS0a5fsMC1lTP4Ku/2q07Mb5JUbXQ846OtGQUv8/48h+uoaUnqmGMb1I0lgx0nRTNHwp0yTjnLSzng+tr+O7zr7GzofuE2zvnaO0dYmgkfhqqyx5jXS5qoecNHWnJSJ+/6gzKIwV8+ufbGE61NI/l20/WsuGrv+HPf7xlbNlva9u59q7fcudje8eWtfQMcedje/nF5qMzJXUNDPONJ17h3udfm1JdR/qi/M9f7siKWxUcPSmqFnq+UKBLRioLF/C1961lb3Mvdz1Ve9xtR+cmffHAkbGJpx/Z3sjWui7+/bcHiaeW/WxjHXc9tZ/P/OLoH4mHtjXyrd+8ylce2U1n//AJ6/r2U7X86IXDfPXRPcTix/9D47WjJ0X1a54vdKQlY12+ejbvO6eau56qPW7XS31n8h4w/cNxDqSuND3cMQAkZ73f39YHwMHUOudgb3NyrPv4K1OnMml1bWvf2GvsmEJ3ULo892rb2HuPNzgc5/MP7OCBl18/P2vXwDAvH+4C1OWST3SkJaN98erVzIwU8Bf3vUx/atTLAy/X8/S+ViDZf17fOcjbllcCsKOhC4BDRwZYUz0DgK11qWUdAyysCAOMhd2hIwMsrowQ8NmUJto4eKSfi5dVAPC7/Ucm3cY5x56mnpMeSz8cS4x9mqjrGOC2n2zhlZZeBofjfOi7L/GubzzD9vqu1z3n9wfa+cmLh/nUz7axLbWfHf3DrPvKE3wjddWtTormDwW6ZLSycAH/dN06Xmvv54sP7SKRcHzqZ9v4yL9vxDnHkf5hBkfivH1FFUVBPzvqexiJJ2jsGuTtK6ooCQXGgu7QkQEuWFzB7Bmhsdb4wSP9rJpTwpnVpa8L9MHhOLsau/nUz7by2V9sA5KB29A5yHkLyllaFRn7QzHR06+0ceU3n+P7vzsIJAN+NNydc2N/mJxzfP93B8c+QXzmF9tY/YX/ora1l4e2NfKr7U38wT8//7q6HpkwPn9H/dGrap+vbQdgV+PrPzmohZ4/dKQl4120tJJPXLaMX2yu587H940t31bfTV2qa2VRRYSFFWEOdwzQ2DVIwiWXLZtdzIG2fvqjMdr7oiyoCHPG3BnUtvYRTzjqOgZYWBFh1ewS9rcd7X75yiO7eO+3nuehbY3cv6WBlp4h6joHkq9bGWFJVTEH2/vZ19zLVd98bqzrpq5jgMd2NgPJ/vaReILr7nmBm7+3kVg8wY9eOMSZX3yMp/a1sqeply8+tIt3/uMzNHUP8uDWRqKxBP/w2CtjLfHhWIIfv3gIgMWVEZ7Z18Yzr7SN1bmjoZulVRFWzSnhhQPJTwz7mntf9++nk6L5Q4EuWeEv3rmcDYtm8p2n948t+3+7W8b6z2tmFlFdVkRj1yCHjiRDfsHMMNVlRTSMW7awIsz88jB1HQM0dQ8yEncsqggzf2YR7X1R2nqj3HjvS/z0pToA4glHPOG4f0v9WB/8osoIiysjHOoY4MGtDexu6uHPfrSZRMLxtr9/ivs2Jp/b3jfM5kOdvPhaB0/ta+PHLx7msV0tAHzu/h1sG9d98uMXjk6WvflwJ1vrurhwSbJr59c7m1kwM8y7V89mX0svH773Jb71m1d5aFsjOxu6WVtdygVLKth4sINYPMErLb1UFof4i9RtiUsKA2k/HpKZFOiSFQJ+H9+8fh1l4SCVxSEWVoQ5eKR/LGTnl4eZlwrv0S6MRZURasrDNHUPUju6rCJCTXkRPUOxsa6MxZUR5s9M9q3fv6WeZ8e1gBdWhFlSFWFbXRevtCRfY2llMYsqIgzHEjy5N9mXv7e593V96jddvAiAn286erLyxdeOjJ3YbO4Z4vFdzcwoDFBSGOD7vz8IwHVvmU9bb5SWnijvPnP22HPfvqKK8xaWjz3+xhOv8Jmfb6O5Z4hzF5azet4MhkYS1HcOsq+5l5VzivnU5SvY+7dXEAkp0PPFlALdzK4ws31mVmtmd0yy/lIz6zazramvL6S/VMl3c0uL+OHN5/Ot69axYGaylb23pZcFM8NEQgHmlRXRPTjCiwc6qCwOMaskRHV5ESNxx5N7Wijw+1g+u3gsvB/d0YQZnFldOrbske2NY+/30O0X88gn3sryWcXUtvaxo6GLBTPDlIaDLKpMbr+3uZfFlREAfr452TL/wc0b+B9XnUEo4OP+LclAv3BJBU/ubaW5Z4irz54HwFP72jh7fhnr5pfROxSjMOjjT94yf+z9L189m2vXJbf9q3ev4F1nzOZHHz1/bP3ovc7/8Ox5LK1K1rCvpZdXWvpYOXsGZqYTonnmhIFuZn7gLuBKYDVwvZmtnmTT55xz61JfX0lznSIArK0p5aJllSyYmewv39vUw6o5JQDMKysE4PHdzaytTgZaTVkRkOy2WDW3hFDAT015ctlju1pYVlVMcSjA/PJkQO9sSL7ey39zOWfVlFFSGGRpVTGHjgyw5VAXa2tKAcZCHOCGDQsAeHBrI4VBHxctrSDg97FsVjGQ7Pp524rKsQt9rh8X2hctrRx7rZsvXsya6lLWVpfy9fevpaY8zNf/6Cy2f+ndlIUL8PmMty6v5Cd/ej5/+a7lrJpTwnVvmU9ZuIAllcn3+uXLDQyOxHnLoqOteckfU/kstgGodc4dADCz+4BrgN3TWZjI8cyfGaZzYITOgRHee1ayFVudCu+EgzXVyeCtToV3NJZgbWrZaHgDnFVTBkBlccHYsrNryiiPHH28tKqYWMLR3DM09hqzSwq5ZEUV1WWFfOjChdzz3AHaeqOcOa+UgD/ZTppVEmIXcPtly5ib+mNTUhjg/CUV1JQXUd85yE0XL6K1J0q4IMAn37WcoN/Hw59469h7Fwb9b2hlX7S0Mnmi+B3LsdSy8kgBMyMF/Dp1Qvb8VP+75JepBHo1UDfucT1w/iTbXWhm24BG4NPOuV0TNzCzW4BbABYsWHDy1YqkLJh5NJTPGGuhF40tOzfV31w9ftmC5LKycHBs2SUrkuPXzYxLVlTxSnMvN1608HXvNdrSBnjLopkA+HzGD27eMLZ8RmGAtt4ot122dGzZF64+k2vP6eKaddVEY3FuvngxN791EX6f8eBtFxPw+SgM+llQEeaOK1ed9L+B32eve7yoIkxH/zCr5pQwc9wfJMkfUwl0m2TZxCsmtgALnXN9ZnYV8Etg+Rue5Nw9wD0A69ev1wwGcsrGB/qGxcmQnVtayGfes5J5ZYVcuqIKgEgowB1XriJc4Od951QDyfD+8cfOx4ALlx5tyY4P6PGWzy5mbmkhf3xezetOTI73zevOYX9bH+9YdfRE5uLUaBhIDh38wtVHeyorikOnsNfH9wdnzWMk7vj0e1am/bUlO9iJrmYzswuBLznn3pN6/DkA59zfHec5B4H1zrn2Y22zfv16t2nTplOpWYSReIKv/3ovH1w/n5WpFvp0cs5hNlnbRuT0MrPNzrn1k62bSgt9I7DczBYDDcB1wA0T3mAO0OKcc2a2geTJ1smvixZJg6Dfx9/8wWTn5qeHwlyywQkD3TkXM7PbgccAP3Cvc26Xmd2aWn838AHgz8wsBgwC1zlNCikiclqdsMtluqjLRUTk5B2vy0VXioqI5AgFuohIjlCgi4jkCAW6iEiOUKCLiOQIBbqISI7wbNiimbUBhyYsrgSOeXVpltG+ZJ5c2Q/QvmSq07EvC51zVZOt8CzQJ2Nmm441vjLbaF8yT67sB2hfMpXX+6IuFxGRHKFAFxHJEZkW6Pd4XUAaaV8yT67sB2hfMpWn+5JRfegiInLqMq2FLiIip0iBLiKSIzIm0M3sCjPbZ2a1ZnaH1/W8GWZ20Mx2mNlWM8uaewSb2b1m1mpmO8ctm2lmT5jZq6nvWTGd/DH25Utm1pA6LltT0yVmPDObb2ZPmdkeM9tlZp9MLc+qY3Oc/ci642JmhWb2kpltS+3Ll1PLPT0mGdGHbmZ+4BXgcpKTUG8ErnfO7fa0sFM0lSn4MpGZXQL0AT9wzq1JLft7oMM59/XUH9py59xfe1nnVBxjX74E9Dnn/sHL2k6Wmc0F5jrntphZCbAZuBb4CFl0bI6zHx8ky46LJaewiqTmUQ4CzwOfBN6Ph8ckU1roG4Ba59wB59wwcB9wjcc15R3n3LNAx4TF1wDfT/38fZK/gBnvGPuSlZxzTc65Lamfe4E9QDVZdmyOsx9ZxyX1pR4GU18Oj49JpgR6NVA37nE9WXqgUxzwuJltNrNbvC7mTZrtnGuC5C8kMMvjet6s281se6pLJqO7KCZjZouAc4AXyeJjM2E/IAuPi5n5zWwr0Ao84Zzz/JhkSqBPNgOv931Bp+5i59y5wJXAbamP/+K97wBLgXVAE/CPnlZzksysGLgf+EvnXI/X9ZyqSfYjK4+Lcy7unFsH1AAbzGyNxyVlTKDXA/PHPa4BGj2q5U1zzjWmvrcCD5DsUspWLam+z9E+0FaP6zllzrmW1C9hAvhXsui4pPpp7wd+7Jz7z9TirDs2k+1HNh8XAOdcF/A0cAUeH5NMCfSNwHIzW2xmBcB1wEMe13RKzCySOuGDmUWAdwM7j/+sjPYQ8OHUzx8GHvSwljdl9Bct5X1kyXFJnYD7LrDHOfeNcauy6tgcaz+y8biYWZWZlaV+LgLeBezF42OSEaNcAFJDlf4J8AP3Oue+6m1Fp8bMlpBslQMEgJ9ky76Y2U+BS0neArQF+CLwS+A/gAXAYeCPnXMZf7LxGPtyKcmP9Q44CHx8tL8zk5nZW4HngB1AIrX48yT7n7Pm2BxnP64ny46LmZ1F8qSnn2TD+D+cc18xswo8PCYZE+giIvLmZEqXi4iIvEkKdBGRHKFAFxHJEQp0EZEcoUAXEckRCnQRkRyhQBcRyRH/Hw77fsui1q2eAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "show_well_head(model_path=model_path, name=name, wel_coords=wel_coords)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fcc708c6-c749-4c9d-be80-3d9352aee812",
   "metadata": {},
   "source": [
    "There is some oscillation in the water level.\n",
    "These can be reduced by adding more time steps in the TDIS package."
   ]
  }
 ],
 "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.10.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}