{ "cells": [ { "cell_type": "markdown", "id": "94a4daa2", "metadata": {}, "source": [ "# Example 6: Temporal Resolution Demonstration and Benchmark\n", "Model simulations depend on the temporal resolution used. A model with monthly time steps of tracer input behaves differently than a model with yearly time steps. Here we demonstrate the differences and benchmark the use of the time step size parameter for the model (`dt`) to convert units." ] }, { "cell_type": "code", "execution_count": 1, "id": "b25444c9", "metadata": {}, "outputs": [], "source": [ "import PyTracerLab.model as ism\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "from datetime import datetime" ] }, { "cell_type": "markdown", "id": "ac5cff6c", "metadata": {}, "source": [ "## 1. Load (Synthetic) Observation Data\n", "See Example 4 on how this data is generated." ] }, { "cell_type": "code", "execution_count": 10, "id": "7456a1bb", "metadata": {}, "outputs": [], "source": [ "# load input series on monthly time steps\n", "# this would be the tracer concentration in precipitation or recharge in a\n", "# practical problem\n", "file_name = \"benchmark_input_monthly.csv\"\n", "data_m = np.genfromtxt(\n", " file_name,\n", " delimiter=\",\",\n", " dtype=[\"" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "### plot input series, output series, and observations\n", "\n", "# get observation timesteps\n", "timesteps_m = [t.year + t.month / 12.0 for t in timestamps_m]\n", "timesteps_y = [t.year for t in timestamps_y]\n", "\n", "# create figure\n", "fig, ax = plt.subplots(1, 1)\n", "# plot input series monthly\n", "ax.plot(\n", " timesteps_m,\n", " input_series_m,\n", " label=\"Tracer Input (m)\",\n", " c=\"grey\"\n", ")\n", "# plot input series yearly\n", "ax.plot(\n", " timesteps_y,\n", " input_series_y,\n", " label=\"Tracer Input (y)\",\n", " c=\"black\"\n", ")\n", "ax.set_xlabel(\"Time\")\n", "ax.set_ylabel(\"Tracer input flux (per month)\")\n", "ax.legend()\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "68d49bd6", "metadata": {}, "source": [ "## 2. Simulation" ] }, { "cell_type": "code", "execution_count": 21, "id": "60800290", "metadata": {}, "outputs": [], "source": [ "t_half = 12.3 * 12.0 # monthly half life\n", "lambda_ = np.log(2.0) / t_half\n", "\n", "### define model monthly\n", "# time step is 1 month\n", "m_m = ism.Model(\n", " dt=1.0,\n", " lambda_=lambda_,\n", " input_series=input_series_m,\n", " target_series=None,\n", " steady_state_input=0.,\n", " n_warmup_half_lives=10\n", ")\n", "\n", "# add an exponential unit\n", "# define the initial model parameters\n", "em_mtt_init = 12 * 1 # 1 year\n", "m_m.add_unit(\n", " ism.EMUnit(mtt=em_mtt_init),\n", " fraction=1.,\n", " prefix=\"em\"\n", ")\n", "\n", "### define model yearly\n", "# time step is 1 month\n", "m_y = ism.Model(\n", " dt=12., # 1 year is 12 months, so we have to adjust the time step\n", " lambda_=lambda_, # because we change the time step, we leave the lambda unchanged\n", " input_series=input_series_y,\n", " target_series=None,\n", " steady_state_input=0.,\n", " n_warmup_half_lives=10\n", ")\n", "\n", "# add an exponential unit\n", "# define the initial model parameters\n", "em_mtt_init = 12 * 1 # 1 year but in months again\n", "m_y.add_unit(\n", " ism.EMUnit(mtt=em_mtt_init),\n", " fraction=1.,\n", " prefix=\"em\"\n", ")" ] }, { "cell_type": "code", "execution_count": 22, "id": "fd5c4fba", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "test\n", "test\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAG2CAYAAACXuTmvAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjUsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvWftoOwAAAAlwSFlzAAAPYQAAD2EBqD+naQAAaKhJREFUeJzt3QeYU9X2Pv4tgtShd4Z2Lx3pCgKKgEhHUFSKV6qoCAgC0kSKSBUQEQRFKTZABTtFQKo0pfdyLyC99yYl/+dd/9/O9yRkZlJOknNy3s/zHGcmkzmTCeOclbXWXvs+l8vlUkREREQOkizaD4CIiIgo0hgAERERkeMwACIiIiLHYQBEREREjsMAiIiIiByHARARERE5DgMgIiIichwGQEREROQ4DICIiIjIcRgAERERkeNYIgCaOHGiKlCggEqVKpWqVKmSWr9+fYL3nTt3rnrooYdUxowZVdq0aVXZsmXVF1984XGfNm3aqPvuu8/jqFu3bgR+EiIiIrKD5NF+ALNnz1bdu3dXkydPluBn3Lhxqk6dOmrPnj0qe/bs99w/c+bM6q233lLFihVTDzzwgPrll19U27Zt5b74Og0Bz7Rp09wfp0yZMmI/ExEREVnbfdHeDBVBz8MPP6wmTJggH9+9e1flzZtXdenSRfXp08evc5QvX141aNBADRkyxJ0BunDhgvrhhx/C+tiJiIjInqKaAfrnn3/Uhg0bVN++fd23JUuWTNWqVUutWbMmya9H7Pb7779LtmjkyJEen1u2bJlkhTJlyqRq1qyp3n33XZUlSxaf57l586YcGoKwc+fOyf1RPiMiIiLrQ1xw+fJllTt3boknLBsAnTlzRt25c0flyJHD43Z8vHv37gS/7uLFiypPnjwStNx///3qo48+Uk8++aRH+euZZ55RBQsWVP/9739Vv379VL169SSowv29DR8+XA0ePNjkn46IiIii4fDhwyo+Pt7aPUDBiIuLU5s3b1ZXrlxRS5YskR6if/3rX6p69ery+ebNm7vvW6pUKVW6dGn173//W7JCTzzxxD3nQwYK5zAGWPny5ZMnMH369BH6qYiIiCgUly5dkjYaxAlJiWoAlDVrVsnInDx50uN2fJwzZ84Evw5prUKFCsn7WAW2a9cuyeLoAMgbgiN8r/379/sMgNAg7atJGsEPAyAiIiJ78ad9JarL4LGKq0KFCpLFMfbf4OPKlSv7fR58jbGHx9uRI0fU2bNnVa5cuUJ+zERERGR/US+BofTUunVrme1TsWJFWQZ/9epVWdoOrVq1kn4fZHgAb3FflLQQ9MybN0/mAE2aNEk+j7IY+nmaNm0qWST0APXq1UsyRsZl8kRERORcUQ+AmjVrpk6fPq0GDBigTpw4ISWtBQsWuBuj//77b49ObgRHr732mmR1UqdOLfOAvvzySzkPoKS2detWNWPGDFkKj07w2rVryxJ5zgIiIiIiS8wBsmoTVYYMGaQZmj1AROQEWJF769ataD8MokSlSJHC52ruYK7fUc8AERFR9OA1MLLvyJgT2QG2wkKLS6hz+hgAERE5mA5+MDg2TZo0HP5Klg7Wr127pk6dOiUfh7qwiQEQEZGDy146+EloUj6RlaD3FxAE4fc2sXKYLXaDJyKiyNM9P8j8ENmF/n0NtWeNARARkcOx7EVO/H1lAERERESOwwCIiIiIQoIdHIoXLy59ZcH6559/VIECBdRff/2lIoEBEBER2ar8kdgxaNAgZXXnzp1T3bp1U/nz55ctoTCwt127djL4N1D4mX/44YewPM4CBQrI7gz+wI4L/fv3D6kpGc9Fz549Ve/evVUkMACiJGHZIRGRFRw/ftx94OKMYXfG23ABNS6bvn37dlQeJ7IZCQU/jzzyiFq8eLGaPHmybNI9a9Ysefvwww+r//3vf8puVq1aJdtOYQuqUL3wwgtyvh07dqhwYwBEiVq/fr1M1ezXr1+0HwoRkQzA0wf+NiEDoj/evXu3iouLU/Pnz5eNtrH9kb44N27cWLZYSpcunQQaCECMsLckMg958+aVr8P+kZ999pn789u3b1f16tWTr8d5XnzxRXXmzBn356tXr646d+4smZ2sWbMmuPfkW2+9pY4dOybfH+fLly+fqlatmlq4cKFMOe7UqVOiGRhsF6WzXPg8PP300/I86I/xedzv448/lp8Hq6aef/55mY5sfLzdunXzOHeTJk1UmzZt3J8/dOiQeuONN9zZtYQggHvyySdVqlSp3LfpxzB16lT5GfG8YRsrlMhGjRol/15Yxj506FCPc2XKlElVrVpVzhluDIAoUUuXLpVXUHPnzo32QyGiCEDWBNmLSB9m7srUp08fNWLECLVr1y5VunRp2SS7fv360qeyadMmVbduXdWoUSOPkhM23p45c6YaP368fB2CB1y0AbOSatasqcqVKyf9Kdiv8uTJkxJUGGEPSpRx/vjjD8nueLt7965c2JHlQADgPd8GAQICIWSJ/PHnn3/K22nTpkn2S38MyCh988036ueff5bHi58b5/fX3LlzVXx8vHrnnXfc2bWErFy5UjYp94bAE8Eovj+eWwSUDRo0kL08ly9frkaOHClls3Xr1nl8HTZGxznDjYMQKVH4BYY9e/bIqwe84iKi2IXZKsOHD4/49+3bt68ED2bARRsZCS1z5syqTJky7o+xOfb333+vfvrpJ8na7N27V4KFRYsWqVq1asl9/vWvf7nvP2HCBAl+hg0b5r4NmQ1kV/C1RYoUkdsKFy4s2Y2EYONvBFNoFvYFtyMQRPCCICAp2bJl89gawujGjRvq888/V3ny5JGPP/zwQwk+xowZc899fcFzhn4eZNSSuj8yRehj8hXw4XnCOUqUKKFq1Kgh15J58+bJJudFixaVIAgvtCtVquT+OpwL5ww3ZoAoUcZ69IYNG6L6WIiI/OGdjUAGCL1BCDAQLCCzgyyPzgBt3rxZLvaPP/64z/Nt2bJFLtL4On0UK1bM40UioOzmj0jsQY6ykw5+oHLlyhKQIAAx2/Xr1z3KXxpKcgh+NJQOEQgh+DHepre2MGbDItF7ygwQ+R0AIfWLNDARxS70oSAbE43va5a0adN6fIzgB9md0aNHS28PLrDPPvusu1FZb6+QEARQKJkhW+HNuB+V9/f1lbFBAIbgyxfcjl4bPEZAoOAdLIU6/VhLZuK50fN0/vz5JP9N8bP5ug2BmRFKgDq7FU4MgChB+J/BWCM31peJKDbhgmRWKcoq0JOD5l40C+uA5uDBg+7PlypVSi7C6EvRJTCj8uXLqzlz5khGI3ny5CEFHegb+uqrr6RMZywtIYvy0UcfSfM0yk+AIMDYe3Pp0iV14MABj3MioPA1ewd/u9FsrUtTa9eudZedfJ37zp070uiNMpWG3wN/5vqgPLhz505lFjwOnDPcWAKjBOF/IOMvPwMgIrIj9OagqRelLpSzWrZs6ZF1QGDTunVrmcWDmToIMpYtWyZ9QYCVWchKtGjRQv4OouyFZuW2bdsGPPgPfUQIfNCjhAbhw4cPqxUrVkjggxedEydOdN8XGfcvvvhCGoK3bdsmj9F7zg4eO5q7T5w44ZGFQUkK98fPi69//fXXJfjSQVfNmjXVr7/+KgdWz3Xs2FH6k7zPjcd29OhRjxVv3vDYsdrOLHi8tWvXVuHGAIiSLH+h0Q/QlIYmPiIiOxk7dqwsr65SpYqUsnDBRlbHaNKkSVIWw0op9Pd06NBBXb16VT6HLAqySAh2cGFGxghLyFHOMvaz+CNLliySjUGm5ZVXXlH//ve/JTDBWwRXxuZrlCLRl9SwYUNpYMYyddzPCE3NKO/h77Qxa4Iy2jPPPCOr3/CYsRoOGSatXbt2EiBh9Ru+B76vMfsDyFIhU4bvmVhJCqvaMLfHjP6iNWvWyIIb/FuE232uSHRj2QzSjFjthH8EDNlyKizjxKsC/MHASgfdvY/ZFURkf1gphGxHwYIFfTaxkj1hBg8yWch4Rcqbb74p106MDwhFs2bNZMVeYrPnEvu9DeT6zQwQJZkBwisDvaqCZTAiIvI14BFbe3g3NAcCTenIrmH4YiQwAKIkAyCkPzE5FSK1SR0REdlHxowZJWsTaEnQCE3XGIyY1Ko8szAAogTp+RbIAOkACBkgVk2JiKxdAotk+cuuGACRTwhyjCUw7OmC1QdYaYAVAURERHbGAIh8Onv2rDSTARrNsJleyZIl5WOWwYiIyO4YAJFPOvuDUeq6y95YBiMiIrIzBkDkk7H8pXElGBERxQoGQJRoA7Rx6JZxJRgboYmIyM4YAJHfGSDMZ8AyRYxbN26SSkREZDcMgCjJGUAagh9M6ASWwYiISMN+ZMWLF/d7bzRsnhofH+/ebiQaGABRkjOAjDgQkYiivVt9Ygdm4FgdNlbFXmKYnIwXlthrDHtzYQPqQOFnxrYX4VCgQAE1btw4v+7bq1cvGWLovVlrQkqUKKEeeeQR2actWhgA0T1u3rypjhw5kmgAxAwQEUXD8ePH3QcuztjvyXhbz5493fdFr+Lt27ej8jixrUNCwQ8u/IsXL5b9Fvfv369mzZolb/H31Y7tBatWrZIXzU2bNg3o69q2bSub0Ebr34gBEN0Du77jD0e6dOnu2QFYrwTbsGGD36lOIiKz5MyZ031g00tkQPTHu3fvVnFxcWr+/PmqQoUKKmXKlO6Lc+PGjVWOHDnk7xoCDQQg3i/8evfuLbuq4+uwm/pnn33m/vz27dtlI2h8Pc7z4osvqjNnzrg/X716ddW5c2fJ7GTNmlV2nE9oz6xjx47J98f58uXLp6pVq6YWLlyoUqRIoTp16pRoBgZDaXWWC5+Hp59+Wp4H/TE+j/thY1L8PJjjhh3nsUGo8fF269bN49zYbb5Nmzbuz+NagH25dHYtIQjgnnzySffIFOwgjy0xvCsF+FmM+4XhaxAQLl++XEUDAyBKtPzl/UuPGm/atGmlbos/NkQUW/DiB/9/R/owc2Vpnz591IgRI9SuXbtU6dKl1ZUrV1T9+vWlT2XTpk2qbt26qlGjRh4lp1atWqmZM2eq8ePHy9cheECwAxcuXFA1a9ZU5cqVk4v6ggUL1MmTJyWoMJoxY4aUtP744w/J7njDhR/BwgsvvCABmxH2v3rttdckEEJQ4A+diZ82bZpkv4yZeWSUvvnmG/Xzzz/L48XPjfP7a+7cudKj884777izawlZuXKl+8UxIBCrVauWPC4jfIwAS+8XhucKgRq+PhqSR+W7ku1WgGmo75YvX15+YfGHQE+HJqLYcO3aNfeFP5IQpODFlRlw0UZ2QcucObN7AQcMGTJEff/99+qnn36SrM3evXslWFi0aJFcuL3//k2YMEGCn2HDhrlvmzp1qmRX8LVFihSR2woXLqxGjRqV4OM6ffq0BFN4IekLbkcgiOClYsWKSf6cOkOPjUi9A6obN26ozz//XIbZwocffqgaNGigxowZc899fcFzhr/3yKgldX9kitDHZPTSSy+pV199VXp8kFHbuHGj2rZtm/rxxx897oevw9dHAzNA5NcKMCMORCQiKzNmI3Rwhd4gBBgIFhDgIcujM0DYOBQX+8cff9zn+bZs2aKWLl0qX6ePYsWKeWTMAWU3f0RijhpKazr4gcqVK0sGas+ePaZ/r+vXr7vLX8ZyGp5TBJowffp0VaNGDXeZzpj5QtAdDcwAkd8rwDQ2QhPFLvSLIGCIxvc1i3cmCcEPsjujR4+W3h5cdJ999ll3ozI+TgyeD5TMRo4cec/ncuXKleD39ZWxQQCG4MsX3I62AzxGQKnIO1i6deuWMkMyE8+NnifMhzNCeQtlRZS9nnnmGfX111+rDz744J6vRbkvoRfb4cYAiAIqgRkDILwqwh8Q/KITUWzABdisUpRVoCcHvSdoFtYBDRp1jUNekR1BM64ugRmh7D9nzhzJXiRPnjykoAN9Q1999ZWU6YylJWRRPvroI2meRvlJB0zG3htsUH3gwAGPc6Jx2teCFGS30GytS1Nr166V71+0aFGf575z5440eiNLo+Fvuz+LXVAexFwfbyiDPfjgg/JzYaUXAiFv+J4IRqOBJTDygFcESZXAcDtexWDVBH55iYisDL05aOpFqQsv3Fq2bOleiQQIbFq3bi2zeDBTB0HGsmXLpC8IsDILmYoWLVpI5htZcjQrYxl3oKth0UeEwAc9SlitdvjwYbVixQoJfJCBmThxovu+aLz+4osvpOcS/TN4jN5zdvDY0dx94sQJjywMSlK4P35efP3rr78uwZcOumrWrKl+/fVXObCgpWPHjtKf5H1uPLajR496rHjzhseO1XbeUHLEkn+srsNz551pQxCKc/sKOiOBARB5OHXqlKzIwKtALFf0BZ9jHxAR2QUacTNlyqSqVKkipSxcsJHVMcI8GmQisFIK/T0dOnRwTylGFgVZJAQ7tWvXlowRlpDjhaBe0eSvLFmySDYGmZZXXnlFXlAiMMFb/D01Zt779u0rfUkNGzaUBmb01Xi/MEVTM8p7aMhGJkZDGQ0ZF6x+w2PGajhkYrR27dpJgIQyFb4Hvq8x+wPIUiFIwff0HolihFVtO3bs8Nlf1L59e6kU4Pt5w6o7PLaErjXhdp+Lu1reA2lGzJfAzAQM2XKSNWvWyB8JNNAl1pnfr18/NXz4cElxTpkyJaKPkYjMgZVCyHYULFjwniZWsi/MAUImCxmvSHnzzTfl2onxAUZYcfftt9+qrVu3etyOoAiZOfQGVa1a1bTf20Cu35bIACHlh1QbfpBKlSqp9evXJ3hfpDGRfUDkjTo1ZgggRWiEmG7AgAHSnIaUG9Jr+/bti8BPEpu7wPvCRmgiIjIOeDQOOUSfFVokMEKgS5cuylePEl5IBxr8mCnqAdDs2bNV9+7d1cCBA2VOAGY1ID2JUowvaA7DE41MBSJK1GBxoB6rYQ4DhllhENW6deskUMI5ETVSaA3Qmi6B4Rc8WksYiYjIGjJmzCgBjS4JYr4SxgJgorSv8hdKdCgBRpUryipWrOjq1KmT++M7d+64cufO7Ro+fLjf5yhXrpyrf//+8v7du3ddOXPmdL333nvuz1+4cMGVMmVK18yZM/0638WLF1EWlLdO07p1a/nZhw4dmuj98DznyJFD7rt69eqIPT4iMs/169ddO3fulLdEsfB7G8j1O6oZINQAsaeUsQMc0SM+RoYnKSh1ofsdjVfYSwVQF0Q3vPGcqAeitObPOZ3O3xIYGqFZBiMiIruK6hwgLKtDVz02ljPCx4ntM4XmJky4xDJsLAlEZ7see47gR5/D+5z6c95wHhzGJiqn8rcEpstgv/zyCwMgIpvjWhhy4u+rLQchYm8SdLejyQoZIPQQ4YKNWmMwsJpp8ODByukwiAuDs8CfyZw6A+S94y8R2QOG6AH6+JKahkxkFbrvVP/+2jIAwvhsZHCwq64RPk5s8zWUyfSocKwCw/hwBDEIgPTX4RzGEeX4GPf1BbMWEEQZM0CYqeA0esIoSoaYmZEU3QiNEiSeM6eNDCCyO/z9RfOqXnSC7ShQ3iayauYHwQ9+X/F76z0U0lYBEMZso0scWRwMeAIsocPH6CD3F75Gl7AwFwBBEM6hAx5cnLEaDJMufcFOtTiczlj+8uePYPbs2WVeEJYzopfLe4gWEVmfftGY0MpbIqtB8OPPjvaWL4Eh84JplMgmVKxYUY0bN06mb2JpO2BKJfp9kOEBvMV9UaJB0DNv3jyZA4QpnoALNyZ0vvvuuzJkCQHR22+/LZM8dZBFviW1BUZCZTAEQCiDMQAish/8zUS2HC9ozNpokyhcUPYKNfNjmQCoWbNm6vTp0zK4EE3KyNosWLDA3cSMi6tx1DiCI4wqP3LkiNSsMbL8yy+/lPNovXr1kvu9/PLLsrfJo48+KufkpNPQdoFPKADCJoFshCayN1xUzLqwENkBt8LwwalbYWCPHKzqwgBJfwdUodSIkQOY5O29SzEREVEk2W4rDLJvCQw9XIAN8xLbLZiIiMhKGACRu5E8kBlAxma0IkWKyPtcDk9ERHbBAIgE+q+wVxp6AAIdAaCXw7MPiIiI7IIBEHk0QGM330CHS3EgIhER2Q0DIBLBlL807glGRER2wwCIQg6AMLoAowqOHz+ujh49GoZHR0REZC4GQBTQLvC+pE2bVpUsWVLeZxmMiIjsgAEQhZwBApbBiIjIThgAkSkBEFeCERGRnTAAInXlyhV18uTJoEtg3ivBOFyciIisjgEQubewyJw5s4wQD0apUqXUAw88oM6dO8ctMYiIyPIYAFFQW2B4S5kypSpdurS8zzIYERFZHQMgCmoXeF84EJGIiOyCARCF3ACtcSUYERHZBQMgMqUEZlwJtmHDBnXnzh1THhsREVE4MAAi00pgxYsXV2nSpJFVZXv37jXp0REREZmPAZDDIVNz8OBBUwKg5MmTq/Lly8v7LIMREZGVMQByOOzd9c8//8gO8PHx8SGfjwMRiYjIDhgAOZzu/ylQoIC6//77Qz4fV4IREZEdMAByOLNWgHkHQJs3b1a3bt0y5ZxERERmYwDkcKHsAu8LzoNp0jdu3FDbt2835ZxERERmYwDkcGZngJIlS+buA2IZjIiIrIoBkMOZHQABByISEZHVMQByOLNLYMCVYEREZHUMgBzs4sWL6uzZs/J+wYIFTc8AoQfo+vXrpp2XiIjILAyAHOzAgQPyNlu2bCouLs608+bNm1dlz55d3b59W23ZssW08xIREZmFAZCDhaP8Bffddx/LYEREZGkMgBwsHA3QGgciEhGRlTEAcjCzdoH3hSvBiIjIyhgAOZhZu8D7oktgu3fvVpcvXzb9/ERERKFgAORg4SyB5ciRQ5qhXS6X2rhxo+nnJyIiCgUDIIfCCq1Dhw6FrQQGLIMREZFVMQByqMOHD0sQlDJlSpUrV66wfA+uBCMiIqtiAOTw8hcGIGL/rnDgSjAiIrIqBkAOFa4ZQEYVKlRwB1t64jQREZEVMAByqHA2QGuZMmVShQoVkveZBSIiIithAORQkQiAgGUwIiKyIgZADhWJEhhwJRgREVkRAyCHilQGiCvBiIjIihgAOdD58+fVhQsX3KvAwql8+fKyyuzYsWNyEBERWYElAqCJEyeqAgUKqFSpUqlKlSqp9evXJ3jfKVOmqMcee0wabHHUqlXrnvu3adNGdiQ3HnXr1o3AT2Kv8hfm/6RJkyas3ytt2rSqRIkS8j77gIiIyCqiHgDNnj1bde/eXQ0cOFC2TChTpoyqU6eOOnXqlM/7L1u2TLVo0UItXbpUrVmzRrZbqF27tjp69KjH/RDwHD9+3H3MnDkzQj+R9UWq/KWxDEZERFYT9QBo7NixqkOHDqpt27aSKZg8ebJkJaZOnerz/l999ZV67bXXVNmyZVWxYsXUp59+qu7evauWLFnicT9MOM6ZM6f7QLaIohMAsRGaiIisJqoB0D///KM2bNggZSz3A0qWTD5Gdscf165dU7du3VKZM2e+J1OUPXt2VbRoUdWxY8dEB/HdvHlTXbp0yeOIZZFaAeZrKTw2RyUiInJ0AHTmzBl1584d2TncCB+fOHHCr3P07t1b5c6d2yOIQvnr888/l6zQyJEj1fLly1W9evXke/kyfPhwlSFDBveBslosi3QGqHTp0ipFihQShB48eDAi35OIiMjSJbBQjBgxQs2aNUt9//330kCtNW/eXD311FOqVKlSqkmTJuqXX36R8guyQr707dtXXbx40X1go1AnBECRygChHIkgCFgGIyIi5fQAKGvWrOr+++9XJ0+e9LgdH6NvJzGjR4+WAOi3335zX1wTgkwHvtf+/fsTvECnT5/e44hVKDv+/fffEc0AASdCExGRlUQ1AHrggQdkw0xjA7NuaK5cuXKCXzdq1Cg1ZMgQtWDBAvcKo8QcOXJEyi9Y9u10CH7wHKPR3Lv0GE5cCUZERFaS3J87jR8/PuATY1VXXFxckvfDEvjWrVvLBbJixYpq3Lhx6urVq/L10KpVK5UnTx7p0wH09AwYMEB9/fXXMjtI9wqlS5dOjitXrqjBgwerpk2bShYJDb+9evWSTTmxvN7pdAM0sj+YjxTpDBCa3hGAodmdiIjI0gFQt27dVHx8vJSr/IEemoYNG/oVADVr1kydPn1aghoEM1jejsyOzk4gY2G8WE6aNEnKOM8++6zHeTBHaNCgQfIYt27dqmbMmCHTjtEgjTlByBih1OV0kW6A1jDiIHXq1Ory5ctq7969MsKAiIjI0gGQ7t3AsnJ/+BP4GHXu3FkOX7wbl5NaRYSL7MKFCwP6/k4SrQAoefLkqly5cmr16tVSBmMARERE0eRXHQLZFZSX/NWvX7975vKQM2cAGXEgIhER2SoDhAAoEFhWTtYUrQwQcCUYERFZRUidqFiGrncVJ+vDFOZoBkB6JdimTZtkejcREZEtA6Bhw4apc+fOmfdoKOyTt9GEjNVfWEEXaYULF5YZSzdu3FA7duyI+PcnIiIyJQDivk72orM/GCtgnJwdKVjNp7NALIMREVE0cRiLg0Sz/KVxICIREdlqGbwvO3fulDk7ZA/RXAGmcSUYERHZPgOEXdP9HY5I0WeFDJAOgLZt2ya9QMHCNGlM/SYiIgpbAISZPmig9Ve+fPnUoUOHgnpAFNsBEH43sDHt7du31ZYtW4I+z/fff6/Gjh0rU8SJiIjCUgLDUvf58+erDBky+HVSbDx6586dgB8MxX4JDCvQkAXC7xPKYJUqVQrqPNg2BU34Bw4cUNmyZTP9cRIRUWzzuwcIG5aSfaHcdPTo0ahngEAHQKGsBLt586a81ZvhEhERmR4Aod+C7A17qCFjgi1NUIKKJjNWgjEAIiKiUHAZvMP6f1D+QhkqmnQj9K5du2QwY6AQyP3zzz/y/qlTp1huJSKigDEAcggrNEBrOXPmVPHx8RLIYFuMQOngBxD8BNKgT0REBAyAHMIKDdBmlcF0+UtjGYyIiALFAMghrJQBCnUgIgMgIiIKFQMgh7BqABTMSjAGQEREFJWtMLAqbP/+/dKA6r1CrFq1aiE/KDIXem2MTdBWUKFCBXdp7ty5czJsM9AACFPI0QOkZwJFu7mbiIhiOABau3atatmypUx69t4NHhcgrsixnpMnT6pr167JbuyYxGwFCHgQjCEAQhaodu3aAQdAaKY+fvy4zDi6ePGiypgxYxgfMREROboE9uqrr0oD6/bt2+WV+/nz590HPibr0dkf7N32wAMPKKsItgymA6A0adK4p0CzDEZERGENgPbt26eGDRumihcvLq+4sT2G8SDrsdoKsFBXgukACMFcrly55H0GQEREFNYACHs3of+H7MNqDdChrgTTc4BSpkypcuTIIe8zACIiorD2AHXp0kX16NFDLjilSpVSKVKk8Ph86dKlAz0lOTQAKl++vPSNYY8y9PLobI6/GSAEQOgDAgZAREQU1gCoadOm8rZdu3bu23AR06tw2ARtPVYtgWFfMpRSd+7cKX1AjRo1CrgEpgMgNEFfv35dpU6dOqyPmYiIHBoAHThwIDyPhByXAdJlMARAKIMFGgAhA5QqVSrpRbtw4YJkgQoWLBjmR0xERI4MgPLnzx+eR0JhgeXvKC9ZOQCaMWNGQCvBjD1AgNIZAyAiIgr7IESUVMaNGye7eUOJEiVU165dLVdiof/L2CFLEsiwwWisBPN3mKExAwRohMbvIvuAiIgobKvAFi5cKAHP+vXrpeEZx7p161TJkiXVokWLAj0dObj8BWXKlFHJkyeXHd0xXNMf3gEQl8ITEVHYM0B9+vRRb7zxhhoxYsQ9t/fu3Vs9+eSTAT8ICn8DtFUDIPTwIIjeuHGjlMEKFCgQcACkG6FPnz6tbt26dc/KRCIiopAzQCg1tG/f/p7bsSoMzaxkLVbbA8yMeUDeAVBcXJxKmzatlNCw7QcREZHpARC2Hti8efM9t+O27NmzB3o6cngJLJiJ0N4BEPqGcufOLe8fO3YsbI+TiIgcXALr0KGDevnll+XCWqVKFbntjz/+UCNHjlTdu3cPx2OkGJwB5CsDtGHDBnX37l3ZtDUht2/fds+aMu5rhgAI27QwACIiorAEQG+//baUHMaMGaP69u3rvvgMGjRIvf7664GejsIIwYReBWblDBAa6NELdOnSJQliihYtmuQSeGMGCHQGCFOliYiITC+BodyAJugjR47I9F0ceB/L4P1ZwkyRg/k/KBdhlRV2grcqPL5y5cr5VQbT5S80OhszRXny5JG3WE2m70NERGRaAGSETBAOsnb5C8MrEWRYmS6DJTUQ0bv/R0MTdIYMGeR9PfiRiIgoIcn93bRyyZIlKlOmTPJKPbFMD5YzkzXYoQE60JVgxn3AvKEMhowk+oD8WU5PRETO5VcA1LhxY/crbrzPUpc92CkA0ivBNm3aJI3OCWWsvLfB8A6AMKaBjdBERGRKADRw4ED3+2h2JnuwwwowrUiRIip9+vTSCI15UhiOGEgJDNgITUREYesBQjbh7Nmz99yOzSjtkGlwEjtlgNDQXKFChSTLYP4EQPhdxCawREREpgVABw8edM9h8b4wYTVYMCZOnCg9G1gKXalSJdlnLCFTpkxRjz32mPQj4ahVq9Y998dE4AEDBsgeUalTp5b7YHm109gpAPJ3IGJiARB+f7JkySLvswxGRESJ8Xtp0E8//eSxIapecQMIiNAkXbBgQRWo2bNnywDFyZMnS/CDXebr1Kmj9uzZ43Oy9LJly1SLFi1kCCMueBjAWLt2bbVjxw73UuhRo0ap8ePHqxkzZshjwuwinBOlFXyNE1y+fFmdOnXKVgGQPyvBEmuC1lkgZCgRABUqVChMj5SIiGzP5af77rtPjmTJkrnf18cDDzzgKlKkiOvnn392BapixYquTp06uT++c+eOK3fu3K7hw4f79fW3b992xcXFuWbMmCEf371715UzZ07Xe++9577PhQsXXClTpnTNnDnTr3NevHjRhacGb+1qy5Yt8jNkyZLFZRcHDhyQx5wiRQrXjRs3fN5n/vz5rkGDBrkWLVrk8/Nr1qyRz/v7b01ERLEjkOt3skCmCuPIly+fZBb0xzjwqhwZm4YNGwYUfGFFD7Y/QInK2AuCj9esWePXOdDrgR3AM2fOLB9j8vGJEyc8zolsFbJLCZ0Tjx/Nt8bD7qy+C7wvmFeEEhb+Pbds2RJwCQzYCE1ERGHpAUKAkTVrVmUGTO1F+SxHjhwet+NjBDH+6N27t1z0dMCjvy6Qcw4fPlyCJH1YeWpyLO0C7w3jFZIqgyW2DB7Q94XzXLlyJSYCWSIiCo+gxgNfvXpVLV++XP39998eezNBJPcDGzFihJo1a5b0BYXS24M9zYwbueLCafcgyG4N0BoCoAULFiTYCJ1UBghbZKB37OTJk9IHhKX1REREIQdAGFRXv359KT0hEELpCZmcNGnSyIUnkAAImaT7779fLlZG+DhnzpyJfu3o0aMlAFq8eLHHzBj9dTgHsgHGc5YtW9bnuXAxTeiCald2LIEhmE5qT7CkAiDAvzv+vVEGK1asWJgeLREROaoEho1QGzVqpM6fPy9LzNeuXasOHTokM1wQlAQCK3nwdVhBpqGnCB9Xrlw5wa/DKq8hQ4ZIpkAvndaw6gtBkPGcyOisW7cu0XPGGruVwFAKxTgErNQDTHRGGSvQVWCgVwNyKTwREZkWAG3evFn16NFDmpWRvcEFCeUiBCX9+vUL9HRSesJsHyxZx0WvY8eOkllq27atfL5Vq1ZSotKw7B3L2qdOnSqzg9DXg0NfLNH/0a1bN/Xuu+/K0v1t27bJOdAn1KRJE+UECCYwr8lOGSDds4NtMJDBQSCMbGMwGSAdACEDhJlQREREIZfA0GOB4AdQ8kIfUPHixaV5+PDhw4GeTjVr1kydPn1aBhcikEGZCpkd3cSM8+vvB5MmTZJSybPPPnvPdh16m45evXpJEPXyyy/LVOBHH31UzumUGUAYSImVVMiS6GDA6oy9ZNgWAzu6owyGoZe+7pdYAITfHfyeIljC75aveVJERORsAQdA6NHAhalw4cLq8ccfl8AFPUBffPGFevDBB4N6EJ07d5bDFzQ4G+nMRmKQBXrnnXfkcCJd/kKGDFk6O9CZHdBBm3cfELI5/mSAEDDjHPhdQVDOAIiIiEIugQ0bNszdXDx06FDZjgJlK7zS/uSTTwI9HYWBHVeAGTNAGTNm9LkUHlktXdJKqmldr+ILdnsWIiKKbQFlgHDxwatpnenB+ygtkbXYaRd4XxmgdOnSydv9+/dLsz2CbO/7oMSVmPj4eHkbTFmWiIhiX7JAAyDsr8SLirXZPQOEkQo6g2PMAhn7f1Dm9CcAwr5g3BmeiIhCCoDQW4HeH1xUyLrsGAAZszvG7JUxAPKn/8cYROmJ5SyDERFRyD1AGD745ptvqu3btwf6pRQhdiyB6eyObtrWjcvGRuhAAiBgGYyIiExbBYaZOigplClTRpZZYxii0blz5wI9JZkIy/71vwGGQtqFDm6wISoyWHoLi1ACIJTRMLeKGSAiIgo5AHr//feT7L+g6Je/kEHRzcR2ygBhYCUmOGfLlk1+zxC8YD4UpnsHEwDpgYgYrGicJ0VERM4WcADUpk2b8DwScuQWGL4anFG6unHjhmSDMMsHfUANGzb0axsMI/QAYfglzuW9NxwRETlbwC+J0aNx6tSpe25HY7Rdhu7FMjs2QIMxuNG9OwiAjGWwQDNAyCCxD4iIiEwJgBLaWwkXJ39fmVP42HEXeO/gRpeusmTJ4rESLNAACBgAERFRSCWw8ePHu19Vf/rppx79Jdh8c8WKFapYsWL+no7CxO4lML1/GX7P9ABEZIAQePuzD5g3ToQmIqKQAiA0PwMuRJMnT/Yod+GihX2ncDtFl91LYAhucKCJGwFP8uTJZZsVbIobaA8Q6GAKq+MuX76s4uLiwvYzEBFRDAZABw4ckLc1atRQc+fOdb86J+vAXlmHDh2yZQBkzADpzA0al9EHhLIeymD4+QLNAOlgCudCGaxEiRJh+gmIiCime4CWLl3K4MeicIFHORIrn+y24sm7v0f37mBZvC6DBdMDZCyDsQ+IiIiCXgaPC+z06dPVkiVLZDUY5qsY/f7774Gekkwuf2EAop1m3vjq79ErwIx9QLqvKdAAKF++fJJB0tkxIiKigAOgrl27SgDUoEED2RWeQxGtw45bYOigWgfSugSWMWNGlSFDBncma8OGDapx48ZBBUA6mMJARcwEQoaMiIicLeAAaNasWeqbb75R9evXD88jIsc2QHs3OCNwwbYeuO3ixYsyIRpbrwQaAGFbDWSSzp8/L2UwbOhLRETOFnCdBBejQoUKhefRkCMDIF3+SpEihUfpDgEQVhvqHh798wUaAOlzASZLExERBRwA9ejRQ33wwQcJDkSk6LFrCSyh5e0YrQDYF8zYxBzMwE19LvYBERFRUCWwVatWyUqw+fPnq5IlS8qrdiMskafIQ0Bq1ynQCQ04RNkKc3uwEare1NTX/QLJAKGMhu/HqeVERM4WcACE5tSnn346PI+GgoZemUuXLrlXgcVCBggN9ghcdOnr+PHj8hbDEQOlm6rRS4RMkt2yZEREZK6AryTTpk0z+SGQGXSQgLk5aBS2k8Tm+yAAwp5g+BzuhwAmWCiDbdmyRfqAGAARETlbUMNibt++rRYvXqw+/vhj2V5AlxauXLli9uOjGG+A9jUF2jtoQWO0LoPpLFAwdBmMfUBERBRwBggXj7p167r3ZnryySelT2PkyJHyMfcDiw679v8klQFC9idt2rSypxd+93QfUCgBEM6BbTW8+9eIiMg5kgUzCPGhhx6SmSrGUgv6gjAdmqLDrrvAJ5UB0n1AekuMULaz0E3VGLrI3eGJiJwt4ABo5cqVqn///j6XLIfy6pycWwJLao8v7wDIODgxEAimuByeiIiCCoDw6hlbF3jDK2q8uqbosHMJLLEMECBoQfYGGUf87m3dujXo78U+ICIiCioAql27tho3bpzHq2o0Pw8cOJDbY0QxgNClITuXwBLKAGEQYpo0adxZIGxsGmoAhIAdzfxERORMAQdAY8aMUX/88YcqUaKEbCzZsmVLd/kLjdAUechmYBAigoTs2bMru0loDlBCfUDYGT5YuqkawQ9LtkREzhVwABQfHy+zVN566y31xhtvqHLlyqkRI0aoTZs22fLiG2vlLwQLsZYBgooVK8rk8VADIDw/elCk7psiIiLnCXyk7v+bxPvCCy/IQdFn5xVg/mSAAEHL22+/LWMWdu7cqa5evSqZnGDgXNu3b5fnrUaNGkE/biIiclAGaPjw4Wrq1Kn33I7bWAKLDjuvAPM3AwQogeXKlUsa8ZFxDJYOFFECQxmXiIicJ+AACNOfixUrds/tKE9wCGJ02HUX+EAyQNrDDz8cchkMe4KhFwh9U9gWg4iInCfgAOjEiRPyKtzXSp1Qtimg4DklA2QMgEJZCQbsAyIicraAA6C8efPKKjBvuE2v0qHIQRbDzjOA8PiTmgNkhCnkoWaAjNkyBkBERM4UcBN0hw4dVLdu3WQvpZo1a8pt2AKjV69eqkePHuF4jJSI06dPS0Owccqxnejgx98MkA6A9u3bpy5cuKAyZswY1PfFc4Xn7OzZs7LDPMpiRETkHAEHQG+++aZcNF577TX3xStVqlSqd+/eqm/fvuF4jJQIncHAeAJ/Agir9v8gGMHqwqRkzZpVylcHDhxQGzZsUE888URQ3xe/s9hgFQMRkUErX758UOchIiKHlMBwocJqL2Qe1q5dKzOBzp07pwYMGBCeR0iJsnP5y7v/x98ZRmaVwfRzhmCKiIicJeAASEuXLp00pD744IO2zDzECifMAArHSjBjAITnEL1IRETkHAGXwNBvgsnP6Ps5deqUzGQxYlNpZDlpBZjZK8FQNkTgde3atQRXNxIRUWwKOAP00ksvqc8++0w99thjqnPnzqpr164eR6AmTpwoDanoyahUqZJav359gvfdsWOHatq0qbuB1bgpqzZo0CD5nPHwNbcoVti9BBZMBgj9Ovh3/fvvvyUID9b999/v3hyVgTsRkbMEnAGaP3+++vXXX1XVqlVD/uazZ89W3bt3lwGKCH4Q0NSpU0ft2bPH575ieKWOC/1zzz0n+5AlBEMZFy9e7P7Yn+Zau7J7CSyYDFD69OlV0aJF1e7du6UM1qBBg6C/P36fsKIMz6MZv9NERBSjGaBMmTKpzJkzm/LNx44dK8vq27ZtK7vLIxDCjua+ttrQpY/33ntPNW/ePNELJgKenDlzug+sHIpF2MZB72jupAyQmWUwHTgeOnRIRjsQEZEzBBwADRkyRFZ8IRsT6it/LGOuVavW/z2YZMnk4zVr1oR0bryix1BGBAXYsBWlkqQuwpcuXfI47ECvXoqLi5OtHZySATJzJRiCY2SU7ty5w20xiIgcJODa0JgxY6TvJEeOHNKLkyJFCo/Pb9y40a/znDlzRi46OI8RPkZpI1gopU2fPl1KJNiaY/DgwdKvhN2/ESgktMEr7mfn8pe/S8hjLQOEAAgruIL9+fF1hQsXlmB879698j4REcW+gAOgJk2aKCurV6+e+/3SpUtLQIRG12+++Ua1b9/e59dggCN6kTRkgLDlh9VZfQUYAtxVq1ZJE7p3oKsFsg2GUdmyZaWJGU3QGGYYyr+XDoCQOQwlmCIiohgOgAYOHGjKN0bpARewkydPetyOj9G3YxZslVCkSBG1f//+BO+D8osdZxlZfQUYMirLli2TpvaXX3450QxQoM9/6tSpZQYVBnEiCxRKAITJ0vhdxJYYGPDpqwGfiIhiS9CDEPGK+csvv5Rj06ZNAX89XvFXqFBB5glpmCmEjytXrqzMcuXKFQkUYnHGi9VXgGFmFKAUmVBfVbA9QGYORMTvot4dHlkgIiKKfQEHQCg5YBNUXHxef/11ORDIYE8mvHoOBMpOU6ZMUTNmzFC7du1SHTt2lIsmVoVBq1atPPYXw8Vy8+bNcuB9rIDC+8bsTs+ePdXy5culoXX16tXq6aefllf3LVq0ULHG6iUw40anyAaZ2QNk5kow0L0/DICIiJwh4ACoS5cu6vLlyzKUEHuA4UCDMV7hIxgKRLNmzdTo0aNlVRl6OhDMLFiwwN0vgtVbyB5ox44dU+XKlZMDt+Nr8T6GM2roB0Gwgybo559/XlZHYc+ybNmyqViCXpVYCIBCyQDplWAIgELdykIHQPidu379ekjnIiKiGOwBQoCCIYPFixd334YZPpjoXLt27YAfAKZJ4/AF/SNGWHWW1IVu1qxZygmwdQMu1BgdoKcZW43O7gCCNQQ73pmeUDJApUqVksDpwoULkgUMZQUX5luhLw2rE1EyRX8RERHFroAzQOjT8V76DrjNe18wCn8DdL58+Xz+e1gtA4QVYb62mwglA4SfG5lDs8pgaJYHlsGIiGJfwAEQ+n+w5xfKURp6cbA1BfqAKDKs3gANerKy3orEVxkslAyQmQMRvfuAGMwTEcW2gAOgCRMmSL8PylG4+OLAChrc9uGHH4bnUdI9rN7/YwxujJkV7xJmKBkgM1eCAZbS43GgtKi3GCEiotiUPJiLBKY9ow9IT2xGP5BxSwsKP6vPADIGN8isoEcHIwmQOcyTJ4+7LIYjlAyQDoDwO4lzYcVfsPC1hQoVkgZ/BGt2GIZJREQRnAOESblPPvmkrAjDweAn8uxQAtMBEDa41Y/TWAYzNkkHmwHCar+0adPK3nQYpWBWGSyhVWtEROTQAAhL3cePH++zNNatWzezHhfFQAnMuM2FLoMZAwv9efQIYTVbsFkbzKEysw8IAT4mkp8/fz7k8xERkTUFfNWZM2eOqlq16j23V6lSRX333XdmPS5KBIZFYhm8nQIgnVnB48aWE6FsgxHOgYjIVqG/DczIKBERUYwEQGfPnlUZMmS45/b06dPLDBUKvwMHDrhn1+CwQwCEMpXuqdFZoGA3Qg3nSjDA5q2ge9yIiCj2BBwAoUkUwxC9zZ8/39LZiFhih/IXVnt5Bzjo1zEGFmZngLAxqnH2UKgB0OHDh2XqORERxZ6AV4Fh/y5Mbsa+X5gJBNjAdMyYMWrcuHHheIxkwxVgWJGlZ+noAAerBbF6EBksNC2blQHC84BMGHp2tm7d6s4IBQvZzPj4eNlWBcGaDrCIiMjBGaB27dpJsPPZZ5+pGjVqyIEd4SdNmqQ6dOgQnkdJtl0BBnpSdebMmVXOnDklO4TAwqwMEJqWjfuCmYFlMCKi2BbU0hvs2o5Xx1gpgwGIuCBj53aKDDuUwBJa4YV942Dnzp0hT4EO10BE0HvdIVvFzVGJiGJPcGuP/x/ssJ4uXTrzHg3FTAksoeBGB0AILPRqMCsGQMhW5ciRQ7JVe/bsMeWcREQUIwEQRR76avQqMDuUwLzLW1myZFHZs2eXn2P79u0+7xMMXQLDFGf0F5mBZTAiotjFAMhmsEcVgguUltCoa1WJNTjrLBC2xkjoPoHC9hroL0JgtWnTJmVmGQzbeJixuoyIiKyDAZBN+3/y58/v3mXdrgGQZkYGCI3QZpfBkKlCKQwr2rA3GBEROTQAunXrlnriiSd4MYgiO6wASyoAQu9Y1qxZ3R+bkQECs1eCIajSZTA0bRMRkUMDICxnxpwVih47NEBDUjN+jFkgMzJAYHYGCB588EH39Grj5q1EROSwEth//vMfmQFE0WGXDFBSS9x1f004AiAEKxcuXDDlnOgrQuP27du32QxNRBRDAm4iwYVg6tSpMtEXu3BjjyejsWPHmvn4yIYzgPzJAGGJOUphmCieMWNGU74nymrYyPTgwYNq48aN7knloZbBkAVavny5rForU6aMKY+ViIhsFgDhIlC+fHmPTS2NFwsKr1gpgeF3pWXLlurcuXPSbGwW9AEhAEIZzIwACEqVKiUBEJ77q1ev3hP0ExGRAwKgpUuXhueRUJIwdfvMmTMxEQABMj9mZX+MZbDvvvvO1D4glMBy5cqljh8/Ls3Q3BuMiMjBy+AxG2XhwoXubQIwMZfCSw9ARKkHG3ZaWUKDEMNNBydmrQQzZoFAD28kIiKHBUBnz56VpfBFihRR9evXl1fF0L59e9WjR49wPEayWfkLzNrpPVC6PHvo0CHpLzJLyZIl5e3ff//t3sKDiIgcFAC98cYbshweF4I0adK4b2/WrJlasGCB2Y+PbLgCLJoBUIYMGVTRokXlfTPLYMi4ocEamAUiInJgAPTbb7+pkSNH3rMNQ+HCheVVN4WPXVaARTMACmcZTM8E2rZtm6nnJSIiGwRAWAVjzPxoWM0T6X4Pp2EJLLCJ0KtXrzb1vJhdlCxZMnXy5ElTy2tERGSDAOixxx5Tn3/+ucdyZmxAOWrUKFWjRg2zHx8ZsATmH/SoAZr0165da9p5EfgXKlRI3t+yZYtp5yUiIhsEQAh0PvnkE1WvXj25yPXq1UtKAytWrJDSGIUHBlBivo1dMkB6EnQ0soL4fWzTpo2837lzZ9nM1Cxly5Z1B0AI/ImIyCEBEC4uGID46KOPqsaNG0tJ7JlnnlGbNm2yRWbCro4cOSJBEDIquXPnVlaGkQjRzADBiBEjpHF5w4YNMrncLFj9iEzQlStXZBQEERE5ZBCiXmnz1ltvmf9oKMnyV8GCBdX999+vrAyBmhatAAhbbQwePFhWLfbt21c1bdpUZc6cOeTz4rkvXbq0lNYQ9CMgIiIiB2SApk2bpr799tt7bsdtM2bMMOtxkY0boI27pmNkQrR06tRJdp3H7KoBAwaYdt5y5crJW2RCkQElIiIHBEDDhw+XScTesJ/TsGHDzHpcFCNL4KO5PxyCrw8//FDenzRpkmmNy/hdRxkSPUBbt2415ZxERGTxAAgDEFGG8ZY/f375HIUHV4AFBxuiPvfccxKsdOnSxbQtW3QWCGUwbgNDROSAAAivfn296sWra2waSeFhpxKYlQIgGD16tDQur1y5Us2cOdOUc2IxQPLkyWUe0LFjx0w5JxERWTgAatGihXr99ddlV3gsL8bx+++/q65du6rmzZuH51ESM0AhyJcvn+rXr5+837NnT3X58uWQz5kqVSoZjKizQEREFOMB0JAhQ1SlSpVk2Fzq1KnlqF27tpQa2AMUHufPn5cDfJUfrcZqARBgo14Ej9i899133zV1JhD2Brt165Yp5yQiIosGQLiozZ49W+3evVt99dVXau7cuVKewawVK13wYjH7g6XdadOmVXZZBWal3wdkbMaNGyfvv//++2rPnj0hnxPBaMaMGeXn5f5gREQxHgBpmH+C5tKGDRtKAzSFj53KX8YMkNX2hsPvaoMGDSRbg5JtqM3LWOGm9x3DzvNshiYiivEACFOJP/roI9WnTx/VvXt3jyNQEydOVAUKFJBX6CitrV+/PsH77tixQwba4f64+OhX9KGc0w7stATeGABFcwZQQvA7g8wU9gn76aefTFkNhuGIJ06cUEePHjXlMRIRkQUDoCVLlqiiRYvKXJUxY8ZIMzSGI6IEtnnz5oDOhVIagqaBAweqjRs3qjJlyqg6deqoU6dO+bz/tWvXJAjANgc5c+Y05Zx2YKcVYFbtAdKwmSn6gaBbt27q+vXrIZ0Pq8uwIkxngYiIKEYDIGwrgJU06HlAhmXOnDnq8OHD6vHHH5eSWCDGjh2rOnTooNq2bSsTeydPniwXlIT2bnr44YfVe++9J6vNEiqvBHpOO2AJzFzYxiU+Pl42l8XvU6jwe6kzlJwMTUQUowHQrl27VKtWreR9zEHBK+h06dKpd955J6Dd4HGRxEaVtWrV+r8HkyyZfLxmzZpAH1ZI50QT66VLlzwOK2EGyFxoJMdsID3ZHIFQKPLkySOToTESAllHIiKKwQAIFw99gcuVK5f74gxnzpzx+zy4Ly4YWNlkhI/RTxGMYM+JiyA2eNVH3rx5lVWgYVdP2GYAZJ7nn39eVa9eXd24ccNdEjMjC4QAHFOniYgoxgKgRx55RK1atUrer1+/vlw8hg4dqtq1ayefsyOU9S5evOg+UNKzCgQ/uKCi3IiA0w7sEAChiR77hKGBGaMcFi1aFNL5SpYsKTOx8Puzb98+0x4nERFZJABCjw1WVsHgwYNlICIaj7Hq6rPPPvP7PNhQFRefkydPetyOjxNqcA7XOdGrkj59eo/DiuWvaG4sGmsBEKB5uXPnzvI+ppvrxx0MrHjT+4PZfdUhEZETBBQAobyEJfDYWkCXw9BkjL3B0AwdyDwgXBwrVKggq8o0ZDrwceXKlQN5WGE9Z7TZbQm8nQIgGDRokMqWLZsM9hw/fnxI58JMIASp+DfzDsKJiMjGARCyK9j2Qm/LECosV58yZYqaMWOGNFd37NhRVtFgBReg2RrlKeOFFUvtceB9zF3B+/v37/f7nHZjtxVgxknQVl0FZoRJzrp5HxlNbJURrEyZMrn3Bwu2kZ+IiCxaAkPZQF+UQ9WsWTNZjTNgwADZVwnBzIIFC9xNzOh/MV6QsOs2ygw4cDu+Fu+/9NJLfp/Tbuy2AsxuGSBo3bq1lHWvXLmievXqFdK5qlSpIm8xJsJqqwmJiOj/3OcKcH4/gglkZbApKspN3ntTWal/Jli4cGE1GBpao/3zIMBDEPfzzz/LVg52gE1xsXoNfTXIitgBhhgiCML/Dmjyr1q1atDnmj59ujp06JAEQ08++aSpj5OIiMy5fvudAcKcH5SSsPJry5Yt6qmnnpJhcrjA4UApwS4XO7vAxdhuJTA8Zr0zul0yQHoZe/v27eV9NEaj3y3ULBCWxOtyIBERWUtyf++I/ohXX31Vtr6gyDh79qy7jIJVdnZgXEllpwBIZ66+++47ybh98skn0j8WjMKFC8uKRMylQhCkAyIiIrJhAKQrZdjygiJDZ38waRgzZuwUAGE1FCaF2wlWg6G026VLF9W/f38ZlpglS5aAz4OfHasOUbZcu3atlNawgICIiGzaBG2XOTSxwu5L4O34+4IsZ6lSpdS5c+dkz7BglS5dWvrjLl++rLZv327qYyQioggHQEWKFFGZM2dO9CDzcAVY5CFrNWHCBHkfZbBg9/bCefTAUDRVc3sMIiJrCahGgT4gdFdTZNitAToWAiCoVq2aatGihZo5c6Y0RCOAwaa6wTRWr169WnqBdu7cKSMkiIjIhgFQ8+bNVfbs2cP3aMj2JTA7DUFMzHvvvad++uknGWj45ZdfylDOQGH/NuyPt2zZMrVixQrZL8yOZUEioljk98ta/uGOPJbAogeN52+//ba8j+GImCkRDJTBEAidPn1askBERGSzACjAeYlkQiYF+65ZsQR2+/btmA+AoFu3btL3hn29MAcrGAh+dC/Q8uXL+f8REZHdAiA0cbL8FTkHDx6UiyVWEmF5tpUmJg8fPjzBbEYsBUAo433wwQfyPjZKDTaDgzIYzsUsEBGRdQTe2UkR7/+xUvkR+7MhGP7tt998ZoJiKQCCunXrqsaNG8vPiq09gsng6F4gQC8Qs0BERNHHAMiirLoCTG9zgZ6Yv/76K+YDIBg7dqxkcJYsWaLmzJkT1DlQBsM5Tp06xSwQEZEFMACyKKs2QBu3uli5cuU9e13FYgCEf4PevXvL+z169FDXrl0L+ByY5K2zQL///ntIe40REVHoGABZlFWXwBu3ukAggGXisR4AAQKgfPnySQlwxIgRQZ0De4KhpwtTprFHGBERRQ8DIIuyaglMBzjly5eXtwiArl69GvMBUJo0adT7778v748aNcr97xMIPCd6Lz2sCONO8URE0cMAyILQJGvVDJDuASpbtqzKlSuXBDwohcV6AARPP/20qlWrlgQub7zxRlDnQOCIDVaRPcOUaCIiig4GQBaERllkVVBmyp8/v7ISY4CDYADQDH3+/PmYmgTtC/49sBwe+3xhSvT8+fMDPgd2ha9Zs6Y7e4bNUomIKPIYAFmQzv7kzZvXcoGEMQBCdgoHGnoXLVp0z+djUfHixVXXrl3lfbwNpoyFc8THx0s2DaUwIiKKPAZAFmTVFWCY/6Nn/6RIkULe1qlTRzIju3btUgcOHIj5AAgGDBigcubMqfbt2+fuCwoEni+dPcNu88j4ERFRZDEAsiCrNkDr/h9jgIPp4A899JC8v2DBAndGJJYDoPTp00sjNLz77rvuLUsCgdJmsWLFpN8LpTQORyQiiiwGQBZk1QZo4wwg9MFo1atXl2nHyGRcv3495gMg+M9//iPL2tGr9eabbwZ1DmTP8Dxi2xMORyQiiiwGQBZk1RKYzgAhuDFuz4El4jVq1PC4b6wHQPj5J0yYIG9nzZoVVC9PxowZ1aOPPirvY2sRY4BJREThxQDIgqxaAkusvwdlMONmuVZr3g6HcuXKqVdeeUXe79Kli8+90ZKCLBICoUuXLnmMEyAiovBiAGQxKCEdO3bMkhmgxAKgZMmSSUkHUA7Dcm8nQA9Q5syZ1bZt29SkSZMC/no0k+vnDcviz549G4ZHSURE3hgAWQxWUulGW1xYrVgC0yvAvCFga968uRxOgaGGQ4cOlffffvvtoFZ0FS1aVBUqVEjGCaCRnA3RREThxwDIwuUvY5+NFfizxB0Xc6sNbwy3Dh06SDns4sWLql+/fgF/Pf6d69atK1mz/fv3SzaJiIjCiwGQxVi1ARqcMOMnGAhc0BANU6dOVevXrw8qk1StWjV5H1kg4/5qRERkPgZAFmPVJfDGACihEpiToZm5VatWUr7q3LmzDI0MVNWqVVWOHDmkDwxBEBERhQ8DIIux6gow72XwdK+RI0equLg49eeff6rp06cHlUl66qmnpCS2fft2tWfPnrA8TiIiYgBkOXYogTED5Bu2xxg4cKC836dPH3XhwoWAz5E7d25VuXJlef/XX39VN27cMP1xEhERAyBLQdlErwKzcgDEDFDCXn/9ddns9PTp0+5gKFCYrI0VgNgpHgMSiYjIfAyALOTEiRPyih+lkHz58imrYQCUNGTHxo8fL+9PnDgxqBVdOAdKYbBp0ybZaJaIiMzFAMiC5S8EP1YsM7EHyD/Y6b1p06Yy1wcTooOZ64NRAmiKhp9//lmyQUREZB4GQBZi5RVg/gxCpP8zZswYlTp1atkjbPbs2UGdA/ur5cqVS1aF/fjjjxyQSERkIgZAFmLlFWDAElhgGZy+ffvK+z179lRXrlwJ+BwohT7zzDOyYzyyg+vWrQvDIyUiciYGQBZi5RVgwAAoMG+++aYqWLCgOnr0qHu7jEBlzZpV1a5dW95fvHix9IkREVHoGABZiNVLYAyAAoNNYceNG+cuie3duzeo8zz00EOqSJEi0lP07bffcmk8EZEJGABZMANk1RIYe4AC16hRI9nnC89dt27dgurjwWDExo0bqwwZMqhz586xH4iIyAQMgCwCPSJ6J3FmgGIHgpcPPvhAgsb58+erX375JajzpEmTRj333HPSF7R79261evVq0x8rEZGTWCIAwryUAgUKSMmgUqVKSW4miTJAsWLF5P6lSpVS8+bN8/h8mzZt5MJjPPAq3Mr0AEQMwMuYMaOyIgZAwUH5qnv37vI+skDBlrDy5Mnj/j1esmSJOnjwoKmPk4jISaIeAGGJMC4OmJq7ceNGVaZMGVWnTh13NsQbXvm2aNFCtW/fXobENWnSRA7snWSEC8Xx48fdx8yZM5WVWb0BGiUXboURvP79+8s2F+jzGj16dNDnqVChgipdurT8e3z33Xfq4sWLpj5OIiKniHoANHbsWNWhQwfVtm1bVaJECTV58mRJ90+dOtXn/VFOQHCDFTbYcmDIkCGqfPnyasKECR73S5kypezNpI9MmTIpK7N6A/Tt27fd7zMDFLh06dK5A59hw4apv//+O6jzIJvZoEEDlT17dnX16lU1a9Ysd2BKREQ2CYDwh3vDhg0yOdf9gJIlk4/XrFnj82twu/H+gIyR9/2XLVsmF4miRYuqjh07qrNnzyors8sMIGAGKDjNmzdX1apVk8GGPXr0CPo8CECRBcULBSyLnzt3LpuiiYjsFACdOXNGlvbmyJHD43Z8nNC8E9ye1P2RIfr888+lT2LkyJEyjbdevXryvXy5efOmunTpkscRaVYvgekVYBjKhyCVgsvefPjhh/L8oXyF389goU8MARWaovfs2SMzgoiIyH8xeSXDhQGbSaJBGv1BWHnz559/SlbIl+HDh8sSY33kzZs34o/Z6iUwNkCbA/07r732mryPfcJ0YBkM/J5iebzujUNPHBER2SAAwpRbvII9efKkx+34GH07vuD2QO6vgwp8r/379/v8PLYsQDOpPg4fPqwiCZkpvaLH6iUwBkChe+edd+T3Ebu8e/euBQpBPspqetPUYIctEhE5TVQDIFxMsarFWAq4e/eufFy5cmWfX4PbvUsHixYtSvD+cOTIEekBwsaSvqBhOn369B5HJGGrBAQYKC/Fx8crK+IQRPOgIR9ZR8Dqx1C3t6hevbp7ZRhGRBw6dMikR0pEFLuiXgLDEvgpU6aoGTNmyCtiNCxjdQtWhUGrVq3cm0pC165d1YIFC2RrAQyEGzRokPrrr79U586d3QMFsUJs7dq1klVBsIQyQaFChaRZ2op0+QuzkJARsyJmgMzVrl079fDDD6vLly+rPn36hNxbhJIv5g1htR5GPnDPMCIiiwdAzZo1k+XBAwYMUGXLllWbN2+WAEc3OmO5MOb4aFWqVFFff/21+uSTT2RmEJpJf/jhB/Xggw/K5xFAbN261X1BwLwgZJlWrlwpmR4rsvoKMGAAZC40QqMhGhD8hzrZGb/3zz77rMqXL5809X/55ZeybQYREfl2n4vrZ++BVWBohkY/UCTKYW+99ZbMhkH266OPPlJWhHEFaCbHWAE0mZM5EKBj5hVmWWECeqgZQEyZnj59uvTF4Xe3devWMl2ciMgJLgVw/Y56BoisvwIM2AMUHnoFIqagf/rppyGfD9vD/Oc//5Ema/whQHaJmSAionsxALIAq+8CDyyBhQeGdWJVmM4EmhGsYOo0Mj86CEJGiEEQEZEnBkAWYIcMEPcBCx/MBUIPG1YqogyKRn4zgyA0WiMIwuBRIiL6/zEAijLUKfU2HXYIgJgBMh/GH+iG6G+++Ub6rNDoH2p7ng6CsmXLJkHQtGnTZOQCERExALJM9gcXqbi4OGX1HiAGQOGBWT5oMkcQfOzYMfXCCy+oxx57LOTpzjoIwk70165dk54gXXIlInIyBkBRZofyF7AEFn7Y5X3Hjh1q6NChstHpH3/8oR566KGQN/NNmzatzNPC7xgCWWSXtm3bZupjJyKyGwZAUWaHGUDADFBkYBVXv379ZMgnZmRhMvrkyZNV4cKFZUQCBh0GAzOwWrZsKb1GOCd2kF+xYgV3kScix2IAFGVW3wVeYw9QZGGj01mzZskGvtjv6/z586pTp04y1BOBSzAwY+iZZ55RlSpVko+XLl2qvv/++6CDKiIiO2MAFGV2K4ExAIqsxx9/XGYEYdNU7CGGKee4rUWLFrLHXTDbZtStW1fKbZhGjVIYVoiZsfKMiMhOGABFmd1KYOwBis4qMWR/sNP7K6+8IkEMskNYLYYJ4pj+HCj0FmFgIkpuWBmGrWWw7QwRkVMwAIoilB70zt3MAFFSMNMH/UDYlqRq1aqyqgvDE0uWLKl+/vnngPt5ChYsqF566SX3rCCsEFuzZg37gojIERgARdHhw4clCEKDKpYpWxkDIOsoV66cbO6LDU9z5colWURs/lu/fn21Z8+egM6VJUsW1aFDB3dz9G+//SaziILJKhER2QkDIAuUv/BKHP0YVsYAyFpQBsOsIAQ8vXv3ltLkggULpGG6V69esgWGv/BviuZoBFBolMYKNGSaDh48GNafgYgomqx91Y1xdlkBdufOHckOAHuArAXDM0eMGCHzgxDAoFfrvffek/6gL774wv3v5k9A9fDDD6t27dpJszUmlKMktmjRIq4SI6KYxAAoiuy2AgyYAbImzAn69ddfZZp0oUKF1IkTJ2T44aOPPio9Q/5CKRaN1iizwerVq2WXepyPiCiWMACKIrutAEOZDiUSsi4sb9++fbsaPny4TIBGUzMyOy+//LI6ffq0X+dATxp6ijCIEROpT548KavEFi9e7P5dICKyOwZAUWSXEhj7f+wFAUyfPn2kPwjTn7Gqa8qUKapIkSKy6aq/Ja1ixYrJNhwlSpSQc2BrDvQGHThwIOw/AxFRuDEAiiK7lcAYANlLnjx51FdffSUrxsqWLasuXLigXn/9dSlvYQq0v5upPvfcc5INQr/RuXPn1Oeff67mzJkTUKM1EZHVMACKElxIcEGyQwDEIYj2hj6gv/76S02aNEllzpxZSmQ1a9ZUzz//vN/DD5ENwjBGDFAEnAPTqRFcsUmaiOyIAVCUsz85c+aUPgsrYwbI/tC79eqrr8o06ddee036ub799lsJbIYMGeLX3B+U1tBjhH6ifPnySWD8+++/q4kTJ8oWHRygSER2wgAoSuxS/gIGQLEDgw8RsGBl2GOPPaauX7+uBgwYIH0+P/zwg19BDIYvtmnTRmYHoSyGTCY2Vf3444/Vvn37GAgRkS0wAIpyA7TVV4ABA6DYg56g5cuXq5kzZ0qvEBqbn376adkoFYMQ/ZkbhKGLXbp0kXIaskNYLfb111/L5qr4/WYgRERWxgAoSuyUAWIPUGxCENO8eXMJePr16ycBLrbCQGDTs2dPv5qc8TuBTBKaqytXriylNvQVYZuOzz77jBkhIrIsBkBRYpcZQMAMUGzDSq+hQ4fKNOlGjRpJU/OYMWNk2TyyOf5Mk0YfW+3atVXXrl1VpUqVZAd77DKPjBCWzm/evFkmihMRWQUDoCixywwgYwDEDFBswwTpn376Sc2bN08mS6Ok1bZtW1WlShX1559/+nUO9AShjIZACBkh/M6cOnVK/fjjj2rcuHGyauzq1ath/1mIiJLCAChKAQV2grdbAMQMkDPUq1dPlrmPHDlSskPr1q2TrM5LL70kwYw/8HXICL3xxhvqiSeekMDoypUrsmrs/fffV3PnzpVSGctjRBQtDICiAH/4UVZInTq1LIO3Sw8QAyDnwL81dpXHNOkXX3xRAhX09KAshkyOv1ti4Hccc4iQEWrSpInsNYZS2LZt29S0adNkNhH2G0NwREQUSQyAolz+QiOq1bEE5lwIWDD5edWqVTJBGrvEI6uDVWRLlizx+zxoji5Tpozq0KGDHPh69AlhfzLsOD927FhZkbZz507uN0ZEEcEAKArstAIMmAGiqlWrSh8QZv1glhAClVq1aqlnn31Wsjn+NEobg6rGjRurHj16yGBFLMNHhglDGjGccfTo0TJXCCvI2DhNROGSPGxnpphYAQbsASKdxcEUaOwNNnDgQBmoiD3BcKDHp0KFCrJVBg7sQF+wYMFEM5ypUqVy3x+ZoC1btkjvEbJMmCyNA/OF0JCNidV4y99BIjILA6AosNMKMGAAREaZMmVS48ePl1IWdp3HxqqXL19Wy5Ytk8N4Px0M6UAnPj7eZ1CULVs2ySihYfrIkSMSCCHLhN4gvI8DAVj+/PklEMKBTBQRUbAYAEWBXUtg7AEiIwxM/PXXX2VuEIIVbLiKA6UyZHPOnz8v/T04tBw5cnhkifAWt2kIjvLmzSsHltMjGMKgRhzYQBj/7+BYuHChBFj4fwhHgQIFLL+nHhFZy30urkO9BybgZsiQQVLx6dOnN/XceLpxTryy3bVrl6T2re6DDz6Q/Z7at28vr+CJknLz5k3J2uiACG/xsa+eHvxOGbNEOLBrvff/N2fPnpW+IByHDh26p+8IKyqxSSuyRDjSpk0b9p+TiOx7/WYAFOEACL0O2bNnl/exESX6IKzuvffeU9euXVMdO3Z0P3aiQOH3HROhjZkiZHZ8/QlCVseYJSpfvrzH/4sIsBAE6YwQ/r/yhiAKwRWySXiLMhvKaEQUuwK5frMEFqXyF1a+2CH4AfYAkRkwEwjToXFo6B3atGmTR6Zo//797sDmm2++cZfGihYt6lE+w1J6zCXS58F8LQRFODCwESUzHGimBiy7R5YIq9Cwoz3eZ1BE5FwMgCLMTrvAA8oM6PEA9gCR2bB6rFq1anJo6B3asGGDR6YIwY3uBcJGq5AsWTJVsmRJj0wRmqixcgzZJvQP6QP7kiFrpD/WcA4EQehDwltkOHHgFaQdZnQRUfAYAEWYXRuggRkgigQ0N2NFGA4NGR0ERTpLhLcnTpyQGUR6qrQO0kuXLu0RFGESNbI/yAYdO3ZMDnwtjhs3bsieZziMcP+sWbPKgdVmOFBSw4FMFhHZHwOgCLNzAISLAlE0ICuDPcpwaMjqGLNEeItGaQRKODC0EVBqxhRrY5M1gitkf9AvgEAIAZY+zpw5I1lPHSR5QwCEIC1jxoweB7JGeMsXCkT2wCtahNmtBGbs/2FJgKwEfXQ4MFUa0EyN/h9jQIQDQc6aNWvk0PC7jGBFZ3X0gUwPbkeQg4Af90MZGOfG/wsIjFBew4FMki8IuNB8iQMlPu8DG8VihRp7j4iiiwFQhNktA8QGaLILBCuYB4QDW3QAghc0VRszRRs3bpRVjeg1wqFflPhLBzcIZDB7CD1HKL0hoMH/JwiefB3eAQ9uQyCkD5zLeODz+i0OfB++CCGKsQAII/Wx1BrpZmyY+OGHH6qKFSsmeH/sF/T222+rgwcPykTYkSNHqvr167s/j1drGNU/ZcoUmV+DfYyw6zTuG03oN0Da3o4ZIDZAkx2hzIWVYjhatmwpt2EWEcpcepWY8UAJLaHbkUkCvNXvBwKZIR0wJRQk6UMHPvgaXXpG8IOvxW36c/pj/RYBGN7XBz7Wt+Et/j9mNpfIIgHQ7NmzVffu3dXkyZNVpUqV1Lhx41SdOnXUnj17fM6cWb16tWrRooUaPny4atiwofr6669VkyZN5FXdgw8+KPcZNWqUjOqfMWOG7EeEYAnnxLTaaC49R8CG4AyvHNFcaQfcCJViDTIxWPVlnEDt7/8LeEHlT7Bk/Bhfo18A4QiUziohiMFj1wcCI+PHgdzuHSgZgyhfhw64EJghiMK59PkSepvQgaAUB1G0RX0QIoIerNaYMGGCO2WNwWVdunSRfYa8NWvWTF29elX98ssv7tseeeQRmQmCIAo/DuZ8YKfpnj17yucxEAl/7KZPn66aN28etUGI8+bNk92vsUoFWwVgGrReYm5VKA/guca/Sbt27aL9cIhsBxknBEG+AqbEgih8jRXn1CJ4SSi48RX8INvk6zAGQ4kd3vczfqzPr2/zfpvY1/k6rz4fPvZ+vN636Y/1c6JvA1+fN37OeN/E3ur7+vrY+3tgthVeWHtn94wfG7/WKKGMoK+vTew+Sd3ufRsCarMntttmECLKK1it0bdvX/dt+AfFCg1jw6IRbkfGyAjZnR9++EHeP3DggJTSjEto8WQg0MLX+gqAMB8EhxZMejuYXeB//PFH6U+wA5bAiIKDi6peSh9o4IQ/4joowt8l/M3UB/5mBfo+MlDeb42fR5YL7+OtPry3L8GLVBzGFaIUfdg7D8kAO3n00Udldle0RDUAQh0e/3N5p6LxMQae+YLgxtf99XJV/Tax+3hDOW3w4MEq3JDxQSpZN0DrV0xWh6C0RIkS0X4YRI6Cvw96dVo0IdgJNujCgSyWDqRwIOuNt8bbjLcndujgC/fFW+Ntxvd9fZzYgceI++Ot/lhn3/T73h/r933dltTnvW9L6uOk3oJexZjQ572zia4kPvb3Pv5+rRVfWFv/6hsByEAZs0p4pYWSj9lQ0uvdu7c72+RPOY6IKNovgHQvEFEsiWonGuqVeJXjPYUVH2OfHl9we2L3128DOSea//TSVn2EC2qg/ENCRETk4AAIKxEqVKiglixZ4r4N6Ud8bNww0Qi3G+8PixYtct8fq74Q6Bjvg4zOunXrEjwnEREROUvUS2AoPbVu3VrG02P2D5bBY5VX27Zt5fOtWrWSaa/o04GuXbuqxx9/XI0ZM0ZWVM2aNUsGnH3yySfuDEu3bt3Uu+++K3N/9DJ4rAzDcnkiIiKiqAdAWNZ++vRpNWDAAGlSxnL2BQsWuJuYsQu0cWZElSpVZPZP//79Vb9+/STIwQowPQMIevXqJUHUyy+/LEtJ0WmOc7L0RERERJaYA2RF4ZoDRERERNa4fnMcJxERETkOAyAiIiJyHAZARERE5DgMgIiIiMhxGAARERGR4zAAIiIiIsdhAERERESOwwCIiIiIHIcBEBERETlO1LfCsCI9HBsTJYmIiMge9HXbn00uGAD5cPnyZXmbN2/eaD8UIiIiCuI6ji0xEsO9wHy4e/euOnbsmIqLi5Pd5c2OThFYHT58mPuMhRGf58jg8xwZfJ4jg8+z/Z9nhDQIfnLnzu2xkbovzAD5gCctPj4+rN8D/+j8Hyz8+DxHBp/nyODzHBl8nu39PCeV+dHYBE1ERESOwwCIiIiIHIcBUISlTJlSDRw4UN5S+PB5jgw+z5HB5zky+Dw763lmEzQRERE5DjNARERE5DgMgIiIiMhxGAARERGR4zAAIiIiIsdhABSgFStWqEaNGsmUSUyJ/uGHHzw+f/LkSdWmTRv5fJo0aVTdunXVvn377jnPmjVrVM2aNVXatGllEFS1atXU9evX3Z8/d+6ceuGFF+RzGTNmVO3bt1dXrlxRThKJ5/rgwYPy3BYsWFClTp1a/fvf/5bVCf/8849yikj9Tms3b95UZcuWle+1efNm5RSRfJ5//fVXValSJfmdzpQpk2rSpIlyikg9z3v37lWNGzdWWbNmlc8/+uijaunSpcopVoT4PONvL77O1/Htt9+67/f333+rBg0ayDmyZ8+u3nzzTXX79m1TfgYGQAG6evWqKlOmjJo4ceI9n8OCOvyh+d///qd+/PFHtWnTJpU/f35Vq1Yt+Trj/1j4Zahdu7Zav369+vPPP1Xnzp09xnYj+NmxY4datGiR+uWXX+SX7eWXX1ZOEonnevfu3bL1yccffyzP9/vvv68mT56s+vXrp5wiUr/TWq9eveSPotNE6nmeM2eOevHFF1Xbtm3Vli1b1B9//KFatmypnCJSz3PDhg3lQvz777+rDRs2yPfEbSdOnFBOcDXE5xlbYRw/ftzjGDx4sEqXLp2qV6+e3OfOnTsS/OAF6erVq9WMGTPU9OnT1YABA8z5IbAMnoKDp+/77793f7xnzx65bfv27e7b7ty548qWLZtrypQp7tsqVark6t+/f4Ln3blzp5znzz//dN82f/5813333ec6evSoy4nC9Vz7MmrUKFfBggVdThTu53nevHmuYsWKuXbs2CHn3bRpk8uJwvU837p1y5UnTx7Xp59+GsZHbx/hep5Pnz4t51mxYoX7tkuXLsltixYtcjmNCvJ59la2bFlXu3btPP5eJEuWzHXixAn3bZMmTXKlT5/edfPmzZAfNzNAJkJqH1KlSuW+Da8YMOxp1apV8vGpU6fUunXrJJVXpUoVlSNHDvX444+7P69ffaDs9dBDD7lvQ+SMc+Frybzn2peLFy+qzJkzh/kncN7zjJR4hw4d1BdffCHpbDL/ed64caM6evSofG25cuVUrly55NX09u3bo/BTxe7znCVLFlW0aFH1+eefS0YDmSBkkfE1FSpUUE5304/n2RuyaCiJoyXBeC0sVaqU/BtoderUkc1UkbEPFQMgExUrVkzly5dP9e3bV50/f17SdiNHjlRHjhyR9B4gJQiDBg2Si8GCBQtU+fLl1RNPPOGujyKFiv+RjJInTy4XZaekVyP1XHvbv3+/+vDDD9Urr7wS0Z8n1p9nvEhEP8Crr77qEdiTuc+z8T79+/eX8jl6gKpXry59hU5n1vOMPpXFixdLaScuLk4u9GPHjpX74vl2umJ+PM/ePvvsM1W8eHEJOjVc74zBD+iPzbgWMgAyUYoUKdTcuXOlOQ7BCl7loikOr8B07Rj9JoALLGr0eJWGvhO8mpg6dWqUfwJnP9d45Yy6/3PPPSd/+Mi85xlB5eXLl+UPIoXvedb3eeutt1TTpk0lGzFt2rR7GkudyqznGQF9p06d5IXqypUrpU8IPS9oCk7oAu8kKfx4no3QXP711197ZH8iIXlEv5sD4A8O0ngooyDqzZYtm6zG0K96kZKGEiVKeHwdIl90u0POnDklDWuEFCteweFzZN5zrR07dkzVqFFDXn188sknEfwpnPE8o1EU6WzvvX9wDjT8o7nR6cx4nn3dB8/5v/71r3t+553KrN9nZNeQ3cAKMPjoo49k0Qp+l/v06aOcrkISz7PRd999p65du6ZatWrlcTuudwguvUvp+nOhYgYoTDJkyCD/4EiZ/vXXX7JcEgoUKCArYPbs2eNxf0TK6JKHypUrqwsXLkhNVMP/cHhlgl8gMu+51pkflAj0q2Vfr1AotOd5/PjxsiIJfxBxzJs3T26fPXu2Gjp0aBR+mth8nvE7jIDHeJ9bt27JkmPj7zyF9jzjYg3efyvwsc4gUeLPs3f566mnnpL7GeFauG3bNo+EAIJMBJ3eAWpQQm6jdpjLly/LyhUcePrGjh0r7x86dEg+/80337iWLl3q+u9//+v64YcfXPnz53c988wzHud4//33pYv922+/de3bt09WG6RKlcq1f/9+933q1q3rKleunGvdunWuVatWuQoXLuxq0aKFy0ki8VwfOXLEVahQIdcTTzwh7x8/ftx9OEWkfqeNDhw44LhVYJF6nrt27SorwRYuXOjavXu3q3379q7s2bO7zp0753KCSDzPWAWWJUsW+brNmzfLqqeePXu6UqRIIR87wWUTnmfA84sVzljp7O327duuBx980FW7dm15XhcsWCAryfr27WvKz8AAKED4B8U/tvfRunVr+fwHH3zgio+Pl/8R8uXLJ//j+FquN3z4cLlfmjRpXJUrV3atXLnS4/Nnz56VgCddunTyP2Lbtm3lF85JIvFcT5s2zef3cNJrg0j9Tjs9AIrU8/zPP/+4evToIUFPXFycq1atWh7LkWNdpJ5njCnBhTlz5szyPD/yyCOybNsplpr0PCOYyZs3ryyT9+XgwYOuevXquVKnTu3KmjWr/G5j3IMZ7sN/Qs8jEREREdkHmx2IiIjIcRgAERERkeMwACIiIiLHYQBEREREjsMAiIiIiByHARARERE5DgMgIiIichwGQEQUc7DzPDanJCJKCDdDJSJbwc7miRk4cKD64IMPZMduIqKEMAAiIls5fvy4+31spjpgwACPjSvTpUsnBxFRYlgCIyJbyZkzp/vATtPICBlvQ/DjXQKrXr266tKli+rWrZvKlCmTypEjh5oyZYq6evWqatu2rYqLi1OFChVS8+fP9/he27dvV/Xq1ZNz4mtefPFFdebMmSj81ERkNgZAROQIM2bMUFmzZlXr16+XYKhjx47queeeU1WqVFEbN25UtWvXlgDn2rVrcv8LFy6omjVrqnLlyqm//vpLLViwQJ08eVI9//zz0f5RiMgEDICIyBHKlCmj+vfvrwoXLqz69u2rUqVKJQFRhw4d5DaU0s6ePau2bt0q958wYYIEP8OGDVPFihWT96dOnaqWLl2q9u7dG+0fh4hCxB4gInKE0qVLu9+///77VZYsWVSpUqXct6HEBadOnZK3W7ZskWDHVz/Rf//7X1WkSJGIPG4iCg8GQETkCClSpPD4GL1Dxtv06rK7d+/K2ytXrqhGjRqpkSNH3nOuXLlyhf3xElF4MQAiIvKhfPnyas6cOapAgQIqeXL+qSSKNewBIiLyoVOnTurcuXOqRYsW6s8//5Sy18KFC2XV2J07d6L98IgoRAyAiIh8yJ07t/rjjz8k2MEKMfQLYRl9xowZVbJk/NNJZHf3uTgulYiIiByGL2OIiIjIcRgAERERkeMwACIiIiLHYQBEREREjsMAiIiIiByHARARERE5DgMgIiIichwGQEREROQ4DICIiIjIcRgAERERkeMwACIiIiLHYQBEREREymn+P5n5v/xZdch6AAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# simulate with both models\n", "sim_m = m_m.simulate()\n", "sim_y = m_y.simulate()\n", "\n", "# plot output series\n", "fig, ax = plt.subplots(1, 1)\n", "ax.plot(\n", " timesteps_m,\n", " sim_m,\n", " label=\"Tracer Output (m)\",\n", " c=\"grey\"\n", ")\n", "ax.plot(\n", " timesteps_y,\n", " sim_y,\n", " label=\"Tracer Output (y)\",\n", " c=\"black\"\n", ")\n", "ax.set_xlabel(\"Time\")\n", "ax.set_ylabel(\"Tracer concentration [-]\")\n", "ax.legend()\n", "plt.show()" ] }, { "cell_type": "code", "execution_count": null, "id": "ba549955", "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "PyTracerLab_app", "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.12.11" } }, "nbformat": 4, "nbformat_minor": 5 }