examples/Notebooks/flopy3_mf6_B_complex-model.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# FloPy\n",
    "\n",
    "## Creating a Complex MODFLOW 6 Model with Flopy\n",
    "\n",
    "The purpose of this notebook is to demonstrate the Flopy capabilities for building a more complex MODFLOW 6 model from scratch.  This notebook will demonstrate the capabilities by replicating the advgw_tidal model that is distributed with MODFLOW 6."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Setup the Notebook Environment"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:40.085630Z",
     "iopub.status.busy": "2022-07-30T20:16:40.085282Z",
     "iopub.status.idle": "2022-07-30T20:16:40.949475Z",
     "shell.execute_reply": "2022-07-30T20:16:40.948676Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "3.10.5 (v3.10.5:f377153967, Jun  6 2022, 12:36:10) [Clang 13.0.0 (clang-1300.0.29.30)]\n",
      "numpy version: 1.23.1\n",
      "matplotlib version: 3.5.2\n",
      "flopy version: 3.3.6\n"
     ]
    }
   ],
   "source": [
    "import sys\n",
    "import os\n",
    "import platform\n",
    "from tempfile import TemporaryDirectory\n",
    "\n",
    "import numpy as np\n",
    "import matplotlib as mpl\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# run installed version of flopy or add local path\n",
    "try:\n",
    "    import flopy\n",
    "except:\n",
    "    fpth = os.path.abspath(os.path.join(\"..\", \"..\"))\n",
    "    sys.path.append(fpth)\n",
    "    import flopy\n",
    "\n",
    "print(sys.version)\n",
    "print(\"numpy version: {}\".format(np.__version__))\n",
    "print(\"matplotlib version: {}\".format(mpl.__version__))\n",
    "print(\"flopy version: {}\".format(flopy.__version__))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:40.987071Z",
     "iopub.status.busy": "2022-07-30T20:16:40.986715Z",
     "iopub.status.idle": "2022-07-30T20:16:40.991328Z",
     "shell.execute_reply": "2022-07-30T20:16:40.990444Z"
    }
   },
   "outputs": [],
   "source": [
    "# For this example, we will set up a temporary workspace.\n",
    "# Model input files and output files will reside here.\n",
    "temp_dir = TemporaryDirectory()\n",
    "model_name = \"advgw_tidal\"\n",
    "workspace = os.path.join(temp_dir.name, model_name)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:40.995165Z",
     "iopub.status.busy": "2022-07-30T20:16:40.994835Z",
     "iopub.status.idle": "2022-07-30T20:16:40.999076Z",
     "shell.execute_reply": "2022-07-30T20:16:40.998341Z"
    }
   },
   "outputs": [],
   "source": [
    "data_pth = os.path.join(\n",
    "    \"..\", \"data\", \"mf6\", \"create_tests\", \"test005_advgw_tidal\"\n",
    ")\n",
    "assert os.path.isdir(data_pth)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:41.002634Z",
     "iopub.status.busy": "2022-07-30T20:16:41.002350Z",
     "iopub.status.idle": "2022-07-30T20:16:41.056455Z",
     "shell.execute_reply": "2022-07-30T20:16:41.055706Z"
    }
   },
   "outputs": [],
   "source": [
    "# create simulation\n",
    "sim = flopy.mf6.MFSimulation(\n",
    "    sim_name=model_name, version=\"mf6\", exe_name=\"mf6\", sim_ws=workspace\n",
    ")\n",
    "\n",
    "# create tdis package\n",
    "tdis_rc = [(1.0, 1, 1.0), (10.0, 120, 1.0), (10.0, 120, 1.0), (10.0, 120, 1.0)]\n",
    "tdis = flopy.mf6.ModflowTdis(\n",
    "    sim, pname=\"tdis\", time_units=\"DAYS\", nper=4, perioddata=tdis_rc\n",
    ")\n",
    "\n",
    "# create gwf model\n",
    "gwf = flopy.mf6.ModflowGwf(\n",
    "    sim, modelname=model_name, model_nam_file=\"{}.nam\".format(model_name)\n",
    ")\n",
    "gwf.name_file.save_flows = True\n",
    "\n",
    "# create iterative model solution and register the gwf model with it\n",
    "ims = flopy.mf6.ModflowIms(\n",
    "    sim,\n",
    "    pname=\"ims\",\n",
    "    print_option=\"SUMMARY\",\n",
    "    complexity=\"SIMPLE\",\n",
    "    outer_dvclose=0.0001,\n",
    "    outer_maximum=500,\n",
    "    under_relaxation=\"NONE\",\n",
    "    inner_maximum=100,\n",
    "    inner_dvclose=0.0001,\n",
    "    rcloserecord=0.001,\n",
    "    linear_acceleration=\"CG\",\n",
    "    scaling_method=\"NONE\",\n",
    "    reordering_method=\"NONE\",\n",
    "    relaxation_factor=0.97,\n",
    ")\n",
    "sim.register_ims_package(ims, [gwf.name])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:41.060267Z",
     "iopub.status.busy": "2022-07-30T20:16:41.059967Z",
     "iopub.status.idle": "2022-07-30T20:16:41.084606Z",
     "shell.execute_reply": "2022-07-30T20:16:41.083792Z"
    }
   },
   "outputs": [],
   "source": [
    "# discretization package\n",
    "nlay = 3\n",
    "nrow = 15\n",
    "ncol = 10\n",
    "botlay2 = {\"factor\": 1.0, \"data\": [-100 for x in range(150)]}\n",
    "dis = flopy.mf6.ModflowGwfdis(\n",
    "    gwf,\n",
    "    pname=\"dis\",\n",
    "    nlay=nlay,\n",
    "    nrow=nrow,\n",
    "    ncol=ncol,\n",
    "    delr=500.0,\n",
    "    delc=500.0,\n",
    "    top=50.0,\n",
    "    botm=[5.0, -10.0, botlay2],\n",
    "    filename=\"{}.dis\".format(model_name),\n",
    ")\n",
    "\n",
    "# initial conditions\n",
    "ic = flopy.mf6.ModflowGwfic(\n",
    "    gwf, pname=\"ic\", strt=50.0, filename=\"{}.ic\".format(model_name)\n",
    ")\n",
    "\n",
    "# node property flow\n",
    "npf = flopy.mf6.ModflowGwfnpf(\n",
    "    gwf,\n",
    "    pname=\"npf\",\n",
    "    save_flows=True,\n",
    "    icelltype=[1, 0, 0],\n",
    "    k=[5.0, 0.1, 4.0],\n",
    "    k33=[0.5, 0.005, 0.1],\n",
    ")\n",
    "\n",
    "# output control\n",
    "oc = flopy.mf6.ModflowGwfoc(\n",
    "    gwf,\n",
    "    pname=\"oc\",\n",
    "    budget_filerecord=\"{}.cbb\".format(model_name),\n",
    "    head_filerecord=\"{}.hds\".format(model_name),\n",
    "    headprintrecord=[(\"COLUMNS\", 10, \"WIDTH\", 15, \"DIGITS\", 6, \"GENERAL\")],\n",
    "    saverecord=[(\"HEAD\", \"ALL\"), (\"BUDGET\", \"ALL\")],\n",
    "    printrecord=[(\"HEAD\", \"FIRST\"), (\"HEAD\", \"LAST\"), (\"BUDGET\", \"LAST\")],\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:41.087919Z",
     "iopub.status.busy": "2022-07-30T20:16:41.087642Z",
     "iopub.status.idle": "2022-07-30T20:16:41.094825Z",
     "shell.execute_reply": "2022-07-30T20:16:41.093963Z"
    }
   },
   "outputs": [],
   "source": [
    "# storage package\n",
    "sy = flopy.mf6.ModflowGwfsto.sy.empty(gwf, layered=True)\n",
    "for layer in range(0, 3):\n",
    "    sy[layer][\"data\"] = 0.2\n",
    "\n",
    "ss = flopy.mf6.ModflowGwfsto.ss.empty(\n",
    "    gwf, layered=True, default_value=0.000001\n",
    ")\n",
    "\n",
    "sto = flopy.mf6.ModflowGwfsto(\n",
    "    gwf,\n",
    "    pname=\"sto\",\n",
    "    save_flows=True,\n",
    "    iconvert=1,\n",
    "    ss=ss,\n",
    "    sy=sy,\n",
    "    steady_state={0: True},\n",
    "    transient={1: True},\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:41.098706Z",
     "iopub.status.busy": "2022-07-30T20:16:41.098390Z",
     "iopub.status.idle": "2022-07-30T20:16:41.175458Z",
     "shell.execute_reply": "2022-07-30T20:16:41.174721Z"
    }
   },
   "outputs": [],
   "source": [
    "# well package\n",
    "# test empty with aux vars, bound names, and time series\n",
    "period_two = flopy.mf6.ModflowGwfwel.stress_period_data.empty(\n",
    "    gwf,\n",
    "    maxbound=3,\n",
    "    aux_vars=[\"var1\", \"var2\", \"var3\"],\n",
    "    boundnames=True,\n",
    "    timeseries=True,\n",
    ")\n",
    "period_two[0][0] = ((0, 11, 2), -50.0, -1, -2, -3, None)\n",
    "period_two[0][1] = ((2, 4, 7), \"well_1_rate\", 1, 2, 3, \"well_1\")\n",
    "period_two[0][2] = ((2, 3, 2), \"well_2_rate\", 4, 5, 6, \"well_2\")\n",
    "period_three = flopy.mf6.ModflowGwfwel.stress_period_data.empty(\n",
    "    gwf,\n",
    "    maxbound=2,\n",
    "    aux_vars=[\"var1\", \"var2\", \"var3\"],\n",
    "    boundnames=True,\n",
    "    timeseries=True,\n",
    ")\n",
    "period_three[0][0] = ((2, 3, 2), \"well_2_rate\", 1, 2, 3, \"well_2\")\n",
    "period_three[0][1] = ((2, 4, 7), \"well_1_rate\", 4, 5, 6, \"well_1\")\n",
    "period_four = flopy.mf6.ModflowGwfwel.stress_period_data.empty(\n",
    "    gwf,\n",
    "    maxbound=5,\n",
    "    aux_vars=[\"var1\", \"var2\", \"var3\"],\n",
    "    boundnames=True,\n",
    "    timeseries=True,\n",
    ")\n",
    "period_four[0][0] = ((2, 4, 7), \"well_1_rate\", 1, 2, 3, \"well_1\")\n",
    "period_four[0][1] = ((2, 3, 2), \"well_2_rate\", 4, 5, 6, \"well_2\")\n",
    "period_four[0][2] = ((0, 11, 2), -10.0, 7, 8, 9, None)\n",
    "period_four[0][3] = ((0, 2, 4), -20.0, 17, 18, 19, None)\n",
    "period_four[0][4] = ((0, 13, 5), -40.0, 27, 28, 29, None)\n",
    "stress_period_data = {}\n",
    "stress_period_data[1] = period_two[0]\n",
    "stress_period_data[2] = period_three[0]\n",
    "stress_period_data[3] = period_four[0]\n",
    "wel = flopy.mf6.ModflowGwfwel(\n",
    "    gwf,\n",
    "    pname=\"wel\",\n",
    "    print_input=True,\n",
    "    print_flows=True,\n",
    "    auxiliary=[(\"var1\", \"var2\", \"var3\")],\n",
    "    maxbound=5,\n",
    "    stress_period_data=stress_period_data,\n",
    "    boundnames=True,\n",
    "    save_flows=True,\n",
    ")\n",
    "\n",
    "# well ts package\n",
    "ts_data = [\n",
    "    (0.0, 0.0, 0.0, 0.0),\n",
    "    (1.0, -200.0, 0.0, -100.0),\n",
    "    (11.0, -1800.0, -500.0, -200.0),\n",
    "    (21.0, -200.0, -400.0, -300.0),\n",
    "    (31.0, 0.0, -600.0, -400.0),\n",
    "]\n",
    "wel.ts.initialize(\n",
    "    filename=\"well-rates.ts\",\n",
    "    timeseries=ts_data,\n",
    "    time_series_namerecord=[(\"well_1_rate\", \"well_2_rate\", \"well_3_rate\")],\n",
    "    interpolation_methodrecord=[(\"stepwise\", \"stepwise\", \"stepwise\")],\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:41.178973Z",
     "iopub.status.busy": "2022-07-30T20:16:41.178669Z",
     "iopub.status.idle": "2022-07-30T20:16:41.187897Z",
     "shell.execute_reply": "2022-07-30T20:16:41.186909Z"
    }
   },
   "outputs": [],
   "source": [
    "# Evapotranspiration\n",
    "evt_period = flopy.mf6.ModflowGwfevt.stress_period_data.empty(gwf, 150, nseg=3)\n",
    "for col in range(0, 10):\n",
    "    for row in range(0, 15):\n",
    "        evt_period[0][col * 15 + row] = (\n",
    "            (0, row, col),\n",
    "            50.0,\n",
    "            0.0004,\n",
    "            10.0,\n",
    "            0.2,\n",
    "            0.5,\n",
    "            0.3,\n",
    "            0.1,\n",
    "            None,\n",
    "        )\n",
    "evt = flopy.mf6.ModflowGwfevt(\n",
    "    gwf,\n",
    "    pname=\"evt\",\n",
    "    print_input=True,\n",
    "    print_flows=True,\n",
    "    save_flows=True,\n",
    "    maxbound=150,\n",
    "    nseg=3,\n",
    "    stress_period_data=evt_period,\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:41.191273Z",
     "iopub.status.busy": "2022-07-30T20:16:41.191003Z",
     "iopub.status.idle": "2022-07-30T20:16:41.217218Z",
     "shell.execute_reply": "2022-07-30T20:16:41.216541Z"
    }
   },
   "outputs": [],
   "source": [
    "# General-Head Boundaries\n",
    "ghb_period = {}\n",
    "ghb_period_array = []\n",
    "for layer, cond in zip(range(1, 3), [15.0, 1500.0]):\n",
    "    for row in range(0, 15):\n",
    "        ghb_period_array.append(((layer, row, 9), \"tides\", cond, \"Estuary-L2\"))\n",
    "ghb_period[0] = ghb_period_array\n",
    "ghb = flopy.mf6.ModflowGwfghb(\n",
    "    gwf,\n",
    "    pname=\"ghb\",\n",
    "    print_input=True,\n",
    "    print_flows=True,\n",
    "    save_flows=True,\n",
    "    boundnames=True,\n",
    "    maxbound=30,\n",
    "    stress_period_data=ghb_period,\n",
    ")\n",
    "ts_recarray = []\n",
    "fd = open(os.path.join(data_pth, \"tides.txt\"), \"r\")\n",
    "for line in fd:\n",
    "    line_list = line.strip().split(\",\")\n",
    "    ts_recarray.append((float(line_list[0]), float(line_list[1])))\n",
    "ghb.ts.initialize(\n",
    "    filename=\"tides.ts\",\n",
    "    timeseries=ts_recarray,\n",
    "    time_series_namerecord=\"tides\",\n",
    "    interpolation_methodrecord=\"linear\",\n",
    ")\n",
    "obs_recarray = {\n",
    "    \"ghb_obs.csv\": [\n",
    "        (\"ghb-2-6-10\", \"GHB\", (1, 5, 9)),\n",
    "        (\"ghb-3-6-10\", \"GHB\", (2, 5, 9)),\n",
    "    ],\n",
    "    \"ghb_flows.csv\": [\n",
    "        (\"Estuary2\", \"GHB\", \"Estuary-L2\"),\n",
    "        (\"Estuary3\", \"GHB\", \"Estuary-L3\"),\n",
    "    ],\n",
    "}\n",
    "ghb.obs.initialize(\n",
    "    filename=\"{}.ghb.obs\".format(model_name),\n",
    "    print_input=True,\n",
    "    continuous=obs_recarray,\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:41.221053Z",
     "iopub.status.busy": "2022-07-30T20:16:41.220730Z",
     "iopub.status.idle": "2022-07-30T20:16:41.230965Z",
     "shell.execute_reply": "2022-07-30T20:16:41.230178Z"
    }
   },
   "outputs": [],
   "source": [
    "obs_recarray = {\n",
    "    \"head_obs.csv\": [(\"h1_13_8\", \"HEAD\", (2, 12, 7))],\n",
    "    \"intercell_flow_obs1.csv\": [\n",
    "        (\"ICF1_1.0\", \"FLOW-JA-FACE\", (0, 4, 5), (0, 5, 5))\n",
    "    ],\n",
    "    \"head-hydrographs.csv\": [\n",
    "        (\"h3-13-9\", \"HEAD\", (2, 12, 8)),\n",
    "        (\"h3-12-8\", \"HEAD\", (2, 11, 7)),\n",
    "        (\"h1-4-3\", \"HEAD\", (0, 3, 2)),\n",
    "        (\"h1-12-3\", \"HEAD\", (0, 11, 2)),\n",
    "        (\"h1-13-9\", \"HEAD\", (0, 12, 8)),\n",
    "    ],\n",
    "}\n",
    "obs_package = flopy.mf6.ModflowUtlobs(\n",
    "    gwf,\n",
    "    pname=\"head_obs\",\n",
    "    filename=\"{}.obs\".format(model_name),\n",
    "    print_input=True,\n",
    "    continuous=obs_recarray,\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:41.234638Z",
     "iopub.status.busy": "2022-07-30T20:16:41.234303Z",
     "iopub.status.idle": "2022-07-30T20:16:41.264075Z",
     "shell.execute_reply": "2022-07-30T20:16:41.263034Z"
    }
   },
   "outputs": [],
   "source": [
    "# River\n",
    "riv_period = {}\n",
    "riv_period_array = [\n",
    "    ((0, 2, 0), \"river_stage_1\", 1001.0, 35.9, None),\n",
    "    ((0, 3, 1), \"river_stage_1\", 1002.0, 35.8, None),\n",
    "    ((0, 4, 2), \"river_stage_1\", 1003.0, 35.7, None),\n",
    "    ((0, 4, 3), \"river_stage_1\", 1004.0, 35.6, None),\n",
    "    ((0, 5, 4), \"river_stage_1\", 1005.0, 35.5, None),\n",
    "    ((0, 5, 5), \"river_stage_1\", 1006.0, 35.4, \"riv1_c6\"),\n",
    "    ((0, 5, 6), \"river_stage_1\", 1007.0, 35.3, \"riv1_c7\"),\n",
    "    ((0, 4, 7), \"river_stage_1\", 1008.0, 35.2, None),\n",
    "    ((0, 4, 8), \"river_stage_1\", 1009.0, 35.1, None),\n",
    "    ((0, 4, 9), \"river_stage_1\", 1010.0, 35.0, None),\n",
    "    ((0, 9, 0), \"river_stage_2\", 1001.0, 36.9, \"riv2_upper\"),\n",
    "    ((0, 8, 1), \"river_stage_2\", 1002.0, 36.8, \"riv2_upper\"),\n",
    "    ((0, 7, 2), \"river_stage_2\", 1003.0, 36.7, \"riv2_upper\"),\n",
    "    ((0, 6, 3), \"river_stage_2\", 1004.0, 36.6, None),\n",
    "    ((0, 6, 4), \"river_stage_2\", 1005.0, 36.5, None),\n",
    "    ((0, 5, 5), \"river_stage_2\", 1006.0, 36.4, \"riv2_c6\"),\n",
    "    ((0, 5, 6), \"river_stage_2\", 1007.0, 36.3, \"riv2_c7\"),\n",
    "    ((0, 6, 7), \"river_stage_2\", 1008.0, 36.2, None),\n",
    "    ((0, 6, 8), \"river_stage_2\", 1009.0, 36.1),\n",
    "    ((0, 6, 9), \"river_stage_2\", 1010.0, 36.0),\n",
    "]\n",
    "riv_period[0] = riv_period_array\n",
    "riv = flopy.mf6.ModflowGwfriv(\n",
    "    gwf,\n",
    "    pname=\"riv\",\n",
    "    print_input=True,\n",
    "    print_flows=True,\n",
    "    save_flows=\"{}.cbc\".format(model_name),\n",
    "    boundnames=True,\n",
    "    maxbound=20,\n",
    "    stress_period_data=riv_period,\n",
    ")\n",
    "ts_recarray = [\n",
    "    (0.0, 40.0, 41.0),\n",
    "    (1.0, 41.0, 41.5),\n",
    "    (2.0, 43.0, 42.0),\n",
    "    (3.0, 45.0, 42.8),\n",
    "    (4.0, 44.0, 43.0),\n",
    "    (6.0, 43.0, 43.1),\n",
    "    (9.0, 42.0, 42.4),\n",
    "    (11.0, 41.0, 41.5),\n",
    "    (31.0, 40.0, 41.0),\n",
    "]\n",
    "riv.ts.initialize(\n",
    "    filename=\"river_stages.ts\",\n",
    "    timeseries=ts_recarray,\n",
    "    time_series_namerecord=[(\"river_stage_1\", \"river_stage_2\")],\n",
    "    interpolation_methodrecord=[(\"linear\", \"stepwise\")],\n",
    ")\n",
    "obs_recarray = {\n",
    "    \"riv_obs.csv\": [\n",
    "        (\"rv1-3-1\", \"RIV\", (0, 2, 0)),\n",
    "        (\"rv1-4-2\", \"RIV\", (0, 3, 1)),\n",
    "        (\"rv1-5-3\", \"RIV\", (0, 4, 2)),\n",
    "        (\"rv1-5-4\", \"RIV\", (0, 4, 3)),\n",
    "        (\"rv1-6-5\", \"RIV\", (0, 5, 4)),\n",
    "        (\"rv1-c6\", \"RIV\", \"riv1_c6\"),\n",
    "        (\"rv1-c7\", \"RIV\", \"riv1_c7\"),\n",
    "        (\"rv2-upper\", \"RIV\", \"riv2_upper\"),\n",
    "        (\"rv-2-7-4\", \"RIV\", (0, 6, 3)),\n",
    "        (\"rv2-8-5\", \"RIV\", (0, 6, 4)),\n",
    "        (\n",
    "            \"rv-2-9-6\",\n",
    "            \"RIV\",\n",
    "            (\n",
    "                0,\n",
    "                5,\n",
    "                5,\n",
    "            ),\n",
    "        ),\n",
    "    ],\n",
    "    \"riv_flowsA.csv\": [\n",
    "        (\"riv1-3-1\", \"RIV\", (0, 2, 0)),\n",
    "        (\"riv1-4-2\", \"RIV\", (0, 3, 1)),\n",
    "        (\"riv1-5-3\", \"RIV\", (0, 4, 2)),\n",
    "    ],\n",
    "    \"riv_flowsB.csv\": [\n",
    "        (\"riv2-10-1\", \"RIV\", (0, 9, 0)),\n",
    "        (\"riv-2-9-2\", \"RIV\", (0, 8, 1)),\n",
    "        (\"riv2-8-3\", \"RIV\", (0, 7, 2)),\n",
    "    ],\n",
    "}\n",
    "riv.obs.initialize(\n",
    "    filename=\"{}.riv.obs\".format(model_name),\n",
    "    print_input=True,\n",
    "    continuous=obs_recarray,\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:41.267410Z",
     "iopub.status.busy": "2022-07-30T20:16:41.267142Z",
     "iopub.status.idle": "2022-07-30T20:16:41.290530Z",
     "shell.execute_reply": "2022-07-30T20:16:41.289900Z"
    }
   },
   "outputs": [],
   "source": [
    "# First recharge package\n",
    "rch1_period = {}\n",
    "rch1_period_array = []\n",
    "col_range = {0: 3, 1: 4, 2: 5}\n",
    "for row in range(0, 15):\n",
    "    if row in col_range:\n",
    "        col_max = col_range[row]\n",
    "    else:\n",
    "        col_max = 6\n",
    "    for col in range(0, col_max):\n",
    "        if (\n",
    "            (row == 3 and col == 5)\n",
    "            or (row == 2 and col == 4)\n",
    "            or (row == 1 and col == 3)\n",
    "            or (row == 0 and col == 2)\n",
    "        ):\n",
    "            mult = 0.5\n",
    "        else:\n",
    "            mult = 1.0\n",
    "        if row == 0 and col == 0:\n",
    "            bnd = \"rch-1-1\"\n",
    "        elif row == 0 and col == 1:\n",
    "            bnd = \"rch-1-2\"\n",
    "        elif row == 1 and col == 2:\n",
    "            bnd = \"rch-2-3\"\n",
    "        else:\n",
    "            bnd = None\n",
    "        rch1_period_array.append(((0, row, col), \"rch_1\", mult, bnd))\n",
    "rch1_period[0] = rch1_period_array\n",
    "rch1 = flopy.mf6.ModflowGwfrch(\n",
    "    gwf,\n",
    "    filename=\"{}_1.rch\".format(model_name),\n",
    "    pname=\"rch_1\",\n",
    "    fixed_cell=True,\n",
    "    auxiliary=\"MULTIPLIER\",\n",
    "    auxmultname=\"MULTIPLIER\",\n",
    "    print_input=True,\n",
    "    print_flows=True,\n",
    "    save_flows=True,\n",
    "    boundnames=True,\n",
    "    maxbound=84,\n",
    "    stress_period_data=rch1_period,\n",
    ")\n",
    "ts_data = [\n",
    "    (0.0, 0.0015),\n",
    "    (1.0, 0.0010),\n",
    "    (11.0, 0.0015),\n",
    "    (21.0, 0.0025),\n",
    "    (31.0, 0.0015),\n",
    "]\n",
    "rch1.ts.initialize(\n",
    "    filename=\"recharge_rates_1.ts\",\n",
    "    timeseries=ts_data,\n",
    "    time_series_namerecord=\"rch_1\",\n",
    "    interpolation_methodrecord=\"stepwise\",\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:41.294282Z",
     "iopub.status.busy": "2022-07-30T20:16:41.293970Z",
     "iopub.status.idle": "2022-07-30T20:16:41.315595Z",
     "shell.execute_reply": "2022-07-30T20:16:41.314923Z"
    }
   },
   "outputs": [],
   "source": [
    "# Second recharge package\n",
    "rch2_period = {}\n",
    "rch2_period_array = [\n",
    "    ((0, 0, 2), \"rch_2\", 0.5),\n",
    "    ((0, 0, 3), \"rch_2\", 1.0),\n",
    "    ((0, 0, 4), \"rch_2\", 1.0),\n",
    "    ((0, 0, 5), \"rch_2\", 1.0),\n",
    "    ((0, 0, 6), \"rch_2\", 1.0),\n",
    "    ((0, 0, 7), \"rch_2\", 1.0),\n",
    "    ((0, 0, 8), \"rch_2\", 1.0),\n",
    "    ((0, 0, 9), \"rch_2\", 0.5),\n",
    "    ((0, 1, 3), \"rch_2\", 0.5),\n",
    "    ((0, 1, 4), \"rch_2\", 1.0),\n",
    "    ((0, 1, 5), \"rch_2\", 1.0),\n",
    "    ((0, 1, 6), \"rch_2\", 1.0),\n",
    "    ((0, 1, 7), \"rch_2\", 1.0),\n",
    "    ((0, 1, 8), \"rch_2\", 0.5),\n",
    "    ((0, 2, 4), \"rch_2\", 0.5),\n",
    "    ((0, 2, 5), \"rch_2\", 1.0),\n",
    "    ((0, 2, 6), \"rch_2\", 1.0),\n",
    "    ((0, 2, 7), \"rch_2\", 0.5),\n",
    "    ((0, 3, 5), \"rch_2\", 0.5),\n",
    "    ((0, 3, 6), \"rch_2\", 0.5),\n",
    "]\n",
    "rch2_period[0] = rch2_period_array\n",
    "rch2 = flopy.mf6.ModflowGwfrch(\n",
    "    gwf,\n",
    "    filename=\"{}_2.rch\".format(model_name),\n",
    "    pname=\"rch_2\",\n",
    "    fixed_cell=True,\n",
    "    auxiliary=\"MULTIPLIER\",\n",
    "    auxmultname=\"MULTIPLIER\",\n",
    "    print_input=True,\n",
    "    print_flows=True,\n",
    "    save_flows=True,\n",
    "    maxbound=20,\n",
    "    stress_period_data=rch2_period,\n",
    ")\n",
    "ts_data = [\n",
    "    (0.0, 0.0016),\n",
    "    (1.0, 0.0018),\n",
    "    (11.0, 0.0019),\n",
    "    (21.0, 0.0016),\n",
    "    (31.0, 0.0018),\n",
    "]\n",
    "rch2.ts.initialize(\n",
    "    filename=\"recharge_rates_2.ts\",\n",
    "    timeseries=ts_data,\n",
    "    time_series_namerecord=\"rch_2\",\n",
    "    interpolation_methodrecord=\"linear\",\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:41.319132Z",
     "iopub.status.busy": "2022-07-30T20:16:41.318863Z",
     "iopub.status.idle": "2022-07-30T20:16:41.339998Z",
     "shell.execute_reply": "2022-07-30T20:16:41.339277Z"
    }
   },
   "outputs": [],
   "source": [
    "# Third recharge package\n",
    "rch3_period = {}\n",
    "rch3_period_array = []\n",
    "col_range = {0: 9, 1: 8, 2: 7}\n",
    "for row in range(0, 15):\n",
    "    if row in col_range:\n",
    "        col_min = col_range[row]\n",
    "    else:\n",
    "        col_min = 6\n",
    "    for col in range(col_min, 10):\n",
    "        if (\n",
    "            (row == 0 and col == 9)\n",
    "            or (row == 1 and col == 8)\n",
    "            or (row == 2 and col == 7)\n",
    "            or (row == 3 and col == 6)\n",
    "        ):\n",
    "            mult = 0.5\n",
    "        else:\n",
    "            mult = 1.0\n",
    "        rch3_period_array.append(((0, row, col), \"rch_3\", mult))\n",
    "rch3_period[0] = rch3_period_array\n",
    "rch3 = flopy.mf6.ModflowGwfrch(\n",
    "    gwf,\n",
    "    filename=\"{}_3.rch\".format(model_name),\n",
    "    pname=\"rch_3\",\n",
    "    fixed_cell=True,\n",
    "    auxiliary=\"MULTIPLIER\",\n",
    "    auxmultname=\"MULTIPLIER\",\n",
    "    print_input=True,\n",
    "    print_flows=True,\n",
    "    save_flows=True,\n",
    "    maxbound=54,\n",
    "    stress_period_data=rch3_period,\n",
    ")\n",
    "ts_data = [\n",
    "    (0.0, 0.0017),\n",
    "    (1.0, 0.0020),\n",
    "    (11.0, 0.0017),\n",
    "    (21.0, 0.0018),\n",
    "    (31.0, 0.0020),\n",
    "]\n",
    "rch3.ts.initialize(\n",
    "    filename=\"recharge_rates_3.ts\",\n",
    "    timeseries=ts_data,\n",
    "    time_series_namerecord=\"rch_3\",\n",
    "    interpolation_methodrecord=\"linear\",\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Create the MODFLOW 6 Input Files and Run the Model\n",
    "\n",
    "Once all the flopy objects are created, it is very easy to create all of the input files and run the model."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:41.344176Z",
     "iopub.status.busy": "2022-07-30T20:16:41.343882Z",
     "iopub.status.idle": "2022-07-30T20:16:41.407652Z",
     "shell.execute_reply": "2022-07-30T20:16:41.406718Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "writing simulation...\n",
      "  writing simulation name file...\n",
      "  writing simulation tdis package...\n",
      "  writing ims package ims...\n",
      "  writing model advgw_tidal...\n",
      "    writing model name file...\n",
      "    writing package dis...\n",
      "    writing package ic...\n",
      "    writing package npf...\n",
      "    writing package oc...\n",
      "    writing package sto...\n",
      "    writing package wel...\n",
      "    writing package ts_0...\n",
      "    writing package evt...\n",
      "    writing package ghb...\n",
      "    writing package ts_1...\n",
      "    writing package obs_0...\n",
      "    writing package head_obs...\n",
      "    writing package riv...\n",
      "    writing package ts_2...\n",
      "    writing package obs_1...\n",
      "    writing package rch_1...\n",
      "    writing package ts_3...\n",
      "    writing package rch_2...\n",
      "    writing package ts_4...\n",
      "    writing package rch_3...\n",
      "    writing package ts_5...\n"
     ]
    }
   ],
   "source": [
    "# write simulation to new location\n",
    "sim.write_simulation()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:41.411323Z",
     "iopub.status.busy": "2022-07-30T20:16:41.411028Z",
     "iopub.status.idle": "2022-07-30T20:16:41.415333Z",
     "shell.execute_reply": "2022-07-30T20:16:41.414398Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "['advgw_tidal.nam', 'advgw_tidal.ic', 'recharge_rates_1.ts', 'advgw_tidal.sto', 'tides.ts', 'advgw_tidal.ims', 'well-rates.ts', 'advgw_tidal.ghb', 'advgw_tidal.obs', 'advgw_tidal.riv.obs', 'advgw_tidal.dis', 'advgw_tidal_1.rch', 'advgw_tidal_3.rch', 'advgw_tidal_2.rch', 'advgw_tidal.oc', 'river_stages.ts', 'advgw_tidal.wel', 'advgw_tidal.ghb.obs', 'advgw_tidal.npf', 'advgw_tidal.tdis', 'recharge_rates_2.ts', 'mfsim.nam', 'advgw_tidal.riv', 'recharge_rates_3.ts', 'advgw_tidal.evt']\n"
     ]
    }
   ],
   "source": [
    "# Print a list of the files that were created\n",
    "# in workspace\n",
    "print(os.listdir(workspace))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Run the Simulation\n",
    "\n",
    "We can also run the simulation from the notebook, but only if the MODFLOW 6 executable is available.  The executable can be made available by putting the executable in a folder that is listed in the system path variable.  Another option is to just put a copy of the executable in the simulation folder, though this should generally be avoided.  A final option is to provide a full path to the executable when the simulation is constructed.  This would be done by specifying exe_name with the full path."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:41.419162Z",
     "iopub.status.busy": "2022-07-30T20:16:41.418848Z",
     "iopub.status.idle": "2022-07-30T20:16:42.921601Z",
     "shell.execute_reply": "2022-07-30T20:16:42.920396Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "FloPy is using the following executable to run the model: /Users/wes/dev/usgs-executables/mf6\n",
      "                                   MODFLOW 6\n",
      "                U.S. GEOLOGICAL SURVEY MODULAR HYDROLOGIC MODEL\n",
      "                            VERSION 6.3.0 03/04/2022\n",
      "\n",
      "   MODFLOW 6 compiled Mar 07 2022 13:50:09 with Intel(R) Fortran Intel(R) 64\n",
      "   Compiler Classic for applications running on Intel(R) 64, Version 2021.5.0\n",
      "                             Build 20211109_000000\n",
      "\n",
      "This software has been approved for release by the U.S. Geological \n",
      "Survey (USGS). Although the software has been subjected to rigorous \n",
      "review, the USGS reserves the right to update the software as needed \n",
      "pursuant to further analysis and review. No warranty, expressed or \n",
      "implied, is made by the USGS or the U.S. Government as to the \n",
      "functionality of the software and related material nor shall the \n",
      "fact of release constitute any such warranty. Furthermore, the \n",
      "software is released on condition that neither the USGS nor the U.S. \n",
      "Government shall be held liable for any damages resulting from its \n",
      "authorized or unauthorized use. Also refer to the USGS Water \n",
      "Resources Software User Rights Notice for complete use, copyright, \n",
      "and distribution information.\n",
      "\n",
      " \n",
      " Run start date and time (yyyy/mm/dd hh:mm:ss): 2022/07/30 16:16:41\n",
      " \n",
      " Writing simulation list file: mfsim.lst\n",
      " Using Simulation name file: mfsim.nam\n",
      " \n",
      "    Solving:  Stress period:     1    Time step:     1\n",
      "    Solving:  Stress period:     2    Time step:     1\n",
      "    Solving:  Stress period:     2    Time step:     2\n",
      "    Solving:  Stress period:     2    Time step:     3\n",
      "    Solving:  Stress period:     2    Time step:     4\n",
      "    Solving:  Stress period:     2    Time step:     5\n",
      "    Solving:  Stress period:     2    Time step:     6\n",
      "    Solving:  Stress period:     2    Time step:     7\n",
      "    Solving:  Stress period:     2    Time step:     8\n",
      "    Solving:  Stress period:     2    Time step:     9\n",
      "    Solving:  Stress period:     2    Time step:    10\n",
      "    Solving:  Stress period:     2    Time step:    11\n",
      "    Solving:  Stress period:     2    Time step:    12\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     2    Time step:    13\n",
      "    Solving:  Stress period:     2    Time step:    14\n",
      "    Solving:  Stress period:     2    Time step:    15\n",
      "    Solving:  Stress period:     2    Time step:    16\n",
      "    Solving:  Stress period:     2    Time step:    17\n",
      "    Solving:  Stress period:     2    Time step:    18\n",
      "    Solving:  Stress period:     2    Time step:    19\n",
      "    Solving:  Stress period:     2    Time step:    20\n",
      "    Solving:  Stress period:     2    Time step:    21\n",
      "    Solving:  Stress period:     2    Time step:    22\n",
      "    Solving:  Stress period:     2    Time step:    23\n",
      "    Solving:  Stress period:     2    Time step:    24\n",
      "    Solving:  Stress period:     2    Time step:    25\n",
      "    Solving:  Stress period:     2    Time step:    26\n",
      "    Solving:  Stress period:     2    Time step:    27\n",
      "    Solving:  Stress period:     2    Time step:    28\n",
      "    Solving:  Stress period:     2    Time step:    29\n",
      "    Solving:  Stress period:     2    Time step:    30\n",
      "    Solving:  Stress period:     2    Time step:    31\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     2    Time step:    32\n",
      "    Solving:  Stress period:     2    Time step:    33\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     2    Time step:    34\n",
      "    Solving:  Stress period:     2    Time step:    35\n",
      "    Solving:  Stress period:     2    Time step:    36\n",
      "    Solving:  Stress period:     2    Time step:    37\n",
      "    Solving:  Stress period:     2    Time step:    38\n",
      "    Solving:  Stress period:     2    Time step:    39\n",
      "    Solving:  Stress period:     2    Time step:    40\n",
      "    Solving:  Stress period:     2    Time step:    41\n",
      "    Solving:  Stress period:     2    Time step:    42\n",
      "    Solving:  Stress period:     2    Time step:    43\n",
      "    Solving:  Stress period:     2    Time step:    44\n",
      "    Solving:  Stress period:     2    Time step:    45\n",
      "    Solving:  Stress period:     2    Time step:    46\n",
      "    Solving:  Stress period:     2    Time step:    47\n",
      "    Solving:  Stress period:     2    Time step:    48\n",
      "    Solving:  Stress period:     2    Time step:    49\n",
      "    Solving:  Stress period:     2    Time step:    50\n",
      "    Solving:  Stress period:     2    Time step:    51\n",
      "    Solving:  Stress period:     2    Time step:    52\n",
      "    Solving:  Stress period:     2    Time step:    53\n",
      "    Solving:  Stress period:     2    Time step:    54\n",
      "    Solving:  Stress period:     2    Time step:    55\n",
      "    Solving:  Stress period:     2    Time step:    56\n",
      "    Solving:  Stress period:     2    Time step:    57\n",
      "    Solving:  Stress period:     2    Time step:    58\n",
      "    Solving:  Stress period:     2    Time step:    59\n",
      "    Solving:  Stress period:     2    Time step:    60\n",
      "    Solving:  Stress period:     2    Time step:    61\n",
      "    Solving:  Stress period:     2    Time step:    62\n",
      "    Solving:  Stress period:     2    Time step:    63\n",
      "    Solving:  Stress period:     2    Time step:    64\n",
      "    Solving:  Stress period:     2    Time step:    65\n",
      "    Solving:  Stress period:     2    Time step:    66\n",
      "    Solving:  Stress period:     2    Time step:    67\n",
      "    Solving:  Stress period:     2    Time step:    68\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     2    Time step:    69\n",
      "    Solving:  Stress period:     2    Time step:    70\n",
      "    Solving:  Stress period:     2    Time step:    71\n",
      "    Solving:  Stress period:     2    Time step:    72\n",
      "    Solving:  Stress period:     2    Time step:    73\n",
      "    Solving:  Stress period:     2    Time step:    74\n",
      "    Solving:  Stress period:     2    Time step:    75\n",
      "    Solving:  Stress period:     2    Time step:    76\n",
      "    Solving:  Stress period:     2    Time step:    77\n",
      "    Solving:  Stress period:     2    Time step:    78\n",
      "    Solving:  Stress period:     2    Time step:    79\n",
      "    Solving:  Stress period:     2    Time step:    80\n",
      "    Solving:  Stress period:     2    Time step:    81\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     2    Time step:    82\n",
      "    Solving:  Stress period:     2    Time step:    83\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     2    Time step:    84\n",
      "    Solving:  Stress period:     2    Time step:    85\n",
      "    Solving:  Stress period:     2    Time step:    86\n",
      "    Solving:  Stress period:     2    Time step:    87\n",
      "    Solving:  Stress period:     2    Time step:    88\n",
      "    Solving:  Stress period:     2    Time step:    89\n",
      "    Solving:  Stress period:     2    Time step:    90\n",
      "    Solving:  Stress period:     2    Time step:    91\n",
      "    Solving:  Stress period:     2    Time step:    92\n",
      "    Solving:  Stress period:     2    Time step:    93\n",
      "    Solving:  Stress period:     2    Time step:    94\n",
      "    Solving:  Stress period:     2    Time step:    95\n",
      "    Solving:  Stress period:     2    Time step:    96\n",
      "    Solving:  Stress period:     2    Time step:    97\n",
      "    Solving:  Stress period:     2    Time step:    98\n",
      "    Solving:  Stress period:     2    Time step:    99\n",
      "    Solving:  Stress period:     2    Time step:   100\n",
      "    Solving:  Stress period:     2    Time step:   101\n",
      "    Solving:  Stress period:     2    Time step:   102\n",
      "    Solving:  Stress period:     2    Time step:   103\n",
      "    Solving:  Stress period:     2    Time step:   104\n",
      "    Solving:  Stress period:     2    Time step:   105\n",
      "    Solving:  Stress period:     2    Time step:   106\n",
      "    Solving:  Stress period:     2    Time step:   107\n",
      "    Solving:  Stress period:     2    Time step:   108\n",
      "    Solving:  Stress period:     2    Time step:   109\n",
      "    Solving:  Stress period:     2    Time step:   110\n",
      "    Solving:  Stress period:     2    Time step:   111\n",
      "    Solving:  Stress period:     2    Time step:   112\n",
      "    Solving:  Stress period:     2    Time step:   113\n",
      "    Solving:  Stress period:     2    Time step:   114\n",
      "    Solving:  Stress period:     2    Time step:   115\n",
      "    Solving:  Stress period:     2    Time step:   116\n",
      "    Solving:  Stress period:     2    Time step:   117\n",
      "    Solving:  Stress period:     2    Time step:   118\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     2    Time step:   119\n",
      "    Solving:  Stress period:     2    Time step:   120\n",
      "    Solving:  Stress period:     3    Time step:     1\n",
      "    Solving:  Stress period:     3    Time step:     2\n",
      "    Solving:  Stress period:     3    Time step:     3\n",
      "    Solving:  Stress period:     3    Time step:     4\n",
      "    Solving:  Stress period:     3    Time step:     5\n",
      "    Solving:  Stress period:     3    Time step:     6\n",
      "    Solving:  Stress period:     3    Time step:     7\n",
      "    Solving:  Stress period:     3    Time step:     8\n",
      "    Solving:  Stress period:     3    Time step:     9\n",
      "    Solving:  Stress period:     3    Time step:    10\n",
      "    Solving:  Stress period:     3    Time step:    11\n",
      "    Solving:  Stress period:     3    Time step:    12\n",
      "    Solving:  Stress period:     3    Time step:    13\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     3    Time step:    14\n",
      "    Solving:  Stress period:     3    Time step:    15\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     3    Time step:    16\n",
      "    Solving:  Stress period:     3    Time step:    17\n",
      "    Solving:  Stress period:     3    Time step:    18\n",
      "    Solving:  Stress period:     3    Time step:    19\n",
      "    Solving:  Stress period:     3    Time step:    20\n",
      "    Solving:  Stress period:     3    Time step:    21\n",
      "    Solving:  Stress period:     3    Time step:    22\n",
      "    Solving:  Stress period:     3    Time step:    23\n",
      "    Solving:  Stress period:     3    Time step:    24\n",
      "    Solving:  Stress period:     3    Time step:    25\n",
      "    Solving:  Stress period:     3    Time step:    26\n",
      "    Solving:  Stress period:     3    Time step:    27\n",
      "    Solving:  Stress period:     3    Time step:    28\n",
      "    Solving:  Stress period:     3    Time step:    29\n",
      "    Solving:  Stress period:     3    Time step:    30\n",
      "    Solving:  Stress period:     3    Time step:    31\n",
      "    Solving:  Stress period:     3    Time step:    32\n",
      "    Solving:  Stress period:     3    Time step:    33\n",
      "    Solving:  Stress period:     3    Time step:    34\n",
      "    Solving:  Stress period:     3    Time step:    35\n",
      "    Solving:  Stress period:     3    Time step:    36\n",
      "    Solving:  Stress period:     3    Time step:    37\n",
      "    Solving:  Stress period:     3    Time step:    38\n",
      "    Solving:  Stress period:     3    Time step:    39\n",
      "    Solving:  Stress period:     3    Time step:    40\n",
      "    Solving:  Stress period:     3    Time step:    41\n",
      "    Solving:  Stress period:     3    Time step:    42\n",
      "    Solving:  Stress period:     3    Time step:    43\n",
      "    Solving:  Stress period:     3    Time step:    44\n",
      "    Solving:  Stress period:     3    Time step:    45\n",
      "    Solving:  Stress period:     3    Time step:    46\n",
      "    Solving:  Stress period:     3    Time step:    47\n",
      "    Solving:  Stress period:     3    Time step:    48\n",
      "    Solving:  Stress period:     3    Time step:    49\n",
      "    Solving:  Stress period:     3    Time step:    50\n",
      "    Solving:  Stress period:     3    Time step:    51\n",
      "    Solving:  Stress period:     3    Time step:    52\n",
      "    Solving:  Stress period:     3    Time step:    53\n",
      "    Solving:  Stress period:     3    Time step:    54\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     3    Time step:    55\n",
      "    Solving:  Stress period:     3    Time step:    56\n",
      "    Solving:  Stress period:     3    Time step:    57\n",
      "    Solving:  Stress period:     3    Time step:    58\n",
      "    Solving:  Stress period:     3    Time step:    59\n",
      "    Solving:  Stress period:     3    Time step:    60\n",
      "    Solving:  Stress period:     3    Time step:    61\n",
      "    Solving:  Stress period:     3    Time step:    62\n",
      "    Solving:  Stress period:     3    Time step:    63\n",
      "    Solving:  Stress period:     3    Time step:    64\n",
      "    Solving:  Stress period:     3    Time step:    65\n",
      "    Solving:  Stress period:     3    Time step:    66\n",
      "    Solving:  Stress period:     3    Time step:    67\n",
      "    Solving:  Stress period:     3    Time step:    68\n",
      "    Solving:  Stress period:     3    Time step:    69\n",
      "    Solving:  Stress period:     3    Time step:    70\n",
      "    Solving:  Stress period:     3    Time step:    71\n",
      "    Solving:  Stress period:     3    Time step:    72\n",
      "    Solving:  Stress period:     3    Time step:    73\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     3    Time step:    74\n",
      "    Solving:  Stress period:     3    Time step:    75\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     3    Time step:    76\n",
      "    Solving:  Stress period:     3    Time step:    77\n",
      "    Solving:  Stress period:     3    Time step:    78\n",
      "    Solving:  Stress period:     3    Time step:    79\n",
      "    Solving:  Stress period:     3    Time step:    80\n",
      "    Solving:  Stress period:     3    Time step:    81\n",
      "    Solving:  Stress period:     3    Time step:    82\n",
      "    Solving:  Stress period:     3    Time step:    83\n",
      "    Solving:  Stress period:     3    Time step:    84\n",
      "    Solving:  Stress period:     3    Time step:    85\n",
      "    Solving:  Stress period:     3    Time step:    86\n",
      "    Solving:  Stress period:     3    Time step:    87\n",
      "    Solving:  Stress period:     3    Time step:    88\n",
      "    Solving:  Stress period:     3    Time step:    89\n",
      "    Solving:  Stress period:     3    Time step:    90\n",
      "    Solving:  Stress period:     3    Time step:    91\n",
      "    Solving:  Stress period:     3    Time step:    92\n",
      "    Solving:  Stress period:     3    Time step:    93\n",
      "    Solving:  Stress period:     3    Time step:    94\n",
      "    Solving:  Stress period:     3    Time step:    95\n",
      "    Solving:  Stress period:     3    Time step:    96\n",
      "    Solving:  Stress period:     3    Time step:    97\n",
      "    Solving:  Stress period:     3    Time step:    98\n",
      "    Solving:  Stress period:     3    Time step:    99\n",
      "    Solving:  Stress period:     3    Time step:   100\n",
      "    Solving:  Stress period:     3    Time step:   101\n",
      "    Solving:  Stress period:     3    Time step:   102\n",
      "    Solving:  Stress period:     3    Time step:   103\n",
      "    Solving:  Stress period:     3    Time step:   104\n",
      "    Solving:  Stress period:     3    Time step:   105\n",
      "    Solving:  Stress period:     3    Time step:   106\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     3    Time step:   107\n",
      "    Solving:  Stress period:     3    Time step:   108\n",
      "    Solving:  Stress period:     3    Time step:   109\n",
      "    Solving:  Stress period:     3    Time step:   110\n",
      "    Solving:  Stress period:     3    Time step:   111\n",
      "    Solving:  Stress period:     3    Time step:   112\n",
      "    Solving:  Stress period:     3    Time step:   113\n",
      "    Solving:  Stress period:     3    Time step:   114\n",
      "    Solving:  Stress period:     3    Time step:   115\n",
      "    Solving:  Stress period:     3    Time step:   116\n",
      "    Solving:  Stress period:     3    Time step:   117\n",
      "    Solving:  Stress period:     3    Time step:   118\n",
      "    Solving:  Stress period:     3    Time step:   119\n",
      "    Solving:  Stress period:     3    Time step:   120\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     4    Time step:     1\n",
      "    Solving:  Stress period:     4    Time step:     2\n",
      "    Solving:  Stress period:     4    Time step:     3\n",
      "    Solving:  Stress period:     4    Time step:     4\n",
      "    Solving:  Stress period:     4    Time step:     5\n",
      "    Solving:  Stress period:     4    Time step:     6\n",
      "    Solving:  Stress period:     4    Time step:     7\n",
      "    Solving:  Stress period:     4    Time step:     8\n",
      "    Solving:  Stress period:     4    Time step:     9\n",
      "    Solving:  Stress period:     4    Time step:    10\n",
      "    Solving:  Stress period:     4    Time step:    11\n",
      "    Solving:  Stress period:     4    Time step:    12\n",
      "    Solving:  Stress period:     4    Time step:    13\n",
      "    Solving:  Stress period:     4    Time step:    14\n",
      "    Solving:  Stress period:     4    Time step:    15\n",
      "    Solving:  Stress period:     4    Time step:    16\n",
      "    Solving:  Stress period:     4    Time step:    17\n",
      "    Solving:  Stress period:     4    Time step:    18\n",
      "    Solving:  Stress period:     4    Time step:    19\n",
      "    Solving:  Stress period:     4    Time step:    20\n",
      "    Solving:  Stress period:     4    Time step:    21\n",
      "    Solving:  Stress period:     4    Time step:    22\n",
      "    Solving:  Stress period:     4    Time step:    23\n",
      "    Solving:  Stress period:     4    Time step:    24\n",
      "    Solving:  Stress period:     4    Time step:    25\n",
      "    Solving:  Stress period:     4    Time step:    26\n",
      "    Solving:  Stress period:     4    Time step:    27\n",
      "    Solving:  Stress period:     4    Time step:    28\n",
      "    Solving:  Stress period:     4    Time step:    29\n",
      "    Solving:  Stress period:     4    Time step:    30\n",
      "    Solving:  Stress period:     4    Time step:    31\n",
      "    Solving:  Stress period:     4    Time step:    32\n",
      "    Solving:  Stress period:     4    Time step:    33\n",
      "    Solving:  Stress period:     4    Time step:    34\n",
      "    Solving:  Stress period:     4    Time step:    35\n",
      "    Solving:  Stress period:     4    Time step:    36\n",
      "    Solving:  Stress period:     4    Time step:    37\n",
      "    Solving:  Stress period:     4    Time step:    38\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     4    Time step:    39\n",
      "    Solving:  Stress period:     4    Time step:    40\n",
      "    Solving:  Stress period:     4    Time step:    41\n",
      "    Solving:  Stress period:     4    Time step:    42\n",
      "    Solving:  Stress period:     4    Time step:    43\n",
      "    Solving:  Stress period:     4    Time step:    44\n",
      "    Solving:  Stress period:     4    Time step:    45\n",
      "    Solving:  Stress period:     4    Time step:    46\n",
      "    Solving:  Stress period:     4    Time step:    47\n",
      "    Solving:  Stress period:     4    Time step:    48\n",
      "    Solving:  Stress period:     4    Time step:    49\n",
      "    Solving:  Stress period:     4    Time step:    50\n",
      "    Solving:  Stress period:     4    Time step:    51\n",
      "    Solving:  Stress period:     4    Time step:    52\n",
      "    Solving:  Stress period:     4    Time step:    53\n",
      "    Solving:  Stress period:     4    Time step:    54\n",
      "    Solving:  Stress period:     4    Time step:    55\n",
      "    Solving:  Stress period:     4    Time step:    56\n",
      "    Solving:  Stress period:     4    Time step:    57\n",
      "    Solving:  Stress period:     4    Time step:    58\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     4    Time step:    59\n",
      "    Solving:  Stress period:     4    Time step:    60\n",
      "    Solving:  Stress period:     4    Time step:    61\n",
      "    Solving:  Stress period:     4    Time step:    62\n",
      "    Solving:  Stress period:     4    Time step:    63\n",
      "    Solving:  Stress period:     4    Time step:    64\n",
      "    Solving:  Stress period:     4    Time step:    65\n",
      "    Solving:  Stress period:     4    Time step:    66\n",
      "    Solving:  Stress period:     4    Time step:    67\n",
      "    Solving:  Stress period:     4    Time step:    68\n",
      "    Solving:  Stress period:     4    Time step:    69\n",
      "    Solving:  Stress period:     4    Time step:    70\n",
      "    Solving:  Stress period:     4    Time step:    71\n",
      "    Solving:  Stress period:     4    Time step:    72\n",
      "    Solving:  Stress period:     4    Time step:    73\n",
      "    Solving:  Stress period:     4    Time step:    74\n",
      "    Solving:  Stress period:     4    Time step:    75\n",
      "    Solving:  Stress period:     4    Time step:    76\n",
      "    Solving:  Stress period:     4    Time step:    77\n",
      "    Solving:  Stress period:     4    Time step:    78\n",
      "    Solving:  Stress period:     4    Time step:    79\n",
      "    Solving:  Stress period:     4    Time step:    80\n",
      "    Solving:  Stress period:     4    Time step:    81\n",
      "    Solving:  Stress period:     4    Time step:    82\n",
      "    Solving:  Stress period:     4    Time step:    83\n",
      "    Solving:  Stress period:     4    Time step:    84\n",
      "    Solving:  Stress period:     4    Time step:    85\n",
      "    Solving:  Stress period:     4    Time step:    86\n",
      "    Solving:  Stress period:     4    Time step:    87\n",
      "    Solving:  Stress period:     4    Time step:    88\n",
      "    Solving:  Stress period:     4    Time step:    89\n",
      "    Solving:  Stress period:     4    Time step:    90\n",
      "    Solving:  Stress period:     4    Time step:    91\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     4    Time step:    92\n",
      "    Solving:  Stress period:     4    Time step:    93\n",
      "    Solving:  Stress period:     4    Time step:    94\n",
      "    Solving:  Stress period:     4    Time step:    95\n",
      "    Solving:  Stress period:     4    Time step:    96\n",
      "    Solving:  Stress period:     4    Time step:    97\n",
      "    Solving:  Stress period:     4    Time step:    98\n",
      "    Solving:  Stress period:     4    Time step:    99\n",
      "    Solving:  Stress period:     4    Time step:   100\n",
      "    Solving:  Stress period:     4    Time step:   101\n",
      "    Solving:  Stress period:     4    Time step:   102\n",
      "    Solving:  Stress period:     4    Time step:   103\n",
      "    Solving:  Stress period:     4    Time step:   104\n",
      "    Solving:  Stress period:     4    Time step:   105\n",
      "    Solving:  Stress period:     4    Time step:   106\n",
      "    Solving:  Stress period:     4    Time step:   107\n",
      "    Solving:  Stress period:     4    Time step:   108\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    Solving:  Stress period:     4    Time step:   109\n",
      "    Solving:  Stress period:     4    Time step:   110\n",
      "    Solving:  Stress period:     4    Time step:   111\n",
      "    Solving:  Stress period:     4    Time step:   112\n",
      "    Solving:  Stress period:     4    Time step:   113\n",
      "    Solving:  Stress period:     4    Time step:   114\n",
      "    Solving:  Stress period:     4    Time step:   115\n",
      "    Solving:  Stress period:     4    Time step:   116\n",
      "    Solving:  Stress period:     4    Time step:   117\n",
      "    Solving:  Stress period:     4    Time step:   118\n",
      "    Solving:  Stress period:     4    Time step:   119\n",
      "    Solving:  Stress period:     4    Time step:   120\n",
      " \n",
      " Run end date and time (yyyy/mm/dd hh:mm:ss): 2022/07/30 16:16:42\n",
      " Elapsed run time:  1.483 Seconds\n",
      " \n",
      " Normal termination of simulation.\n",
      "\n",
      "Success is:  True\n"
     ]
    }
   ],
   "source": [
    "# Run the simulation\n",
    "success, buff = sim.run_simulation()\n",
    "print(\"\\nSuccess is: \", success)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Post-Process Head Results\n",
    "\n",
    "First, we get the simulated head data using the `.output.head()` method and the `get_data` function, by specifying, in this case, the step number and period number for which we want to retrieve data. A three-dimensional array is returned of size `nlay, nrow, ncol`. FloPy plotting methods are used to make contours of the head in a specific layer (in this case, layer 1). FloPy plotting methods are also used to plot the model grid and the location of GHB cells in the model domain."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:42.926413Z",
     "iopub.status.busy": "2022-07-30T20:16:42.926108Z",
     "iopub.status.idle": "2022-07-30T20:16:42.980025Z",
     "shell.execute_reply": "2022-07-30T20:16:42.979365Z"
    },
    "tags": []
   },
   "outputs": [],
   "source": [
    "# Retrieve the head data using the .output() method\n",
    "h = gwf.output.head().get_data(kstpkper=(0, 0))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:42.983852Z",
     "iopub.status.busy": "2022-07-30T20:16:42.983564Z",
     "iopub.status.idle": "2022-07-30T20:16:43.166162Z",
     "shell.execute_reply": "2022-07-30T20:16:43.165105Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAaUAAAI/CAYAAADXx8YWAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAAsTAAALEwEAmpwYAACDhUlEQVR4nO3dd3hcV5n48e+ZolHv1ZaL3FsctzjV6YmTkA6kUaIAG350tkLYZdkFLxsWdllYFpYAQaGE9N6d3osd915lW9ao9xlNu+f3x4wcEVxkWbrnju77eR4/kkbjue9czdx3zjnvOUdprRFCCCGcwGM6ACGEEGKAJCUhhBCOIUlJCCGEY0hSEkII4RiSlIQQQjiGJCUhhBCO4TMdwNHk5OTo8vJy02EcEovFAPD7/YYjSZJ4js5p8YDzYpJ4js6J8SSAmEPiAfDHYgQbGlq11mUj8XiOTkrl5eXs2bPHdBiH1NXVAVBbW2s0jgESz9E5LR5wXkwSz9E5MZ5XgDqHxANQW1dH3S231I/U40n3nRBCCMeQpCSEEMIxJCkJIYRwDElKQgghHEOSkhBCCMeQpCSEEMIxJCkJIYRwDElKQgghHEOSkhBCCMeQpCSEEMIxJCkJIYRwDElKQgghHEOSkhBCCMeQpCSEEMIxJCkJIYRwDElKQgghHEOSkhDHqTN6kCcO/AvB8FbToYya/kSM29bezfPBDaZDOaJtba18/olH2NraYuT44WiMb/zqCZ5/f7ttx9y7M0TrM62gtW3HtJskJSGOk6Xj7Op9na7YQdOhjJqAx8drzVvZ3HnAdChHFE0kWLl7F/u7uowcP9Pv4+X1u9m4N2jbMSMRi/59/QSau207pt0kKQlxnLK8BQCE42YuhnZQSlGemU9zv3OfY0lWFgCt4ZCR4yulqCjMpamj17ZjFpf6AcjfZl8itJskJSGOU8CbCyj6E2P30ypAeaCAJkcnpWwA2kJmkhJARVEeTR09th2vqMQPCvK2Ndp2TLtJUhLiOHmUl0xPLuExnpQqsgoc3VIK+HzkZQRoM9RSAigvsrel5PN58BX6yN8qLSUhxCA+TyZxHTEdxqgqDxTQHOnG0pbpUI6oJDvbaEupsiiPls5eLMu+woOM0gzytktSEkIM4lFeLB03HcaoqsgqIKEt2iP2tQSOV1l2Nk195uIrL8wlblm09fTZdkx/mZ9AWy8ZLfZ1G9pJkpIQw+BRPiydMB3GqCrPTBZ0NEec2005saCQ+q5OY8evLMoDsLULL6M0A4D8MTquJElJiGHw4sNijLeUUkmpKdxpNpCjqCksormvj75o1MjxK4oHkpJ9rRZ/iR+tIG+MVuBJUhJiGDzKS2KMd9+VZ+YD0Nzv3JbS5MIiAPZ2dhg5fkVhLmBvUvL4PYQmlkhLSQjxgWT33dhOSoX+HDI8PkeXhdcUFgKwt7PTyPELc7PIzcygvsne43fPrCJvjFbgSVISYhjckJTSYQLtpFRLaY+hlpJSiinjStjV2GrrcXtmVZHZ2kNGm3OLUIZLkpIQw+B1QaEDOH8CbbbfT2VOrrHuO4CpVSXsPtiGtnE9uu6ZlcDYHFeSpCTEMHhcUOgAzp9AC1BTVGSspQQwdVwJnX39tHXbN1+qZ0YyKY3FcaVjJiWl1Eyl1NpB/7qVUl9XShUrpVYqpXakvhal7q+UUj9VSu1USq1XSi0a9Fg3p+6/Qyl182g+MSFGkxsKHSA9JtBOLiwy21IaVwrArsY2246ZyAnQN7HYnS0lrfU2rfUCrfUCYDEQAh4Gvgm8oLWeDryQ+hngUmB66t+twC8AlFLFwHeAU4GlwHcGEpkQ6cYNY0qQHhNoawqL6Ojvp7M/bOT4U6tKANh10OZxpRmV7mwpfcgFwC6tdT1wFXBX6va7gKtT318F/E4nvQ0UKqWqgOXASq11u9a6A1gJXHKiT0AIE9ySlNJhAu1kwxV4JfnZFOZk2tpSAuieVUVmUzf+DvtWk7CDOp7BOaXUncD7WuufKaU6tdaFqdsV0KG1LlRKPQHcrrV+PfW7F4BvAOcCmVrrFanbvw2EtdY/OtLxqqur9YoVK4b3zEZBMJhsKldWVhqOJEniObrRjCc07XUSOW3krbvKMTENx7HiafH1c3/pPpZ3VDE1kmc8nsNpScT5RV8XV2fmMj8jYCSeBzb3YGnNdXPzR/T4h4unCwhWVtLf0E/rk62UXlpK5oTMUT3u0VQGg9x+222rtdZLRuLxhtxSUkplAFcC93/4dzqZ2Uak9EQpdatSapVSapWd1SxCHBftQXucO84yUnITPgD6vM5tFRZ5vCigzTJXDVmS5aU9nLC1Am9guaFoq5nVLEaL7zjueynJVlJT6ucmpVSV1rox1T3XnLq9AZgw6P9Vp25rINlaGnz7yx8+iNb6DuAOgJqaGl1bW3scIY6uuro6AJwSk8RzdKMZz/ONbezt6zrux063c6S15o8rv8Pk+bOonXWp8XiO5Pd1v6agspLaSy43Ek/2K+vYcM+LXHrVx6gsHr0WZV1dHduAulQ8Zzz/vzTlVLCh9mOjdsxjqU2do5FyPGNKNwJ/GvTzY8BABd3NwKODbv90qgrvNKBLa90IPAtcrJQqShU4XJy6TYi045YxpXSYQAvJYoe9HWbLwsHeCjyAnpljr9hhSElJKZUDXAQ8NOjm24GLlFI7gAtTPwM8BewGdgK/Ar4IoLVuB74HvJf6993UbUKkHa/yuaIkHJw/gRaSxQ57Oztt7T4bzFQFXvfMKrIau/B1m6k8HA1D6r7TWvcBJR+6rY1kNd6H76uBLx3hce4E7jz+MIVwFjdsXTGgIquAte17TYdxVDVFRfTGorSGQ5Rl59h+/MLcLErzs9l10OaW0owKAPK2B+lYUmPrsUeLrOggxDC4pfsO0mQCbUFqDTyjXXiltnff9U1OTtzN3j92Op0kKQkxDAPLDLmhQjQtJtAWmd3CAmBKVQm7G9ts3Ro9UpZPIuAjZ5+9yXA0SVISYhi8yguAxdjvwkuHCbTj8vLxezxG18CbNq6E/micg202jr95FKGJJWTvk5aSEK7mUcnhWDd04aXDDrQ+j4cJBQXGVnUAmDJQgWfzuFJoQjHZ0lISwt3clJTSYQdaSJaFG10tvMpMWXjfxBKyGjtRsbHRapekJMQweEh137mgAq/Qn4NfeR3dfQfJ1cLrOzuxDI3z5WYFqCzOM9JSUpYmq8FcQh5JkpSEGAav8gO4Yk8lpRR5/ix6Y86eC1NTWEQkESfY22MshqlVJba3lEKTki20sVLsIElJiGHwpAod3DKBNteXSU+833QYRzXZ8NbokFzZYU+wnXjCvvL50IRkUhor40qSlIQYBjeNKQHk+TPpdXhSqikcKAvvNBbD1KoSYvEEB1rsiyGel0mkKGfMzFWSpCTEMLgtKeX6suiNOTspVeTmkunzGZ9ACwYq8CaOnQo8SUpCDMNAoUPCBYUOALn+AD1xZ48peZQyvjV6TVUxSsFOm9fAG0tzlSQpCTEM3oGWkgsKHQDyfFn0xiOmwzimyQWF7O0yl5SyMvyMLy1gt93FDhNLCHT04etxdmt2KCQpCTEMruu+82c6vvoOksUO+7q6iFvm1umbNq6UnQbKwgGy96d/F54kJSGGYSApJXTMcCT2yPdlEbHi9Cecvcvp1OJi4pbFjnZzF+dp40rZ19xBZ699SXxgYda87UHbjjlaJCkJMQx+TyYAccv5XVojoSIrudRQo4OXGgI4e9JkPErx9I7txmK4cNF0Epbm8bc323bM0IRi+iYWU7lyk23HHC2SlIQYBp8KABDX7khK47OS3UMHw85eNaAsO4czJ0zksW1bjK3gPqO6jJOnVPHga+vtWzFcKYLLT6JozT4CQWdvyHgskpSEGIaBllLMSv+B5aEYl5WcA3Qw5PwKrytnzmZfdxdrg+a2Cf/Y2Sezr7mT97bvt+2YjcvnAVD53EbbjjkaJCkJMQyHWkou6b4rCeQR8Pgc31ICuHjKNDK8Xh7bvtVYDBcumk5BTiYPvLrOtmP2jyui4+QJVD27AdJ4ny9JSkIMw6GWknZHS0kpxbisIhrCzm8p5QUCnD95Ck/u2EbCUBVewO/jytPn8vK6XbR02rc5YnD5PHL3tJK3vcm2Y440SUpCDIPvUKGDO5ISwLjsYg6GnN9SArhs+gxaQyFWNx40FsNHl51EwtI88qZ93WlN58/B8nupfHaDbcccaZKUhBgGr/LhwUfMJYUOAJNzytjb10I47uyycIBzJ08hw+vlmV07jMUwsbyIU2dN5KHXN9i2QGs8P4vWM6ZRuXITKm5urtaJkKQkxDD5PAFXtZROL5tB1IrzbttO06EcU25GBssmTuLZnTuMVeEBfOzs+TR19PLGpj22HbNx+UkE2nopWm3fMUeSJCUhhsnvyXRN9R3AoqIacn2ZvNq8xXQoQ3LJtBk09vawvtnc+MrZ86dQWpDDA6+ut+2YrWdMI1qQxaS7307LggdJSkIMk09luqbQAcDn8XJG6Qxeb9lKQju/a+jCmqn4PB6e2WluIq3f6+WaM+fx5ua9NLTaM39IZ/jYU7uMkvf2UPKm81u1HyZJSYhh8nsCrikJH7CsfDYd0T42dR4wHcoxFWRmcnr1BJ4x3IV3zZknoVA89Lp9xQcHPrqYvonFzPif51Hx9FrJXpKSEMPktu47SI4reZWH19KkC2/51OnUd3Wytc3erSQGqyzOY9lJNTz65iZiNiUI7fOy48sXkrOvjeqH37flmCNFkpIQw+RTAdcsMzQg35/FwqLJvNqSHknpoqnTUMCzO81V4QFce9ZJtPeEeGPTXtuO2XrmdNqWTKbmN6/i63b+Cu8DJCkJMUxubCkBnF0+mz29zezvc/42CWXZOZwyrtpoaTjAaXMmUZSbxVPv2pjMlWLHVy/C39vPlN++Zt9xT5AkJSGGyacyXVUSPmBZ+WwAXkuT1tIl06azva2V3R3mVqPwe71ccspMXl2/m56Qfa+Z3mkVNFy+gOoHVqXNdumSlIQYJr8n03XddwDjs4uZmluRNuNKF0+dBsAzhrvwLls6m2g8wfPv2xvH7lvPwQr4mP6zF2w97nBJUhJimHyegCu77yDZWlrbUU9XNGQ6lGMal5fPyRWVxrvw5kyqYFJ5EU/a2YUHRItz2fPpMyl7fTtFq5w/oVaZLJU8lurqar1ixQrTYRwSDCZ3daysrDQcSZLEc3SjHU//hDVEqrZQ8O5NjonpeA03niZ/mAdL9nNhZyUz+vONx3Msb0TCvBAJ8dXcQgo9XmPxvNsQ5u0D/dQuyCc/MPQ4BsfTBQSPMx4d1wTvD+Lxeyi/thzlUcd97COpDAa5/bbbVmutl4zE4/lG4kFGS8KyeKW+3nQYh1QCBYGA6TCEU1he8FhoLJTLOh3KY5lkJ7zsCfSOaFIaLbP9GbwQCbE1FuW0QJaxOGaWZPD2gX62tUY5Zbx9cSifomBpAe0vtNO3rY/c2bm2Hft4OTopxYA600EMUgvMrKyktrbWcCRJdXV1ABLPEYx2PO+33c9rLeu56VPXE/DmOCKm43Ui8TRufJiVjev5xFWfxO8ZmUvJaJ6fF+7+HW1+P7Ufv9FoPOt+dC/Bvn5uvvnTKHV8LZa6ujq2AXXDiUdrlrTcRWRjBw//840kckbmA3Zt6hyNFHd9vBNiBB3avsJFSw0Ntqx8Nn2JCKvbnT9OAcmJtO83HqS5z779jQ7nslNnsyfYztb9zfYeWCm2f/ViAh19TP7dG/Ye+zhIUhJimPye5CfNmMuWGhpwSslUAh5/2lThXTJtOhp4bpfZ9eAuWjQDv8/LU+/avzNu95xxNC6fx8R73yGzsdP24w+FJCUhhsmnUrvPurQCL9Pr59TSabzavMXo2nJDNb24hClFRcar8ApyMjlrXg3PvLfVtn2WBtv5hfNBKab9/EXbjz0UkpSEGCa/y7vvILm6Q1N/Fzt6Gk2HckxKKS6ZOoN3DuynI2x22Z3Lls6irTvEu9v22X7sSHk+9Z84neLVe/F3Oq+kX5KSEMM0kJTc2lICOLNsJgqVRnssTSehNSt3m+3CWzavhrysAE++Y+a87f3kGbx5zxeIFWYbOf7RSFISYph8Kjmm5LbtKwYrCeQxr3ACrzbbPz4yHHPLyqnOz+dZw+NKGX4fFy2ewUtrdxLqt397eSvTTzzfXGn80UhSEmKYpPsu6ezy2WztbqCp355N7E6EUorlU6fzxr56eiJmP0xctnQW/dE4L63bZTQOp5GkJMQw+aT7DoBl5bMAeD1NWkuXTJtO1Erw4t7dRuNYMHU840ryjXXhOZUkJSGGya8GWkru7b4DqMkpZ1xWEW+1mtt2/HgsrBxHeU6O8QVaPR7FpafM4t2t+2jpMjt3ykkkKQkxTH5Psk++P9FjOBKzlFIsLp7C2o69WNr+Eufj5VGKi6dM45X6PYRiMaOxXHbqbCytuf+V9UbjcBJJSkIMk8+TQaF/PC39MiawsGgy3bEwe3ptXqVgmC6ZNoP+eJxX6/cajaOmspiLF8/gD8+vJtju7g83AyQpCXECyjNn0NS/zXQYxi0srgFgTcdes4EM0dLx1RRlZvLMLvNdjl+9ZhmW1vzPo6+bDsURJCkJcQIqMmfQG28hFO8wHYpR47KKKA/ksyZN1sHzeTxcNGUaL+7ZTSQeNxrLuJJ8PnHBIp5+dysb9waNxuIEkpSEOAHlWTMAaO43O2humlKKBcWTWdtRnxZLDkGyC683GuWN/favqvBhn7lkKSX52fzn/S+nzfkbLZKUhDgBZYHkVttN/ea7gUxbUDSZlkg3DeF206EMyenVE8jNyOBZw2vhAeRkZvDFK85g3e5GVr7v7teSJCUhTkDAm0NRxgSawukxR2c0LSxKjiu91LTJcCRDE/D5uKBmKit376Q3av+qCh925RlzmVFdxn898Cr7mjtNh2OMJCUhTtCE7EXsD60hZpld5NO0KbnlnFY6nf/bvpLNXQdMhzMktScvpKu/n/9++03ToeD1ePiXT19MNBbnlh/ew+Z6d44vSVIS4gRNy1tGXEeo71tlOhSjlFJ8d/71lAby+MaaP9IZ7TMd0jGdXFnFDfPmc9e699nSYr6cfdaEcu78u+vJCvj5qx8/wBub9poOyXaSlIQ4QeOzTyLTm8/OntdMh2JcYUY2/77wJjqiffzTuntJpMFk2n84YxmFmZl8+6XnsRxQZDC5spjf/v31TCwr5K9//ihPvL3ZdEi2kqQkxAnyKC9Tc89kT+/bxC3zYxOmzSmo5u/nXMm7bTv55Y7nTYdzTAWZmdx21jm8H2zkvk0bTIcDQFlBLr/6m4+zaPp4/vmuZ6l79j3XVOVJUhJiBEzLW0bUCrE/9L7pUBzhquolXFW9hLrdL/NKk/M/6V8zaw6njq/mB2+8RlvIGRvf5WYF+J8vX8MlS2by00de50cuKReXpCTECJiQs5AMT7Z04Q3yd7OvYFb+eP5lw/3s62s1Hc5RKaX47rkX0heLcvsbr5oO5xC/z8uKWy7lkxcs4k8vreXpnX1Y1thOTJKUhBgBXuWnJvd0dve8iaUTpsNxhIDXzw8W3oRXefjGmj8Sjju7a3N6SQmfW7iEB7ds4t0G51QPejyKv/nYOfz1R89mZ3uM3Vt78YTH7sr0kpSEGCHT8pbRb/VwILTOdCiOUZVVxPfmX8/u3mZu3/yo6XCO6StLT2N8Xj7//NLzxBLO+nDxqQsXc/HUbEK9cWb85334O8fmdheSlIQYIZNyluBTmdKF9yGnl83gs1PP4+mDaxw/vpTl9/Mv55zP9vY27ly72nQ4f2FWaYCamblktHYx44f3EGhsMx3SiJOkJMQI8XsymZx7Crt63kCnQSm0nW6Zei7T8ir5weZH6Y45e5LxBVOmctGUqfz0nbfotJzVWgLIK/Cz428+joonmPGf95Gz66DpkEaUJCUhRtC0vGWEEu00hp3dIrCb3+Pj2/M+Ske0j59sfcp0OMf0z+ecD8Cz/c6cAByeWMH2v7ueeHYm037yAPnrxs6eXpKUhBhB1dkLAGiJ7DQbiAPNLhjPp2qW8XjDat5qcfaio+Pz8vnqqaezLR5jW8yZBRrRskJ2/P31hMeVMuWXj1PyujPmWJ0oSUpCjKBsbxFelUF3rMl0KI702annU5NTxvc3PUxvvN90OEf1mQWLKfN4eaa/z/i26UcSz8tm59c/RvecSUz84/NUPvEWpPlcJklKQowgpRR5/nJJSkcQ8Pr5p5M+SnN/Nz/b9ozpcI7K7/VyWWYOXdriZ+++bTqcI7IyM9j9hStpO30OVU++zYS7X4BE+o5pDikpKaUKlVIPKKW2KqW2KKVOV0oVK6VWKqV2pL4Wpe6rlFI/VUrtVEqtV0otGvQ4N6fuv0MpdfNoPSkhTMr3VdAjSemITiqcyI2Tz+Sh/e+yqs3ZYyGTfH5O9gf49ZpV7GhzcKWb18u+T11M8JKllL6+gZo7HkdFndm6O5ahtpR+AjyjtZ4FnAxsAb4JvKC1ng68kPoZ4FJgeurfrcAvAJRSxcB3gFOBpcB3BhKZEGNJvr9CWkrH8P+mX0h1djH/tvFhx0+qvTCQTY4/g288/yyd/Q6uHFSKxqvOZP/151GwYTfTfvIg3j5nd5EezjGTklKqADgb+A2A1jqqte4ErgLuSt3tLuDq1PdXAb/TSW8DhUqpKmA5sFJr3a617gBWApeM4HMRwhHy/BWEE53ErPS7INgl05vBP827loZwO7/YsdJ0OEeV4/Hwb+dfyKaWJq685w+sb3L2Pket5y5g7+c+Qva+Zmb86F78bd2mQzou6lgL/CmlFgB3AJtJtpJWA18DGrTWhan7KKBDa12olHoCuF1r/Xrqdy8A3wDOBTK11itSt38bCGutf3SkY1eNG6cv+fznT+T5jahKoCAQoLKy0nQoAASDyTeHxHN4puKJlu4mPO1Nctdegbe/wBExHYnpeF7Na2JjdhfXtE+gKpZlPJ4PGxxPQyLGA6FeerXF8swcFvsDJC999sbTBQSHcH56u2Ps2d6HxwNTZuWRle0dlZgqg0Fuv+221VrrJSPxeL4h3mcR8BWt9TtKqZ/wQVcdAFprrZQakZIPpdStJLv9GDduHOdMmjQSDzsigsEgXZEI2+rrTYcCfJAkhbN4IrkAWIG+v0hK4s+d1lvG3sw+XioIcl2rc97rhzPe6+evcgp4pL+Xp/r72J+I8ZHMXDJsTkxDlZvvZ9qcPPZs7WHn5h5qZuSQm+83HdYxDSUpHQAOaK3fSf38AMmk1KSUqtJaN6a65wa2bWwAJgz6/9Wp2xpItpYG3/7yhw+mtb6DZMuMmpoaXVtbO9TnMurq6urYVl9PnelAUmqBmZWVOOUc1dXVAbg+np5YM3fueo7Tzp3PSUWXOyKmI3FCPPNbd/DVVb8lcvp4Kt/yGI9nsMOdn/+nNb9Y9Q4/fvtN+v2K/73sCqYVl9gWzzag7jjOj7+9m6n/8zA7dnRRX3sOnYtnjGhMtalzNFKOOaaktQ4C+5VSM1M3XUCyK+8xYKCC7mZgYLXFx4BPp6rwTgO6tNaNwLPAxUqpolSBw8Wp24QYU3J8JXjw0h2XYoehOK10OleMX8wf9r5Gs8/543AepfjSKadx19UfpT0c4up7/8jj27eaDuuIYsX57Pi76wlNqmDyb56k9KU1pkM6qqFW330F+KNSaj2wAPg+cDtwkVJqB3Bh6meAp4DdwE7gV8AXAbTW7cD3gPdS/76buk2IMcWjvOT6y6Qs/Dh8fdZlFGfk8mJBkATpMfnzzAmTePzGTzG7tIyvPfMk33n5BSLxuOmwDiuRk8nOr36UrvlTmXDfy1Q98rpjJ9kOKSlprddqrZdoredrra/WWndordu01hdoradrrS8cSDCpqrsvaa2naq1P0lqvGvQ4d2qtp6X+/Xa0npQQpklZ+PHJ82fxzblX0+6PsjrXwfOBPqQyN4+7r72Ozy1czO/Xr+WGB++loduZ1W46w8eeWy+nddlJVD77HhPvehYs502ylRUdhBgFeTKB9rgtK5/F9HAe7+e0s6On0XQ4Q+b3evnWsnP5+WVXsqujnSvu+T0v791jOqzD83jYf+MFHLzidKzMDHBgkYYkJSFGQZ6/gt54GwmdnrPqTVnWXU7A8vK9DQ8Sd+C2EUdzybTpPHbDJ6nKzeMzjz3Ef771OgkHtkRQiqbLTuPA9edJUhLCLfL9FYCmJ9ZiOpS0kqm9nN1dztbug/xx7+umwzlukwuLePC6G7luzjz+9713+PQjD9IaCpkO6/AcmJBAkpIQoyKZlJAuvGGYGsnj/Iq5/GrnC+zpbT72f3CYTJ+f2y9czg8uXM77jQe54k+/572DB0yHlTYkKQkxCvJSSUmKHYbn7+ZcSabXz4qND5FI0118Pz5nHg9ddyOZPh83PXgfv35/FcdaQUdIUhJiVOT5ywAlSWmYSgN5/O3sy9nQuY/76t8yHc6wzS4r59EbPslFU6bx/ddf4R+ef5ZoIr3GyuwmSUmIUeBVfnJ9JdJ9dwIuqVrAmWUz+fn25zgQSp8y8Q/LDwT438uu4Gunns6DWzZxy6MP0tXv/EnCpkhSEmKU5Pkr6JFVHYZNKcVtc6/G5/GwYsNDWGnajQfJ5/K1U8/gvy6+lFUHG/jY/X9iX1en6bAcSZKSEKNEdqA9ceWZBXx91kd4v2MPD+1/13Q4J+zqWXP4/TUfpy0c4qP33c37jQdNh+Q4kpSEGCX5/kp6Yy1YWsYQTsSV4xeztGQaP9v2DI3hDtPhnLCl46t54OM3kpsR4KaH7uPJ7dtMh+QokpSEGCX5vgosEvTF03c8xAmUUnxr7jVo4PsbHxkTFWxTiop58LobOam8gq888wS/WPXOmHheI0GSkhCjJM9fDkhZ+EgYl13El2cs5522HTzR8L7pcEZEcVY2f7jm41w5cxY/fPN1bnvhOWJSmSdJSYjRMjCBtiuWPuu4OdlHJ57KgqLJ/HjrkzSExsYGAwGfjx9ffBlfWXoa923eyKcfeYCNze7+ECNJSYhRUpAxjlxfKWvbH5RxpRHgUR6+Pe+jKOCrq35La6THdEgjQinFX592Jv950aVsaW3hynv+wGcfe5i1QXd+mJGkJMQo8So/y8r/Hy2RXWzofNx0OGPChJwSfrz4Zloi3XxtVR09sbDpkEbMNbPn8GrtX/G3p5/JmuBBrr3vbm5+5AFWHWwwHZqtJCkJMYqm553DhOyFvNnyW0Lx9K8cc4L5RZP4wYJPsKe3mb9Z/Tv6E1HTIY2Y/ECAL51yGq/W/hXfOHMZm1uaue6Be/jEQ/fx9oH9riiGkKQkxChSSnFuxVeIWxFeb/m16XDGjNPLZvDdk69jfec+vrnmbmLWBzu+9sUj7Otrpb6vNe22vxiQm5HB5xcv5ZXav+Ifl53LzvZ2bnroPu4KddMejzp219iR4DMdgBBjXXFgIguLP8rq9nvJyV2Or7fMdEhjwoWVJ9E7t5/vb3qYG1//KRpNa6SH8KCWk095mZBdzOTccibnlDE1r4JlZbPJ8mUYjHzosv1+PrtwMZ84aT73bdrIf776MvtCPVTd/yc6l55GeFKNY7egGC5JSkLYYGnpJ9nW/SKhmnfJ3XCp6XDGjKsnnIJGs7JxPcUZuZQE8ijNzKMkIw8LTX1vC3v6mtnV08SrzVtIaIs8XyZXT1jKxyeeRmVWoemnMCSZPj+fPnkh0ffXcn8swpa+Xiofe5hIeQWdp5xGaMrUMZOcJCkJYYMMTxZnl/8/nop/j2jFDtPhjCnXTFjKNROWHvN+MSvOxs793LfvLf645zXu3vs651XM5YZJZ3JS4QRUGlzUfUoxLiOT5z71KXK3bqbwvXeoePJRIqVlyeQ0bXraJydJSkLYZFre2fg6K+mfsJa+eAc5viLTIbmK3+NjYXENC4traAx3cP++t3l0/3s8H9zA3IJqbph0JhdUzsPn8ZoO9di8XnrnnkTv7LnkbttCwXvvUPH040SLS+hcehp902aAJz1LBiQpCWETpRSZe0+hd/6TvNH8Ky4e9w+mQ3KtqqwivjrzUj439XyeOriGe/a+ybfX38tPtz3NxyeextUTllKYkW06zGPzeOidPZfembPJ2bGNwvfeofyZJ4kVvkl/ZRWJ3FziObkkcnNJ5KS+z8lxdMKSpCSEjbz9BQQaZ7PF8xzzCi9jXPY80yG5WrYvwMcmnsa1E5byVusO7tn7Bj/f8Ry/2fUil41bSMAboTgRMB3msXk89M2cTd+MWWTv3EH++rVkHdiPN9SHsv58yw8NJLJz6J09l44zl5mJ9ygkKQlhs0DDPPyTmnmp6X+4cfLP8ag06C4a4zzKw5llMzmzbCa7eoLcW/8WTx1cQ6QszoRINjNatnFa6XQ8yrktDACUIjR9BqHpM5I/a40nHMLX24u3rxdfXx/e3h58fX3ECgrMxnoEkpSEsJmy/Jxd/gWeOvhd1nc8xoLia0yHJAaZmlfJt+ZdwxdnXMy3H/0/NmZ38vXVdzEpp5TrJ53BR8YtSpuScpTCys4hmp0DVJiOZkgcnvaFGJum5S1jQvYi3mqto09WenCkwowcFveV8MmWKXx3/nVkewP8x+bHuPzl2/mfbc8QDHeaDnFMkqQkhAHJlR6+TNyK8Ebzr0yHI47Ci+KScQuoO/2L/PrUz7O0dBp/3PMa17z6I7619k+s79jniuV/7CLdd0IYUhyYyKLij7Gq/R4pekgDSinmF01iftGk9C8pdzBpKQlh0NLST5DrK+Olpv+R7S3SyEBJ+ePnfoN/mHMlPbF+vr3+Xq5+5YfU7XqZzmjIdIhpS5KSEAb5PVmcXf4FWiO7WN/xmOlwxHEaKCm/b9nX+fHim6nJLefnO57jipd/wL9vfJjdve7esG84lJP7Qqurq/WKFStMh3FIMBikKxIhaDqQlEqgIBCgsrLSdChA8vwAEs9RHC4mjSY060Xiua3krbsSTyzLaDwmjYV42nwRNmR3si2rm4TSTIhkc3JfEROi2ShObAmgYDBIFxB0yPkBqAwGuf2221ZrrZeMxOPJmNJxKggEmOmQF8RAktxWX286FCCZJAFHxVMQcP7ER4Uic+8Seuc/Sf/ENWTvOsN0SEaNhde0Aqb5/HSUJ2isCLG/OERGWFF20EtBq2fYyenQlcch52c0ODop+f1+amtrTYdxSF1dHYBjYqqrq2NbfT11pgNJqU19rTMYw2C1wMzKSsf8veDor6E3mnNY5bmHixZ/hfHZJxmPx4Qx9ZqOKzjoQzV6KSq2qKpKEJ0aZ9d4aDzoo63Vg9bHl5xOKJ5RUjvCjydjSkI4xEDRw8tS9DCmaK1ob/OyaaOfHdt8JOKKmilxTjo5Sll5AqWcO4RigiQlIRzig6KH3VL0MCYpOju9bN7kZ/tWH7GYYnJNnPknRymvkOQ0QJKSEA4yLW8ZE7MX81brb2WlhzFL0dXlZcsmP9u2+IlEFJMmx5m/IEpFZRyPx93JSZKSEA7ywUoPUd5ovsN0OGJUKbq7PWzd4mfrZj/9YcXESQnmL4hSWRXH59Mk1/R2F0cXOgjhRkWBCSwq+Tir2v7E3MLLbCt6EKYoenoU27ZmkJtrMW58nAkTE0yYmMCyIBqFWEwRiyoaowp/FIpTPw/8zrLSe7fZwSQpCeFAS0tuYmvX87zc9D/cOPkXsr2FS/T2eti+LYOcHIvcXAt/BmRkaPx+TVa2prPQwvLC1A/9v3icZJKKpb5GFbFYKmlFFdGYIhYFTnCelB0kKQnhQANFD7K9hTv19Xno6/vL0ZVaIOHV3O/Xf5awMjLAn6HJ8Gsy8y38/r/cXFZriMVIJStFV4eHlhbnfdiRpCSEQ03LW8bEnGTRw/T8c8jxFZsOSTiAN6HoTyj6+492L43Pl0paGZoMfyppZSSTWSBDkxFw5niVJCUhHCpZ9PAV/rD7c7zefAfLx33TdEgibSjicYjHFaTZ2rBSfSeEgxVlVLOo5ONs7X6ehtAG0+EIMeokKQnhcEtLbkptb/FTWelBjHmSlIRwOL8ni3MqvkhbZA/rOh41HY4Qo0qSkhBpYGruWUzKWcLbrXX0xdtNhyPEqJGkJEQaUEpxTmqlh9dlpQcxhklSEiJNFGVUs7jkOil6EGOaJCUh0sgpJTeS5yvn+cYfEbeipsMRYsRJUhIijfg9WSwuuZ7OWAN37rqJVW33Ek2k2UQUIY5CJs8KkWbmF15Jnr+CdR0P80bLr1jVdg8Liq9hQdE1ZHrzTIcnxAmRpCREmlFKMSX3NKbknkYwvJX32u7mndbfsab9AeYXXsnC4o+S7SsyHaYQwyJJSYg0Vpk1iyuqv0tL/25Wtd3NqvZ7WdvxMPMKP8Li4uvI9ZeaDlGI4yJJSYgxoCxzCpeO/ydOi9zMe+33sK7jETZ0Ps7sguUsKb6egowq0yEKMSSSlIQYQ4oCE7i46u85reRTrGq/l81dz7Cp8ylm5V/AkpIbKQ5MNB2iEEclSUmIMSg/o5LzK7/G0pJP8H77/WzofIIt3c8zPe9sTim5ibLMD28TJ4QzSFISYgzL9ZdydsUXWFJyI2s7HmJdxyPs6HmFmtzTWVpyk+nwhPgLkpSEcIFsXyFnlH2GRcUfZ13Ho6xtf5B7e7+Cb1YlgYaT0FqjlPO3yhZjn0yeFcJFMr15nFr6SW6Zdjdnld1KIruTvrkreWDfX1Pf+x5aO3M3UuEekpSEcKEMTxaLS64jb83VZO45he5YE48cuI176r/Erp430NoyHeKf8XksppV1UpobNh2KGGXSfSeEiyntI9A0k09d8l22dD3PqrY/8UTDdygJTOaUkpuYnncOHuW1Py6lqS7qZUZ5JzMqOqkp7SbDZ7H+QAl3vjHH9niEfSQpCSHwKj/zCi9lTsHFbO9+mffa7uaZg9/nbf9dLCm5kVkFF+JVo3e50FrTEd1PpGIbk6t382/lXWRnxAE42JlNY1c2k0p6WV1fNmoxCGeQpCSEOMSjvMwquICZ+eexq/cN3mu9m+eDP+Kd1t+zpOR65hRcgs+TMSLH6o21sj/0Pvv71rA/tIbeeCvUQHZfgPcPlLC9qZDtTYX0Rvx8/cJ1tPZmsr5BVqgY64aUlJRSe4EeIAHEtdZLlFLFwL3AZGAvcJ3WukMlS3h+AlwGhIBarfX7qce5Gfin1MOu0FrfNXJPRQgxUpTyMC1vGVNzz6K+7z3ebfsDLzX9lHfb/sii4o9zUuFH8Huyjusx+xM9HAitO5SEOqL7AMj05jMhewETchbx7tN7eH17G/fwQSXglNIuJpf0cP/qqWgtFYJj3fG0lM7TWrcO+vmbwAta69uVUt9M/fwN4FJgeurfqcAvgFNTSew7wBJAA6uVUo9prTtG4HkIIUaBUorJuUuZlHMKB0LreLftj7zW/H+813Y3C4s+yslFVxHw5h72/8atCAfDGw8loeb+HWgsfCqT8dnzmVdwKRNyFlIamIJSyZqr1ZE64M+3ez9/1gF6+328u6dilJ+tcIIT6b67Cjg39f1dwMskk9JVwO90srb0baVUoVKqKnXflVrrdgCl1ErgEuBPJxCDEMIGSikm5CxgQs4CGsObebf1j7zV+ltWt9/HgqKrWVB0LQFvLs39O9jX9z77Q2toDG8koWN48FKZNZulJZ9gQs4iKrNm4VX+IR23Ir+PeePbeXrjRGIJ+wsuhP3UUOYlKKX2AB0kWzi/1FrfoZTq1FoXpn6vgA6tdaFS6gngdq3166nfvUAyWZ0LZGqtV6Ru/zYQ1lr/6EjHra6u1itWrDiR5zeigsEgAJWVlYYjSQoGg3RFIgRNB5IycFacFE9BIOCYvxc48zUEw4snkd1O//iNxIv3geUDrcAXA8DTV4SvqxJfdyW+7nKUNbQk9OHX9IRTtlM0sYVNjy8lER3aY4wkJ76mwTnxQDKm2//lX1ZrrZeMxOMNtaV0lta6QSlVDqxUSm0d/EuttVZKjcisO6XUrcCtAOPGjRuJhxzTCgIBZjrsAuekeLoiEbbV15sO5ZCBM+OUmE44ni2TCOSVUTazAYDepkJ6mwuIRwaKIeLAweOKZ+A1bflD9Ex+g4zmaZxVNW148Z0gR76mARwSDwDBkU2RQ0pKWuuG1NdmpdTDwFKgSSlVpbVuTHXPNafu3gBMGPTfq1O3NfBBd9/A7S8f5lh3AHcA1NTU6Nra2uN4OqOrrq4OAKfEJPEcXV1dHdvq66kzHcggtamvdQZjGKw29bXuRB6kJxtWTT/hWCAZz8zKSmpra3m9+Ve83w43nfFPFJxr5gOqI1/TQJ1D4gGoTZ2jkXLMFR2UUjlKqbyB74GLgY3AY8DNqbvdDDya+v4x4NMq6TSgS2vdCDwLXKyUKlJKFaUe59kRfTZCiDEhkuhjQ+cTTMs7i4IM6TFxk6G0lCqAh1OLNfqAu7XWzyil3gPuU0p9FqgHrkvd/ymS5eA7SZaE3wKgtW5XSn0PeC91v+8OFD0IIcRgmzqfImr1saj4umPfWYwpx0xKWuvdwMmHub0NuOAwt2vgS0d4rDuBO48/TCGEW2iVYE3HQ1Rnn0xl1izT4Ywd2iKrp4n89l3kt+2mt3ACjVPPMx3VX5AVHYQQjhIrqSccb+GCyr82HUp605rMvhby2neT37aLvPY9+GMhAPqzS+gpmmw2viOQpCSEcBBNpGozJYHJTMo5xXQwaScj3El+6w7y23eT176bjEgPAJHMQjrLZ9FTPJXukinEMgsMR3pkkpSEEI6RV9mBldPJouJbZdPB41QU3MCUdffh0QmigTx6iqfQXTyFnpKpRLKKIE3OpyQlIYRjlM86gIpkMzPfeWMdTlbS8D41Gx6kt3Aie+ddTX9OedokoQ+TTf6EEI4woaiHvIouAsGhL0MkoGzfO0zZ8ADdJVPYvuQW+nMr0jYhgbSUhBAOcf6sAySiXjKazazekI4q9rzOxG1P0Vk2k50LbkJ70z+ZS0tJCGFcSU6Yk6tbad1VhUqMzH5NY5rWVO18kYnbnqK9ch47F35iTCQkkKQkhHCAc2c2YGlFy3ZZveFYNJrq7c9SvfN5WsctZNf869GesdPpNXaeiRAiLeVkxDi1ponV9eUE+gOmw3E0jSaY10LVnh00T1hK/ZwrQY2ttsXYejZCiLRz1rSDZPgsXtw23nQojpbQFi/lN9Ge00lw8lnUz7lqzCUkkKQkhDDI702wbPpBNjYU09SdYzocx4pbCb6z7j62ZndT1lvM/pmXpnWF3dFI950Qwpilk5vIzYzz4rZq06E4VtSK8621f+LV5i2c3lNKd1/xmE1IIC0lIYQhSmnOndnA3rY8drfkmw7HkfoTUf529e95tXkLfzf7Chb2FZsOadRJUhJCGFGeF6Isr58NB0qAsfvJ/0Q837iBd9p2cH7FPK6bdLrpcGwhSUkIYURzTzb72nNZNv0gGb6E6XAc6YLKk1hWNosXmzbyw82PYaFNhzTqJCkJIYzQWvHQ+1MpzI5y0ez9psNxpCxfBv+x6JN8YvJZ3L/vbZ4saiChxnYCl6QkhDBmb1s+7+0t57yZByjNDZsOx5G8ysPXZl3Gt+ZeQ0NGiD3F+wmExu6m3ZKUhBBGPb6uhrjl4eoFu02H4mhXTziFKzqqiXvjzH7r5+R21JsOaVRIUhJCGNXdn8GzmyYyb3w7eZVjtwUwEsZHs6lpm0jCn8XMd39NycE1pkMacZKUhBDGvbpjHM3dWYxfuBs9xsdMTlQgkcHm075Ab9Ekpqy/n/HbnwNtmQ5rxEhSEkIYl7A8PLxmCpn5YaKVW02H43iJjGy2L6mlpXoJ43a/zNR19+BJRE2HNSJkRQchhCNsCRbT1VAM1Rvoi7eR4ysxHZKjaY+PvXOvIZxTxoRtz5AR7mTnwk8Sy0zvicjSUhJCOEbD2imgLN5o/rXpUNKDUjTVLGPnok+S1dvMnLd/QVb3QdNRnRBJSkIIx4j2ZhFonM2W7pU0hjebDidtdJbPZsuptwIw561fMOvtX1K97RkKmrfijYYMR3d8pPtOCOEogYZ5+CYFebnpZ9ww6WeoMbg9w2gI549j82lfoKL+DfLa91Kx9w2q9rya/F1uOT2Fk+gtmkxP0SSiWUWOXdRVkpIQwlGU5eesslt5tvHf2dT1DPMKLzMdUtqIZeZzYOalAKhEjJyuA+R11JPbsZfi4HrKD7wHQDSQT8uEUzg47QKT4R6WJCUhhOPMzD+fDZ2P82bLb5iedzYBb67pkNKO9vrpLa6ht7gmdYNFVm8TuR315HXsJeF15i6/0i4WQjiOUopzK75Mf6KHt1t/ZzqcsUF5COdV0TLxNHaffANNNWeZjuiwJCkJIRypLHMa8wovY13HI7RF9poOR9hEkpIQwrFOL72FDE8OLzf9L1qP/W0bhCQlIYSDZfkKOL2slgOhNezsec10OMIGysmfPqqrq/WKFStMh3FIMBgEoLKy0nAkSRLP0QWDQboiEYKmAxlk4Mw4JSYnxlMQCPzZa0hj0XvSU2hfjLx1V6As++qzHPmaBoIOiQegMhjk9ttuW621XjISjyfVd2JEdUUibKt3xpL6znnb/rmCQICZDrmoDFx0hxOPVgms7E4S2e0ksjvxRHLwt07FEx9+VdfAB4kPv4ZyQhOYfv4GtmW/SdOmScN+/OM1cFYc95p2SDyjwdFJye/3U1tbazqMQ+rq6gAcE5MT49lWX0+d6UBSalNf6wzG8GG1JBOAk/5mcOzXUCTRS0tkFy39u2jp30lLZCftkXoskit6+1QmUd1PrGYD0/POYV7hRxiXNQ91nBM0j/gaaink0/tKmTfrAHfsqaAjlHlcjztctQNx2XK0Y6tNfa0zGMOH1Y7w4zk6KQkh7NcXb6O5f+eh5NPSv5OuWOOh32d7iynLnMbk3FMpD0yjLHMaBf4q2qL1bOx4gi3dK9na/TwlGZOZV/gRZhdcNCLzjB5bO4V549q5asFu6t6cc8KPJ5xJkpIQAoDO6EGePriC5v7th24r8I+jLHM6cwoupTwzmYByfMWH/f+lgRrOrfwKZ5Z/ju3dL7Oh8wleaf5f3mj5NTPyz+OkwsupyJx53K2nQ/GFA6zcMoGPnFTP9PIOdjQXDetxhLNJUhJCsLv3bZ49+O8oPCwr+zwVWTMpDUwZVgvH78libuGlzC28lKb+7WzseIJt3S+yuesZygLTOKnocmbmnU+GN/u4H/ulrdWcWtPEtYt288NnF2JpKSAea+QvKoSLaSz6q9fx+IF/osBfxY2Tf86iko8zPnv+iHS5VWTO4IKqv+Gz0+7lvIqvorF4Mfjf/HrX9bwY/G9a+ncd1+PFLQ+PrJlCVUGIs6Y1Hvs/iLQjLSUhXKo/0U1o1svECw8yu2A551d8FZ9ndNZDC3hzmF90JScVXkGwfwsbOp5gc9dzbOh8gsrMWZxUeAXT888Z0mNtPFjM1sZCLp1Xz+p9ZfRFMkYlZmGGtJSEcKHm/h38ae8XiecHydy9lIsq/27UEtJgSimqsuZw8bh/4HPT7uHs8i8SsUKsDP6Q3+y8gfCkVQTyjrX/j+KhNVPJ8FlcftLYLY12K2kpCeEym7ue48Xgf5PlLSBn88X4ekuHXXxwIjK9+SwsvpYFRdfQEF7Pho4n2F7xCrMvs/hycz5v7qpi3YFSEtZffnZu7snmle3jOHdmA2/uqmR/R57t8YvRIS0lIVwibkV5MfgTVjb+B1VZc7hx8s/x9ZaaDgulFNXZJ3Pp+H8k7/1rObhuMgVZUT59+jb+9Yp3uGL+Hkpzw3/x/57dNJHefj/XLtqFwrkr04jjIy0lIVygJ9bCUw3/SrB/K4uLr+eMss/gUV7TYf0FTzyT5l0TuGtrNdMrOjlzaiPnzjzABbMPsC1YyBu7qtjYUIylPUTiPp5YX8NNp25n8aRmVtVXmA5fjABJSkKMcQf61vLUwRXEdYTLxv0z0/PPNh3SMWkU25uK2N5URH5mhFOnNHHGlCCfOXMLXWE/b++u5O3dlby3t5wzpzVyxcl72NBQQiQul7R0J39BIcYorTXvt9/PGy2/pjCjmsvHf4figH3rxo2U7v4AKzdP5PktE5hd2cEZUxu5aM5+Lpqzny2NxWwNFrJ87n4unrOfx9fXmA5XnCBJSkKMQdFEiOeDP2JHz6tMy1vGRZV/P6zJqk6itWJzYzGbG4spzO7n9ClBTpvSxNxx7QCcM6OBd/ZU0NyT3s/T7SQpCTHGtEf28WTDv9ARPcBZZbeyqPjjRqrrRlNnKJOnN07m2U0TmTeunTOmNTKrspO549pp3iZJKZ1JUhJiDNnZ8xorG3+IV2VwzYT/YELOAtMhjSpLe1jfUMr6hlJyA1HCMbmkpTv5CwoxBlg6wZstd7K6/V4qMmfxkfH/TJ6/3HRYtuqVlR3GBElKQqSxSKKXrd3Ps6HjCdqiezmp8HLOLv8iPo9coEV6kqQkRJp6q+Uu3m+/j7iOUJ45g0vH/RMz8s81HZYQJ0SSkhBpqL5vNe+2/Z6puWextPQTlGdONx2SECNCkpIQaSZqhXkp+FMK/eO5ZNy3pKtOjCmy9p0QacTSCZ5q+B7dsUYuqPxrSUhizJGkJESa0FrzUtNPqe97l/Mqv0r1GC/3Fu4kSUmINLGq/R42dj7JkpIbOanwctPhCDEqJCkJkQa2db/Imy2/YUb+eZxReovpcIQYNZKUhHC4htB6Vjb+kPFZ87mo8u9RSt62YuySV7cQDtYe2cfjB75Dvr+Ky6v/VQobxJgnSUkIh+qLd/DogW/hUV6uqv43Mr2y5bc4cT6fprwiQVFxwnQohyXzlIRwoJgV5rED/0go3snHJv4nBRlVpkMSaczj0RQVW5SUJMgv0CgFra0eOtqdt/uwJCUhHMbSCZ4++H1a+ndy+fh/pSJrpumQRBpSSlNQYFFcalFYaOH1QiQCjQe9tLd5CIed2VEmSUkIB9Fa80rTz9nT+xbnVnyFKXmnmw5JpBVNXp6muCRBcYmFzwexWLJV1N7qpbdXAc7eW0uSkhAOsqb9AdZ3Psqi4o9zctFVpsMRaUGTla0pKUl2z2UEIJGAzg4Pba0eurs9aO3sRDTYkJOSUsoLrAIatNaXK6VqgHuAEmA18CmtdVQpFQB+BywG2oDrtdZ7U49xG/BZIAF8VWv97Eg+GSHS2Y7uV3mt5ZdMzzubs8r+ynQ4wsGU0gQyNUVFFiUlFlnZGsuC7i4P+/d76OzwYFnpk4gGO56W0teALUB+6ucfAD/WWt+jlPo/ksnmF6mvHVrraUqpG1L3u14pNQe4AZgLjAOeV0rN0Fo7swRECBsdDG3i2cZ/pyprLhdXfVPmIgm83mTiyQwkvwYCmj2ZmlhAszgDBna47+lW7N3jo6PdQzyenolosCElJaVUNfAR4N+Av1FKKeB84KbUXe4C/oVkUroq9T3AA8DPUve/CrhHax0B9iildgJLgbdG5JkIkaY6ogd4vOHb5PnKuWL8d2UukmtoMjI4lHAyU18DqSTk+9DVORYFHVFk93jY26+IRBQ93R6i0fRPRIMNtaX038A/AAMTJUqATq11PPXzAWB86vvxwH4ArXVcKdWVuv944O1Bjzn4/xxWLBajrq5uiCGOvmAwCOCYmJwYTyVQazqQlMrU11qTQXxIJcnzNPA3s3z99M17Bu2N4dt4Hve+9rCt8chr6OhO5DVkeTRxP8T8mniGJpaRbOVEMzXRAMQCGj2oQaws8EcUGRHwt3nI6FdkRBQZ/Qp/ROG11KF4qk/saY2oymPf5bgcMykppS4HmrXWq5VS547w8Q93vFuBWwHGjRs32odLe12RCNvq602HASRfnAWBADMrR/plOjwDF1ynxAMfxASgVZzQzJexMsLkbL4Qb0Qmx4LzX0MxLELeOCFPgj5vnD5PnJAnTt+g20KeOBGP9RePl2F5KEhkUBH3kx/yU5Dwk5/IID/hJzfhwzO4Mi6Q+pf/wU3BYJAuAIecHwAGvaZHwlBaSmcCVyqlLgMySZ6inwCFSilfqrVUDTSk7t8ATAAOKKV8QAHJgoeB2wcM/j+HaK3vAO4AqKmp0bW1tcN4WqNj4NOkU2Kqq6tjW309daYDSakl+eZ10vkB5/y94IOYPn3zp5L7IvW28ZHx/8y0OcuMxuOUc+TEeBJo/OdM4/GGVbT0d9MT7/+L+/mVl9JAHqWZRcwI5FESyKMskE9pZj6lgTzKAnmUBvIpyMg+4Xi2AXUOOT8AtSPcyj5mUtJa3wbcBpBqKf2d1voTSqn7gY+RrMC7GXg09V8eS/38Vur3L2qttVLqMeBupdR/kSx0mA68O6LPRog0sa7jEXb1vs7Z5V9gWp6ZhCSOrdEf5pX8Jtq37WBB0WQWF0+hNJCfTDKphFMayKfAn4VSY2tsx5QTmaf0DeAepdQKYA3wm9TtvwF+nypkaCdZcYfWepNS6j5gMxAHviSVd8KNtEqwuv0+qrMXsrD4o6bDEYfRGQ3xv9uf4dGS/eQkfPxg4Sc4t3yOJB4bHFdS0lq/DLyc+n43yeq5D9+nH/j4Ef7/v5Gs4BPCtWJluwnH27i46h9MhyI+RGvNkwfX8NOtT9ET7+fkviKW9pZwXsVc06G5hqzoIISNNBaRqs2UZ85gQvYi0+GIQXb1NPEfmx9lTcdeTiqcyDfnXsUbD8r8frtJUhLCRrGSfVhZPSwp/mvpCnKI/kSU3+x8iT/sfY0cX4Bvzb2GK6sX41Ee3jAdnAtJUhLCJlprIuM24QnnMy3vLNPhCOC15q38aPNjNPZ3cvn4RXxl5iUUZeSaDsvVJCkJYZP6vvewcjrI2nU6aqEsI2RSMNzJf255gleaN1OTW87/Lf0rFhXXmA7LFr6uTnK2byWRm0/v7Dmmw/kLkpSEsMmqtntQkWz8rZNNh+JacSvBPfVvcsfO59EavjRjOTdNPhO/Z2xfCj2hEDk7tpG7bQuZwUYAuueeJElJCLdqDG+mIbyezMbFKO283T7dYF1HPbdveoRdvU2cVTaLv5t9BeOyi0yHNWpUNErOrp3kbN9C1r56lNZES0ppP2MZfTNmEc/PP/aDGCBJSQgbvNf2JzI9eWQ0TzMdiut0RkP8bPszPHZgFRWZBfzHwk9yTvnssVlokkiQXb+HnG1byd6zC088Tjwvj65Fp9A7cxax0jLTER6TJCUhRllrZA97et/i1NJPs9X6y/XQxOiwtMWTDWv46ban6Y3388nJy/jctPPJ9gVMhzaytCbzYAM527aQs2M73kg/icxMemfPpXfmbCJV4z7Y5yINSFISYpStbrsXv8rk5KKr2Iq9q4C71eauA/x465Os66hnfuFEvjH3KqbnVZkOa+RoTUZrCznbt5K7bSu+3h4sn4/Q1Gn0zphNeOIk8KZnN7EkJSFGUXc0yLbuF1lQdA1Z3gLT4YxpcSvBS02buLf+TdZ37iPfn8U/zruWK8YvwjNGNk3ssBLUxyKM/+NdZLS3oT0ewhMn037mMkI1U9EZ6b8XlyQlIUbR6vb7UXhYWPwx06GMWR3RXh7Zv4oH971Nc6Sb6uxi/mbWR7h8/GJy/ZmmwzthbaEQT+7YxuPbt7K6txMAq7iI1nMvoG/aDKzsE1t53GkkKQkxSvriHWzqeppZBReS53f+AHO62d59kHvr3+LZxnVErTinlkznm3Ov5oyyGWnfMuqLRlm5exePbtvC6/v2ktCaGSWlnB/Ips+fwT0fv9F0iKNGkpIQo2Rdx8MkdIwlxdebDmVM2d/Xxvc3Pczq9t1kev1cMX4xH590GlNyK0yHdsJaQn3811tv8Oi2LfTH44zLy+Nzi5Zw1czZzCoto66ujldMBznKJCkJMQq01mzrfpFJOUsoCkw49n8QQ7KycT3f3/gwXo+Hr828lCurl5DnzzId1gmLWxZ3b1jHf771Bv3xGB+dPZerZ81hybjxeNKocm4kSFISYhS0RvbQHQtySsnY7WaxUyQR47+3PsWD+99hXsEEvr/gRiqzCk2HNSLebzzIP7/0PJtbWzhzwkT+9dwLmFJUbDosYyQpCTEKdve+CShqck83HUra29fXyrfW/ontPY18smYZX5x+MT5PepY7D9YWCvEfb77G/Zs3UpmTy88uvZxLp80Ym5N6j4MkJSFGwe7eN6nMnEWOz72feEfCswfX8e+bHsbv8fFfiz7NWeWzTIc0Il7fV8+Xn36cUCzGrYtP4SunnEbOGCjnHgmSlIQYYT2xZpr7t3Nm2edMh5K2+hMx/mvLEzxy4D1OLpzEigU3UJE5NuZ5bWtr5YtPPsa4/Hx+dunlTCsuMR2So0hSEmKE7e59C4ApuWcYjiQ9dXijfOatX7CzN0jtlHO4ddqFY6K7DqClr4/PPvYQ2Rl+7rzyGsblOXNRVJMkKQkxwnb1vEFhRjXFgYmmQ0k72zK7eSW/idxIFj9ZXMvpZTNMhzRiwrEYtz7xCB3hMPd87AZJSEcgSUmIERRJ9NIQWicrOByn/kSUH25+nBcKg4yLZnHH+V+hfIx01wFYWvO3zz3N+qYgv/jIlZxUnv5zqkaLJCUhRtDevnexSEjX3XHY3dvEt9b+iT29LSzuLeaU3pIxlZC01nzv1Zd4ZtcOvnXWOVw8dbrpkBwtvdfiEMJhdvW8Qba3iMqssVElNtqeOLCa2rd+Tme0j58uqeXU3lI8jJ2SaK01P3zzde5at4bPLFjMZxcuNh2S40lLSYgRErei1Pe9x/S8c/CosTEwP1rC8Sg/2PwoTx1cw+LiKXxv/nWUZuazhTdMhzaifr7qHf5v9bvcOG8+/7jsHNfPQRoKSUpCjJADobVErRBT8s40HYqjxaw4n3vnl+zsCfK5qefz2Wnn403zBVQP5841q/nPt97gmllz+N55F0pCGiJJSkKMkN29b+JXmUzMXmg6FEd79uA6dvQ0suLk67m46mTT4YyKP21cz4rXXubSaTP4wYXLnbd+naXxRGNYmc6bsCtJSYgRoLXF7t63mJizBJ9njG23PYIsbfG7Pa8yI6+Kiyrnmw5nVDyydTP/9OJKzp1cw4+XX4bP46xWYFZ9ExPueZFoaQF7P3uZ6XD+giQlIUZAU/92+uJtTJWuu6N6pWkze/ta+LeTbxiT3VlP79zO3618htOqJ/Lzy64gw0Fbknv7+ql67A1KX1tPPC+blvMWmA7psCQpCTECdve+icJDTe6ppkNxLK01dbtfoTq7mPMq5poOZ8S9tHc3X3/mSRZWVnHH5VeR6fObDinJ0hS/vYlxD7+Or6+flnMX0njF6VhZzmzRS1ISYgTs6nmD8dknkemVWfpH8l77LrZ0N3Db3KvHzLJBA97cv48vPPkYM0tK+c2V1zpmcdWs/c1U3/Miubsb6Z06jl03nE+42tm7ICuttekYjqi6ulqvWLHCdBiHBINBACorKw1HkhQMBumKRAiaDiSlEigIBBx1fmD0/16JQA+9Cx8lc+8SAsGjz09y4msI7InnsaIDtPsifKqlBu8Rpkim4/nZH4/xh1A3hR4vN2fnkz2KY0jBYJAuIHiM85OIWzQe6KetKYLPp6iamEVRacaodJlWBoPcftttq7XWS0bi8aSldJy6IhG21debDgP4IAnMdNgb2G3ixfsB8HdUG47EuZp8/RwIhDi9u/SICSkdHUzEuTvUQ57Hw6dGOSENhdaajtYojfvCxOOakooAVdWZeH3pc84dnZT8fj+1tbWmwzikrq6ObfX11JkOJKUWmFlZ6ZhzVFdXB+C6eO6v/zpRawqfuPHLjolpqOyK5xtr/kheWybfv/ZL5PiOPJaRTudnW1srNz54L6V5edz7settWWC1rq6ObUDdYeLJOtCS6qrrpHdKFQduOJ/whPJRj6k2dY5GiqOTkhBOF4p30hjezCklnzAdimPt7W3m5abN3DL13KMmpHSyu6OdTz18PwGvjz9c83GjK357whGqHn+TspfXEc/JpP5TF9N+2hzwpGd1oyQlIU7Ant630VhMzZMFWI/kd3teJcPj4/pJY2Nr+APdXXzq4fvRWvP7az/GpMJCY7Fk7z7IlF8+jq8nROuy+TReeSaJnExj8YwESUpCnIDdvW+Q6yujLDDNdCiOFAx38vTBtXxs4qkUZeSaDueENfX28smH7qcvFuNP115ndNfYjJZOpvziMaysDLZ98SbCk8bGdhjpM/olhMPErH729b3P1NwzxuRE0JFw997XAfjE5GWGIzlxbaEQn3r4ftrCIequ+iizy0Z/vOZIvH39TP3fR1Bas+vL14yZhASSlIQYtn19q4nrCFOk6+6wOqN9PHLgPS6pWkBlVqHpcE5IV38/n37kAQ70dPObK69lQWWVsVgsS1NzxxNktHax+/NXECkvMhbLaJDuOyGGaVfvm2R4chifPTYXFT1R99a/SSQR59NTzjYdygmJaM0tjz7ErvZ2fnXF1Swdb670X2tNw94QeS2d7K1dTt/0sTcNQZKSEMNg6QR7et+iJvdUvEreRh/WF49w/763OadiNjW55rq5TlRMa/4U6qahN8HPL7uSZZMmG41ndWOE9pYojZedSsepc4zGMlqk+06IYWgMb6I/0c3UXFmA9XAe2f8u3bEwN9ecYzqUYYvE49wb6mFfIs5/Lb+Mi6aaLWZ5/v3tvLk/TGGJn+DlY6OS8XDkI54Qw7Cr5028ys+knFNMh+I4USvO3XvfYEnxFOYWTjAdzrDEEgm++swT7E7EuCIzhytmmN3efsOeRr5d9wyVuV7KpuTAGC6skZaSEMdJa83u3jeozl5AhjfbdDiO81TDGloi3dROOdd0KMOSsCz+duXTrNy9i0sys1mYYXbez8G2Lv76F49Rmp/D5TNy8aTppNihkqQkxHFqi+6lK9YoXXeHkdAWv9/zKrPzx3NKyVTT4Rw3S2u+9eJKnti+jX84YxlLM7KMxtMTjvC1/32UaDzBT750Ndn+sX/JHvvPUIgRtrvnTQBqck8zHInzvBjcyP5QGzdPOSft5m5prfnuKy9y/+aNfGXpafy/JUuNxhNPWHzjV09S39TBD2+9nClV5ibq2kmSkhDHQWvNjp5XqcicRa6/1HQ4jqK15q7drzApp5RzK9KrMkxrzX+8+Rq/W7+Wzy5czNdPNTv3TGvNf9z7Em9vqee2G8/n1FkTjcZjJ0lKQhyHg+ENtEZ2MbfgEtOhOM7brTvY3tPIp2vOwaPS69Lys/fe5per3+Omk07mW2eZb+X98cU1PPDaem6+eAnXnHWS0Vjsll6vHCEMW9P+IJmePGYVXGg6FMe5a/crlGcWcMm49JpM/Ov3V/Hjt9/k2llz+O65FxhPSC+v28WPH3yF8xdM4ytXnWU0FhMkKQkxRF3Rg+zqfZOTiq7A70nvlZhH2vqOfbzfsYdPTD4Lvyd9Zprcu3E933/9FS6bNoPbL1yOx3BC2ra/mW/d+RRzJlbwvVsuGfOVdocjSUmIIdrT+w6gmVf4EdOhOM6D+98h35/F1dXpM2+rLRTiu6++xFkTJvFfyy/DZ3jXWID/fug1cjIz+PEXriIrw286HCPM/xWESBNtkT1keQvI86XvsjmjIW4leLt1O6eVTifLl2E6nCH71ZpVRBIJvnPOeWR4vabDYfuBFt7Zuo+bzl9IaUGO6XCMkaQkxBC1RfdSEphsfMzBaVa176Yj2scFlekzIN8WCvH7dWu4fMZMphrcE2mwP774PpkZPq49a77pUIySpCTEEGitaYvspSRjsulQHOfZg+vI9WVyRukM06EM2a/XrKI/HucrpzhjrllLVy9Pv7uVK0+fS0Ga7xx7oiQpCTEEPfFmolaIkkCN6VAcpT8R4+WmTZxfOY+ANz3GQNpCIX6/fi1XzJzlmFbS/a+sJ2FZ3HT+QtOhGCdJSYghaIvsBaAkMNloHE7zestW+hIRllelTxn4r9esIhyL8WWHtJL6o3EeeHUd58yfysQxtmHfcEhSEmIIJCkd3rMH11EayGNRcXq0IHsiEX6/fi2Xz5jFNIe0kp58ZzOdff184vxFpkNxBElKQgxBW2QPub4yAt5c06E4RncszJst27iocj7eNFnB4eX6PYRiMT41f4HpUIDk1uZ/fHENsyaUs2j6eNPhOEJ6vJKEMKwtsldaSR/yUtMmYjrB8jRaweG5XTsozc5mUdU406EA8ObmvewNtvPJCxdJVWeKJCUhjsHSCdqj9ZKUPuTZg2uZkF3C7Pz0+IQficd5ee8eLpoyzfjKDQP+8PxqygtzuWhR+lQujjZJSkIcQ1eskYSOSeXdIC393axu38PyqpPT5hP+G/v30ReLcfEUs9uaD9h+oIV3t+3nhnMX4PeZn7zrFJKUhDiG1v7dAJRkTDIciXOsDK5Ho1k+boHpUIbsuV07yM3I4PQJztgG4o8vDEyWtX/Ssb8zRM7uFtuPOxSSlIQ4hoPhjfhUgNLMKaZDcYxnDq5jdv54JuWkx55SCcvihT27OG/yFEcsKdTS1cvT723lqjPmkW9gsuzku97g1NpfkdHea/uxj0WSkhDH0BBaR2XWHLwqPSaHjrb63ha2djekVYHD6saDtIXDLJ/qjK67gcmyN55n/2TZQHM31Q+vIrj8JKLFzqsmPWZSUkplKqXeVUqtU0ptUkr9a+r2GqXUO0qpnUqpe5VSGanbA6mfd6Z+P3nQY92Wun2bUmr5qD0rIUZIJNFLS2Q31dnuXo9ssGcb16FQXFSZPufk2V07yPB6OWeS+XHBcDTGA6+u49z5U5lYXmj78WvqXkdZmt2fWWb7sYdiKC2lCHC+1vpkYAFwiVLqNOAHwI+11tOADuCzqft/FuhI3f7j1P1QSs0BbgDmApcAP1dKmW9HC3EUDeGNgGZ8VvosNjqatNY827iOxcU1lGXmmw5nSLTWPLdrJ2dNmEROhvlVzJ98Z0tysuwF9k+WzTrQzrjH19Jw9SL6qwptP/5QHDMp6aSBjkd/6p8GzgceSN1+F3B16vurUj+T+v0FKlmecxVwj9Y6orXeA+wElo7EkxBitBwMrcer/FRmzTYdiiNs6W5gf6gtrQocNrc009DTzcUO6LqzLM3dL7zPnEkVLJxmfyn9lN+8ivZ52HOzc3e0VVrrY98p2aJZDUwD/hf4IfB2qjWEUmoC8LTWep5SaiNwidb6QOp3u4BTgX9J/Z8/pG7/Ter/PPDh4w2orq7WK1asOIGnN7KCwSBdkQhB04GkVAIFgQCVlZWmQwGS5wcYU/H0znsaLA+5m0emtzndz9Ebec1syO6itnkKmXrkOzpG4/y81B/i9WiYv8ktIuc4N/Ib6Xj2dMR4fHsvy6fmMLP0+FttwWCQLiA4jHhi7TGaHmgi9+RcCk8tPO7/fySVwSC333bbaq31kpF4vCHtW6y1TgALlFKFwMPArJE4+OEopW4FbgWorKrilfr60TrUcRtIAjMddkERo0N7YiRy2gk0zDUdiiNYaHZk9jApkj0qCWm0bItHmej1HXdCGg1rgv3kZiimFdtfNNO1qgvlV+SdnGf7sY/HkJLSAK11p1LqJeB0oFAp5dNax4FqoCF1twZgAnBAKeUDCoC2QbcPGPx/Bh/jDuAOgNJx43TdcT2d0VULzKyspLa21nAkSXV1dQASzxGcaDz1fat4ZP+9XHr6J5mUMyIfAtP6HL3btpPQezu4dfGVo7ah30ifnz2dHXz3d3fy7bPPo3bB8Y/hjGQ82/Y389N3/shXrzmL2ouHt218XV0d24C644wnb8tBTr3jTnZ97mxWfubsYR37SGpT52ikDKX6rizVQkIplQVcBGwBXgI+lrrbzcCjqe8fS/1M6vcv6mQf4WPADanqvBpgOvDuCD0PIUZcQ2g9Cg9VWdJSguSK4DneAGeWjVpHyYh7btcOAEes4nD3i2vICvi59kz7i2am3vEy0YIs9l1/qu3HPl5DaSlVAXelxpU8wH1a6yeUUpuBe5RSK4A1wG9S9/8N8Hul1E6gnWTFHVrrTUqp+4DNQBz4UqpbUAhHagitpyJzBhmeLNOhGBdJxHipaRPnVswlM0028wN4btdO5pWVMz7fbKXgwGTZjy2bb/tk2cI19ZS+s5vtX76ARE7A1mMPxzGTktZ6PfAXM7y01rs5TPWc1rof+PgRHuvfgH87/jCFsFfcitDUv40FRdeYDsUR3mzdTm+8P60mzDb19rIm2MjfnHam6VD4w/PvY1maG+3eWVZrpt7xMv2leRz46Mh0QY828yN/QjhQMLyFhI4xXibNAsmuu+KMHJYUp89SSyt37wQwXgre3h3i/lfXcckpM5lQVmjrsUve2U3Ruv3sqT0LK5AeLVxJSkIcxoHwekAxLmue6VCM643383rLVi6snI/Pkz5Vd8/t2snkwiKmG95h9nfPryIaS/C5S20ez9Gaqb98iXBVAQevWGDvsU+AJCUhDqMhtIGywFTZaRZ4uWkzUSueVl13Xf39vN2wn4unTjO6tUZ7d4j7XlnH8lNmMrmy2NZjl7+8lfxtQXZ/9my0P30+TEhSEuJDEjpGMLxZuu5Snj24lvFZxcwrmHDsOzvES3t3E7cslhuuuvv9C6vNtJISFlN+9Qq9k0tpXJ5eS2RJUhLiQ5rC24nriCQloDXSw3ttu1g+Ln028wN4eOtmKnNyObmyylgMHT0h7n15LctPmUmNza2kyuc2kru3ld2fOwe86XWZT69ohbBBQ3gdgCzCCrzQuAELzfKq9Om629TcxGv76vnk/AVGtz3/3fOricTitreSVCzBlN+8SveMSprPTZ85ZQMkKQnxIQ2h9ZRkTCbLV2A6FOOebVzHjLwqanLLTYcyZL98/z1y/Rl8cr65RDrQSrpkySzbW0njnlhL9sFOdt16LnjSp3U7QJKSEINYOsHB8CbGZUsr6UCojY1d+9OqwKG+s5Ondmznpvknkx+wf0fXAb9PtZI+e6m9GyF4IjFqfvs6nfOraTt9qq3HHimSlIQYpC1ST8wKM16SEs82JrsxL06jzfx+vWYVPuXhM8NY526k9IYj3PvKOpYvnsmUKnvL0asfWk1maw87P38epNEY4GCSlIQYpD/RBUCOz+zcFtO01jx7cB0LiyZTkVVoOpwhaQn1cf/mjVw7ew7lOeZK+V/fuIdwJMZ159rbwvT2RZj8+zdpO6WGzoWTbD32SJKkJMQgUSsEQIYn23AkZu3oaWRvX0taFTjct2kjsUSCzy0yu5zOS2t3UpKfzfyacbYed+J975LRGWLX58+z9bgjTZKSEIPErDAAfpcvwvrMwXX4lJfzR2mLitHwWv1e5pZXMKXI3sKCwSKxOG9s2ss586fisbHIwNcdZtLdb9N89ky659ibDEeaJCUhBpGWElja4rnGdZxeOp3CjPQ4D6FYjDXBg5wxYaLRON7Zuo9QJMZ5C+ydtDv5D2/hDUXY9Vfn2Hrc0SBJSYhBpKUEW7sP0hzp5qKq9ClwWBtsJGZZnD7e7KoTb27aS3bAz9KZ9sWh4gmqH1pF0wVz6JuaPqX7RyJJSYhB4joGgE9lGI7EnA0d+wBYVFxjOJKh29HeCsDssjKzcTS0MqO6DL/PvrXm8rY24gtFaT53tm3HHE2SlIQYJKGjgEK5+K2xoWsf5ZkFlGemz+Thne3t5GUEKMvOMRaD1ppdB1uZOs7eys2iNckPER0LzHZdjhT3vvOEOAxLx/Eqf1qt8zbSNnUeYF5Btekwjsuu9namFRcb/bu1dvXRHYowdVyprcctWltP7+RSYsXmEvJIkqQkxCAJHcOr0mMztNHQEe2lIdzO3ML0WRFca82O9jamFpurugPYcTDZhTjNxpaSilsUrttP58Kx0UoCSUpC/Jm4jro6KW3sPACQVttUBHt7aQuHmFNqdpB/18E2AFtbSrk7gvhCUToWpO9k2Q+TpCTEIH2xNlev5rCpaz9e5WF2wXjToQzZ+uYgAPMrKo3GsfNgK6X52RTl2le5eWg8SVpKQoxN3fEg+f70L6sdro2d+5maW0GmN32qD9c3BfF5PMwxXHm3s6HVyHhS34RioqV5th53NElSEiJFa01PrJl8v9lP3KZY2mJT137mpdF4EiST0oziEjJ95rpdY/EEOw+2MaPaxsSYSI0njZGquwGSlIRIiVg9RK0Qef4K06EYUd/XSl88klbjSVprNjQ3Gd1hFmBHQwuxeIK5k+177eTuasbf009HGi++ejiSlIRI6Y41AZDv0qS0sXM/QFpV3u3t6qQ7EmF+udm/2aa9ydfO3En2tbLH4ngSSFIS4pAPkpI7u+82du0n15fJpBx7x0VOxPomZxQ5bKpvojAnk3El+bYds2htPaFxhUQq0meS81BIUhIipTuWvMC5uaU0t6Aaj0qfy8L6piCZPh/TS8wm0s31QeZOrrRv8q6lKVyzb8yNJ4EkJSEO6Yk14/dkEfCMnUqmoQrHo+zqCaZV1x0kk9LcsnJ8HnOXslB/lN2N7cyZZN+HmZw9LWR0h8fceBJIUhLikO5YkHx/hSuXGNrS3YCFTqsih7hlsamlmZMMd91t2d+MpTVzJ9s5nlQPjL3xJJCkJMQh3bEm8n0uHU86VOSQPmve7WhrpT8eZ3652b/Z5vqBIgf7WkpFa/bRX5FPf1Whbce0iyQlIVJ6Yk2uLQff1LWf8VnFFGXkmg5lyD4ocjBdeReksiiPknx7FkTVWlO4dl9yVfAx2KqXpCQEEEn0ErF6XV3kkHaTZpubyMsIMLmwyGgcm+ubbJ2fFO+ME+joG5PjSQBKa206hiOqGjdOX/L5z5sO45BKoCAQoLLSGV08wWDyk6LEc3jHE08iu4Pe+U+SvX0Z/vbRe7M78RyF/BbPzO3lrO4y5ofMXuCP5/zc0dtJlvLwqZzRK8M+VjzhmMWv3u/ijAlZLBmXOWpxDI5nx16LfestKq6vwF9gfvHgymCQ22+7bbXWeslIPJ5vJB5ktHg9Hs6Z5JxPAwMvUDH2WIFeAFRkbOxJczw6shMAVMRG/6I6UuJa02wlOD3D7Bp9zX2pc5dj306zfR0aT5YHX76jL9/D5uhn5ff7qa2tNR3GIXV1dQCOiUniObrjiWdt+8O80vwKN13zObJ9hY6IyQ51dXXsy21F0c/Xb/gcAa/ZT95DPT9rGg9i3f8nPnnRRSyfOt1YPL966m3Y9hZfv/XT5GUFRi2OwfF0h5ponj+B5265adSPNxS1qXM0UmRMSQiSlXc+FSDLO7Zmxw9FlzdGRWaB8YR0PNYNFDkYrrzbtLeJyRVFtiQkAMvSxDpi9E4bu2OfkpSEwN1zlLq9McZnm9219XitbwpSlp1DZa65akGt9aGVHOzS3RkHC3qnjt3tVSQpCcFAObgzig/slpZJqTnI/AqzHyKaO3tp7Q7ZOj+psz0GQO80SUpCjGnd8SZXbu4XVRZhb4LqrPRJSt2RCLs7OhyxCCvAHBtXBu9oj4GCvknps2ju8ZKkJFwvmgjRn+h25erg3d7kJ+90ailtbE4mA/PjSckdb2dOsG9jv862GP4iP9pvX7Wf3SQpCdfrjrt3H6UubxSA6uwSw5EM3cBKDic5YCWHaeNLCfjtK2LubI/hL06fgpThkKQkXK8ntY+SG5cYGmgpVadRS2l9U5AJ+QUUZ2Ubi8GyNJv3Ndu6kkN3d5hQXwJ/iSQlIcY0N2/u1+2LEbA85PmzTIcyZOubg5xseDxpf0snveGIrTvN7tmZfJ1KS0mIMa47FsSr/GR7C02HYrsub4z8RPpc5FpCfRzs6XFAkUOyC9HOPZR272wGkJaSEGPdwOrgKo12XB0p3d4YBXGzS/UcD8dsf763iUy/jylV9o3F7d7VTCDTgydrbL9Ox/azE2IIumNNrixyiFsJetKspbS+KYhHKeaWmS3f31QfZNbEcnxe+y6he3Y2U1jsH/MTvCUpCddLbu7nvqTU1N+FVlCQVkmpiWlFxeQYXIg1lkiwbX+zrSs5JBIWe3c3UzTGx5NAkpJwuZjVTzjR6crVHA6E2gDIj6fHhU5rzYamoPHtz3cfbCMSS9i6kkNjQweRSJzCMT6eBJKUhMv1xNw7R+lAqB0gbbrvGnq6ae8Pmx9PMrCSw0CRQ6G0lIQY27pdnJQaQu14tSLXcvQONoc4psihPkh+doAJZfatKL97VxMej6KwSJKSEGNaTyz5CdSNE2cbwu3kJfwo0mPgfENzE36Ph1klZtd921LfzJxJ9i4Gu3dXC+MnFOP1pcff6kRIUhKuFkp0AJDtM7sNuAmN4Q7yE+nRSgJo7O2hMjePgM9szE0dPYwvtXffrdaWHioq3bHXlyQl4WrheCcBTy5elT4X55HSGO4kL03GkwBaQyFKs80tLQQQT1h09oUpzrM3jo72XoqKc2w9pimSlISrhRNdrtxtNhSP0BULpVdS6uujLNvshbmzN4zWUJpvXxxaazo6+igskqQkxJgXTnSR5XNfUmoMdwKQl0bdd62hECWGW0pt3X0AFOfbF0eoL0IsmpCWkhBu4NaWUrC/E4DcNGkpxRIJOvrDxltKbT0hAEpsTEod7clEKC0lIVwgnOgi04VJqTGcLPBIl5ZSeziMBuNjSm1dyQRRYmP33UBSkpaSEGOc1ppw3KUtpXAnPuUlJ03mKLWGkhfmshyHtJRsLHT4ICnl2nZMkyQpCdeKWiEs4q7csiIY7qQyqyBt5ii1hJLJwHRLqb07RGaGj+xM+9be6+iQlpIQrhBOdAGQ6cZCh/5OKjMLTYcxZIdaSqbHlLpDtlbeQbIcXCkoKDCbkO0iSUm4VjjRCeDK7rvGcAdVWekzYbgllZRKjSelPlsr7wA62/soKMjG63PH5dodz1KIwwgnugH3JaWoFac10kNlVqHpUIasNRQix+8n22+2WrC9J2RrkQMkx5QKXdJ1B5KUhIuF452A+5JSUzjZbVmVVt13IeOtJIDW7pCtRQ6QTEpuGU8CSUrCxfpTY0puS0qN/cly8HRqKQV7e6jIMVt9Fksk6OoL295919EhSUkIVwgnuvAqP35PlulQbBVMreaQTmNKjb09VOXlGY2hs2dgiSH7x5SKitxRDg5DSEpKqQlKqZeUUpuVUpuUUl9L3V6slFqplNqR+lqUul0ppX6qlNqplFqvlFo06LFuTt1/h1Lq5tF7WkIcW3I1h0JbtyBwgmC4Ew+K8sx806EMiaU1Tb29VOaavTAPzFEqtnFMKRyO0t8fo7DIHZV3MLSWUhz4W631HOA04EtKqTnAN4EXtNbTgRdSPwNcCkxP/bsV+AUkkxjwHeBUYCnwnYFEJoQJyYmz6XFhHklN/V0UB3Lxe9Jj4mxbKETMsqjKNdtSaunsBaCswL6k1DmwxJB0331Aa92otX4/9X0PsAUYD1wF3JW6213A1anvrwJ+p5PeBgqVUlXAcmCl1rpda90BrAQuGcknI8TxCCc6yXLhxNmOaB/FGenTHdTY2wPAOMPddy2pJYbKCuw7d25bzQFAaa2HfmelJgOvAvOAfVrrwtTtCujQWhcqpZ4Abtdav5763QvAN4BzgUyt9YrU7d8GwlrrHx3peNXV1XrFihXDeFqjIxhMbsdcWWl2O+YBEs/RHSuengWP4O0tJXvnWY6JyQ4PFu/Drz1c2VHtiHgGO1w8W2MR7gv38lc5BVR57W3dDY7nnQNh3mno50unFOL12NPlu39vmFefa+OSa8opKcsgGAzSBQQd8vcCqAwGuf2221ZrrZeMxOMNudBBKZULPAh8XWvdPfh3OpnZhp7djn6cW5VSq5RSq44nYQpxvCxfBBULmA7DdmFPgizLazqMIeuyLADyldm6rL6YRZZP2ZaQAPpDCQCystxTkzakjx1KKT/JhPRHrfVDqZublFJVWuvGVPdcc+r2BmDCoP9enbqtgWRrafDtL3/4WFrrO4A7AGpqanRtbe1Qn8uoq6urA8ApMUk8R3e0eBI6xs+2/YFFJ53Oqed90hEx2RbDyn9h8dT51M7+iCPiGexw8QRff4WMdWv44i232F6UMjieNT9/hIkZvdTW2vd6+f2dr/Lu66/yV5+/Bb/fS11dHduAOof8vQBqU+dopAyl+k4BvwG2aK3/a9CvHgMGKuhuBh4ddPunU1V4pwFdWutG4FngYqVUUarA4eLUbULYLhx35xyl/kSMUCJKcUb6DJwf7O2hKjfPeJVkS2cfZYX2ju10tPeRl5eJ358+LdsTNZSW0pnAp4ANSqm1qdu+BdwO3KeU+ixQD1yX+t1TwGXATiAE3AKgtW5XSn0PeC91v+9qrdtH4kkIcbzcusRQZzQ1cB5In4HzYG+v8co7gJauXmZPKrf1mJ0uW2IIhpCUUgULR/qIcsFh7q+BLx3hse4E7jyeAIUYDYcWY3XZCuHt0WRZc1EatZQae3o4Zdx4ozHEEgnae0K2Vt4BtLf3umo1B5AVHYRLhQ8tMVRoNhCbHWoppUlSsrSmqa/X+GoO7d0htLZ3jhJAW0sPpWXmW4l2kqQkXOmDde/cNXm2Pc2SUmuoj7hlUWl64uzAHCUbx5S01rS19lJa5q7XqCQl4UpxKwKAz5NpOBJ7fdBSSo8xpcae1MRZ00nJwGoO3V1hYrEEJaXSUhJizLNIzn3xuOwt0B7tJcPjI9tr33beJ6KxN5kMTHffmVjNobUlmZBLytLjA8RIcdc7UogUrZNJSRmekGm3jkgfRRk5xsurh2pgiSHTi7G2dPbi9SiK8uxbUb6tNfncZUxJCBfQqZaSctlboDPWlzbjSZDcRyng9VGUaXZ7kZauPkryc/B67Hu9DLSUSqX7Toixz9LJ5VvclpTaI71pM54EyTGlqjwHTJzt6jVSeQdQLElJiLFPY6HwGL/Y2a0z2kdhGrWUGnt7qDLcdQfJ7rtym1dzaG3pobAw21WrOYAkJeFSlk7gUe56s0OyJDydlhhqTC0xZFpLVx+lNk+cbWvtocRl40kgSUm41EBLyU3C8SgRK5Y23XcJy0rtOGv2why3NF19/ZQV2pvMW104cRYkKQmX0tpyXeVdui0x1BoKkdDaeDl4XzRZFGP3EkNuXM0BJCkJl3JjSyndlhgaKAc33X3XF0vu62ZnSykWS9DZGZLuOyHcwtIJF06cTdekZLa70URLqT01R8ltqzmAJCXhUhoL5bJCh55YGIA8v9k5P0PVHk7GW5KdbTSO/niypVSYY9+SVD3d/QDkF6TH32okSVISrqTRruu+CyWS6/3l+NJjC/juSPLCXBAwuz5hzEompexM+5Zm6u+PApCVlR7LQY0kd70rhUhJFjq4a45SXzzdklKEDK+XgG8oe5GOnmgimZSyMvy2HTMcjiWPKUlJCHdwZUspHsGDIuCx7+J6IrojEfIyzCfQaEKTHfDj8dj3IaY/nGwpZWalx99qJLnrXSlESjIpua+llOMLpE0LsScSIT9gPinFUknJTv390lISwl1cOE8plIiQnSZddwA9UWckpWhC2zqeBBAOSUtJCFdxY0spFE+vpNTtlJaSBdkBe5PSQEspU1pKQrhDcusKd738++JRcrzmL/JD5ZSkFE1ocuxuKaXGlAI2dxs6gbvelUKkaK3TZmxlpPTF+9Om8g4GkpL57epjCU2W3WNK4RiZWfYWVziFJCXhSm5cZijdxpSS1Xfmu69MtZQybT6mU7jrXSlEijvHlKJkp0n3XSQeJ5KIk+eElpKlyc60u6UUdWWRA0hSEi5l6bjr9lNKtpTS49N3TzQ5puKEMaVYQpNjd6FDOObKcnCQpCRcKmb14/eY/xRuF611ap5SejzngSWGTCclrTUxC9vHlJLdd9JSEsI1YlYYv3LPYpcRK05CW2lT6NATSS6JZDopRRPJr3aPKfX3x8jKlpaSEK4R0/34Pe5JSqHUunfZ3vS40HWnkpLpZYZMLMYKqTElaSkJ4R4xK+yq7ruBFcLTpfqu2yEtpVhqMVYZU7KPJCXhSm4bU/pghfD0eM7dUWckpYEVwu2uvgu7uPpOaa1Nx3BE1dXVesWKFabDOCQYDAJQWVlpOJIkiefojhZP1yl/IqNpBln7FjsmptF00B/ikZIDXNleTXX0g03znPo3211UwPOREN/MKybD4CTn9XubeLkpg2tm5TKhwL4kce9vG5g2K4fFpxf+2e3BYJAuIOiQvxdAZTDI7bfdtlprvWQkHk9aSsJ1NBZ4EyjL7D49dop6klt6+3V6zM3q1xoPYLqtELOS5yvDa99501oTj2l8/vT4W400R78r/X4/tbW1psM4pK6uDsAxMUk8R3ekeKKJEL/YcTenLDydxRde54iYRttzjet4at29XH/VR5mcW248niMZiCc0qZrNO7Zzyy23GI1n649+CS0hbvj4R5lUUWTLMfv7Y9z9qx9w2mmncP0nz/iz39XV1bENqHPI3wugNvU3GynSUhKuE9NhAFdV3w2MKaVLoUNPNEquI5YYSn61c0zpg20rzD9/EyQpCdeJWcmJmW4qdAilWVLqjUYckZRMVN/19yeTUpZLCx0kKQnXiVkubimlyTylvmiMHAckpYHqOztXdOgPp/ZSkgVZhXCHqAuTUigRJcubgSdNdttNtpTMt+pilibDi63bnAzspeTWkvD0eIUKMYLiA913yj3ddz2xMHn+9EnCvbEYORnmL8rRhMZv855GA7vOyuRZIVzCjS2lrliI/HRKStGI8SWGIJmU7CwHB+gPSUtJCFeJ6WRLKcNFhQ49sTD5vvRJSn3RKDl+8xdlE0lpoPsuK9t8UjZBkpJwnYFCB5+LWkrdsTD5GenxfC2tCcfjjih0iBlISqGQVN8J4SoDJeFuail1x8LkpUlLKUqy4s0JhQ7RBPjt7r5LtZSypaUkhDskW0oKn3LPm747FqbAn33sOzpAJLUeZ65bu+9SLaWAbF0hhDvErH78KoBKk/LoExVJxIhYsbSpvhtISo7ovrPMdN9lZvnx2Fz15xTueFcKMUhMh11VedcTS46h5fnTo7syqp3UfWegJDwcdW05OEhSEi6U3ODPPUmpK5WU0qb7joGWktnuq2gsjqUhw2vvccPhqGu3QgdJSsKF3LbBX098oKWUHonYKS2lvtQadCa676SlJISLuK2l1J1qKaXL5NkPCh3MXphNJaX+kLSUhHCVZFJyT0spXZOS6e67UCS53I/dJeHhcExaSkK4SUxHXLXuXV8sOS8r15cezzmW+pptuCQ8PJCUbC50iPTHXFsODpKUhAsldAyvcs+bvt9KbYWQJttWxFMtpYDX7MbY0Xhyhz+vzVfJWDyB329zdYWDSFISrpOw3JWUIolkUgp4zF7khyqOJsPjtXW7iMOJJVJJyeY44rEEPp97L83ufebCtVzXUkrECHj8xi/yQ5UAMnzmWwqxmJmWUjxu4XPA8zdFkpJwnYSO4XFTUrKiZHrT5/nGtSbgNX9RPtR9Z3dLKZ7AJ913QriH5bKWUiQRT6+kBGQYHk8Cg2NKsYS0lIRwk4SO4fWkz0X6RPUn0qullNCaDAe0lEyNKSWk0EEI99DawiLhqpbSwJhSuogDAZ/5lpKpMaVYXAodhHCNhE5WorkpKUWsdOu+c0ZLaaD7zs5pSpalsRJaxpSEcAs3JqVk9116zFEC5xQ6xAwUOsRTx5QxJSFcIu7KpBRLq5ZSAmcUOhwaU7LxKhmPDSQl916a3fvMhSsldHKRTa9Kn5bDiYpYsbSZOAuplpIDWgrRWAIFeOxsKSUsAPz+9Pl7jTRJSsJVPui+c8+bPpqIk5FGhQ7JlpIDklI8Yf/E2YGWkt+9l2b3PnPhSpaOA+7qvovqBBke8xf5oYqjja97BxBPJIwsMQQypnRUSqk7lVLNSqmNg24rVkqtVErtSH0tSt2ulFI/VUrtVEqtV0otGvR/bk7df4dS6ubReTpCHJ0bCx2iiZgjxmiGKq6d01KyeYHwQ8UVkpSOrg645EO3fRN4QWs9HXgh9TPApcD01L9bgV9AMokB3wFOBZYC3xlIZELYyZVJyUrgT6cxJbQj5ikZ6b6LD4wpSVI6Iq31q0D7h26+Crgr9f1dwNWDbv+dTnobKFRKVQHLgZVa63atdQewkr9MdEKMuoTlrkIHrTUxK05GGiWlhENaSrGYdN+ZoHRq75Kj3kmpycATWut5qZ87tdaFqe8V0KG1LlRKPQHcrrV+PfW7F4BvAOcCmVrrFanbvw2EtdY/Otpxq6ur9YoVK4b51EZeMBgEoLKy0nAkSRLP0R0unljBQUKzXyRn43J8vWWOiGk0JdD8snIHS3tKWNJXYjyeYwkGg/w6y8dpgSwuzMwxGsuT23tp6Y1wybiYbeenrTnKM480c+4lJYyf+Jc7BQeDQbqAoEP+XgCVwSC333bbaq31kpF4vBNunOpkVjt2ZhsipdStSqlVSqlVQ0mYQhwXT/KTqNLuqPGxVPI95NXpsW2FBiyl8GE+3oS2dzUHgISV/Ht57D6wgwy3Td+klKrSWjemuueaU7c3ABMG3a86dVsDydbS4NtfPtwDa63vAO4AqKmp0bW1tcMMceTV1dUB4JSYJJ6jO1w8O7pf4amDr3DVFddQmjnFETGNps5oH7968d8489TTuW7S6cbjOZZf/fa30NPOqUsWU7vkVKOxvPvTB+nf30BlZYlt52ft6r2sfOwPXPaRSzl54aS/+H1dXR3bgDqH/L0AalOvoZEy3I+LjwEDFXQ3A48Ouv3TqSq804AurXUj8CxwsVKqKFXgcHHqNiFslRgoCfe4Y0wpaiWfb7qMKcVTnS6ZPvOFKNG4gTGlVPWd38VjSsd8pSql/kSylVOqlDpAsoruduA+pdRngXrgutTdnwIuA3YCIeAWAK11u1Lqe8B7qft9V2v94eIJIUbdBys6mL/o2WEgKfnTZJ5SPNVj74y17yzbu+8Gqu+8Ll5m6JhJSWt94xF+dcFh7quBLx3hce4E7jyu6IQYYW4rCU/flpL5eGPxBF6bs9LAdhlSEi6ES1h6YDsCd7zpo4n0SkqxVHGTE1Z0iMbjeO1uKUlJuCQl4S6aZPeIB3e86aNW8iKXLkkpnvrqiAVZ4wl8NreUDo0pSUtJCHfQIzd7IS3E0q37bqCl5JjuO3uPOZCU3Dym5N5nLtzp0Nw3d8wD+aDQwfxFfigGxpSc0H1npNAhJssMSVISrjLQUlLKHS/9gaTkhIv8UByqvnNAS8lESXhMxpQkKQl3OZSUDMdhl7QrCXdU9V3cWPedT1pKQriFO7vvAmmyyZ9T5ilprY1OnpWWkhAu8UFLyV1JKd1aSqa7G+OWhdaY23lWCh2EcAt3tpTSpfpuYJ6S6e67Qysr2D2mFE/g83lQNh/XSSQpCVcZWHneLW/6dEtKH8xTMhtvdKA02/aWkuXq8SSQpCRcR1pKTvZB953ZC3M0nlq418Amf24eTwJJSsJlPkhJ7khKMSuBB4U3TUrg41rjAbwes/HGTLWUEpKU0uOVKsQIGVhmyC1JKZqI4/f40qa7Mg6O2OAvmio4MNFS8vvdfVl297MX7uOyFR1iOp42lXcACa3xOeBPE0uk1ki0OZZYLIFXWkpCuMehlJQmLYcTldAWvjTpugOwcMZFybLMJCWttau3Qgdn/P2FEKMkYVl40i4pmb8oJ6yBKk17j2slJCmlz6tViBHhrlXCE+i0KXKAZFJywiXZ0mYmWVta4zFc5GGau5+9cCGNMy579khYVlolJY12xEXJMtVSsjQel3QtH4kT/v5C2MotlXcAFlZaXeQsbf84zuHjGKjStJe2NB67t7t1GElKwlXctslfQqdXSynZfWf+ojwwpmR3grQs7ZoinCNJn1erEOK4pVtS0jjjomSu+86SQgfTAQhhNyd8EreLpdNrjMJyzJiSme47KXSQpCRcxp3dd+kzGdMp1XcJPdB9Z3P1XUKjpKUkhItol1Xfacv2pXJOhNbOaMlqy8yHF5k8K0lJuJCbBpItnX7zlJxwTf6gpWTvca2EjCmlz6tViBHgxu679FrRwSljSoYKHaSl5Ii/vxBilCSTUvpc5JzSuZowVOigZfKsJCXhLhoLRfoM/J8oKy0XZDV/UbYsM4UOCUsKHdLn1SrECLB0Iq26s05U2nXfaWdclEx132lLSsLd/eyF6yRbSu552Vta40mj56vRtieCw/lgQVabj2vJmFL6vFqFGAFaW6g0ajmcqIS28KXRJ2+nrOiQGFj7TgodbKe0dm41UnV1tV6xYoXpMA4JBoMAVFZWGo4kSeI5usPFE655h1jRfvLf/5hjYhpN95XUk5vwcVnneEfEcyw/62wl39J8urjMaBybWyI8vzvEZeMi5PrtOz9PPdhEdo6Xcy8pPezvg8EgXUDQIX8vgMpgkNtvu2211nrJSDyeEz6UCGEb7Zj6LntYaEdMRh0qp6zoMPBZ3fYxJQ0uasgfls90AEfj9/upra01HcYhdXV1AI6JSeI5usPFs7KxhX19XcZitPscPf3aj6nJraD26pscEc+x/OxnPybb4zMez4OvreeFPS9QUV5ObobHtnjeeP6XVE8sobb28C35uro6tgF1Dvl7AdSmXkMjxeU5WbiN1u4qdEhonXbLDDlhSOVQ9Z2B47ppxZHDcc+7UwgGNr1zz8veIv0WZHXCPKWBZYZMlIR7nZCVDXLPu1MIXNhSstJrRQfnLDOUrL6zOz8kLEsmz5oOQAg7aRJ40qjlcKISpOGCrKaDYHD3nb0JQss8JUf8/YWwjeXCeUrplpSccElOGNt5VpJS+rxahRgBmoTr1r5Lp+67ZKGD+Xgtbab7TnaelaQkXMbSVlotu3OiEtrCl0bdlU4ZU0oYqr7TljOWWTLJCX9/IWyjcVn3nSXdd8NhbD8lWZBVkpJwF7dV38V1Am8aXeQcUxJuaD8ly9J4vOafv0np82oVYgRYLmspxdOu+84ZFyUrtdme3RNZtZbJs074+wthG60TeFxS6KC1To0ppcfbfGC7CCcUn5lardtKWFJ9ZzoAIezkpjGlge0X0mVMaaDLzAnRJizLyMoKUn3njL+/ELZx0yZ/A0nJ50mPluFAvE4YUzJVcKBT3YZu5o53pxApVpptD34i0q2lFDdUhn04CcvMQraWpWWZIdMBCGEnN7WU4mmWlA5NWDUcBxgcU5IVHZy9n5IQI606+2QCnjzTYdjCoxTnVcxlQnaJ6VCGxKM8zPJlUOyA7sbp40vpO2kKELT1uGcsm8HkKWZ33TVNkpJwlWXlnzcdgm1yfZn8YOEnTIcxZLkZGVyX7YwPDNeedRLXnnXSoU0Q7fKP37vW1uM5kRNaykIIIQQgSUkIIYSDSFISQgjhGJKUhBBCOIYkJSGEEI4hSUkIIYRjSFISQgjhGJKUhBBCOIYkJSGEEI4hSUkIIYRjSFISQgjhGJKUhBBCOIYkJSGEEI5he1JSSl2ilNqmlNqplPqm3ccXQgjhXLYmJaWUF/hf4FJgDnCjUmqOnTEIIYRwLrtbSkuBnVrr3VrrKHAPcJXNMQghhHAopbW272BKfQy4RGv9udTPnwJO1Vp/+XD3r66u1itWrLAtvmMJBpO7UFZWVhqOJEniOTqnxQPOi0niOTonxtMFBB0SD0BlMMjtt922Wmu9ZCQez3E7zyqlbgVuTf0YueWWWzaajCcNlAKtpoNwMDk/xybn6Ojk/BzbzJF6ILuTUgMwYdDP1anbDtFa3wHcAaCUWjVS2XesknN0dHJ+jk3O0dHJ+Tk2pdSqkXosu8eU3gOmK6VqlFIZwA3AYzbHIIQQwqFsbSlpreNKqS8DzwJe4E6t9SY7YxBCCOFcto8paa2fAp4a4t3vGM1Yxgg5R0cn5+fY5BwdnZyfYxuxc2Rr9Z0QQghxNLLMkBBCCMdwbFJy63JESqk7lVLNSqmNg24rVkqtVErtSH0tSt2ulFI/TZ2j9UqpRYP+z82p++9QSt1s4rmMBqXUBKXUS0qpzUqpTUqpr6Vul3OUopTKVEq9q5RalzpH/5q6vUYp9U7qXNybKjZCKRVI/bwz9fvJgx7rttTt25RSyw09pVGhlPIqpdYopZ5I/SznZxCl1F6l1Aal1NqB6jpb3mdaa8f9I1kEsQuYAmQA64A5puOy6bmfDSwCNg667T+Ab6a+/ybwg9T3lwFPAwo4DXgndXsxsDv1tSj1fZHp5zZC56cKWJT6Pg/YTnLJKjlHH5wjBeSmvvcD76Se+33ADanb/w/4Qur7LwL/l/r+BuDe1PdzUu+9AFCTek96TT+/ETxPfwPcDTyR+lnOz5+fn71A6YduG/X3mVNbSq5djkhr/SrQ/qGbrwLuSn1/F3D1oNt/p5PeBgqVUlXAcmCl1rpda90BrAQuGfXgbaC1btRav5/6vgfYAoxHztEhqefam/rRn/qngfOBB1K3f/gcDZy7B4ALlFIqdfs9WuuI1noPsJPkezPtKaWqgY8Av079rJDzMxSj/j5zalIaD+wf9POB1G1uVaG1bkx9HwQqUt8f6Ty54vylulEWkmwJyDkaJNU1tRZoJnkh2AV0aq3jqbsMfr6HzkXq911ACWP7HP038A+Alfq5BDk/H6aB55RSq1VypR2w4X3muGWGxNFprbVSyvUlk0qpXOBB4Ota6+7kB9ckOUegtU4AC5RShcDDwCyzETmHUupyoFlrvVopda7hcJzsLK11g1KqHFiplNo6+Jej9T5zakvpmMsRuUxTqilM6mtz6vYjnacxff6UUn6SCemPWuuHUjfLOToMrXUn8BJwOskulYEPooOf76Fzkfp9AdDG2D1HZwJXKqX2khwaOB/4CXJ+/ozWuiH1tZnkB5ul2PA+c2pSkuWI/txjwEDVys3Ao4Nu/3Sq8uU0oCvVtH4WuFgpVZSqjrk4dVvaS/Xl/wbYorX+r0G/knOUopQqS7WQUEplAReRHHt7CfhY6m4fPkcD5+5jwIs6OUr9GHBDqvqsBpgOvGvLkxhFWuvbtNbVWuvJJK8tL2qtP4Gcn0OUUjlKqbyB70m+PzZix/vMdIXHUSo/LiNZWbUL+EfT8dj4vP8ENAIxkv2vnyXZf/0CsAN4HihO3VeR3DRxF7ABWDLocT5DcuB1J3CL6ec1gufnLJJ93euBtal/l8k5+rNzNB9YkzpHG4F/Tt0+heRFcydwPxBI3Z6Z+nln6vdTBj3WP6bO3TbgUtPPbRTO1bl8UH0n5+eD5zWFZGXhOmDTwDXYjveZrOgghBDCMZzafSeEEMKFJCkJIYRwDElKQgghHEOSkhBCCMeQpCSEEMIxJCkJIYRwDElKQgghHEOSkhBCCMf4/5gEGGgcWWgaAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 720x720 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig = plt.figure(figsize=(10, 10))\n",
    "ax = fig.add_subplot(1, 1, 1, aspect=\"equal\")\n",
    "\n",
    "# Next we create an instance of the ModelMap class\n",
    "modelmap = flopy.plot.PlotMapView(model=gwf, ax=ax)\n",
    "\n",
    "ghb_quadmesh = modelmap.plot_bc(name=\"ghb\", plotAll=True)\n",
    "riv_quadmesh = modelmap.plot_bc(name=\"riv\", plotAll=True)\n",
    "linecollection = modelmap.plot_grid()\n",
    "contours = modelmap.contour_array(h[0])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Post-Process Flows\n",
    "\n",
    "MODFLOW 6 writes a binary grid file, which contains information about the model grid.  MODFLOW 6 also writes a binary budget file, which contains flow information.  Both of these files can be read using FloPy methods.  The `MfGrdFile` class in FloPy can be used to read the binary grid file, which contains the cell connectivity (`ia` and `ja`).  The `output.budget()` method in FloPy can be used to read the binary budget file written by MODFLOW 6."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:43.169859Z",
     "iopub.status.busy": "2022-07-30T20:16:43.169582Z",
     "iopub.status.idle": "2022-07-30T20:16:43.175629Z",
     "shell.execute_reply": "2022-07-30T20:16:43.174966Z"
    }
   },
   "outputs": [],
   "source": [
    "fname = os.path.join(workspace, \"{}.dis.grb\".format(model_name))\n",
    "bgf = flopy.mf6.utils.MfGrdFile(fname)\n",
    "ia, ja = bgf.ia, bgf.ja"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:43.179112Z",
     "iopub.status.busy": "2022-07-30T20:16:43.178835Z",
     "iopub.status.idle": "2022-07-30T20:16:44.077225Z",
     "shell.execute_reply": "2022-07-30T20:16:44.076167Z"
    }
   },
   "outputs": [],
   "source": [
    "flowja = gwf.output.budget().get_data(text=\"FLOW-JA-FACE\")[0].squeeze()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:44.081916Z",
     "iopub.status.busy": "2022-07-30T20:16:44.081599Z",
     "iopub.status.idle": "2022-07-30T20:16:44.089488Z",
     "shell.execute_reply": "2022-07-30T20:16:44.088401Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Printing flows for cell 412\n",
      "Cell 412 flow with cell 262 is 251.46262091207623\n",
      "Cell 412 flow with cell 402 is 0.7176346498656017\n",
      "Cell 412 flow with cell 411 is 439.8629968543804\n",
      "Cell 412 flow with cell 413 is -693.4212447574185\n",
      "Cell 412 flow with cell 422 is 1.3779378787739915\n"
     ]
    }
   ],
   "source": [
    "# By having the ia and ja arrays and the flow-ja-face we can look at\n",
    "# the flows for any cell and process them in the follow manner. Note\n",
    "# layer, row, column locations are zero-based.\n",
    "k = 2\n",
    "i = 11\n",
    "j = 2\n",
    "cell_nodes = gwf.modelgrid.get_node([(k, i, j)])\n",
    "\n",
    "for celln in cell_nodes:\n",
    "    print(\"Printing flows for cell {}\".format(celln))\n",
    "    for ipos in range(ia[celln] + 1, ia[celln + 1]):\n",
    "        cellm = ja[ipos]\n",
    "        print(\n",
    "            \"Cell {} flow with cell {} is {}\".format(\n",
    "                celln, cellm, flowja[ipos]\n",
    "            )\n",
    "        )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Post-Process Head Observations\n",
    "\n",
    "MODFLOW 6 observations can be read using the `output.obs()` method in FloPy."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:44.094139Z",
     "iopub.status.busy": "2022-07-30T20:16:44.093826Z",
     "iopub.status.idle": "2022-07-30T20:16:44.259819Z",
     "shell.execute_reply": "2022-07-30T20:16:44.258920Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXAAAAD4CAYAAAD1jb0+AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAAsTAAALEwEAmpwYAABLgklEQVR4nO2deXyUxf3H35MbcgK5SMJ9hfsUb+p9tfVorfVoa61n1V+1rS3VXmqtVVuPtmrVeqDWCwRFAUVAQFEEwhGOHCTcCblDstlcm92d3x/z7LO7kJBAgrj6ffPKi+zsPDPfmfnOZ77P7LMTpbVGEARBCD3CjrcBgiAIwtEhAi4IghCiiIALgiCEKCLggiAIIYoIuCAIQogS8WVWlpycrAcPHvxlVikIghDyrF+/vlprnXJw+pcq4IMHDyYnJ+fLrFIQBCHkUUrtaS9dtlAEQRBCFBFwQRCEEEUEXBAEIUQRARcEQQhRRMAFQRBCFBFwQRCEEEUEXBAEIUT5Up8DP1qa1q+ntajYn6CU9T+HpgHKfj8oQzvXqw7ebu/9zspq5/qO6urk+k7b0mFZh7n+KG05Jtd3qazuXt+JLRZBxynroDdo9w3d+e8dlslBRzd3p6wjtK9b5Ry2rC60VWvQ2n+d1uZt+1rdTh7f+9qfz/d+4LUH5dFad1i+/V5g+b7rDy6fg94PKj/wva6UbxISL76YqB7+ImNICLhj4UIOvP7G8TZDEAThqOk1eXKPC7j6Mv+gw7Rp0/TRfBPT09CAt7nZvLDN7STq6DAS6e717UUcXS/r2FxP0Ep/XG05gsiwu23pPBI8zPVHckdzUJ5O7woOKaeDMnuwrC6V04V2HnKX0pP2+X7sfCr4TuvgPNb76uDr7NcqoKh28nRWvlLB7x1SvrGhw/KV8jexg/Lbu+s7GpRS67XW0w5OD4kIPDw+nvD4+ONthiAIwleKkBDwNWVr2H5g+2HzqINDm/by9MBq2BP1dKmMTvL0SFu6UEZXbNW0ExEHvLbfPyifpuP8XclzpGUe7nr/fx3bcEj7jrTMTmw/XLlH2x+d2noU5R6uz46ozCMYi0P640j6OCCPry87Gl+7bE1QXvv6g8ak3fcCygkchz+f/GempR8SRHeLkBDwJXuW8FbhW8fbDOEbiCLgNjvwtX1vTrvph+Q/aLE84nKV6rCO7pR5ONsPV+7B19p5u1pmJ/m6Uu4heQ8zFoHl+V6HhYWZcgK3Wnz/lD/94DoDr7F/Dyj34H7yvR8XFUdPExJ74M3uZtq8bR2+3xNt6EoZB0cUR5Onp+rptIxO6ulSW7pQRk8JVnsT/Kgndw+LrSAcb0J6D7xXRC960et4myEIgvCVQr7IIwiCEKKIgAuCIIQoIuCCIAghigi4IAhCiCICLgiCEKJ0WcCVUuFKqY1KqQXW6yFKqTVKqWKl1FtKqahjZ6YgCIJwMEcSgd8B5Ae8fhh4XGs9HDgAXN+ThgmCIAiHp0sCrpTKAr4NPG+9VsBZwNtWlpeBS4+BfYIgCEIHdDUCfwL4LeC1XvcD6rTWbut1CZDZ3oVKqZuUUjlKqZyqqqru2CoIgiAE0KmAK6W+A1RqrdcfTQVa6+e01tO01tNSUlKOpghBEAShHbryVfpTgYuVUhcBMUAC8E8gSSkVYUXhWUDpsTNTEARBOJhOI3Ct9d1a6yyt9WDgSuBjrfU1wHLgcivbtcD8Y2alIAiCcAjdeQ58JvArpVQxZk/8hZ4xSRAEQegKR3QaodZ6BbDC+n0nML3nTRIEQRC6gnwTUxAEIUQRARcEQQhRRMAFQRBCFBFwQRCEEEUEXBAEIUQRARcEQQhRRMAFQRBCFBFwQRCEEEUEXBAEIUQRARcEQQhRRMAFQRBCFBFwQRCEEEUEXBAEIUQRARcEQQhRRMAFQRBCFBFwQRCEEEUEXBAEIUQRARcEQQhRRMAFQRBCFBFwQRCEEEUEXBAEIUQRARcEQQhROhVwpVSMUmqtUipXKbVNKXWflT5LKbVLKbXJ+pl0zK0VBEEQbCK6kKcVOEtr7VRKRQKrlFIfWO/9Rmv99rEzTxAEQeiITgVca60Bp/Uy0vrRx9IoQRAEoXO6tAeulApXSm0CKoElWus11lt/VUptVko9rpSK7uDam5RSOUqpnKqqqp6xWhAEQeiagGutPVrrSUAWMF0pNQ64G8gGTgD6AjM7uPY5rfU0rfW0lJSUnrFaEARBOLKnULTWdcBy4AKtdZk2tAIvAdOPgX2CIAhCB3TlKZQUpVSS9Xsv4FygQCnV30pTwKXA1mNnpiAIgnAwXXkKpT/wslIqHCP4s7XWC5RSHyulUgAFbAJuOXZmCoIgCAfTladQNgOT20k/65hYJAiCIHQJ+SamIAhCiCICLgiCEKKIgAuCIIQoIuCCIAghigi4IAhCiCICLgiCEKKIgAuCIIQoIuCCIAghSle+iSkIgnBMaGtro6SkhJaWluNtyleCmJgYsrKyiIyM7FJ+EXBBEI4bJSUlxMfHM3jwYMyxSt9ctNbU1NRQUlLCkCFDunSNbKEIgnDcaGlpoV+/ft948QZQStGvX78juhsRARcE4bgi4u3nSPtCBFwQBCFEEQEXBOEbTVxcXNDrWbNmcfvttwPwzDPPMH78eCZNmsRpp51GXl4eADU1NZx55pnExcXZeX1ccMEFTJw4kbFjx3LLLbfg8XgAyM3N5eSTT2b8+PF897vfxeFwdNt2EXBBEIQOuPrqq9myZQubNm3it7/9Lb/61a8A87TIX/7yF/7xj38ccs3s2bPJzc1l69atVFVVMWfOHABuuOEGHnroIbZs2cJll13G3//+927bJwIuCILQAQkJCfbvjY2N9h51bGwsp512GjExMR1e43a7cblc9jXbt29nxowZAJx77rnMnTu32/bJY4SCIHwluO/9beTt7/62QiBjMhL483fHHjZPc3MzkyZNsl/X1tZy8cUX26+feuopHnvsMVwuFx9//HGX6j3//PNZu3YtF154IZdffjkAY8eOZf78+Vx66aXMmTOHffv2HXmDDkIicEEQvtH06tWLTZs22T/3339/0Pu33XYbO3bs4OGHH+aBBx7oUpmLFy+mrKyM1tZWW/RffPFFnn76aaZOnUpDQwNRUVHdtl0icEEQvhJ0Fikfb6688kp+/vOfdzl/TEwMl1xyCfPnz+fcc88lOzubjz76CDDbKQsXLuy2TRKBC4IgdEBRUZH9+8KFCxkxYsRh8zudTsrKygCzB75w4UKys7MBqKysBMDr9fLAAw9wyy3d/zvwEoELgiB0wJNPPsnSpUuJjIykT58+vPzyy/Z7gwcPxuFw4HK5ePfdd/noo4/o168fF198Ma2trXi9Xs4880xbqN944w2eeuopAL73ve9x3XXXdds+pbXudiFdZdq0aTonJ+dLq08QhK82+fn5jB49+nib8ZWivT5RSq3XWk87OK9soQiCIIQonQq4UipGKbVWKZWrlNqmlLrPSh+ilFqjlCpWSr2llOr+R6qCIAhCl+lKBN4KnKW1nghMAi5QSp0EPAw8rrUeDhwArj9mVgqCIAiH0KmAa4PTehlp/WjgLOBtK/1l4NJjYaAgCILQPl3aA1dKhSulNgGVwBJgB1CntXZbWUqAzA6uvUkplaOUyqmqquoBkwVBEAToooBrrT1a60lAFjAdyO5qBVrr57TW07TW01JSUo7OSkEQBOEQjugpFK11HbAcOBlIUkr5niPPAkp71jRBEIRjT08eJ9vU1MS3v/1tsrOzGTt2LL/73e/s9/bu3cuZZ57J5MmTmTBhAosWLeq27V15CiVFKZVk/d4LOBfIxwj55Va2a4H53bZGEAThK8TRHCd71113UVBQwMaNG/nss8/44IMPAHjggQe44oor2LhxI2+++Sa33nprt+3rSgTeH1iulNoMrAOWaK0XADOBXymlioF+wAvdtkYQBOErxJEeJ9u7d2/OPPNMAKKiopgyZQolJSWA+XNpvj/iUF9fT0ZGRrft6/Sr9FrrzcDkdtJ3YvbDBUEQus8Hv4PyLT1bZvp4uPChw2Y5FsfJAtTV1fH+++9zxx13AHDvvfdy3nnn8e9//5vGxkaWLl16ZG1pB/kmpiAI32iOxXGybrebq666il/84hcMHToUMGeh/PSnP6WkpIRFixbx4x//GK/X2y3b5TArQRC+GnQSKR9vjuQ42ZtuuokRI0Zw55132mkvvPACH374IQAnn3wyLS0tVFdXk5qaetQ2SQQuCILQAUd6nCzAH/7wB+rr63niiSeC0gcOHMiyZcsAc2BVS0sL3X20WiJwQRCEDjjS42QTEhL461//SnZ2NlOmTAHg9ttv54YbbuDRRx/lxhtv5PHHH0cpxaxZs+wPRY8WOU5WEITjhhwneyhynKwgCMI3ABFwQRCEEEUEXBAEIUQRARcEQQhRRMAFQRBCFBFwQRCEEEUEXBCEbzSHO072k08+YcqUKURERPD2228fcq3D4SArKyvoSNkLLriAiRMnMnbsWG655RY8Hs8xs10EXBAEoQMGDhzIrFmzuPrqq9t9/49//CMzZswISps9eza5ubls3bqVqqoq5syZc8zsEwEXBEHogMGDBzNhwgTCwg6VyvXr11NRUcF5550XlO47gtbtduNyubr9bcvDIV+lFwThK8HDax+moLagR8vM7pvNzOkzD5uns+Nk28Pr9fLrX/+a//3vf+0eC3v++eezdu1aLrzwQi6//PJ2SugZJAIXBOEbTWfHybbH008/zUUXXURWVla77y9evJiysjJaW1uP6AzxI0UicEEQvhJ0Fil/lVi9ejWffvopTz/9NE6nE5fLRVxcHA895D8SNyYmhksuuYT58+dz7rnnHhM7RMAFQRCOkNdee83+fdasWeTk5PDQQw/hdDppaGigf//+uN1uFi5cyOmnn37M7JAtFEEQhA5Yt24dWVlZzJkzh5tvvpmxY8ceNn9jYyMXX3wxEyZMYNKkSaSmpnLLLbccM/vkOFlBEI4bcpzsochxsoIgCN8ARMAFQRBCFBFwQRCEEKVTAVdKDVBKLVdK5Smltiml7rDS71VKlSqlNlk/Fx17cwVBEAQfXXmM0A38Wmu9QSkVD6xXSi2x3ntca/2PY2eeIAiC0BGdCrjWugwos35vUErlA5nH2jBBEATh8BzRHrhSajAwGVhjJd2ulNqslHpRKdWng2tuUkrlKKVyqqqqumetIAhCD3O0x8lecMEFJCUl8Z3vfCco/ZprrmHUqFGMGzeOn/3sZ7S1tQEwf/58+/nwadOmsWrVqm7b3mUBV0rFAXOBO7XWDuA/wDBgEiZCf7S967TWz2mtp2mtp6WkpHTbYEEQhC+Lwx0n+5vf/IZXX331kPRrrrmGgoICtmzZQnNzM88//zwAZ599Nrm5uWzatIkXX3yRG264odv2dUnAlVKRGPF+TWs9D0BrXaG19mitvcB/gendtkYQBOErxOGOkz377LOJj48/JP2iiy5CKYVSiunTp1NSUgKYSN93tGxjY2OPHDPb6R64MrW8AORrrR8LSO9v7Y8DXAZs7bY1giB8Yyl/8EFa83v2ONno0dmk33PPYfMczXGyXaGtrY1XX32Vf/7zn3baO++8w913301lZSULFy7sdh1deQrlVODHwBal1CYr7R7gKqXUJEADu4Gbu22NIAjCl4zvOFkfvsOpusutt97KjBkzgg6zuuyyy7jsssv45JNP+OMf/9juWeJHQleeQlkFtBfrL+pWzYIgCAF0FimHEvfddx9VVVU8++yz7b4/Y8YMdu7cSXV1NcnJyUddj3wTUxAEoQd5/vnnWbx4MW+88UbQ3nlxcTG+wwM3bNhAa2sr/fr161ZdIuCCIAgdcLjjZE8//XR+8IMfsGzZMrKysli8eDEAt9xyCxUVFZx88slMmjTJ/gs/c+fOZdy4cUyaNInbbruNt956q9sfZMpxsoIgHDfkONlDkeNkBUEQvgGIgAuCIIQoIuCCIAghigi4IAhCiCICLgiCEKKIgAuCIIQoIuCCIHyj6enjZK+//nomTpzIhAkTuPzyy3E6nQDs2bOHs88+mwkTJnDGGWfYh1x1BxFwQRCEDjia42Qff/xxcnNz2bx5MwMHDuTJJ58E4K677uInP/kJmzdv5k9/+hN33313t+0TARcEQeiAozlONiEhAQCtNc3Nzfa3LfPy8jjrrLMAOPPMM5k/f3637evKaYSCIAjHnE9nb6d6n7NHy0weEMfpV4w8bJ5jcZzsddddx6JFixgzZgyPPmr+1s3EiROZN28ed9xxB++88w4NDQ3U1NR06zwUicAFQfhG4ztO1vfjO7ukO7z00kvs37+f0aNH89ZbbwHwj3/8g5UrVzJ58mRWrlxJZmYm4eHh3apHInBBEL4SdBYphxrh4eFceeWVPPLII1x33XVkZGQwb948AJxOJ3PnziUpKalbdUgELgiC0ENorSkuLrZ/f++998jOzgaguroar9cLwN/+9jd+9rOfdbs+EXBBEIQOONLjZLXWXHvttYwfP57x48dTVlbGn/70JwBWrFjBqFGjGDlyJBUVFfz+97/vtn1ynKwgCMcNOU72UOQ4WUEQhG8AIfEhZktjG22tnuNtRo/RzT/C8RXj69OYr9W4hEhbvF6Nx+093mZ8KYSFKVRYzw5MSAj4mvd2snVl6fE2QxCEHmbaNX2oLmno9p8WCwUSU3sT3evwknukW9ohIeAjp6eTMvDQbzyFJF/eRw7HnC/z8xOh64TSsHh7OfCEN5OU2OdrL+IRkYffsdZaU1NTQ0xMTNfL7CyDUmoA8AqQhpGf57TW/1RK9QXeAgYDu4ErtNYHulzzEdB/WCL9hyUei6IFQTiOtLWlUlJSwp7SYyIdIUdMTAxZWVldzt+VCNwN/FprvUEpFQ+sV0otAX4KLNNaP6SU+h3wO2DmUdjcc7S1QGQXVi+3CyKierZurxfaOS/hENqaITy6a3l7Gq3B3Xr8+qgrtLXA+lkw7vsQl/Ll1691aGyGH0873a0QEd15Pq8Hwg7/TcPIyEiGDBnSfh1eN0TFHpoeFtFpubhbAXV8fFhraKmHXknHvKpOVURrXaa13mD93gDkA5nAJcDLVraXgUuPkY0HGwSbXocP7zZi6GP1U/DQQChe5k+r2QGzfwJVhf608q3w8GBY8CsjumBEY/7tsHWeP5/HDbO+A8+dAZX5/vQtb8PHD4Bjvz9t27vwyBBY9bg/rakWFvwSCj/wpzXXwVMnwsvftRzMYvmDsPKR4HvfeTfDExOC21O0FD76A9Tu8qft+gSenA5vX+9Pa3XCJ/+ANc8FtKcN3rgK/jUJnFX+9Lz58M4tpg98rHnW9OXmOf60qkKYdxNUbAtI2w5PnQSvXua/3uOGnBfNj69/tYZVT8Dj44KvP7AbPrwHnJX+fO/+HD6cCW9fZwQAwNUE798Bu1f5r21rgRcvMOPbVOtP3zwbNr5m7PCR+yY8+y0zTj6cVcamslx/WlMtPHu6qSuQL54xZQTy7q3wwnlwYI8/bcfHpsxAvyxYCC9dBJ//25/W6jT9s+tTf5qrEebeCC+cbxZPu8zlh/rGuufhsdHBttfsgA9mQmONP61sMzx5Aiz5k38sXE3w8V+NXT7crfDubfD4eKgu8qfv+RwW3hXsL9s/Mr6R82Jwv719Pez+zJ/WWGNsfP2H0Npg0rSGpfcZ//IGPJTw3i+MD/vmiscNz58DD2bCuhcC7Fltynzn5mDbF/8e1v7Xn+abu0+fCI6ygH57wcz7QF9f+Qj8fUTw9fs3waLfwN4v/GlVhfDMafD0Kf7x8biN7nx4d7CvL7gTHh0F+zf6rz9G+1pHFAYqpQYDk4E1QJrW2tc75ZgtlmOHb8Dz3zOT/IunIeclk1ZZAIvvMSv2+3eYQW1rhlcvNQI15zp/p3/0e2hrhJwXoNBy4s//DRtfNaLhE8ztH8LuT80gLH/QpNXsgLnXwyd/96e5Go34t9SZNJ/DvHmNcfK3f+YXmGX3Qf0+2LPKTHSAbe/Ayodh+V9hzTMmrSIPNr8JdXuMI3m9xv73bje2vnmNcQitTburC2Hr21BiPWO/5I/w8V/gg9/4F5+1/4XtHxix/Mj6AkHpemNf7htmkoOxddn94G6G+bea11rD+3fC5reMaLbUm7xfPA1V+Ua4tlhin/u6WbgW/BKKPjJp+9bA0j+bti++x6Q1HzCL4xdPmbyOMti3FrbNg0Gnmr7fvthqz59MVP6/y82i4atn72ozvp9ZfVlVCPNuNHZ/+g+T1lgD7/0flG0y/eJpM+159TJj07yb/ZNvwZ1QvsXU5RO49S+bBeXdn/sXkD2rYdNrpl0Lfmm1pw7euNqU+c4tJk1rWHov7PnM9GlDhb+eBb804tZsbR189k/YMhv2fQFrnzVpO1caH17+VxO0gFkwFv8eGsrgzR8ZEfG0Gd9d84zJ75srKx+G6u2m7B2WX6/7L3zyCLx5tV9wt70Lm/4H9XvNNWCE/KULTf4Vlq+3OmHuDaa+D++Gun0mfe71xv/e+hE0lJu0Tx8FZ4WZRz6/LlgIqx6DD35rRBzM3N3wsvHh+bebwCjnRSjfDGgz1zxu05fv/R801Rhf843FB7+F1U/CorvMggXm+pK1xr4PrU2Bss2w8Fdm3i/8tX/MPvsnNFaawKix2tSz6C5Y+xy89gP/4rP8QeMblduM74EZp8X3mHmQ965JK15m/Ed7jR9oDfUl8K/JwQFID9FlAVdKxQFzgTu11o7A97T5NKvdJUYpdZNSKkcplVNVVdVels5Z+YhxJq1h9dPQZzAMOs1MXLfLCIsKh0ueMiJRsMAIY91emH6T6fTtH5oJvnMFnP1niM8wk7PFYcoZeQHEpvgji/UvQXx/OOX/THkHdpt84dEw6tvGiZpqIX8BuBpM3V63cZD9G2Hv5zD5R2YhWfuc+X/zbJh4NQw7Cza8YibaiocgbRwMPMU4utaw5j8QGQsXPAy1O2DXSiPoDWUw9numPbtXwf4NxqnOuQ969TH2NdXCpjcg+zsQFW8icTAinTEFTrwFts41Qr76KXOLOv4K096mWhNpupxw2XPgcRmbi5dZ7fkxtDrMXUir0/TBxKuN/V/8x0y0VY9D8iiIS/NP3E2vmfac+XvT//s3GhubD5h+L1gAj4+BOddCTBJc+Tr0TjZtri8xfTr+B6A9xk6tje0ZU0w7N7xqoqpVT0BEL5O+7gXjG7mvm3acc6/xh63zTH9WbIGhZ5gFKH++aXvBQjjx59BniGmP12N8L2u68ZdP/m7as/pJ4yun/MIsXgf2mDR3Cww90wQZjv1mganeDqf90tiQ8wLU7jT9P+BEE0ise94sILlvmGuHnmkWaU8bfP4v049p48xdn9djhNzdChc8ZAR3xzKziJXlwphLjfDtWmnu0goWwql3GNs//7fpo9VPwYCTTD/7For1s6DvUNOeLW8bX1/7HIRHweiLja/Wl0LhImith0v/Y2zY+D8o3WD6YNr10FxrAiFXkylz4lUw+HQzPl6PaUPyKMicavpCayPeYZHwo3lGLP81xSy0w86CH75mfH77hybYqCmCC/8OsammHY01xo/GXGr837dor58FmdPgxJuhYJER5k8eMXnGfd/4bfMB45cuJ3z/BTN2OS+aeVWyDiZcaXw9900zVwoWwEm3Gds/+6fx9TXPmnFMHWvmsdbG33r1hW8/ClUF5i5mwyumTxMHHJ3+HYYuCbhSKhIj3q9prX37DBVKqf7W+/2Byvau1Vo/p7WeprWelpJylHuaCZkm2ln4axOhnHQrnHanWeHz3zMiM+wsmPBDSBxoJu+aZ4yznP83iEs3jr/lbVBhMOlqmHwNFC810Z3LCd+aadILPzARVvEyI1gn3WquWXqfcZbJP4Kz/mAGfMPLZnCSBhohG36ucdaVfzeCdf6DMGSGqbtwkalnwhWmXEeJiaSqCoyoTr3WDHLxMhMRjbkEpv7UCMWqx4ww9p8Elz5tHGTts8ZRI3vDtOtgyk+Ms370BxM9n3mPsTX/PXOrXr7Z9M/Ua81Cs/JhyHsPJl1jFimPywjL1rchfQJM/KFx1jX/gdX/NpPm249B2nhz57NljmnPtOvgpJ+bRWXhL41AnfV7mH4j7FxuJt62d2HMxWZCRfQyjr/2WSOMV70J131g+quhzNjdKwnGX27GYsmfzMQ4+08w6kKzWO/5HGqK4YTrTZnNtaY9m98yfXbm701Ute0d00cDToJT7oB+w01E+emjZoG48g3jL+tnmcjf64ZJV5l+2/2pEVBHiemfydeYiLhim7mzGP8DExwoZSZ0zkvGvu88ZqKvDa+YcqMTYMZvYOSFxi9XPmL2cK94xfhszkuwa4VZXCZeZcp0VpjFqHgpTPuZub5hv+nPLXNgyOlGMHsnG6Fd94JZdC57BmISjbBunm3mzvSb4MSbjKgvvtuUfeY9xl/yF/gX56k/NX2pFKx42JQx7vtw7v2mX3JfN2KWOMC0fdiZJs/n/zJjd869Jqja9Ibx9bZG41vTrrPuJO8yi+Vpd8IJN5jxK/rIzI3R34XhZ8NtX8Do75hxuuxZs7jHZxhhXT8LImKMX06+xtydLfkjeFrhW781glv4gbmrqNwGE68088zbZu6CChaZ16feYa7ZPNvYmjnV+NqI80xfrnjI+Pp3nzDzbe1/zaLkdVt9dIvx8fm3mXaddKvx/+pCM14Fi2DsZTDucjP2q580geKIc6HPoKPTv8PQqYAr82zPC0C+1vqxgLfeA661fr8W6P7p5B0x4YdGjHNeMBNu6k+N8ydkmts3R4mZzGFhcPJt5pa1LNcMbHiEmZRFi83kHTID4tNNx0fEmIguazpkTjGTAg0vXWAcecpPICHDDMi2ecYZTv0FpI0x5Sy912yHTL/J1D39JnCWm62Z0+40k2nyj8xAL7zL2D74NBM1po414hjf30yU0d81Yv3a983KP+EK82HjKf9n9rlrd8Lpv4LIXtbke9+IxLjvm3qmXW+EY9NrJi1tLEz5sRHml79jJtn4H0DKKDMx1j1v2n/izdB/AvSfCCv+ZgR3/A9Mv59xj1lUdq4wky4iCk66xUSvC+40bcg6wThrbIqxJ3kkZH8Xpl5n7lb+e5Zpzwk3Gjsn/MBEmwd2wxm/M/086BS4cRn8bLGxB+Dk2/2LythLzSI59TpzCz3rImP76O+aCC9zmlnkVBiccrvxjb7D4J2bjFBM/akZn2nXm+hq1yfGN6J6W8K8wgQH/SeZxWvyjyA60YxvvxFGmCdeZWx940pj17jLIWmAEamcF6Cp2ox/36Ew4nzTl5vfMr4bFQsn32ry5L5hrolPN+LsKDXbOb2TzSI34jzjJ8sfMAJwwg2m/l59zP5t7Q5zxxQRZXyxeKkR4Ok3Gt+Y9CPTZyseNL6WmGXaH9HLCGHmVOO7J9xgfP1/3zMR8KRrTN5RFxmxbmsyfdR3iOnjjx8w0f6Ua01fnvILY/u2d2D6DRCTYPyydgd8+DsjvINOMdFxSrapOzbF3EGO/Z65s3j9CrP1OPWnZsz7DoXvPw8//wziUs3cnXqtqXfjq0aAYxKNDT5fH/s94+uTrjbj8taPjN+NuRRSs83c3fiqyT/tZ8bPB5xotl4qthjhB6MHjVVmPp98q+nLE282wrzsfrOtlzLSBFZxaebusP8kc4cy9lIzv1673ARPU681vnXybWYxa6wyd2HHAq31YX+A08xIsxnYZP1cBPQDlgFFwFKgb2dlTZ06VR819aVar3xE6z1f+NNy39L6zwlaP3eW1l6vSfN6tf74r1qveNifr7le68fHa/1Af63Lt/nTt83X+v07tW6o9KfN/z9T5qKZ/rS6fVov+4upz8euT7X+S6rWDw3SurnOn75lrimjrcW8dru0fuVSU2bBIn++8m1av/kjrauL/WnrXjT5/vcDrT1uk9baqPUHd2u98XV/Gw/s1fovaVrfn6z1/k3BNr13h9Z1Jf60935hyvz4wYD+qNN6zs+0LlriT9u3Tut7k7T+z6latzr96cv+ovWqJ0w7fO15fLwpc9t8f769a7Re9Fut96z2py3+vcn30rf9aa2Nxqal9+lO2fqOGQdXk3nt9Wr94kWmzOV/C7Z99k+1znvfn7bmOZPv4SGmTp/ty/+m9dwbtW5rNWmOcjOGf07Qumip//rij7V+/lytKwsPbc+cn/nHwlGm9SuXaZ3zkj/fni9MvnuT/Nd7vVqvflrruTcZX/bZ89TJJu/Kv/uvry7WetZ3gv1l/Ssm3zMztHa3Wde3af3Wj7Vecq/WHo9Jc1aZfH9O0Lp4mf/6TW9ovfCuYH+b/VOTb94t/rTGGq0X/6H99jw0SOuWBn/6uhfMWPrs8Xi0fuF8rf+cGDwWe9dq/c6tWlfk+9M++5cp8+WL/X3ZHs5qMw5/zdD6wJ6AMtdoPe9mY6+PRTNNmav+6U9rqND67euDdWPXKpPvyen+uev1muu+eNbvG65mrZ+YYPLmvee/vrrY9GX51oC+sObu3Bv9aR631h/8TuvNczpuXxcBcnQ7mhr6h1nVl5rVsnffw+dzlJlIMGXU4fO1tZjb1b5DO6+7rdn8dFZ3W4u5fcyYfPh8Wps9xYxJnT8m1eqE8MiuPc5VkWfa3VmZZbkm0u3V5/D5GqtNtNtZu8F88Nu7b+dldpXGGtOXg087fD6v10SDsSmdP87VXGe2mIbMOHw+d6uJqEZe2PmjmKXrjQ911m5Pm3nqIWOyiTgPR9ESsx+e0P/w+cpyTcSd0sn52l6PuaPpndz5Y60N5SbCTRp4+HzNB8zdVWe+7vWaPfykQZ0/Dqm1qbszX/d6zOcrGVM6b09Vodl26uwxQ6/X1N2VR29rdpi7mK7MySOko8OsQl/ABUEQvubIaYSCIAhfM0TABUEQQhQRcEEQhBBFBFwQBCFEEQEXBEEIUUTABUEQQhQRcEEQhBBFBFwQBCFEEQEXBEEIUUTABUEQQhQRcEEQhBBFBFwQBCFEEQEXBEEIUUTABUEQQhQRcEEQhBBFBFwQBCFEEQEXBEEIUUTABUEQQhQRcEEQhBBFBFwQBCFEEQEXBEEIUUTABUEQQpROBVwp9aJSqlIptTUg7V6lVKlSapP1c9GxNVMQBEE4mK5E4LOAC9pJf1xrPcn6WdSzZgmCIAid0amAa60/AWq/BFsEQRCEI6A7e+C3K6U2W1ssfTrKpJS6SSmVo5TKqaqq6kZ1giAIQiBHK+D/AYYBk4Ay4NGOMmqtn9NaT9NaT0tJSTnK6gRBEISDOSoB11pXaK09Wmsv8F9ges+aJQiCIHTGUQm4Uqp/wMvLgK0d5RUEQRCODRGdZVBKvQGcASQrpUqAPwNnKKUmARrYDdx87EwUBEEQ2qNTAddaX9VO8gvHwBZBEAThCJBvYgqCIIQoIuCCIAghigi4IAhCiCICLgiCEKKIgAuCIIQoIuCCIAghyjdSwEsONOH16k7z5Zc5aGx192jdWndeL0BZfTP7aps6zdfm8VLf1NZds4LweHWX6tZak7uvrstt6mrdK7dX0ebxdpq3qKKBJlfn4+Nye/F0YbyPhIaWNuqbO+/3xlY3pXXNPVo3dM2PtNZsKanvUnktbZ7umnQIe2uaulRuyYGmHvfhNo8Xdxd8qLbRRWF5Q6f5tNbHpI+6S8gJuNaaz4urgya416v517Ii7p63+RDBPVgI1u6q5bSHlzM7Z19QmW+t20tlQ4ud5mhp48J/fsp1s9YFXb8kr4LiSmdQ2j3vbOE3c3I7Ffu6JhfnP/EJj3xYEJReVNHAnprGoLRbX9vAtS+uDZqoG/Ye4P3c/UH5HlyUz8T7P6KqodVOa3K52V4R7JRVDa3c+942nli6PSjd6z3UMf/w7lZOf2R5UH/UN7XxyurdQQvfu5tKueSpz3gvwCaPVzN3fUlQmY2tbq6ftY7X1uxpv2MCuOHldVz74lreXl8SlP55cTXOgP6tbXTx7X+t4omlRUH5Fm4uo+igtl/61Gf88NnVQX3paGk7pN07q5zc884WluVXBKV7vDrIj7TWzHhkOec9vjIoX36Zg0+Lgg9su+PNTZz60Mc0tPgFqtrZyifbg/NtLa3n5ldzWL/nAJ3xnxU7OOlvy4L8raXNw2fF1UH5Fm4p47tPruLzHf50l9vL85/uDLq2rL6Z7D9+yLMrdwRdX+NsPSTQWZZfwY2v5Bzir3VNrqD+LatvZsbfl/PAwrygfPM3lQaNj9aa0x5eztmPBffl7a9v4Nezc4PSZn22i9tf38DemsMHF26Plwue+IQbXskJSi+udB7Sv/e9v43vPf0ZzS6/L+yocvL6mr1B+eZtKCX7jx+ybb9/QdRaB807ML7y2JLt3PxqTpeCkO4ScgK+bvcBrn5+DT//3wY7La/MwWNLtvPG2n22E2utuWtOLiN+/wEb9/oH7d8fmwm/vLDSTvusuIaZc7fwp3e32Wmrd9QARvB9wtHkcnPjKzmc89hK21ldbi+vr9nLnPUl/PfTnfb1d83J5eS/LWNJnl8M7l+Qx/YKJ8994s9XXt/CuY9/wg+eWW2n7a1pYuPeOnZWNwYtFve9n8ev5+QGRX5zcozQ/X2xf1H498fFnPf4JyzeVm6nzc7Zx6zPd/PE0iI7aq1tdHHO4ys57eGPaXUbB3a0tPHG2r1BfQDw27m5/Gn+Nj4LEIM1O80pw2+t8y+Gy/Ir+PWcXO573z9xN+6tY1lBJb9/Zyu5++oAIzgXP7mK7//nc/ZbEWqlo4XlhUbYNgRMtA+3lnP182v49zK/WC/eVo7L42XRljJ7LA40uvjFmxu5f4G/7romF3llDnL2HODDrf7+uO6ldWT/8cOg6PgP727l9TV7g8ZxS0k94+9dzMy3N9tpa3fVcqCpjQpHK3VNLsD424X//JQfv7AWl9s/cZdai0HggnTPvC385MW1dl8AvLuxlMXbKrjupbW2aG4trefcx1Zy2+sb7DYWVTTw8IcFVDhag66//fWNXPP8GvL2O+y0Bbll5v/NZf6+3FbOAwvzgwRq7S4zjn/7oIAKh1m0a5ytTH1gKZc89ZndnsZWN9e/nMOSvAoWbvGX+fLnu5l0/xLmbii1015ZbRbrD7b4+3xfbRN3vrWJRxYX2mk7q81CUO1sJb/M2F7X5GLB5jLmbigJCkSeX7WLBZvLeOjDfDtt1me7OPlvy4J84+XVe9hR1ciKwirbr1vdHs55bCXf/8/ndv+2uj0sy6+k0eUJ0oObX13PPe9sse0BWGS19zdz/H7w4dZyTvjrUh79yN+epfkV/GtZEYu3VdiR/YFGF/9eVnTIotcThJyAr91lRGVpvj8SLgi4BdpmOfDumiZ70nxaZESnsdVt/567r96eFC9+tgswC4GPVUV+ofJFZD5HB/jEej/QwTZZE8rr1SzYvJ+y+hbeXu8Xt017zftur7ZvGf9pOV5lQ6sdBQROjo+sBWB/XTO5++pwub22M7W0eewo8tMAez+3FrH7A0Q0sI98v7+3qZSdVY1UO13klzUEXQvYi2GTy83HBcbBVxb6I8cvdpqx+HxHDbWNRsjWWH00d32JHXVuDYhafH2UX+Zgc0k96/ccYIVVpm/s4qMj2BggTn9dlBd0LfgnVMmBZvu6ZQWVeLyaVcXVtjAHlvOJFR03udx2JPb08mLACLCvnPyyBts3nli6nSaXh3kbS+20wDsO3yK3Ya+/nnW7ay3b/JHiu5v81/jqeXTJ9kPSHC1u9lrbV0vyKiiqdLJwcxkVDhPprd3t98Ecqw17a5rshcIXbTe7PKzYbsZs8dZyewvJF1C8uW6v3Z6NAbb7Fk5f/2wprbfFLdDXfVszbR4vDy4ygroiQATfs9pb2+SyF7nZOfvQGj4tqrL9NjBImG9dsyLAx95ca+ZPS5uHsvoW2yYf72wspay+hfntjAnA1lKHVbd/Ad1e2WDnc7a6Uco/59wer60rvkAG/AtNXpnD9mtfMPPvj4vtIG9bwAK62eqjLaX1PLpkO6UHen4rLeQEfN3uAyhlft9SWgdAYbmDqIgwhiTH2h0YeKvji1R2VplBOH1EMuWOFnbXNNHm8dpivbe2yd77/WJnDaePSCYiTNkraeAtqs9JtlrONGVgkl33ntomWtq8lh0mrdXtYXdNI1MHmaPTN5UYmzZb/4PZIvGVObhfb7LT48mxJqxvEenTO9KOrDeX1OP2ak4Y3Iey+hbqm9to83gptCZaaV0zByxhLSx3kJ0eH2RT4IK1yap7R0AfrbYEesOeOto8ZrL7hLy83vTfGaPMEcG+ye3rI5fHa9eztbSezKRe9I2NssclcBL6BN733tUnDqS40kl9cxvVzlb21RrHzy0xC5jWmk376jhndBrgF5sleeUkxESgNSyx+mjj3jrCFIzPTLQXrs0lgXUbGyscrdQ3tzEiNY765jZ7AQjsI5+wFlc6mTggibjoCD63/OCjgLsdn5h+Yd2hnDKsH0UVDXi95pbbV/baXTV4vdpaPOqZkJUYND6BfeTrm6IKJ72jwhmZFmcL+IaAO0zfoppf7qClzcs5o1OpaXSxr9b4+orCShJiIthR1WiL0sZ9dUzMSkQp/+IeuPBts+zwzYMJWYl2H+6saqTVitBzLV9ucpl9/zNGpaC1vx+W5FWQ1DuSljav7SdrdtWSnhDDxKxEux3LCytJjotmQlYiBeUOu889Xs3YjAT21TZT39yGy+21A4+dVU47ANpR5WT64L4A9vzJCVj4Pi82fbRxbx1Kwbmj02yNCFyIl+UbX69rcrGr2j93fb6+pdTvGwWWnxSWOxiSHEtS70h7bvsi+dH9E+hpQkbAtdZ4vJoNew7ww2kDiIoIs28XC8obGJEax/jMRPIsR8/b7yAiTHHxxAw2WR+07agyK+v3pmQCpmN3Vzfi8ni5acZQwEyGNo+XXdWNjMtMZEhyLNsrzHVrdx9g+pC+ZKfH24O4dX898dERXDS+P1UNrVQ2tNgDd/HEDEoOGGfbVd2IV8PlU7NQykTjbo+XokonV00fSJjyR687qpwMTYljZFo8RVY0kFtST0p8NGeMSrUnkq+eK6YNAIxj5ZeZiXvVdJOWV+bA5fays6qRs7JTSewV6e+jMgenDU8mNT7ajm53VDrpnxjDSUP7sa+2mYaWNnsS3TRjKDurG6lvaqPIimIumZQBmFv7+qY2CsobuObEgbY9YARpfGYiYzMSbHHK3VdPclw0Jw3ty9K8Cj4vrmbbfgeD+vXmlOHJ9hj62vrjkwbR0uZlS2k91U4XDS1uTh3ej/iYCNuWTfvqOHt0Gn16R1JojdnmkjpGpsUzeWAS28tNZO0T/O9NyaSgzIHb47Xb+L0pWbbNBxpdlNW3cNlk4y++SHVHVSOj0uIY0z/Bvq6gvIHR/RM4aWhfO58vsLhofH+aXB5K65rt6Pyq6QNoafOy70ATJQeacbS4uWxyJuFhirwyc3e4uaSeC8am233h843hqXFMHdSHTXsPoLUmv9xBZLji8qlZrNlVi9aaXdZCfOG4/vZY7KxqpKHFzeVTB9hj1ubxkr/fwYlD+zG4X6w9Zhv2HGBiViIjUuPsRaygvIHeUeFcNL4/pXXN1ARse1w+NYt9tc1UNbTagdJ3J2SglBF+l9tEtt+bnEWYMv7s66NxmYlMHtiHzSV1tHm8fL6jhlOH92NUWrw993x+4Buf/DIH2ysacHm8XDwxA6+GwooGWt0e9tQ0cuLQvgzs29teVArLGzgrO5VB/XqzxrqL31ndSFafXozPTKTkQDNNLre93Xrj6UMorWvG2eq2F9IrpmXZ/WAWDwcXT8yw/cVXT3Z6POMzE+025pc56J8YQ5/YKHqakBDwfywu5KJ/rWJXdSMNrW6mDupDdnp8UKeNSo9nbEYC++tbONBo9j2Hp8YxfUhfahpdlBxoZkeVk/AwxVnZaShlohlfxHHhuHTClBGwfbVNuL2aocmxlog24PFqCsuNEI1I8wt4QZmZuGMz/NHT5pJ6oiPCbHEzzmYccdKAJAb17c32igZ21zTicns5YXAfhiTHUlDmwOvV7K5pZGhyLCNS42zH8jnG8NQ4yupbaGhpY0dVI31jo2zBKyxvsCOJa04cZNlTz44qJ26vtvto234HbR4v28udjM1IYNKAJDui2lHlZFhKHCNS46zXjRSWN9hiC1Bc1cAOa2E5ZVgy8dERFFU6ybfE7JwxaSTERFBYbp4S2VXdyJiMBMZkJJhJ5/aypbSOCVmJjEyLp7KhlaufX8OaXbWMy0xkVJq5UyiqbLAFwrdIFZQ72GktxEMtO4sqnNQ1uahwtJKdHs+I1Hj7g7LiSicj0+IZlR5PoyWiufvqGJIcy6nDkml1m8XaJxCXTDKis22/w677kkkZ9IoMJ7ekjvomc1cwLCWOEWlxbK9worW2x2dkWjzFlSZtZ5UZR9+dT1Flgx1YXDrJLArbK5y2QE4akMSI1Di27XdQ7mih2tnKycP6MbhfbztPUYWT4SlxZKcn4GhxU9nQSmF5A8NS4pg8MImGFjf761vYWe0kIkxx9uhUq26nvdh8Z2J/u+69tU24PF5GpcUzMi2OwnLj61tK65k8sE/QoltY3sCItHgmZBpfzytzkFdmFqnvW8K6uaTODpTGZiaQ1acXRZUN7Kw2PjhxQCKD+sVSbM2p3TVNDEuJZcqgPrS0eRn358VUNbRy0tB+jEyLp9rZyoFGFwXWYvidCf3teeYT5yutYGXb/np2Vzfh1TA81QRA2yucuNxedlQ5GZUez5j+CRRZc3FnlZOhyWYcAXZUGj9IT4jhBCuCL6posPOflZ1GXLTxa58fnzsmjX6xUeTtd9DkcrOntolRloAXVTTYdwnHIvqGEBHwXlHh5Jc5WG3tOY3PMtFcXpmD2kYXlQ1m4gaK6NZSB2MyEuyO217RwI4qJwP79iaxVyQD+vRme2UDheUNhIcpxmQkMKhfLEWVTjuCGJZqBndvbROF5Q20tHnJTo9nVJoR1sZWN0WVToanxTEmw9STt9/BlpJ6xmYkMCEryU4rrmggTMEQa1EorGiwb/9GWRO/qNLJ/vpmWtq8Rpwsx9pe4WR7RYM1yYwYFFc6TaSeHEtGYgzxlmPllTlI6h3J2IwE+ifGBEWx2ekJjM1IoKC8gYIyE734+mh3TSMtbR52VDUyLCWWEVY92ysaKKxosIURjIjsrG4kLjqC1PhohqfFsb2iwb6NHNM/gVHp8RRVONlRafpyZFocYzMSafNoNpfUUVzpZHxmImdlp9rjXNvo4qxRqaQlRBMfE2HqthaPcZkJxEaF23UD9gJbXOlfiEelxzM8LY6iSnNLXVrXzLCUOFtEC8oa2LbfwdiMBMZmWmNWZvooLSGajKReDE2OJW+/wxbNcZmJjLTauKPaTGbfIlff3MaOKifljhZGpcczIjUOZ6ubCkcrO6sbGRrUl07yyxwMS4ljrCWC2ysa7D3XEWnxjLEE07egGl9PZNt+s/da7mhheFqc7QeF1lgGjo8v2h7YtzdJvaPISIyhuNJJYXkDEWGKcRmJDOjby7THqntYahyj0o0f5Jc5aHJ5GGstumX1LeZxu4oGstPiybbmVEGZWZBGBvj/9gozf5SCwf1iGZEab9dt/MAEIUUVTkoPNONyexmaEsspw/oxNCXW3o45aWi/AP9voKC8gZFpcaQlxJAaH822/fVsKaknqXckJw0xd2J5+x323ZgR8Dh2W4tzm0dbfRRn+/oua3yGpwYEDFYw6OvfogonxVVOEntFkhwXxaj0ePvzG4CJWUmMsbSoqMKJ1pCdbvrI7dXklznYUeVkdP94jgUhIeC+vcE31+0jKiKMYSnm9rWuqc3+4GSUJU4AH+WVU+1sZfKAJEZaTlBoraTDUmIBIyjFVgQ+JDmW6Ihw41iVTnb6JmmymShaw/xc8wn76P4J9oT8fEeNvW+a2CuSAX17sbW0nq3765mQlURKfDQp8dFs2+8gr6yBoSlxxESGMyo9nl3Vjazfc4Co8DCGpxox2FPTSIEl6oGOtSy/gla31xYIsETUipaVUozNTCC3pI68/Q7G9E8waZYYFJQ3EBmuGJoSy9iMRFxur/3h6qQBSYxKj8erzf61s9XN0JQ4BvbtTVREmB1tjEqPJzOpFzGRYfYiNzQlFqUUo9LiKbAWjz69I0mNj2ZkWjwF5Q77TmV4apw9PnNySvBqM65njEplx4MXkRIfbd0dpaKUsqOnLaX1jMkw7RmeGkdxpWl3dEQYmUm9GJ4aR02jy/5MIjs9gZGWsK7ZVYPWMCw1ltH9EwhT5oPM0rpmxmYkMiwljqiIMLuPRqUb+8ZmmK24jXvryEzqRXKcaU9hudMWvEBh9j3lMSo93h6zvLJ69tY2MTTZ+EZ6gllMC8obyO4fT1x0BJlJvcgvc1BU0UBmUi/ioiMYm5FIVUMrH+dXmsCivxHRvbVN5Ow2t/e+aBlg9c4ayh0tZPdPsH2juMI/PgDDrTtGX6QeFRFmC6vvM4+hKbFMzErEq/1PFY3NSLSDok+2V1Hb6GJUejx9Y6NISzAimruvjglZSST2iqR/YowdKGX16UVMZDgj0uLYWdVo33n4Fr5d1Y32HcHQlDiS46L5+Ndn8M6tp/C7C7MZ3K83o3yLbrkR8FFpZnzGZCSQZy1y4zMTCbP6Ka/MQUGZCciGpcQxKj0et1ezYMt+2zeGpxlfX72zhiaXh6EpcQzq15vIcGXEttJJdno8A/r2JjoizF7khqeaeTY+M5GtpQ427D1AUm8z58f0T6CwvMG+UxmVnmAHDHPW78Pt1YzPTDq8yB0lISHg4zP9kfXo9Hgiw8MYYznWPOvRpez0ePrEmmjj1S/MI0yTB/YhPiaSzKRebNhTR3GVk3FWWcNT49le2cCq4iomD0gCMKtzdSOF5U76xUaR2DuSadYHFy+t2k2YMkLkK2Ou9ZTLcGvijO2fyAdby2lyeWybfY6Vt7+eMVbkMjItHo9X89qaPUwamGQWD8uxfE9XDE81jtU7Ktx+Zn10/wQG9O1NTGQYX+ysodrpYliqmaRTBvZhc0k9uVb076t7R5WT3H11DEuJIzI8zH7v1S/20C82ioF9e9ti4HskcVxmoj0J3li7l5Y2LxMHJBEWpvyLnBX9A5wwuC91TW28u2k/o63FY1xmIo4WN8sKKogIUwzqF8uQfrH0jgrnLas9vj4KD1PccNoQfnTiQHufcERqHGt31VJQ3mBv3QxPNdtZhRVOhiTHEhambIGZnbOPxF6RpCVE29GTT1iHpcTROyqC7PQEXrMenxuTkUBkeBij0uLNHUGV0550vq24JfkVTBtsxn9UurmdX7m9itiocAb1i7UFc9bnu1EKxmUk2lHjsnzzRIxPRE8Y0peP8soprWu27wrHZ5oPA4ssgfDVDfBWzj5GpsUTExlu+43vqapxmYn0i4smOS6KVz7fDcCJQ/rSJzaK5Lgo8soc7KppZGiKKXNchhGYTfvqbFEclR7PjionW/fXkxofTUJMJJOsefC/NXuICg9jRFqcbeu8jf55Zv5P4L3c/TS0uu054tta3Fzi9/URqfG4PF7e2VhKdv94s3ikxeH2avsDcZ8fgZmzt3xrGEop0hNiSImPZnlhJVUNrXYU67uLzCtz2MHdmIwECsoayC2pY0SqCZR8dySzPttNn96RjEj1bw1+YM2zYcmxtp68snqP2U5Kjyc8zAQRW/fX23elxr4kmts8zN1QwvjMRJQyd+++R1pjIsMY2Lc3Q5JjiQoP439f7LXmSId/971bhISAJ/WOIsx68mT6EDOZs9PjUQpWFVeTZEV9YAbS930Cn7ONTItjaX4FWsNEy0lPG56M1tDS5rX3CUdbtz0fbC2zbwlTE2IYkhyLy+Nl0oAkYiLDyUiMIT0hhg+tJw98juKbfAATB/gdK7/Mwf76Fvt93/9tHs3JQ/tZZRjHmrexlKHJsSTHRRMZHsaZ2alUOFrJTDIrfXiYYsrAPvaEGppsrvN9Qm7K99WdaEcbvolr7gLC8GoTfStlxDUqPIwPt5UTGa5s+2aMTKbJ+mT/WyPN0yZj+ieQs7vWao+p5+Rhpg0ut5dp1t7h5IGmnxdtKWeINUnCwpQtCOkJMaQmxNg23/ytYdx3yTj7tU9MwOyzA4xIi6PCYb4E49ujnDa4D4m9Iimrb+Hs0SZ6nzggiYgwxdvrS1DWtpXPJt/jdD6BGZuRwBc7a3G5vfbeu69dLrfX7tfAaHvigCTCwxQp8dFkp8dT19TGpAHmjqufFZ36HqPz9fvJQ/vZTyb5AobJA5PYW9vEtv0Oe/wD90p9e80+X1y4pYzkuGjb18dlJtJojc94OzCJ452NpUG2nz06DbdXU9PosresThjchzaPZuHmMnvx6BcXTVLvSLSGkelmwe8bG0X/xBj7i0e+9mT3NwEHwLRBZixGWvv3e2ubmDLQ1D3DepKrsqGVc0enB9n6fu5+kuOi6dvBh3tKKSZmJdmPFY5KDx4fU1aSPZ7NbR4+Laq2A6xhqbHER0fQ6vZyyvBkwsLMXWh4mGJ2TgkRYYoJ1licnZ1qb9+cNsL427TBffhiZy3Vzla7j3x+qbV/Z8A3X1YVVzMyzYh/pHVnDdAvNop+cdHttrG7hISAA9x3yThuPH0Ivz5vFACx0REkW50yKi0eZT1beN2pQ4gMV5wzOpWIcNO8MQHCOsnalz51eD877bQRRpy+NSqFqPAwmlwe+xE1MBEcwP+dNQIwjuUbyJT4aNISjB2njkgmPExxdnYqQ5KDI6pAO4amxNmP3/nEb3hqHCnWxDxxqN+2C8cZp79sciZh1ip2qvWhpclrJo9vwiTERHCBdU1g3T7nDw9Tdju+ZdkQGR5GthXdmAg/HDBPfgBEhCkSe0UCcM7oNFvUzxlj+igjqZddz3WnDAb8i1pgu8E8nghwqfVkR0d8f2qW/fs46/rTAtrt67fI8DDOtkTJ90FnbHQE463JdcLgvnZ7fNH05VOz7L4e004fnTCkj734n2SNReAepq+vlTJPfhycdvboNBpa3faiG2hvdESYHYRMCVh0ffvKib0imTqoD1ERYXb/psabaNtnh8/Xf3nOSMBsf/h8fYa10AbaPnlAEslxUST2irR9w7cAQrA/zbDmwh1nj7TTfG1Ijou2hejq6eZJo7SEaAb0NeM/IWDR9bUtNSHGnkvnjzP/D7O2LRpdHs4bm2a3pz0mZvnF2mfHCYP7ktWnFwP79rYj20BR9/lLdEQ4158+BICp1vhER4TbfjQhK5G4aPNHyXxB3MC+vUmNjwnqP/AHEQP79qZfbBRKwVVWHwxJjrMDTF8QAHDX+aYPLxof+CeEexit9Zf2M3XqVN2T/OGdLXrEPYv0Z8VVQektbW7d0ua2X9c4W/WgmQt09h8+CMq3eV+dnr+pNCjt+llr9aCZC/S+2kY7reRAk37l813a6/Xaaf9aul0PmrlAz163N+h6j8cb9NrR7NJXPbdan/TgUl3X5LLTm11u/eHWsqAy/7owTw+auUC/uXZPUFseXVygqxta7LT1e2r1oJkL9KVPrQqqK2d3ja5v9tfh9Xr1GX9frs95dIUur28OylvX6Aqy9Ysd1XrQzAX6l29tDMr3zoYSvWnvgSC7x/zxA33OoyuC8q0qqtIrCyuD0mY88rEeNHOB3lXlDEp3H9RHHbG/rknn7jsQlDZo5gI9aOYCXeNstdMKyx360cUFQe355Zsb9aCZC/R7AePr8Xj1mp01QfkKyx166N0L9fWz1uk2tyfIxr01fh/QWutfvmXKXLKt3E5rdrn1wx/kB9mzorBSD5q5QN/33jY7zev16pdW7dR7qhuDrvW1J9BfvV5vkF9orfUrn+/Sg2Yu0M+sKA5KX7y1TBdVOOzX1Q0tdpmBLMjdrxdu3h+U5stXG2B7XZMryEattX519W494d7F+qVVO4PSnS1tutLh98s2t8cus9nlb8+Oygb93ModQW36y/vb9KCZC/SqouC5ezBrd9XoQTMX6CeWbD9sPq/Xq3/xxgY9aOYCnV9Wb6e73B49e93eIHsWbt6vB81coH8zZ1PQ9c+sKA7y1QONRjeG37MwqK4KR7NubG0LSvv17E367EdX6M8Oak+Fo1k3tbp1dwFydDuaGtIC7vZ4uywGFY5mvbva2Wm+baX1hzhqezS2tumP8ysOmWjdodnl1i+u2hk0mdvD4/Hq/36y4xBRbo8jsW9LSZ2ua3R1mu/djSWHiHV77K52HiLA3aWgzBG0wHVEeX2zfvLjoiBR7oiDF92OaG3z6A+2lHXqcy63Rz+4MO+QBaA91uys0VUBi/PhaHa5u+Tv724sOSSoaY91u2r0h1vLulR3V6lwNOsvdlR3mq/kQJN+7KPCLrUnMJjqjECh7gi3NX+60u/vbSrtkm4cazoScKV78CS5zpg2bZrOycnpPKMgCIJgo5Rar7WednB6yOyBC4IgCMGIgAuCIIQonQq4UupFpVSlUmprQFpfpdQSpVSR9f+xechREARB6JCuROCzgAsOSvsdsExrPQJYZr0WBEEQvkQ6FXCt9SdA7UHJlwAvW7+/DFzas2YJgiAInXG0e+BpWmvfXx0oB9I6yqiUukkplaOUyqmqquoomyAIgnCEdPtDTOsZxQ6fRdRaP6e1nqa1npaSktJRNkEQBOEIOVoBr1BK9Qew/q/sJL8gCILQw0Qc5XXvAdcCD1n/z+/KRevXr69WSh38p8mTger28ocgX5e2fF3aAdKWrypfl7Z8We0Y1F5ip9/EVEq9AZyBMbQC+DPwLjAbGAjsAa7QWh/8QWeXUErltPcNo1Dk69KWr0s7QNryVeXr0pbj3Y5OI3Ct9VUdvHV2D9siCIIgHAHyTUxBEIQQ5asg4M8dbwN6kK9LW74u7QBpy1eVr0tbjms7vtTTCAVBEISe46sQgQuCIAhHgQi4IAhCiHJcBVwpdYFSqlApVayUCtkDsZRSu5VSW5RSm5RSIfUXK75Op0120JZ7lVKl1thsUkpddDxt7ApKqQFKqeVKqTyl1Dal1B1WesiNy2HaEorjEqOUWquUyrXacp+VPkQptcbSsbeUUu3/leZjYdPx2gNXSoUD24FzgRJgHXCV1jrvuBjUDZRSu4FpWuuQ+2KCUmoG4ARe0VqPs9IeAWq11g9ZC2sfrfXM42lnV+igLfcCTq31P46nbUeC9e3m/lrrDUqpeGA95sC4nxJi43KYtlxB6I2LAmK11k6lVCSwCrgD+BUwT2v9plLqGSBXa/2fL8Om4xmBTweKtdY7tdYu4E3MKYfCl8jX6bTJDtoScmity7TWG6zfG4B8IJMQHJfDtCXksP48pdN6GWn9aOAs4G0r/Usdl+Mp4JnAvoDXJYTowGIG8SOl1Hql1E3H25geoMunTYYItyulNltbLF/5bYdAlFKDgcnAGkJ8XA5qC4TguCilwpVSmzDnPy0BdgB1Wmu3leVL1TH5ELNnOE1rPQW4ELjNupX/WtDZaZMhwH+AYcAkoAx49LhacwQopeKAucCdWmtH4HuhNi7ttCUkx0Vr7dFaTwKyMLsI2cfTnuMp4KXAgIDXWVZayKG1LrX+rwTewQxsKPO1OW1Sa11hTTov8F9CZGysPda5wGta63lWckiOS3ttCdVx8aG1rgOWAycDSUop37EkX6qOHU8BXweMsD7BjQKuxJxyGFIopWKtD2dQSsUC5wFbD3/VVx7faZNwBKdNfhXxCZ7FZYTA2Fgflr0A5GutHwt4K+TGpaO2hOi4pCilkqzfe2EewMjHCPnlVrYvdVyO6zcxrUeHngDCgRe11n89bsYcJUqpoZioG8zhYK+HUjuO9WmTXyYdtOUMzG26BnYDNwfsI38lUUqdBnwKbAG8VvI9mL3jkBqXw7TlKkJvXCZgPqQMxwS/s7XW91sa8CbQF9gI/Ehr3fql2CRfpRcEQQhN5ENMQRCEEEUEXBAEIUQRARcEQQhRRMAFQRBCFBFwQRCEEEUEXBAEIUQRARcEQQhR/h9ew0JMteyKMQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "csv = gwf.head_obs.output.obs(f=\"head-hydrographs.csv\").get_data()\n",
    "\n",
    "for name in csv.dtype.names[1:]:\n",
    "    plt.plot(csv[\"totim\"], csv[name], label=name)\n",
    "plt.legend();"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-30T20:16:44.263945Z",
     "iopub.status.busy": "2022-07-30T20:16:44.263597Z",
     "iopub.status.idle": "2022-07-30T20:16:44.274135Z",
     "shell.execute_reply": "2022-07-30T20:16:44.273473Z"
    }
   },
   "outputs": [],
   "source": [
    "try:\n",
    "    # ignore PermissionError on Windows\n",
    "    temp_dir.cleanup()\n",
    "except:\n",
    "    pass"
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "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": 4
}