{ "cells": [ { "cell_type": "markdown", "id": "a97813c2", "metadata": {}, "source": [ "# Hysteresis Model\n", "\n", "## Packages" ] }, { "cell_type": "code", "execution_count": 1, "id": "b9dc3c04", "metadata": {}, "outputs": [], "source": [ "import matplotlib.pyplot as plt\n", "\n", "import pedon as pe" ] }, { "cell_type": "markdown", "id": "ce1c0ccc", "metadata": {}, "source": [ "## Implementation\n", "\n", "There is only one hysteresis model available in `pedon` which is based on the work of Kool and Parker (1987). This model is implemented in the `Kool` class. The parameters for this model are the same as for the `Genuchten` model, but with an additional parameter $\\xi$ which controls the degree of hysteresis by multiplying the (dry) air-entry pressure $\\alpha$. The 'ink-bottle-effect' (pores empty at larger suction than filling up) leads to a shift of the primary wetting curve compared with the primary drying curve. The van Genuchten shape parameters n and m remain similar for both, drying and wetting (Bashir et al., 2015). Also, no air entrapment is assumed - values for the residual and saturated water contents remain similar (Kool & Parker, 1987; Bashir et al., 2015). The following table from Bashir et al. (2015) shows the parameters for the three soil types used in this example.\n", "\n", "| SoilParameter | Loam | SiltySand | SiltyLoam |\n", "|-----------------|------|-----------|-----------|\n", "| Ks (cm/day) | 1.915| 25.92 | 137.3 |\n", "| theta_r | 0.21 | 0.07 | 0.03 |\n", "| theta_s | 0.52 | 0.38 | 0.4 |\n", "| n | 1.704| 2.02 | 1.228 |\n", "| alpha(_d) (1/cm)| 0.015| 0.02 | 0.091 |\n", "| alpha_w (1/cm) | 0.08 | 0.03 | 0.503 |\n", "| $\\xi$ | 5.33 | 1.5 | 5.53 |\n", "\n" ] }, { "cell_type": "code", "execution_count": 2, "id": "07585678", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAs0AAAGeCAYAAABvvTxFAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjksIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvJkbTWQAAAAlwSFlzAAAPYQAAD2EBqD+naQAAm3BJREFUeJzs3Xd4FGXXwOHf7GbTG0lIIJAQutQQIDSRoiICgoAIFhRBsQEWXkXU9wXLZ+8lNgTBAqIoIIIo0pFeQu+E0AME0vvufH9MshADpG2ZTc59Xbl2szs7c5LsyZx95imKqqoqQgghhBBCiKsyODsAIYQQQggh9E6KZiGEEEIIIUohRbMQQgghhBClkKJZCCGEEEKIUkjRLIQQQgghRCmkaBZCCCGEEKIUUjQLIYQQQghRCimahRBCCCGEKIUUzUIIIYQQQpRCimYhhBBCCCFKIUWzEEIIIYQQpZCiWQghhBBCiFJI0SyEEEIIIUQppGgWQgghhBCiFFI0CyGEEEIIUQopmoVwgHPnztGvXz98fHxo2rQpS5curdC2n3/+OW3btsVkMvHSSy+VeG2PHj3w9PTE19cXX19f+vTpY48fRwhdskWe5ebmMmrUKCIjI/H396dTp06sW7euwscRoiqy1Tlt9+7ddOvWDX9/f5o3b86KFSuKvba0c56juTk7ACGqgzFjxlCrVi3OnTvH33//zdChQzl48CBBQUHl2rZ27dq89NJLzJw586rH+vrrrxk+fLg9fxwhdMkWeebh4UFUVBRr1qyhbt26/PTTT/Tv35+jR4/i6+tb7uMIURXZItf8/Py4/fbbeeKJJ1i+fDnLly9nyJAh7N+/n+DgYIAynfMcShVCRywWi5qZm2/3L4vFUuaYcnJy1AkTJqjh4eFq3bp11X/++Ue988471WXLlpXp9enp6arJZFKPHz9ufax79+7qtGnTKrztI488ok6ePLnE67t3765+9913ZfzJRHWmt1zTY54VqV27trp58+YKvVYIybUrb7tz5041KCio2PY333yzOmXKlBL7udo5z9GkpVnoSna+meaT/rT7cfa80htv97K9/f/zn/9w5swZDh48yOeff86zzz6L0WikZ8+e3HbbbaxZs+aKr5s4cSITJ07k4MGD+Pr6UrduXetzrVq1Yvfu3SVeU55tr+bpp5/m6aefpk2bNrz33nu0bt26zK8V1Yfeck2veXbw4EEuXLhAo0aNyv1aIUBy7WrbxsbGoqpqse1VVdV1LknRLMQ1XLhwgalTp5KYmIi3tzexsbE888wzLFq0CIDff/+91H1kZGTg7+9f7DF/f3+Sk5Mrte2VvP322zRv3hyj0cgnn3xCnz592LdvH35+fmV6vRDOoNc8y87OZvjw4Tz//PMEBASU+zhC6I2ecq1p06YEBgby/vvvM27cOJYuXcrKlStp0KBBJX5C+5KiWeiKl8nInld6O+Q4ZbF06VLatm1LaGgoAHl5ebRp06ZcA+x8fX1JS0sr9lhaWpq1f2RFt72SDh06WO9PmDCBadOmsX79enr16lXmeEX1oKdc02Oe5efnc+edd9KoUSMmTZpUoeMIAZJrV9vWZDIxb948xo0bx2uvvUb79u0ZNmxYsVZpvZHZM4SuKIqCt7ub3b8URSlTPKdOnaJ27drW76dPn06TJk2s3/fp08c6U8W/v15//XUAGjduTEZGBidPnrS+bteuXbRo0aLE8cqzbVkYDIYSl7+EAH3lmt7yzGKxcN9996EoCjNmzCj2M9g6R0XVJ7l29W1bt27NypUrSU5O5s8//+TIkSPFGn90x8l9qoXQtZ9++kmNjIxUk5OT1QULFqhhYWFq165dy72fIUOGqA8++KCalZWlLliwQA0KClKTk5PLvW1+fr6anZ2tPvTQQ+qLL76oZmdnqwUFBaqqqurFixfVv/76S83JyVFzc3PV999/Xw0LC1NTUlIq/gsQwgH0lmcPPfSQ2q1bNzU7O7vSxxFCT/SWa9u3b1ezs7PVzMxM9e2331Y7duxY7LXXOuc5gxTNQlxDdna2OnToUNXf319t1aqVunPnTrVp06bqk08+Wa79nD17Vu3Tp4/q5eWlNm7cWF2yZIn1uVtvvVV97bXXyrTt5MmTVaDY1zfffGN9Xbt27VRfX1+1Ro0aas+ePdUtW7ZU6ucXwhH0lGdHjx5VAdXT01P18fGxfq1atapMxxFCz/SUa6qqqk899ZQaEBCg+vn5qYMHD1ZPnz5d7DjXOuc5g6Kqcu1WCCGEEEKIa5E+zUIIIYQQQpRCimYhhBBCCCFKIUWzEEIIIYQQpZCiWQjh8qKionjggQes369YsQJFUVixYoXTYnKE6vJzCv2orrkmBEjRXGHTp09HURQ2b97s7FCEqNJ27tzJkCFDqFevHp6entSpU4devXrxySeflGs/M2fO5MMPP6xULEePHmXkyJE0bNgQT09PatWqRbdu3Zg8eXKl9iuEHugl14oK8Tlz5lR4H0LYg6wIKITQrbVr19KzZ08iIyMZPXo0tWrV4vjx46xfv56PPvqIcePGAbB//34Mhmu3AcycOZNdu3bx1FNPVSiWQ4cOERsbi5eXF6NGjSIqKorTp0+zdetW3nrrLV5++eUK7VcIPdBTrgmhV1I0CyF067XXXiMgIIBNmzYRGBhY7LmzZ89a73t4eNg9lg8++ICMjAzi4+OpV6/eVWMRwhXpKdeE0CvpnmFH27Zto0+fPvj7++Pr68tNN93E+vXri21z4cIFnnnmGVq1aoWvry/+/v706dOH7du3F9uu6HLVTz/9xMsvv0ydOnXw8/NjyJAhpKamkpuby1NPPUVoaCi+vr6MHDmS3NxcR/64Qtjc4cOHadGiRYmTOEBoaKj1/r/7Wf5bjx49WLhwIYmJiSiKgqIoREVFkZGRgY+PD08++WSJ15w4cQKj0cgbb7xhjaVu3bolCuZ/xwIwf/58+vXrR3h4OB4eHjRs2JBXX30Vs9lcIq6WLVuyZ88eevbsibe3N3Xq1OHtt9++YjwDBw7Ex8eH0NBQnn76aclxYTN6yrWyOnLkCHfeeSdBQUF4e3vTqVMnFi5cWGybvLw8Jk2aRLt27QgICMDHx4cbbriB5cuXF9vu6NGjKIrCu+++S1xcHA0aNMDb25tbbrmF48ePo6oqr776KnXr1sXLy4vbb7+dCxculCte4fqkpdlOdu/ezQ033IC/vz8TJkzAZDLx5Zdf0qNHD1auXEnHjh0BLennzZvHnXfeSf369UlKSuLLL7+ke/fu7Nmzh/Dw8GL7feONN/Dy8mLixIkcOnSITz75BJPJhMFg4OLFi7z00kusX7+e6dOnU79+fSZNmuSMH18Im6hXrx7r1q1j165dtGzZssL7efHFF0lNTeXEiRN88MEHAPj6+uLr68ugQYOYPXs277//Pkaj0fqaWbNmoaoq9957rzWWv//+m2XLlnHjjTde83jTp0/H19eX8ePH4+vry7Jly5g0aRJpaWm88847xba9ePEit956K4MHD2bo0KHMmTOH5557jlatWtGnTx8AsrOzuemmmzh27BhPPPEE4eHhfPfddyxbtqzCvxMhLqenXCuLpKQkunTpQlZWFk888QTBwcHMmDGDAQMGMGfOHAYNGgRAWloaX3/9NXfffTejR48mPT2dqVOn0rt3bzZu3EibNm2K7feHH34gLy+PcePGceHCBd5++22GDh3KjTfeyIoVK3juuees595nnnmGadOmVfh3JVyQ09YidHHffPONCqibNm264vMDBw5U3d3d1cOHD1sfO3XqlOrn56d269bN+lhOTo5qNpuLvTYhIUH18PBQX3nlFetjy5cvVwG1ZcuWal5envXxu+++W1UURe3Tp0+xfXTu3FmtV69eZX5E4SBnz55V+/btq3p7e6tNmjRR//777wptu23bNrVLly6qn5+fWr9+fXXKlCkVPo5e/PXXX6rRaFSNRqPauXNndcKECeqff/5ZLAdUVVXr1aunjhgxwvp9Ub4sX77c+li/fv2umBN//vmnCqh//PFHscdbt26tdu/e3fr9rl27VC8vLxVQ27Rpoz755JPqvHnz1MzMzBL7zMrKKvHYI488onp7e6s5OTnWx7p3764C6rfffmt9LDc3V61Vq5Z6xx13WB/78MMPVUD96aefrI9lZmaqjRo1KvFziiuzVZ6NHj1arVWrlurn56e2bNlS/e233yp8HD3RU64V7fPnn3++arxPPfWUCqirV6+2Ppaenq7Wr19fjYqKsp5XCwoK1Nzc3GKvvXjxohoWFqaOGjXK+lhCQoIKqDVr1lRTUlKsjz///PMqoEZHR6v5+fnWx++++27V3d29WD4LTXlyoHv37qqHh4d1ufpbb721zPspLRftQYrmCrpW0VxQUKB6e3urQ4cOLfHcI488ohoMBjU1NfWKrzt//rx67tw5tXXr1urAgQOtzxX9E3n77beLvaboZPrvfy5PPfWUajAYiiW50Kc777xTHTVqlJqZmanOnz9fDQoKUpOTk8u9bcuWLdWXX35ZNZvN6pYtW1RfX191z549FTqOnmzcuFEdNGiQ6u3trQLWE9v8+fOt21TmRG42m9Xw8HB1+PDh1sd27typAiU+eOzfv18dPny4GhgYaI3F19dX/eqrr64af1pamnru3Dn1+++/VwE1Pj7e+lz37t1VX19f1WKxFHvNgAED1JiYGOv3t9xyi1q7du0S27399ttSNJeRrfJs79691kJp48aNqr+/v3r+/PkKHUdv9JJrZSmamzRponbo0KHE42+88YYKqDt37rzi8ZOTk9Vz586p/fr1U9u0aWN9rqhofvzxx4u9Zt68eSqgvvPOO8UeLzr3Xt4wJjTlyYHu3bur3333XYX2U1ou2oP0abaDc+fOkZWVRdOmTUs816xZMywWC8ePHwfAYrHwwQcf0LhxYzw8PAgJCaFmzZrs2LGD1NTUEq+PjIws9n1AQAAAERERJR63WCxX3Icon9zcXJ577jnq1KlDREQEa9euZejQoSX6xFVERkYG8+bN4+WXX8bb25sBAwbQqlUr5s+fX+5tjx49yt13343BYKBt27Y0a9aMffv2lfs4ehMbG8uvv/7KxYsX2bhxI88//zzp6ekMGTKEPXv2VHr/BoOBe++9l3nz5pGVlQVol2g9PT258847i23bpEkTvvvuO86fP8+OHTt4/fXXcXNz4+GHH+bvv/+2brd7924GDRpEQEAA/v7+1KxZk+HDhwOUyMm6deuiKEqxx2rUqMHFixet3ycmJtKoUaMS213pf4yrcpU8u+6666yD4RRFIS8vj5MnT5b7OHqkp1wrTWJi4lXPsUXPF5kxYwatW7fG09OT4OBgatasycKFCyt9jgWK5amr0EuuVXY/18pFe5Gi2clef/11xo8fT7du3fj+++/5888/WbJkCS1atMBisZTY/vJ+YGV5XFVVm8Zrd6oKeZn2/yrH7+U///kPhw8f5uDBgzz11FM8++yznDlzhp49e17zdbfddhuBgYFX/HrzzTcBOHjwIL6+vtStW9f6ulatWrF79+4S+ytt23HjxvH9999TUFDAxo0bOXbsGJ06dSr3cfTK3d2d2NhYXn/9dT7//HPy8/P5+eefbbLv+++/3/pPWlVVZs6cyW233WY9Mf6b0WikVatWPP/888ydOxfQTv4AKSkpdO/ene3bt/PKK6+wYMEClixZwltvvQVQIq+dlrs6yzVXyTOAxx9/HC8vL2JjY7nxxhtp1apVuY+jZ3rKtcr6/vvveeCBB2jYsCFTp05l8eLFLFmyhBtvvNFx51jJtWvmwNNPP03NmjXp1asXO3bsKNd+rpaL9iIDAe2gZs2aeHt7s3///hLP7du3D4PBYP3UOmfOHHr27MnUqVOLbZeSkkJISIhD4tWV/Cx4Pbz07SrrhVPg7lPqZhcuXGDq1KkkJibi7e1NbGwszzzzDIsWLSr1tb///nup22RkZODv71/sMX9/f5KTk8u9bZ8+fbj//vt57bXXAJg6dSq1a9cu93FcQfv27QE4ffp0mV/z71bay7Vs2ZKYmBh++OEH6taty7Fjx8q8oMO/Y1mxYgXJycn8+uuvdOvWzbpdQkJCmWP9t3r16rFr1y5UVS32c1zpf0yZ6SjXXCnPAD777DM++eQTVqxYwa5du6x/k6qWZ6CvXLtcvXr1rnqOLXoetHNsgwYN+PXXX4vF5dAFiSTXrpoDb7/9Ns2bN8doNPLJJ5/Qp08f9u3bV+b9XC0X7UVamu3AaDRyyy23MH/+fI4ePWp9PCkpiZkzZ9K1a1frm8FoNJb4pPrzzz/b/RKDKJulS5fStm1b65RLeXl5tGnTxjqrQWX5+vqSlpZW7LG0tDR8fX3Lte2FCxfo168f77zzDrm5uWzdupXnn3+erVu3lvs4erJ8+fIrtuQU/YMvT/cEHx+fa3ZXuu+++/jrr7/48MMPCQ4OLvE3Xr16Nfn5+aXGUtQidXnceXl5fPbZZ2WO9d/69u3LqVOniq2QlpWVxVdffVXhfeqJq+TZ5YxGIzfddBN///239T3gqnkG+sq1sujbty8bN25k3bp11scyMzP56quviIqKonnz5sCV83HDhg3FXled6CnXADp06ICvry9eXl5MmDABPz8/1q9fX679XCkX7UVamitp2rRpLF68uMTjL730EkuWLKFr1648/vjjuLm58eWXX5Kbm1tsDtbbbruNV155hZEjR9KlSxd27tzJDz/8QIMGDRz5Y+iHyVv7tOyI45TBqVOnrK21oE0l1qRJEwAWL17MwYMHefjhh4mMjGTt2rUcPnyYHTt28Mwzz9CnTx9Wr159xf2+8MILvPDCCzRu3JiMjAxOnjxJnTp1ANi1axf3339/iddca9vDhw/j4+PDkCFDAGjdujVdunRh5cqVtG3btlzH0ZNx48aRlZXFoEGDuO6668jLy2Pt2rXMnj2bqKgoRo4cWeZ9tWvXjtmzZzN+/HhiY2Px9fWlf//+1ufvueceJkyYwNy5c3nssccwmUzFXv/WW2+xZcsWBg8eTOvWrQHYunUr3377LUFBQdbVz7p06UKNGjUYMWIETzzxBIqi8N1331Wqu8Xo0aP59NNPuf/++9myZQu1a9fmu+++w9u7bO/jK9JRrrlKnl1JQUEBhw4dqtBr9URPuVbkl19+sbYcX27EiBFMnDiRWbNm0adPH5544gmCgoKYMWMGCQkJ/PLLL9ZVC2+77TZ+/fVXBg0aRL9+/UhISOCLL76gefPmZGRklPO3VEGSa2XOAYPBgKqqFdrP5bloN3YdZliFFc2ecbWv48ePq1u3blV79+6t+vr6qt7e3mrPnj3VtWvXFttPTk6O+p///EetXbu26uXlpV5//fXqunXr1O7du5dpCp6rzeIxefJkFVDPnTtnt99BdfDTTz+pkZGRanJysrpgwQI1LCxM7dq1q6qqqrpu3Tr19ddfV7/77ju1X79+6tatW9Xhw4eXe/TukCFD1AcffFDNyspSFyxYcM2RxlfbNiUlRfX391fnzZunWiwWdffu3WpoaKi6ZMmSCh1HL/744w911KhR6nXXXaf6+vqq7u7uaqNGjdRx48apSUlJ1u3KMqI/IyNDveeee6wzX1xpdH/fvn1VoESeqqqq/vPPP+qYMWPUli1bqgEBAarJZFIjIyPVBx54oMQI+n/++Uft1KmT6uXlpYaHh1un7/p3TN27d1dbtGhR4lgjRowoEV9iYqI6YMAA1dvbWw0JCVGffPJJdfHixVVi9gxXyrMffvhBTU9PV/Pz89WffvpJ9fDwULdt21ah4+iJnnKtaJ9X+yqaZu7w4cPqkCFD1MDAQNXT01Pt0KGD+vvvvxfbl8ViUV9//XW1Xr16qoeHhxoTE6P+/vvvJXKsaPaMf8+SUd5zr97pKdcuXryo/vXXX2pOTo6am5urvv/++2pYWJh1yr9r7acsuWgPUjQLcQ3Z2dnq0KFDVX9/f7VVq1bqzp071aZNm6pPPvmkunfvXvWFF15Q77nnHvWbb75RZ8+erY4cObLcxzh79qzap08f1cvLS23cuHGxQvfWW29VX3vttTJtu3jxYjU6Olr19fVVIyMj1TfffLPMxxGagQMHqg0bNnR2GNWOq+RZamqq2qNHDzUgIED19/dX27Ztq/7yyy9lPo64RHLNOfSUa2fPnlXbtWun+vr6qjVq1FB79uypbtmypUz7KUsu2oOiqq42vYIQ+nDmzBmGDh1Kv379aNKkCe+99x5vvvkmXbt2dXZoogJOnz5NvXr1ePHFFx07SEhck+RZ1SO5pk+Sa6WrlkVzVFQU/v7+GAwGatSoYZO5CUX1k5OTQ0BAAMePH2fnzp089dRT7Ny509lhiXJKSEjgn3/+4euvv2bTpk0cPnyYWrVqOTssUUjyrOqQXNM3ybXSVduiedeuXS4xolkIYV/Tp09n5MiRREZG8t5771kHUwohbEtyTbg6KZqFEEIIIYQohcvN07xq1Sr69+9PeHg4iqIwb968EtvExcURFRWFp6cnHTt2ZOPGjcWeVxSF7t27Exsba13FSwghhBBCiKtxuaI5MzOT6Oho4uLirvh80dyQkydPZuvWrURHR9O7d2/Onj1r3WbNmjVs2bKF3377jddff926bKMQQgghhBBX4tLdMxRFYe7cuQwcOND6WMeOHYmNjeXTTz8FwGKxEBERwbhx45g4cWKJfTz77LO0aNGCBx54oMRzubm55ObmWr+3WCxcuHCB4OBguy/VKIQjqapKeno64eHh1kUBHElyTVQXkmtC2J/d8szuk9rZEaDOnTvX+n1ubq5qNBqLPaaqqnr//ferAwYMUFVVm3Q9LS1NVVVVTU9PV9u2batu3LjxivsvWiBEvuSrunwdP37cLrlaGsk1+apuX5Jr8iVf9v+ydZ5VqZbmU6dOUadOHdauXUvnzp2t202YMIGVK1eyYcMGjhw5wqBBgwAwm82MHj2aJ5988or7//cn8tTUVCIjIzlw4ABBQUH2+8HKID8/n+XLl9OzZ8+rLkHq8DiykzB9dQM5qhu/3Licoe3qOi8WvfxOXCSO9PR06tevT0pKCgEBAQ6K8BJn5NqJ/Vupv+AOkvHH/5ntpW6vl7+rXuMB/cWkt3gALly4QJMmTapFrinbf8RtyfNYorpjHjLd+vjPW07w+h8H6NwgiM/uaVPidXr7u+ktHtBfTHqLx1555mazPbmIBg0asH176SdIAA8PDzw8PIiLiyMuLg6z2QxAUFAQwcHB9gyzVPn5+Xh7exMcHOz0wswaR2Y+eCh4qAY8vP0d/jvS5e/EReIoet5Zl2edkWuZgQH4eyjkYSjTMfTyd9VrPKC/mPQWz+WqRa5FtQQPBfJOw2X77tTMiGHZCQ6lWAgKCirxu9Db301v8YD+YtJbPEVsnWcuNxDwWkJCQjAajSQlJRV7PCkpqVITqI8ZM4Y9e/awadOmyoZYtVnfnCp5ZotTQxGuyZG5Jt03RXXmkFyrUU+7TUkEy6VzwnW1/DAaFJIz8ziTlmO/4wthY1WqaHZ3d6ddu3YsXbrU+pjFYmHp0qXFumsIe9GqEAMqeQVSNAudK6yapXYWwk7864JiBHMeZJyxPuxpMtI4VFsnYdfJNGdFJ0S5uVzRnJGRQXx8PPHx8YC2LGd8fDzHjh0DYPz48UyZMoUZM2awd+9eHnvsMTIzMxk5cmSFjxkXF0fz5s2JjY21xY9QdRmM2g0WKZpFhUiuCeEYDsk1oxsEFI5tuZhY7KkW4Vo/010nU+13fCFszOWK5s2bNxMTE0NMTAygFckxMTFMmjQJgGHDhvHuu+8yadIk2rRpQ3x8PIsXLyYsLKzCx5TuGWWkaEWzUZHuGaJiHNo9o6iN2WWHQgtRcQ7LtaIuGhePFnu4ZR1/AHafkqJZuA6XGwjYo0cPSpvwY+zYsYwdO9ZBEQkr5dJnsPz8AicGIkTppFuGEA5QIwoSVmn9mi/Tsk5RS7N0zxCuw+Vamp1BLhmX0WUTiOcVSNEsys+RuaYYtLJZGppFdeSwXAssamkuXjQ3q+2PosCZtBzOpede4YVC6I8UzWUg3TPKqLB7BkhLs6gYZ3TPUKRsFtWQ47pnRGm3F44Ue9jXw436IT4A7DyZYt8YhLARKZqF7VzWPaPALEWzEEJUe8GNtNsLh0s81SYiEID4YymOi0eISpCiuQyke0YZGS61NBcU5DsxEOGqHJprMlGzqMYclmvBDbXbzHOQnVLsqZjConnb8eKPC6FXUjSXgXTPKCPl8qJZWppF+UmuCeEYDss1Dz/wLVxc7F+tzTGRNQCIP56CxSLdpIT+SdEsbOfylmbpniF0ThqahXCQoi4aycWL5qa1/PBwM5CeU8CR85lOCEyI8pGiWdjOZX2azdLSLHSuqGiWgYBC2FlRF43kQ8UeNhkNtK6rTT0XL100hAuQorkMpE9zGSkKauGMBNI9Q1SEY/s0a//+pMFZVEcOzTVrS/OhEk8VDQbcduyi/eMQopKkaC4D6WdZdmphv+b8ArOTIxGuyDlTzsnqlaL6ceh57RpF8+X9moXQOymahW0Vtt5J9wyhewZpaRbCIS7v0/yvFX2LWpr3nUknO08aW4S+SdEsbKuwaM6XolnonFL4788gLc1C2FeNKO3ckJcBGUnFnqod4EmYvwdmi8rOk6nOiU+IMpKiWdhW4QwaefkyT7PQN0VamoVwDDf3S8tp/6uLhqIoxERoXTQ2J15wdGRClIsUzWUgAwHLobBPc64soy0qwBm5Jn2aRXXk8Fwr6qJx/mCJp2LrBwGwKUGKZqFvUjSXgQwELIfC1juLuQCzTFYvysmhuaYUDQQUovpx+HntGoMBO0RpRfPmxIty3hC6JkWzsCmlsHuGAZWsPGltFjpmkD7NQjhMzSba7bl9JZ5qVtsPXw830nMK2HcmzcGBCVF2UjQL2yrsnmHEQna+jIQW+mW8bJ5mVZXWLSHsqmYz7fZsyaLZzWigbT2tX7N00RB6JkWzsClFudR6J9MHCT0ruiqiYEGuCAthZ6HXabdpJyCnZGtyx8J+zRuPStEs9EuKZmFb1u4ZFrKkaBY6Zrhs9gzpRymEnXnVAN9a2v1z+0s8HVvYr3ljwkW58iN0S4pmYVuXdc+QolnomdF46aqIRU7SQthfUWvzub0lnmpdNwB3NwPnM3I5mpzl4MCEKBspmstAppwrB+vgKlW6Z4hyc2SuKZfNniFFs6hunHJeu0a/Zk+TkTZ1AwFtFg0h9EiK5jKQKefK4fI+zTIQUJSTI3PNUNSVSFGle4aodpxyXrtGSzNAh6L5mo9K0Sz0SYpmYVsGNwDcsMiUc0LXjIZL//4sZimahbC7a7Q0w6WieX3CReTij9AjKZqFbRlMALgpBWTkStEs9MtwedFskasiQthdzababfopyE4p8XRsVBDuRgOnU3M4l+PY0IQoCymahW0ZtZZmE2bSc6RoFvp1edFsVmWBEyHszisQ/MK1+1eYQcPL3UjbeoEAHEiVtTqF/kjRLGzL6A6AG2bSc/KdHIwQV1c0EBDAYpaWZiEcopR+zdc3DAGkaBb6JEWzsK3C7hkmCkjLlpZmoWPKZS3NFmlpFsIhQptrt0m7r/j09Y21ovlgqiIDdIXuSNEsbKtY9wxpaRZ6dllLs4w6EsIxarXSbs/suuLTresE4OvhRpZZYc/pkisHCuFMUjQL2yoaCEgBadKnWejZZS3NqnTPEMIxrEXzTrjCFR43o4GO9WsAsPawLKkt9EWK5jKQxU3KobBPs0mRlmZRfg7Ntcv6NMtlYFHdOO28FtJEO0/kpUPK0Stu0qVhMABrjyQ7MDAhSidFcxnI4iblYO2eIX2aRfk5NtcuL5qlpVlUL047rxlNEFo4X/OZnVfcpHMDbb7mLYkp5MgiWUJHpGgWtmXtniEtzULnLu+eIVPOCeE4tVprt1cpmhvV9CHApJJbYGFjgnTREPohRbOwLeOloln6NAtdKzblnBTNQjhMKUWzoig0q6F1mVq+/6yjohKiVFI0C9sqLJrd0VYELJBiROiVTDknhHNcPhjwKpoHakXziv3nHBGREGUiRbOwrcu6ZwCkZEsXDaFXl1qaVenTLITj1Gqp3aadhMwrD/ZrGqBiMioknM8k4XymA4MT4uqkaBa2VdjS7GvSWgkuZOY5Mxohrs5gtN61SNEshON4+EFQA+3+mR1X3MTTDdrX06aeW75PumgIfZCiWdhWYUuzn7tWNCdnSNEsdEpRKCj8F2gukP73QjiUtYvGlYtmgO5NtNUBpV+z0AspmoVtFbY0+7lJS7PQPwtaa7M5X96nQjhU7Wjt9lT8VTfp0aQmABuOXCArTz7YCuerlkVzVlYW9erV45lnnnF2KFVPYdHsYy2ac50ZjRDXZC4qms1yQhbCoeq0025PbrnqJg1CvIkI8iLPbGHtIVnoRDhftSyaX3vtNTp16uTsMKqmwu4Z3oVFc7K0NAsdMyta0VxQIANWhXCo2m2025TEqw4GVBSFnk1DAVgmXTSEDlS7ovngwYPs27ePPn36ODuUqqmopdmoDayS7hlCz6zdM6RoFsKxvAIhuLF2/9TWq27W8zqtaF66NwmLLHcvnMyliuZVq1bRv39/wsPDURSFefPmldgmLi6OqKgoPD096dixIxs3biz2/DPPPMMbb7zhoIirIZMXAN4G7XK3tDQLPbMoRUWzdM8QwuHqtNVuT169aO7SMBhfDzeS0nLZfiLFMXEJcRUuVTRnZmYSHR1NXFzcFZ+fPXs248ePZ/LkyWzdupXo6Gh69+7N2bPaZZ358+fTpEkTmjRp4siwqxc3DwC8FK1YTs6QPs1CvyzW2TOkpVkIhytDv2YPNyM9mmoDAv/cneSIqIS4KjdnB1Aeffr0uWa3ivfff5/Ro0czcuRIAL744gsWLlzItGnTmDhxIuvXr+fHH3/k559/JiMjg/z8fPz9/Zk0adIV95ebm0tu7qWiLy0tDYD8/Hzy8517ki06vt7iUBR33ABPRfs+KTXHYTHq9XfiCnE4O1Zn5VpRS3N+Xm6px9HL37WI3uIB/cWkt3jA+bHo6bymhLXGDVBPbaUgL8+6tP2//243X1eT33ec5s9dpxl/UwMURbnaLu1Cz+8jvcSk13hsTVFV1SU7CSmKwty5cxk4cCAAeXl5eHt7M2fOHOtjACNGjCAlJYX58+cXe/306dPZtWsX77777lWP8dJLL/Hyyy+XeHzmzJl4e3vb5OeoasIvbiT26Kec9mpK54uTcTeovN3BjIP/x4lyysrK4p577iE1NRV/f3+HH99ZudYh/hlqq2eJC5lE3YhGdjuOEEWqa65dicGSR7/tj2DAzF8t3ifbPeSK2+UUwAubjZhVheejC6glp19RCnvlmUu1NF/L+fPnMZvNhIWFFXs8LCyMffv2VWifzz//POPHj7d+n5aWRkREBD179iQ4OLhS8VZWfn4+S5YsoVevXphMJt3EoRw0wtFPCQ3wgouQZ1G44cZe+HvZP0a9/k5cIY6i1iZncVauXdj1AuRDVGQkt/Tpe81t9fJ31Ws8oL+Y9BYPQHKyc6dO09t5TUn6GM7s4MamAajNtBy80t9tUepWVh44T07IdfTt0cChMerxfaS3mPQWj73yrMoUzeX1wAMPlLqNh4cHHh4exMXFERcXh9mszQhhMpl08aYA/cRijcPDBwCjJY8a3iYuZuVzPstMsL/jmgZ09ztxgTicHaezcs2iFP4LVC1lPo5e/q5F9BYP6C8mPcXj7Dh0d16r2x7O7MDtTDy0HlLsqctj6tOyNisPnOfvfed4sldTx8eJvt5HRfQWk17isVcMLjUQ8FpCQkIwGo0kJRUfKJCUlEStWrUqte8xY8awZ88eNm3aVKn9VAuFs2dQkEOYvycAp1OznRiQcCWOzjVVkcVNRPWkm/Nanfba7fGN19zs5uZhGBTYeTKVkylyThHOUWWKZnd3d9q1a8fSpUutj1ksFpYuXUrnzp2dGFk146YVyuTnUDtAu38mNceJAQlxdapBK5pVKZqFcI7IwoXGTm2D/KufK0J8PWhfLwiAP3aedkRkQpTgUkVzRkYG8fHxxMfHA5CQkEB8fDzHjh0DYPz48UyZMoUZM2awd+9eHnvsMTIzM62zaVRUXFwczZs3JzY2trI/QtVnbWnOplaAdv9MmhTNomwcnWtqYfcMi1kfI76FcBTdnNeCGoBPTTDnwen4a27ar3VtABbskKJZOIdLFc2bN28mJiaGmJgYQCuSY2JirFPGDRs2jHfffZdJkybRpk0b4uPjWbx4cYnBgeWlm8tYrqBwnubLW5pPp0jRLMrG4blmkGW0RfWkm/OaokBER+3+sfXX3LRPq1oYFNh+PIVjyVkOCE6I4lyqaO7Roweqqpb4mj59unWbsWPHkpiYSG5uLhs2bKBjx47OC7g6crvUpzm8sGg+kSL/3IROGbSWZlncRAgnKuqicXzDNTcL9fOkc0Ntho8FO07ZOyohSnCpotlZdHMZyxWYCvs0o1K/hlaQHD0vRbMoG4d3zzC6A2DOl+XeRfWiq/NaRGHRfGw9lLJ0xIDocAAWbJeiWTieFM1loJvLWK6gaCAgEOmvvb1Op2aTW2B2VkTChTg81wqLZrVAuhCJ6kVX57Xa0dq5I/sCnD94zU17t6iFyaiw70w6B5PSHRSgEBopmoVtGd2hcBqvEA8z3u5GLCqcuChTBAkdumy2FyGEk7i5Q5122v3j1+7XHOjtTrfGNQEZECgcT4rmMtDVZSy9UxRw99Xu5mUSGaQtaiKDNkRZODzX3Apbms3SPUNUL7o7r1kHA167XzNA/8IuGr9vP4VaSncOIWxJiuYy0NVlLFfgoRXN5KZTL1grmo8mZzoxIOEqHJ1rSlFLs3TPENWM7s5rRYMBj60rddObm4fh4WbgyPlMdp5MtXNgQlwiRbOwvcKWZvIyiArWltVOlJZmoUPGwoGrSoG0NAvhVBEdQTHAhcOQdu1uF74ebtzSQlvp99etJx0RnRCAFM3CHqwtzRlEhWhF8+FzGU4MSIgrM7gXFs2WXCdHIkQ15xWoDQgElMTVpW4+pF1dAObFn5SB5sJhpGgWtndZS3OTMO3+wSQpmoX+GEzaYjwG6dMshPNF3QCAIfGfUjft2iiEMH8PUrLyWb7vrL0jEwKQorlMdDdgQu/cL/VpbhzmB2hLaadmyQIS4tocnWvGwmXfDWZpaRbViy7Pa/W7AaAkril1U6NBYXBbrbV5zpYTdg1LiCJSNJeB7gZM6J3HpZZmf0+TdWXAA2dlTk1xbY7ONXcP7b1psORhscgofFF96PK8FtkJFCNKSiJeuedK3fyOwqJ5+f5znEuXD77C/qRoFrbnfqlPM0CTWlpr8/4zUjQLffHw0mZ3cSefzLwCJ0cjRDXn4Qd12gJQM2NvqZs3CvWlTUQgZovK/HgZECjsT4pmYXuXtTQDNC3sonFAVm8SOuNWOHuGB/mk5UjRLITTFXbRCEkvvWiGSwMCf958QuZsFnYnRXMZ6LLvl565a0UyuVqR3KSwaN4nLc2iFI7ONcVda2n2Ipe0bOlzL6oP3Z7XCgcDhmTshTIUwf1bh+PuZmB/Ujq7TqbZOzpRzUnRXAa67PulZ/9qaW4e7g/AnlNp0m9UXJPDc81dmxLRW8mRollUK7o9r0V0RDWY8Mq/oM3ZXIoAbxO3Fs7ZPHPjMXtHJ6o5KZqF7XkGarc52kpNjUN98TIZycgt4Mh5mXpO6Ehh/3sfcqV7hhB64O6NWriktuHIijK95J6OkQDMjz9Jeo58+BX2I0WzsD2vGtpt1gUA3IwGWtbRWpvjj8uSp0JHrEVztrQ0C6ETaoMbAVCOLC3T9h3rB9Gwpg9ZeWbmxZ+yZ2iimpOiWdheUdGcfdH6UHTdQAC2H09xfDxCXE1hVyJvJUdaqITQCUvDwqI58R8oKH0qOUVRuKdjPQBmbjgmAwKF3UjRLGzPO0i7zU6xPhQdEQjAjhMpJTYXwmkK+zT7kkNqtnTPEEIXQluQ4xaAkp8Fx9aV6SV3tK2Dh5uBvafT2CaNM8JOpGgWtlfU0pybCmatEGlTWDTvOZ1GTr7ZSYEJ8S+F3TO8lDwuZmQ5ORghBACKwln/Vtr9Q2XrohHo7c5trcMB+GG9DAgU9iFFcxnodmoevSoaCAjWwYB1a3gR5u9Bvlll67GLV36dqPYcnmtFC/EA6WnS315UH3o/r531K1/RDHBvJ21A4O87TpGSlWePsEQ1J0VzGeh2ah69MrqBhzbwr6hfs6IodGoQDMD6IxecFZnQOYfnmpsHFsUIQEaGFM2i+tD7ee2cf0tUFDi7G9JOl+k1MRGBNKvtT26BhZ82H7dzhKI6kqJZ2IdXoHZ72WDAS0VzshMCEuIKFAWLSevXnCNFsxC6kefmhxoeo31zuGytzYqi8EAXbUDgjLWJFJgt9gpPVFNSNAv78CoaDFiyaI4/liL9moVuqIVXRcxZ0m1ICD1RG/TU7hxcUubX3N6mDkE+7pxMyeavPUl2ikxUV1I0C/uwztV8qVU5KtibMH8P8swWtiZKgSL0QSmc7cWUl0pugXyYE0Iv1Ea3aHcOL4OCsvVR9jQZubdwsZNpaxLsFZqopqRoFvbhG6rdZp6zPqQoCp0LW5tXHzrvjKiEKMHoo70na5BOcoYMHhJCL9TwGPAJhdw0SPynzK8b3qkeJqPC5sSLsjaAsCkpmoV9FBXNGcUvj/Voqj2+fN9ZR0ckxBUphVdFaigZnEnLcXI0QggrxQBNb9Xu7/+jzC8L8/e0Tj/3zT/S2ixsR4pmYR++Ydrtv4rm7k1qYlBg35l0TlyUeXGFDhR2zwhQMjh5MdvJwQghimnaV7vdvwjKsdLfqOvrA/D7jtMkyYdhYSNSNAv7uErRXMPHnXb1tJY9aW0WulA4aLUGGZxMkaJZCF2p3x3cvCD1OCTtKvPLWtUNIDaqBgUWlelrj9ovPlGtSNFcBnqfBF6XrN0zShbGN16nFdTLpGgW/+KUXCtsaa6hZMjVD1FtuMx5zd0bGhbOolGOLhoAD93QAIDv1yWSlpNv68hENSRFcxnofRJ4XfKtpd1mlJzy56ZmWkH9z6Fk+UcminFKrhW2NAeSzgnpniGqCZc6rzXto93uX1Sul/VqFkajUF/ScwuYuUGW1haVJ0WzsI+ilubsi1CQW+ypxqG+NA71Jc9s4c9dZ5wQnBCXKZw9I1hJkz7NQuhRk1sBBU5tg7RTZX6ZwaDwaPeGAExdkyDrA4hKk6JZ2IdXDTC6a/f/1UVDURT6R2sjmxfsKNvyqELYjV9tAMKUi5y4mI1ajsFGQggH8A2FiA7a/b0LyvXSAdHhhAd4ci49l1+2nrBDcKI6kaJZ2IeigF9hF420kyWeLiqa/zl0nuSM3BLPC+EwhUVzkJKBOT+H06ky0l4I3Wl+u3a7e165XubuZmB0N61v85crj8jS2qJSpGgW9hOgrcpEyvEST9UP8aFVnQDMFpVF0kVDOJNXDTB6ABCqpHDwbIaTAxJClFBUNB9bB+nlO2cMi42ghreJYxey5HwjKkWKZmE/gRHabeqVB2AMKGxtnrNFLpkJJ7rsqkgYFziYlO7kgIQQJQTUhbqxgAp7fivXS73d3RhZOG/zZ8sPSRcsUWFSNAv7CSgqmq9cFA9qWweTUWH78RT2nEpzYGBC/Mtl/ZoPSUuzEPrUfKB2u2deuV96f+d6+Lgb2XcmnT93l5zVSYiykKJZ2E9RS/MVumcAhPh6cEtzrYXvx00yHZBwosKW5lrKRemeIYReFXXRSFxb7i4agd7ujOqqtTZ/+PcBLBZpbRblJ0WzsB9rS/OVi2aAuzto/Z7nbj1JVl6BI6ISoqRA7X0YoZzlwJl0OaEKoUeBEVCnPaCWexYNgIe6NsDPw419Z9L5Q/o2iwqodkVzSkoK7du3p02bNrRs2ZIpU6Y4O6SqK+Cyluar9CHr0jCYyCBv0nMLmLut5CwbQjhEjSgAogxnSc8t4GhypnPjEUJcWYuB2m05Z9EACPA28eANl1qbzfLhWJRTtSua/fz8WLVqFfHx8WzYsIHXX3+d5ORkZ4dVNQXUBRTIz4TM81fcxGBQeKBLFABTVydIC59wjiDtRNrYpL1Pd5xIdWY0QoirKerXnPgPpJa/oWVU1/r4e7px8GwGv+8o+0IpQkA1LJqNRiPe3t4A5ObmoqqqjKS1F5PnpX7NyYeuutnQ2Aj8Pd04cj6TpfvOXnU7IeymhlY017IkoWAh/niKc+MRQlxZYATUux5QYefP5X65v6eJhwvnbf5o6UFpbRbl4nJF86pVq+jfvz/h4eEoisK8efNKbBMXF0dUVBSenp507NiRjRs3Fns+JSWF6Oho6taty7PPPktISIiDoq+Gghtpt8kHr7qJr4cb93SsB8CUVUccEZUQxQXUBcWIm5pHGBfZcSLF2REJIa6m9VDtdsdPFXr5A9fXJ9DbxJFzmfy2XboFirJzc3YA5ZWZmUl0dDSjRo1i8ODBJZ6fPXs248eP54svvqBjx458+OGH9O7dm/379xMaGgpAYGAg27dvJykpicGDBzNkyBDCwsJK7Cs3N5fc3Eur1aWladOi5efnk5+fb6efsGyKjq/3OAw1GmJkGeZzB7BcI9Z7O9Rh6pojbDx6gbUHzxIbVcPmsTiKK8bh7Fj1kGtugZEoFxOobzjDllM1ycjKwcNktD6vl79rEb3FA/qLSW/xgPNj0UOulabUv1vjfrgZn0U5u5v8E/EQ1qJc+/cwwEPXR/HukoO899cBbrmuJu5uV29D1PP7SC8x6TUeW1NUF+6boCgKc+fOZeDAgdbHOnbsSGxsLJ9++ikAFouFiIgIxo0bx8SJE0vs4/HHH+fGG29kyJAhJZ576aWXePnll0s8PnPmTGsXD3Ft9c8tofWJ7zgd0JaNDZ665razjxhYm2Sgkb+FcS1kqVNHysrK4p577iE1NRV/f3+HH18Pudbh8AfUTtvG65YRfJXXm7HNzTQOcNl/j0KnJNdsIzbhE8JTNnEwtA976txd7tfnmuH/thlJy1cYHGWme23J9arEXnlWpYrmvLw8vL29mTNnTrFCesSIEaSkpDB//nySkpLw9vbGz8+P1NRUrr/+embNmkWrVq1K7P9Kn8gjIiI4ffo0wcHB9v7xrik/P58lS5bQq1cvTCaTbuNQjqzAbdYQ1JAmFDyy9pr7Op2aw00frCbfrPLtyHZ0blC+37Gr/E70GEdaWhohISFOO5HrIdcMy/8P49oPWRNwG8OT7mFsjwY8eVMj6/N6+bvqNR7QX0x6iwcgOTmZ2rVrV+tcK01Z/m7K/kW4zbkf1bcWBeO2g8F4xe2u5cdNJ/jfb3uo4W1i6dNd8fO88rH0+D7SW0x6i8deeeZy3TOu5fz585jN5hJdLcLCwti3bx8AiYmJPPzww9YBgOPGjbtiwQzg4eGBh4dHicdNJpMu3hSgn1iuGkdYUwCUCwmYDAoYr/6WiwwxcW/Hekxfe5QPlx7mhiZhKIpiu1gczJXicHacusi12tr/gesM2gqWG45evOKx9fJ3LaK3eEB/MekpHmfHoYtcK6NrxnTdreAZiJJxBtOJddCwZ7n3f3fHekxfl8jhc5lMXXuMZ3tfV/F4nERvMeklHnvF4HIDASurQ4cOxMfHs337dnbs2MEjjzxS6mvi4uJo3rw5sbGxDoiwivGvC+6+YMm/5mDAIo/3aIiHm4Gtx1L4a48sdVrdODXXQpsBEJR5GFDZdiyF9Bx99M8TwtZc/rzm5gEtC8c1VXBAoJvRwIRbtUJ56poEzqTm2Co6UUVVqaI5JCQEo9FIUlLxYispKYlatWpVeL9jxoxhz549bNq0qbIhVj8Gw6VBGmd2lbp5qL8nDxVOPv/awr3kFpjtGZ3QGafmWnBjMJgw5KXTJTiTAovKygPnHB+HEA5QJc5rrYdpt3vmQ256hXZxS/Mw2terQU6+hQ//PmDD4ERVVKWKZnd3d9q1a8fSpUutj1ksFpYuXUrnzp0rvF+X/0TubGEttduknWXa/PEejQjz9+DYhSymrkmwY2BCb5yaa27uUEvrojG0tvbBe4lc7RBVVJU4r0V01KY1zc+E3XMrtAtFUXi+r9ba/NPm4xxIqljxLaoHlyuaMzIyiI+PJz4+HoCEhATi4+M5duwYAOPHj2fKlCnMmDGDvXv38thjj5GZmcnIkSMrfMwq8YncmWoVFs1nylY0+3i4MbGP9k/s02WHSEqTS2bVhdNzrW57ADq7a/OFL9t3lnyzzOQiqh6n55otKAq0vV+7v/XbCu+mXb0gbm1RC4sKr/6+RxY8E1flckXz5s2biYmJISYmBtCK5JiYGCZNmgTAsGHDePfdd5k0aRJt2rQhPj6exYsXX3EeZuEgYYUDLcvQPaPI7dF1iIkMJCvPzP8t3GunwIT4l7paq1to6i5CfN1Jzylg7eFkJwclhLiq6LvB4AYnNkHSngrv5oW+zXA3Glh98LxcYRJX5XJFc48ePawzX1z+NX36dOs2Y8eOJTExkdzcXDZs2EDHjh0rdcwqcRnLmcKaAwpknoWMsi2TbTAovDKgJQYFFmw/xV+7z9g3RqELTs+1wpZm5cwObmuuTb81d+sJ58QihB05PddsxTcUmvbR7leitTky2Ns6nub/Fu4lJ1/G04iSXK5odoYqcRnLmdx9Li2nfSq+zC9rVTeAh7s1BOC/83aRmiUzGVR1Ts+1GvXBOwTMudwbqbUwL959RmbREFWO03PNltqO0G53/Aj5Fe/ON6ZnI0L9tPE00/6R8TSiJCmahWMUXvbmxMZyveypmxvToKYPZ9Nz+b+FFb/0JkSZKArUvwGARumbaVjTh5x8C4t2nnZyYEKIq2p4I/jXgeyLsO/3Cu9GxtOI0kjRLBwjorBoPl6+otnTZOTtO1qjKPDzlhMs3St9zYSdNegBgJKwkiHtIgD4abN00RBCtwxGiBmu3d8yvVK7Gtjm0niat/7YV/nYRJUiRXMZVJm+X84UUdiv/OQWsJSvr1j7qCBGXa/1NXvm5+3y6b8K00WuFRbNnNjEkJYBmIwKWxIvsvNkqvNiEsLGdJFrthRzHygGOLoazu2v8G4MBoXJ/bW1BX7ddpJNRy/YKkJRBUjRXAZVqu+Xs9S8Dtz9IC8Dzpa/m8WEW5vSvLY/F7PyeerHeMwWmRKoKtJFrtWI0r4sBdQ8t57bWocDMH3tMefFJISN6SLXbCkwApr21e5vnFKpXbWJCGRYe+0q04tzd8q0k8JKimbhGAYj1G2n3T++odwv93Az8sk9MXi7G1l3JJm45YdsHKAQlyk6+e773XqVY9GuM6TkOjEmIcS1xT6k3W6fBTlpldrVxD7XEeTjzoGkDKb9k2iD4ERVIEWzcJzIwlUZj/5ToZc3rOnLK7drC6V88PcBlu8r2/R1QpRbs/7a7f5FtKrlRYf6QRRYVP4+Jf8yhdCtBj0guLF2RXP7j5XaVQ0fd17o2wyAT1ccJll6BQqkaC6TKtf3y1nqd9duE1aCpWKXu4a0q8s9HSNRVXjix20knM+0YYDC2XSTaxEdwacm5KTC0dU8dVNjANYmKZxOlbOncH26yTVbUhTo8LB2f9MUqOTKfne0rUPH+kHk5FuYk2CQlQKFFM1lUeX6fjlL3fbg7gtZyXB2d4V381L/FrSrV4P0nAJGf7tZ5tCtQnSTawYjXHebdn/3r3RuGEyHqBqYVYXPVx5xbmxC2IBucs3Wou/SzjPnD2gNNJWgKAqvDWqJyaiwJ8XAX3vk6mZ1J0WzcByjCepdr90/sqLCu3F3M/D5vW0J8/fg0NkMHv9hqwzUELbXeqh2u3seSl4mT96kLbTz85aTHExKd2JgQoir8vTXltYG2PBVpXfXKNSPh7pGAfDqon3SSFPNSdEsHKtBYReNw8srtZtQf0+m3N8eL5OR1QfP8+LcnXLpTNhWZGcIaqj1j9w9lw5RQbSqYaHAovLygj3yfhNCrzqM1m73L4Lkw5Xe3ePdGxDsoZKUlsubMndztSZFs3Cshjdqt0fXQG5GpXbVum4gn94Tg0HRFp+IWyGXzYUNKcqlBRO2fQfAwCgL7m4G1hw6z5+7zzgxOCHEVdVsCo1vAVRY/1mld+dpMnJXQ+1q5g8bjrH20PlK71O4Jimay6BKDphwlprXaXPgmnPh8LJK7+6mZmHWGTU+WnaYNWeUSu9TOI/ucq3NPaAYtWkSz+4hxBMeuj4KgJcX7CFNLtUKF6W7XLO1LuO0220/QGZypXfXJEDlrti6ADz36w6y8goqvU/heqRoLoMqO2DCGRQFmvbT7u9fZJNdDu9UjzE9C/ubJhiZu+2UTfYrHE93ueZXC5ppAwKNG7QWq0e6RREZ5M3p1BxeXVD+hXqE0APd5ZqtRd0AtaOhIBs2T7XJLifc0oTwAE+OX8jm7cUVX3VQuC4pmoXjXVe4cMSBxWC2zaf1Z25pyn2dIgGYOHcXC3ectsl+haDLEwAou37BM+8C3u5uvHtnNIoCP285wdK9SU4OUAhRgqJYc5cNX0J+5aeK9PN04407WgMwY91RWWK7GpKiWTheRCfwqgHZF+HYOpvsUlEU/tunKR1rWrCo8OSP21i8S/qcChuo2x4iu6BY8mlw7k8AOtQP4sHClQIn/rqTc+myVKAQutN8IAREQNZ5bZVAG+jepCZ3tquLqsKEOTvIyTfbZL/CNUjRLBzP6Hapi8auOTbbrcGgcFdDC7e1qkWBRWXMzK38tl26aggb6PoUAPXPLYV07cPYM72b0jjUl3PpuYybtZUCmfZQCH0xukGnx7X76z6t8KJa//bf25oT6udBwvlM3vlTumlUJ1I0C+ewzoE7Fwps10pnUOCdO1oyuG0dzBaVJ3/cxs+bj9ts/6KaanwLljqxuKl5GNa8B2gj6j8f3g4fdyPrj1yQk6cQetT2PvAMgORDsO93m+wywMvEm3e0AmDqmgT+kdk0qg0pmsugyo8ydoaoG8C/jrZM8YE/bbprN6OBd4dEW5fbfnbODr5dd9SmxxD2odtcUxQsPf8LgCH+O7igTW/YKNSXd++MBuDLVUf4fYdc2RCuQbe5ZmsefhBbOG/zqncqvbR2kRuvC+PuDto4mv/8tJ3ULJlJpzqQorkMqvwoY2cwGKDVEO3+jtl22L3CawNbMrJwerBJ83fz1uJ9WCyyIIWe6TnX1HrXk+TXGsVSAItfsD7ep1VtHu7WAIDxP21nY4IMDhL6p+dcs7nOY8DkA2d22LSR5n+3NaN+iA9n0nL47/xdNtuv0C8pmoXztB6m3R74E7JsX2goisKk25rzn15NAPh8xWHG/xRPXoH0PRUVs6vu3agGExz4A/ZdmjJxQu+m9GoeRl6BhYdmbOKALLMthH54B0GHh7T7q962WWuzt7sbHwxrg9GgsGD7KebHn7TJfoV+Vapozs/P5/jx4+zfv58LF6R1RZRTWAttHk1LPmz73i6HUBSFcTc15p0hrXEzKMyLP8XI6RtlUQpRIRmedbB0GqN988dzkJcJaF2CPrk7hnb1apCWU8CIaRs5lZLtxEiFEMV0HgduXnByi00W1irSJiKQcTc2AuC/83ZxUvK+Sit30Zyens7nn39O9+7d8ff3JyoqimbNmlGzZk3q1avH6NGjq8flHmEb7R/UbjdPtdnI5iu5s30EUx+IxcfdyD+Hkhn6xTopakSFWLqOh4BISD0Gf79kfdzTZGTqiPY0rOnD6dQc7pmyntOp8h4TQhd8a0L7kdr9lbZrbQYY27MRbSICSc8p4D8/xWOWboBVVrmK5vfff5+oqCi++eYbbr75ZubNm0d8fDwHDhxg3bp1TJ48mYKCAm655RZuvfVWDh48aK+4RVXR6k7wCICLR2366f9KujepyexHOlPTz4N9Z9IZ8OkaNsvk9KK8TN7Q/wPt/sav4NDf1qcCvd359sGORAR5cTQ5i7u+Wi8fzoTQiy5PgNEDjq+Ho6tttls3o4EPhrXBy6TNpPP5ikM227fQl3IVzZs2bWLVqlVs3LiR//3vf/Tu3ZtWrVrRqFEjOnTowKhRo/jmm284c+YMAwcOZPVq270pRRXl7g0x92r3N35l98O1rBPA3Me7cF0tP85n5HH3lPX8uPGY3Y8rqphGN0OHR7T78x6HzGTrU3UCvfjx4c5EBHmRWFg4JyZnOilQIYSVf21oe792f8VbNm1trh/iwyu3twDg/SUHZEBwFVWuonnWrFm0aNGi1O08PDx49NFHGTVqVIUDE9VI+wcBBQ7+CWf32v1wdWt48+vjXejbqhb5ZpWJv+5k0vxd5MviFKI8er0MIU0hIwnmPVase1FR4Vwv2JtjF7K44/O17DqZ6sRghRCAtlCR0R0S18CRFTbd9ZB2dRkcUweLCk/M2saFzDyb7l84n8yeIZwvpBE0u027/89HDjmkt7sbcfe0tc6s8e26RO6buoHkDFkOWZSRyQvumAJuntoHvhWvF3u6TqAXPz/amea1/TmfkcewL9ex6sA5JwUrhAAgoO6lsTRLX7Fpa7OiKLw6sCUNCqehe/bn7ag23L9wvkoVzTk5OWzcuJHff/+d3377rdhXVVJtJoF3pq5Pa7c7f4YUx3SXKJpZ46v7Lq3q1u9j6efsTC6Xa7Wjof/H2v1V78DuecWeDvXzZPYjnejSMJjMPDMjp29i+j8JciIVTudyuWZLN4zX5m0+tdVmqwQW8fFw49N72uLuZmDpvrNMXZNg0/0L56pw0bx48WIiIyPp1KkTAwYMYODAgdavQYMG2TJGp6tWk8A7S512UL87WApg7ScOPfQtLWoxd8z1NKyptQ4M+2o9X606LIWNE7hkrkUPg85jtfvzHoNT8cWe9vM08c3IWOvS7i8t2MMLc3fKfOHCqVwy12zFNxQ6PabdX/Z/YDHbdPfNw/35323NAXhr8T7ij6fYdP/CeSpcNI8bN44777yT06dPY7FYin2ZzbZ9A4pq4obx2u2W6Q5rbS7SJMyP38Z2ZUB0OGaLyuuL9jH62y2yNKoom5tfhgY9IT8Lvr8DzhcfPe/hZuS9O6N5oe91KArM2nicu75aJ3O6CuEsXcaBZyCc26dd4bSx4R0jreNmxs7cKueSKqLCRXNSUhLjx48nLCzMlvGI6qx+d6jfDcx5sPwNhx/ex8ONj+5qw/8NbIm70cDfe5Po98lqtksrgSiN0Q2GztC6a2Sdh+8GQuqJYpsoisLD3RoybUQsfp5ubD2WQt+PVvPX7jPOiVmI6swrUBsUCLD8dSiw7aA9RVF4Y3BrIoO8OXExm6dmb8Mi8ze7vAoXzUOGDGHFihU2DEVUe4oCN72k3d8+C5L2OCEEheGd6vHr412ICPLixMVshnyxli9WHpZ/eOLaPAPg3l8guBGkHofvBkHm+RKb9bwulIXjbiC6bgCp2fk8/N0WXl6wm9wCuUInhEN1eBh8wyAlETZPs/nuA7xMfD68LR5uBpbvP8cny2T+ZldX4aL5008/5ddff+WBBx7gvffe4+OPPy72JUSF1G0HzQYAqjay2Ula1gng93E30KeldnntzT/2cd+0DZxJzXFaTMIF+NaE++aBfx04fwCm94P0ki3JkcHe/PxoFx7qWh+Ab/45ypDP13HobLqDAxaiGnP3gR4Ttfsr34TsizY/RIvwAF4f1AqAD5ceYPn+szY/hnCcChfNs2bN4q+//uKXX37hk08+4YMPPrB+ffjhhzYMUVQ7N00CxQgH/ii22pqjBXiZ+Ozetrx1Ryu8TNry27d+tIo/5XK6uJbACLj/N/AL1/pLftMHUo6X2MzdzcB/b2vO1/e3J9DbxM6TqfT9eA2frzhMgcwZLoRjxNwPNa/TCuZV79rlEHe0q8vwTpGoKjw5axvHkrPschxhfxUuml988UVefvllUlNTOXr0KAkJCdavI0eO2DJGUd2ENIaOj2r3Fz0LBc6bO1lRFIbFRvL7E11pWceflKx8HvluCy/M3Ul2nlxOF1cR0ghGLoLASLhwBKbdetWFe25uHsYfT95Aj6Y1ySuw8NbifdzxxToOJkmrsxB2Z3SDW/5Pu7/hSy1f7eB/tzWnTUQgaTkFPPr9FnLy5fzhiipcNOfl5TFs2DAMBlkfRdhBj4ngW0v7B7bW+d19Gtb05dfHrueRbg0AmLnhGLd9sprdp2SVN3EVQfVh5B8Q3BjSTsDU3pCw6oqb1g7w4psHYnlnSGv8PN3YfjyFfh+v4bMVh2SlSiHsrdHN0PBGsOTD3y/Z5RAebkY+H96WYB939pxO48W5u2RaUxdU4Yp3xIgRzJ4925axCHGJpz/0fk27v+o9uHjUqeGAdjn9+b7N+P7BjoT6eXD4XCaD4tby9eojMkhQXFlAXXjwL4jsDLmp8N1giJ95xU0VReHO9hEsebo7PZvWJM9s4e3F++n38Wo2HEl2cOBCVCOKorU2KwbYMx8S19nlMLUDvPjknhgMCvyy9QTT1x61y3GE/VS4aDabzbz99tt0796dcePGMX78+GJfenX8+HF69OhB8+bNad26NT//bPv5GYWNtLwDom6Agmz4bRxY9NHi1rVxCIuf6sbNzcLIM1v4v4V7GT51A6dTZc5dcQXeQdrgwBaDtJaseY/BHxPBfOV5W2sFeDLtgVjeuzOaIB93DiRlMOyr9Tw9O56z6TIQVQi7CGsBMfdp9/98wW7nmy4NQ3ihbzMAXv19D6sPnrPLcYR9uFX0hTt37iQmJgaAXbt2FXtOUZTKRWVHbm5ufPjhh7Rp04YzZ87Qrl07+vbti4+Pj7NDE/+mKND/I/iiq3ZZe9MU6PiIs6MCIMjHnSn3t2PmxmP83+97WXs4md4frOL1wa3o3ayms8NzmFGjRlXodQMHDmTAgAE2jkbHTJ5wxzRtOrpV78CGz+HMTrhzujbjxr8oisId7epyU7NQ3vlzPzM3HmPutpP8vSeJZ3o3ZXineo7/GYTTSJ45SM8XYdcv2vLa22dCzHC7HObBrvXZdyadOVtOMOaHrcwf25X6IVKDuIIKF83Lly+3ZRwOU7t2bWrXrg1ArVq1CAkJ4cKFC1I061VwQ+j1Cix6BpZMhoY3aYOsdEBRFO7tWI/ODYJ5enY820+kMnbmNm6Prk1nd2dH5xj16lWseAsMDLRtIK7AYIAb/6stgDL3UUhcA1/eAIO+hAbdr/iSQG93XhvUiqHtI/jvvF3sPJnK5N92M3PDMZ67tbGDfwDhLJJnDuIXBt2fgyX/08431/UDrxo2P4yiKLw2qCVHzmWw9VgKD83YxNwx1+PvabL5sYRtVbhodpZVq1bxzjvvsGXLFk6fPs3cuXMZOHBgsW3i4uJ45513OHPmDNHR0XzyySd06NChxL62bNmC2WwmIiLCQdGLCmn/IOz7HY6sgHmPwsjF2ohnnWhQ05c5j3Xhk6UH+XT5IeZvP80qdyMRrS9wfeOqvWLm5MmTnR2C62nWH0KawOzh2lzO396urUzW80UwXvmkGR0RyLwx1zNz4zHe/XM/+5PSGTVjK80CDTRul0HzurY/sQv9kDxzoE6Pwbbv4fx+baXAvu/Y5TAebka+uK8dt3/6D4fPZfLkrG18PSIWo0G/V+pFJfo0v/HGG0ybVnIFnWnTpvHWW29VKqhryczMJDo6mri4uCs+P3v2bMaPH8/kyZPZunUr0dHR9O7dm7Nni08ofuHCBe6//36++uoru8UqbMRggAGfgoc/nNgEKxy/xHZpTEYD429pys+PdiaihhcX8xSGT9vMm3/sI69AH32xhY7UbAoPr4B2DwAqrPkApt5yzemujAaF+zrVY+WzPXiwa31MRoW9KQZui1vLC3N3ci7deVMzClFlGE2XCuVNX8PpHXY7VKifJ1Pub4+nSVsx8O3F++x2LGEbFW6u+/LLL5k5s+Qo8BYtWnDXXXfx3HPPVSqwq+nTpw99+vS56vPvv/8+o0ePZuTIkQB88cUXLFy4kGnTpjFxorbyT25uLgMHDmTixIl06dLlqvvKzc0lN/fSiSgtLQ2A/Px88vOvPIjHUYqOX23i8KmF0vc93OaOhtXvUhDeHrXRzc6J5Rpah/vx68PtGTttJRvOGfhi5WFWHTjLu0Na0TjU16GxlOf3YYvfmcFgKDaeQVEUatasyR133MF7772Hh4fHVV+r51wrYvP3l+IOt76LEtUd48KnUU5tRf2iK5abXsYSc782kv8KfEwKE3s3ZkibUJ6buY4dFwzM3HCM+fEnGdWlHiO71MPPSZd59ZCDl9NbPFD5WCqTZ1BNc628IrpgbHY7hr3zsSx8hvy759otnqah3rw5qCVP/bSDL1cdoVFNbwa2CS/1dU7/Hf2LXuOxNUWt4ESBnp6e7N27l/r16xd7/MiRIzRv3pycHPuP8lYUpVj3jLy8PLy9vZkzZ06xLhsjRowgJSWF+fPno6oq99xzD02bNuWll1665v5feuklXn755RKPz5w5E29vbxv+JKKsWh+fQf3zS8k1+rLiulfJcQ92dkhXtT1ZYfYRA5kFCiZFZUA9CzfUUtHjONmsrCzuueceUlNT8ff3r9A+5s2bV+xkbjabOXjwIB999BHDhw/n7bffvuprq3uueeYl0y7xC0Iy9gNwzrcZ2yNHkelReveew2kw96iR45na797bqHJjHQvdaql4GO0atqiAyuZaZfIMJNfKyjPvAjftfQ43Sy5bIx/meHBXux5v4TEDf5004KaojGlupkHF/g2LQrY4p11JhYvmxo0bM3nyZIYPLz669LvvvmPy5MkOWRXw30XzqVOnqFOnDmvXrqVz587W7SZMmMDKlSvZsGEDa9asoVu3brRu3bpYzK1atSqx/yt9Io+IiOD06dMEBzu3WMvPz2fJkiX06tULk8l5gwccHkdBLsYZfTGc2Y4lvC3m4fPB5OWcWK7i8jgu5liY+OsuVh/S5tm9oVEwbw5uSajftVuDbB1Hab+PtLQ0QkJCbP4PBrRxCAMHDuTChQtX3UbPuVbE7u8vixnD5q8xrHgNJT8L1c0TS/eJWDo8CoaSFwUvj8dodOPPPUl8tOwwh89lAhDkY+LRbg24O7YunibHVM96yUG9xgOQnJxM7dq1bZ5rZckzkFwrD8PajzEufwXVO4Q/Gr5Cjz6D7BaPxaIy9sftLNl7lkAvE7NHd6BBzatPUKCX35Fe47FXnlW4e8bo0aN56qmnyM/P58YbbwRg6dKlTJgwgf/85z82C9DWunbtiqWM8y96eHjg4eFBXFwccXFxmM3aspcmk0kXbwrQTywOi8NkgqEzYEpPDKe2YvhjPAyewuXNt3r6ndTxNvHtgx35dl0iry/ay+pDyQyIW8d7Q6Pp0TTUYXGU9vuwxe8rJibmitNNWiwWUlJSaNWqlfU4W7duLbaNK+RaEfvFZILrx0KzfrDgSZSElRiXvoRx73xt6sXa0deMZ0BMBP2i6/Lb9pN8+PdBEpOzeP2P/Uz95yhjezZiWGwk7m6OWcFVb383PcVT2Tgqk2cguVYu14+DnbNRzu+nxckfMZmG2jWej+9uy91T1hN/PIWHvt/Kr49dT81SGlic/jv6F73EY68YKlw0P/vssyQnJ/P444+Tl5cHaF02nnvuOZ5//nmbBVgeISEhGI1GkpKSij2elJRErVq1KrzfMWPGMGbMGNLS0ggICKhsmKKygurD0G/hu0Gw82eoeR10e8bZUV2VoiiM6BLF9Y2CGTcrnr2n03jgm0083K0Bz9zS1GGFjL3dfvvtVzyZZ2VlsXPnTvr164eXl9c19yG5hvb+vn++NoL/zxfh1Db4sju0Hwk3/k9bLOUqjAaFQTF1ua11OL9uPcHHSw9xMiWb/83fzRcrj/DkTY0Z3LYObsaq8Z6rjmyRZyC5ViZu7jDgY5jWm3oXVlFwdDU0vtFuh/NyNzJ1RHsGf76WxOQsHpqxiVkPd8LbXT+zRVV3Ff5LKIrCW2+9xf/+9z/27t2Ll5cXjRs3LnUQgj25u7vTrl07li5dau2yYbFYWLp0KWPHjnVaXMIO6nfTRjj//jQse1Wbz7nJbc6O6poahfox9/EuvL5oL9+uS+SrVUfYkHCBT+6KITLY9fsSXm2MwAcffECDBg148803HRuQK1MUaHsfNLoZ/vov7JoDm6fB7rlw0yRoO+KaLzcZDQyLjWRgTB1mbzrOp8u04nnCLzv4bMUhnrq5Cf2jw2V6KxckeeZgkZ0wtx2Jces3GBeNh8fXWbsE2kOwrwfTR3Zg8Gf/sP1EKk/M2sYXw9vJB12dKFfRfOzYMSIjI4s95uvrS2xs7BW3P3nyJHXq1Kl4dFeQkZHBoUOHrN8nJCQQHx9PUFAQkZGRjB8/nhEjRtC+fXs6dOjAhx9+SGZmpnU2jYr492UsoRPtR8G5/bDhC/j1YZS7fnR2RKXyNBl55faWdGkYwoQ529l+PIV+H6/m9cGt6B9d+ohpPRs8eHCx7y0WC0ePHmXPnj38+GPZ/jaSa//iXxuGTNWmpvtjApzdo31Q3DId5ZbSiyMPNyP3d45iaPsIvl+fyOcrDnM0OYunZscTt/wQ43s1oXeLWhikeHYZtsgzkFwrD0vP/5G3cx5eFxNg5dtws33nza4f4sPXI2K5Z8p6/t57lpcW7ObV21vqerXl6qJcH11iY2N55JFH2LRp01W3SU1NZcqUKbRs2ZJffvml0gH+2+bNm4mJibEu4T1+/HhiYmKYNGkSAMOGDePdd99l0qRJtGnThvj4eBYvXkxYWMUXmRgzZgx79uy55s8tnKT369BsAJjzMP58PwFZR50dUZnc2rIWi568gXb1apCeW8C4Wdt46bfd5Jtdd05nf39/AgICrF8hISEMHjyYzZs3lzjRX43k2lXUvwEeWQ23vgUeAXB6O24z+tA+4VO4mFDqyz1NRh66oQGrJvTk2d5N8fd04+DZDB77YSu3fbKGpXuTqOCYcOFgtsgzkFwrF09/dtS9X7u/9mM4s8vuh2xXrwYf3RWDosD364/xxUr7T64gSleuluY9e/bw2muv0atXLzw9PWnXrh3h4eF4enpy8eJF9uzZw+7du2nbti1vv/02ffv2tXnAPXr0KPWf+9ixY6U7RnVhMGoDAX+4iHJ0NZ0OvwsXekPYdc6OrFR1a3gz++FOvL/kAJ+tOMz0tUfZdTKVz+5tS6i/p7PDK7fp06c7O4SqzegGnR6FlnfA0pdRt31PnZSNqF90gdiHoPuEa/Z3BvDxcGNMz0YM71SPqWsSmLYmgT2n03hwxmbaRATyn1ua0LVRiLRo6ZjkmXOcCWyHpeltGPb/Dr+Ng4f+1s4/dnRry1pMuq05Ly/Yw1uL9xEe6MntbWx79V6UT7mK5uDgYN5//31ee+01Fi5cyJo1a0hMTCQ7O5uQkBDuvfdeevfuTcuWLe0Vr1PIZSydM3nCXT+gftMPz6SdqN8PgpELIaiBsyMrlZvRwIRbr9MKlp+2sznxIv0+WcNn97YlNuraBZAe1K9fv0IF1lNPPcUTTzxR4nHJtTLwrQm3f0pB+9Fc+HEsYek7YMPnED8TbhgPHR/VcuIaArxMjO/VhJFdovhy1RFmrD1K/PEU7pu6kQ71g/hPryZ0bKCP6ceE7fMMJNcqwtz7DQxHV8OprbD2E+j6lN2POfL6+py4mM3UNQk8+/MOwvw96SS56TQVGgjo5eXFkCFDGDJkiK3j0SUZZewCPAMouGs2OV/ejF/6KZjeHx74XZuJwAXc0qIW88f68uj3WziQlMHdX63nldtbck/HyNJf7EQVbfWKioq64uOSa+UQ2pz1jZ6h33XeuC17GZJ2wt+TYeMUbTaZmOHaksDXUMPHnYl9ruPBrvX5fMVhvt+QyMaECwz7aj03NA5hfK8mxETWcNAPJK7G1nkGkmsV4lcbbn0d5o+B5a9B41sgrLndD/ti32acTs1m0c4zPPztZn55rAtRQa53NbIqkHlMRNXhG8o/jSbS+/QnKMkHYUZ/eGAh1Kjn7MjKpEFNX+aNuZ7nftnJgu2neGHuThLOZzCxTzPdznLQvXt3Z4dQ7akNemjTYO2YDcv+D9JOwO9PwZoPoMdEaDVU69pxDTX9PJjUvzmju9Xn02WHmL3pOKsPnmf1wfP0ah7GM7c0pWktP4f8PKIkyTMdaXMv7F0ABxbDvEfhoaWlfjitLINB4f2hbTibtoHNiRd54JtNzB595QkYhH3JHCZlEBcXR/Pmza86S4jQj1xTIAX3zoXgRpB6HL7pA+cOODusMvN2d+Pju9rw9M1NAJiyOoFHv99CVl6BkyNzDMm1CjIYoc09MG4r3Pom+IRCSiLMeww+6wg750AZFnWqHeDFa4NasfyZHtzZri4GBZbsSeLWj1bx9Ox4jiVnOeCHEY4guVZBiqItNuQZCKe3w+r3HXJYT5ORKfe3p0GIDydTshn93TZypGeNw0nRXAYyytjF+NWCEQsgpAmknYRvboVT8c6OqswUReHJmxvz0V1tcHczsGRPEkO/XMe59NzSX+ziJNcqyeQJnR6DJ+Oh1yvgFQTJh+CXB+HzLrDnNyjDLBkRQd68c2c0fz3dnb6taqGqMHfbSW56fwX/m7eLs2k59v9ZhF1JrlWCXy3o9552f9XbDju/1PBxZ/rIDoT4urP3TDrf7De49IxLrqjCRfNff/3F2bNnbRmLELbjHw4j/4DabSArWeuqkbjW2VGVy+1t6jBrdEeCfdzZdTKNYV+u42RKtrPDEq7A3QeufxKe2gE9/wueAXBuL/x0H3zZDfYvLlPx3CjUl8/ubceCsV25oXEI+WaV79Yn0u2d5bz5xz5Ss/Id8MMIoUMt74Dmt4OlQLuiU+CYRo3IYG+mjojFy2RgX6qByQv2ynSRDlThonngwIHUrl2b8PBw+vbty4svvsiKFStsGJoQleQTorU417sectO0Zbf3/ObsqMqlXb0gfnmsC3UCvThyPpOhX6wj4Xyms8MSrsLDD7o/C0/ugG4TwN0XzuyAWcPg65vh0NIyFc+t6gbw3YMdmTW6EzGRgeTkW/hi5WFueHsZU9ckkFcgrV2imlEU6Pc+eIdoiw4tf81hh46OCOSDoa1RUPl5y0k+WXao9BcJm6hw0Zyens727dt55513aNmyJRs2bODWW2+lZ8+eZGZWrZO69P1yYZ7+MPwXaHIrFOTAT/drUwW50CfzqBAffn60s7Uv251frGPfmTRnh2UXkmt24hUIN76oFc/XPwluXnByM3w/GL7pC0fXlGk3nRsG8+tjXfj6/vY0DfMjLaeAV3/fQ+8PV/H3HlkgxZVIrtmAT4jWvxngn4/hyAqHHfqm60IZUl/7sPr+kgPM2XLCYceuzipcNCcnJ9OyZUvuvfde3n77bf7++2+OHTtGfn4+r776qi1jdDrp++XiTF4w7AeIHQ2o8Nd/YeF4MLvO4LrwQC9+erQzzWr7cz4jl3unbODQ2Qxnh2Vzkmt25hOs9XV+agd0ehyMHnBsLUzvBzMGwPGNpe5CURRubh7Goidv4I3BrQjxdSfhfCYPfbuZB6Zv4WTVajOpsiTXbKTZbdoy96gw91HITHbYobvWUnn4higAJv6yg9UHzzns2NVVhYvmWrVqER4eTp8+fZg4cSKzZs0iOTmZDz/8kGnTptkyRiEqz+gGfd+B3m8ACmyepl2iznGdFtsQXw9+HN2JlnX8Sc7MY/jXGzh+QWYzEBXgGwq3vqENGGz/IBhMkLASpvaC74fAqW2l7sJoULi7QyTLn+nBo90b4m40sPbIBd7ZYWTygj2kZkt/Z1FN9H5dG3ieflpbLdCBV1z+c3NjBkSHU2BReez7rew97TrnNFdU4aL54MGDfPrpp3Ts2JH9+/fz4osv0qJFC7p27UpycjLDhw/nk08+YfXq1baMV4iKUxTo/DgM+167PH3ob5jWGy4ccXZkZRbgbeLbUR1pHOrLmbQc7v16A+czqv6sGsJO/MPhtvdh3BaIuQ8UIxxaAl/1gB/vhaQ9pe7Cz9PExD7X8ff47tzaIgwVhZkbT3DTeyuZH39SumyIqs/dB+74GozusH+h1ijjIAaDwjt3tqZj/SAycgsY+c0mTqfKgHF7qXDR3LBhQwYPHsxLL73E3LlzOXLkCCkpKcyePRtVVVFVlRkzZtC7d29bxitE5TW7DUYuAt8wbQDHVz3g4N/OjqrMgnzc+f6hjkQGeXPsQhajv91MTr5M2CkqoUY9uP1TGLsJWt8FKLDvd22aul9GQ/LhUncRGezNJ3dFM7a5mQYh3pzPyOXJH+O5f9pGmd9ZVH21o+Hml7T7f74AZ/c57NAebka+uq+9tTHlgWmbSMuRKz32UOGiOSgoiBtvvJGnn36a6dOns2XLFk6fPs2iRYto2LAhP/zwA5s3byY9Pd2W8TqFDJioguq0hYdXQN1YyEmFH4bAqnfKtACEHoT5e/LNyFgCvExsO5bCf37ejsXi+i16kmtOFtwQBn8Jj6+HZgMAFXb+BJ/GwoInIfVkqbtoHKDy25gujO/VBHc3A6sPnqf3h6v4bn2itDrriOSaHXR8DBrepA06/+VByHfcfOYB3ia+GRlLTT8P9iel89j3W2RWGzuocNE8bdo0unfvTmJiIq+88gqxsbE0a9aMH374gY8++si6ndFotEmgziQDJqoo/3Btme12IwFVW4L4p/tcpp9zw5q+fDG8HSajwsIdp/lqtet0M7kayTWdCL0Ohn0HD6+ERr1ANcOW6fBxDCx+HjKuPeDIw83AEzc15s+nutG5QTDZ+Wb+N28X90/bKJeOdUJyzQ4MBhj4uTYNXdIurcXZgerW8OabB2Lxdjfyz6FkJv6yQz6o2lil5mmePHkyv/76K0eOHCE1NZV9+/aRlJRE3759bRmjEPbj5gH9P4T+H2v90fb9DlNuhNM7nB1ZmXRuGMxLA1oA8M6f+9lwxHEjt0U1EN4Ghs+BkYu1+c7NubD+M/goGpa/AXnXniqjfogPPzzUkcn9m+NR2Op8ywerWLzrtGPiF8LR/MK0qzUosHkq7PjZoYdvWSeAz+5ti9Gg8Ou2k7y/5IBDj1/V2WwZbT8/P5o0aYKPj4+tdimE47Qboa0g6BcOyQfh65tgXZxLdNe4p0Mkg2LqYLaojJu1jYuZec4OSVQ19TprV2WG/wrhMZCfCSvfhI/bwtbvwHL1PvUGg8LI6+uz6MkbiI4IJD2ngEe/38r//b5HlgAWVVOjm6HbM9r9BU/COccWrj2ahvL6oJYAfLLsELM2HnPo8asymxXNQri8uu3h0dXQpA+Y87RLaz/cAelnnB3ZNSmKwmuDWtKwpg9n03N5acFuZ4ckqiJFgUY3wejlcOd0qBEFGWfgt7Ha0tyHl13z5Q1r+jLn0c480q0BAF+vSeDur9ZzJtVx/T6FcJgez0PUDdoHzJ/uL/WqjK0Ni43kiRsbAfDfebtYsf+sQ49fVUnRLMTlfELg7lnQ7z1w89QKgc+7wP7Fzo7smrzd3Xh/aBuMBoX58afk8rewH0WBFoNgzEZtflrPQK3/5neDMM4ahm/Oqau+1GQ08HzfZnx5Xzv8PNzYnHiRgXH/yNyyouoxGOGOqdosTef2wsL/OHwl2qd7NWFwW+0q5JgftrLnlORZZUnRLMS/KQrEPqQNggprBVnJ2kIoC/8D+fodxBQdEcij3bVWvMm/7SYz13VWPBQuyM0DOo+BJ7ZpqwsaTBiOLKXnvhcxLH/1mi1rvVvUYsG4rjQqnCLrzi/WseqArGYmqhi/MK1wVgywfRZs+86hh1cUhTcHt6Zzg2Ay88yMmi5zOFeWFM1lIFPzVFOh18HopdBpjPb9pq+1OZ3P7HJqWNcy7sbGRAZ5k5SWy+crE5wdTrlJrrkg7yBtdcExG7A0ugWDasa49iP4tAPs+e2qrWtRIT788mgX66IMo6ZvYn586VPaCduQXHOQ+jfAjf/V7i96Fk5vd+jh3d0MfHFfO+sH1FHTN5MhDSoVJkVzGcjUPNWYmwfc+joM/6XwMts+mNIT1n2my0GCniYj/7utOQDT1h4l2cW6i0quubDghpiHzWR9g6dRAyIh7YQ2heOsuyHtyt2FArxNfPtgB25voy0D/PTseBZsv3r3DmE7kmsOdP3T0PgWbf7mH4dDpmNnOQrwMvHNA7GE+Lqz93QaY37YSoEMwq0QKZqFKItGN8Nja6HJrYWDBJ/XFkRJT3J2ZCXc3CyU6xsFk29WWXJSUlw4VlJADAWPrIFuz4LBBAf+gM86QvzMK7Y6e7gZ+WBoG4a2r4tFhadmx/P7DimcRRViMMDgr6BGfUg9BnMeALNjW3sjgryZOiIWT5OBlQfO8b/5u2UO5wqQM6oQZeUTAnf/CH3fLRwkuFQbJHjgT2dHVoyiKDx9cxMANpxTOH5RljAWDmby1i5JP7ISarfRVt2c9xjMHHrFD5oGg9b3cki7upgLW5zXy5zjoirxqgF3zQSTDySsgiWTHB5CdEQgH98Vg6LArI3H+HKV6y+I5WhSNAtRHooCHUZrS3CHtYSs81ohsPAZXQ0SbB8VxPUNg7GoCl+4YN9mUUWEtYCHlsJNk7XFgw7+BV90veL0dAaDwlt3tKZfq9rkm1Ue/X4LCecdO02XEHYV1hwGfaHdXx8H22c7PIRbWtTif/20Lnxv/rGPhTtkpqXykKJZiIoIbaYVA50e177fNAW+6glJ+pkjeWxPbSaNedtPk5IlC54IJzG6wQ3j4ZFVENocMs/Cd4Ph75dLXKI2GhTeGxpNm4hAUrLyeXD6Jhm0JKqW5gPghqKFT56AU/EOD2FU1/o80CUKgP/8HM+uk6kOj8FVSdEsREWZPLVZA+79BXxCtbk4v+qJYdNXDp+P80raRQZSx1slr8DCnC0nnB2OqO5Cm8HoZdB+FKDCmvfhu4GQdaHYZp4mI1Pub0/tAE+OnM/kpd/080FUCJvo+cKlgYGzh0PmeYeH8L/bmtO9SU1y8i08/O1mzqXnOjwGVyRFsxCV1bhwkGDj3mDOxfjXC3Q68h5kOHcFJkVR6FpLGyH9w4ZjMuhDOJ/JC277QFtR0N0Xjq7Wlqz/1zLDNf08+OiuGAwKzNlyQmbUEFWLwQiDp0BQQ0g9Dj8/AOZ8h4ZgNCh8fHcMDWr6cCo1h0e/30JugdmhMbgiKZqFsAXfmnDPbOj7LqqbJ2FpO3D7ujsc+MupYbULUfHxMJJwPpPNiRedGosQVi0GwYNLIDASLhyBr2+GIyuLbdKhfhBje2rLAE+av0u6GImqxStQGxhY9OFx0bMOv0IZ4GXi6/vb4+fpxpbEi/x37i5pXCmFFM1lIJPAizIpHCRYMOpvUj0jUDLPwcw7YdEEyHfOhMkeRuh1XSiASwz4kFyrRsKaw0PLIKIT5KbCD3eWmIlm3E2NaRrmx8WsfN75c7+TAq2aJNd0IPQ6uONrQIEt38D6zx0eQoOavnx6T1sMCvy85QTT/jnq8BhciRTNZSCTwItyqXkdq5pOxhz7iPb9xi+1BVGS9jglnFtbhgHwx67TWCz6bkWQXKtmfGvCiN+gaT8w58KP98Ke+danTUYDr9zeAoCZG4/JgCUbklzTiaZ94JZXtft/vgD7Fzs8hO5NavJC32YAvLZwD/8ccnwfa1chRbMQdmAxuGO55TW4dw741ISze7QluDd86fBLcF0bheDr4UZSWi7bjqc49NhClMrNA4bOgJZ3gCUffh5ZrHDo2CCY/tHhqCp8tPSgEwMVwk46j4W2IwAVfnkQzuxyeAgPdq3PkHbaAkNPzNrGmVQXW07WQaRoFsKeGvcqHCR4i9aS9scEbVlhBy6j6uFmoFuTEADWHJQWBKFDRpM2MKr1MFDN2sCoE5utTz95U2MMCizZkyStzaLqURTo9x7U7wZ5GTDrLoevNqsoCv83sCXNa/uTnJnHmJlbyZeltkuQolkIe/MNhXt+gj5vaws8HPhDW+Dh6BqHhdC1UU0Auewm9MtghNvjtCXrC7K1Ps7JhwFoFOpL/+hwAKasllXMRBVkNMHQbyG4kTajxo/3OHzBLE+Tkc+Ht7UODHxj0T6HHt8VSNEshCMoCnR8RFsQJbgxpJ+CGf1h+eslFniwh66NtJbmrccukimLRQi9MprgzhkQHgPZF+DnEdbC4aGu2mI9f+w8w/kMmVNWVEFeNbQGFs9AOLkZ5j3u8O589YJ9eO/OaACm/ZPAop36H0DuSFI0C+FItVvDIyuhzXBQLbDyLa14TrXv4iORwd6EB3hSYFHZcUIubwsd8/DVpuLyDoYzO+GP5wBoVTeA6LoB5JllsR5RhQU3hGHfg8ENdv8KK950eAi3tKjFI920D6kT5uwgMVmWsy8iRbMQjubuAwPj4I6p4O4Hx9bCFzfA4eV2PWzruoEA7DiRYtfjCFFp/uGXpuLaOgP2LgBgWGwk4BrTJwpRYfVv0BYBAlj5JsTPcngIz/ZuSmxUDTJyC3hqdjwF0r8ZkKJZCOdpNQQeXQW1o7VL0d8PhlXvgsU+/5yiIwIBpKVZuIaGN0LXp7T7C5+B7BR6twjDoMDOk6kcv5Dl1PCEsKu298P1T2n3fxsLh5c59PBuRgMfDGuDn6cb246l8PGyQw49vl5Vy6J50KBB1KhRgyFDhjg7FFHdBTWAUX9BzH1ad41lr8LseyE7xeaHiq4bAMB2aWkWrqL7c9pSwxlnYOkrBPt60KlBMID0tRRV302ToeUQsBTA7Pu17koOVLeGN68NagXAp8sOsvnoBYceX4+qZdH85JNP8u233zo7DCE0Jk+4/VPo/zEYPWD/Im0xlPO2nZO2RR2taD5xMZsMGQwoXIHJC/p/pN3f8g2cO0DvFrUAWHXwnBMDE8IBDAYY+BlE3QB56dqMMinHHRrCgOhwBretg0WFp2bHV/tzR7Usmnv06IGfn5+zwxCiuHYj4ME/ISASLhyBr2+CIytstvsALxMhvu4AJJyTgR3CRdS/QVsxsPBKTJeGWkvzlsSL5BaYnRycEHbm5qENDKzZDNJPa4WzHa5EXsvLA1pQt4YXJy5m887i6j0NncsVzatWraJ///6Eh4ejKArz5s0rsU1cXBxRUVF4enrSsWNHNm7c6PhAhaiI8BgYvQzqdoCcVPj+Dtgy3Wa7bxDiC8CR8xk226cQdnfT/0AxwN7faKQeJdjHnZx8i/TPF9WDVyDc+zP41YZze2H2cChw3LSLfp4m3rqjNQAz1iWyMaH6dtNwuaI5MzOT6Oho4uLirvj87NmzGT9+PJMnT2br1q1ER0fTu3dvzp496+BIhagg35owYgG0ulPry7bgSfjzRZsMEKwf4gPAEWlpFq4ktBk0vx0AZf0XdGwQBMAm6WMpqovACK1wdveDo6sxLhirXX1xkOsbhXBXbAQAz/2yg5z86nmVx83ZAZRXnz596NOnz1Wff//99xk9ejQjR44E4IsvvmDhwoVMmzaNiRMnlutYubm55OZe+jSXlpYGQH5+Pvn5+RWI3naKji9xXKKXWGwThxH6f4ahRkOMq96EdZ9iyUzG3O8Dbf7OCsZRL9gTgENn04s97uzfmZ5zrYhe3l9F9BYP2DcmJfYR3HbPRd35E23bDWcRsPtk6jWPpeffkbNIrpWfbuIJvg5lyHSMPw7DsGcuzUPzyM/v7bDDT7ilEcv3nyXhfCbv/bmPCb2bWJ/Tze+okL3icLmi+Vry8vLYsmULzz//vPUxg8HAzTffzLp168q9vzfeeIOXX365xOPLly/H29u7UrHaypIlS5wdAqCfOEA/sdgmjuZERI4m5tjXGHbM4kziQbZEPYbFYKpQHOeTFcDI7oTTLFp0aYGIrCznTt/lCrlWRC/vryJ6iwfsF1M37/rUyEog8tAPQB+2HCr+PnZ0PBUhuVZ2evq7gX7iiag7irbHvqLx2YXsmDmBhJq9HHbsAbUVpqQZmfpPAiHph6j1r7eMXn5H9sozRVUdvEajDSmKwty5cxk4cCAAp06dok6dOqxdu5bOnTtbt5swYQIrV65kw4YNANx8881s376dzMxMgoKC+Pnnn4ttX+RKn8gjIiI4ffo0wcHB9v3hSpGfn8+SJUvo1asXJlPZC6iqGoeeYrFHHMq+hRjnjUYx52Gp3wPzkBnaIinljGPb8RSGfrWR8ABPVj7TzbptWloaISEhpKam4u/vb5OYy0PPuVZEL+8vvcYD9o/JsOkrjH+9QG6ttjQ9+gxGg8L2/96Ih8nolHgqIjk5mdq1a0uuXYPe/m56iwdAXfUO7qvfQkXBPGQGatO+Djv2I99vY9n+c3RpGMT0Ee1QFEV3vyN75VmVamkuq7///rtM23l4eODh4UFcXBxxcXGYzVofHpPJpIs3BegnFr3EAfqJxaZxtBoI3gHw470YElZgmHMf3POTNiVXOeKoE6QNBDyXkYvR6IbBoFi3cSZXyLUieotJb/GAHWNqNQSW/BePM1tp4XWB3dlBHL2YS8vC6RQdHk8FODsOybWK01M8+Tc8w9E9G4lKXo7bvIe1cTARHRxy7MkDWrDmg1WsPXyBZQeSubVlbetzevkd2SsGlxsIeC0hISEYjUaSkpKKPZ6UlEStWrUqvN8xY8awZ88eNm3aVNkQhai4hj3h/nng7gsJq+Cn+6Egr1y7CPXzACDfrHIxq3yvdQTJNXFNfmHanLXAYO94AE5clJUBK0JyzcUpCjsi7sfSqBcU5MDMYXDeMav21Qv24ZFuDQB49fe91WpQYJUqmt3d3WnXrh1Lly61PmaxWFi6dOkVu18I4XIiOmgtzG5ecPAv+PWhcs2qYTIarHM1J6U5bsoiIWym8S0AdFK11dGOX8h2ZjRCOI2qGDEP+lqbqjT7AvxwB2Q4ZtGfx3s0IjzAk5Mp2Xy77qhDjqkHLlc0Z2RkEB8fT3x8PAAJCQnEx8dz7NgxAMaPH8+UKVOYMWMGe/fu5bHHHiMzM9M6m0ZFxMXF0bx5c2JjY23xIwhROVHXw90zwegOe+bD35PK9fJQP20GjaS0HHtEVymSa6JUDXoA0DhnO+7kc1xamitEcq2KcPfRGlIC68HFozBzKOTZf0pRL3cjT/fSZs/4bMVh0nP0MWuGvblc0bx582ZiYmKIiYkBtCI5JiaGSZO0wmHYsGG8++67TJo0iTZt2hAfH8/ixYsJCwur8DHlMpbQnYY3wu2Fc5Wv/QQ2TinzS4N8tJZm6Z4hXFJYC/AJxd2SQ7RymBMXpaW5IiTXqhDfUBj+K3gFwamtMGcUmO2/3PXgtnVpFOpLSlY+U/9JtPvx9MDliuYePXqgqmqJr+nTp1u3GTt2LImJieTm5rJhwwY6duzovICFsJfWQ6Hnf7X7fzwHiWvL9LIAL22ARGp29WgZEFWMolgHPLUyJHA+Q7oZCUFII7j7R3DzhAOLYdF/wM6ToxkNCs/corU2f7M2kfRqcEpxuaLZGeQyltCtbs9AyyGgmuHnkZCeVOpL/AuL5rRs+7dElJfkmiiT2tEAtDAkkJyhvysmrkByrQqK7Ah3fA0osGU6rH7P7ofs3aIW0XUDyMozs+J01S8pq/5PaANyGUvolqLAgI+hZjPIOAO/PFjqwEB/L22mST22NEuuiTIpLJpbKke5kClFc0VIrlVRzfpDn7e1+8tehfhZdj2coiiM6dkIgH/OKKTn6K8xxpakaBbC1bn7wLDvwOQDR1fDhi+uuXlR94y0ajJwQ1RBhUVzI+UklvxssvKq9olaiHLp+DBc/6R2/7excHiZXQ93c7MwGoT4kG1WmL259BU6XZkUzWUgl7GE7oU0hlte1e4vfRnOHbjqpnru0yy5JsrErxaqZwBGRSVCOStdNCpAcq2Ku+klreuepQBm3w+nd9jtUAaDwkNdowCYvjaRvIKyT4PqaqRoLgO5jCVcQvtR0PAmbaL735++6iAQP8+iPs36K5ol10RZKTWiAIhUzupyJhi9k1yr4gwGGPiZthhQXjr8cCekHLfb4QZE1ybApJKUnssfu07b7TjOJkWzEFWFokD/D7WFTxLXoOz59YqbeZuMAORU4dYAUQ1cVjRn5lafFcmEKDM3Dxj2/aUxLz8MgeyLdjmUh5uBzmHaOWXmhmN2OYYeSNEsRFUSGAnd/gOA8e9JGC0lp+PyLCqa86TQEC7ssqI5O1/6NAtxRV6BMHwO+IXDuX3w471QYJ9pGjuHqhgU2JBwgUNnM+xyDGeTolmIqqbLExAYiZKRRNT5kgNAvNy1tM8pkKJZuLDCojlCWpqFuLaAulrh7OEPif/A/DF2mcM50AN6Nq0JwKyNVbO1WYrmMpABE8KluHlAtwkANE76vcSSqh5uhS3N+forNCTXRJn51gKgppIis2dUgORaNRPWQptlyeAGO3+GlW/b5TDD2tcF4NetJ8g3V70ugFI0l4EMmBAuJ/ou1Br18ShIx7BlarGnirpnZOuwe4bkmigzH61FK5h0snT4XtY7ybVqqEEP6Fe44MmK12HnHJsf4oZGwYT4unMxK59/Dp23+f6dTYpmIaoiownz9eMBMGz6GsyXWuK83GUgoKgCfIIBCFLSpGgWoqzaPQCdx2r35z0Ox237ocnNaKBvq9oALNhe9WbRkKJZiCpKbTGYXDc/lPRTsH+R9XFPNy3t8wosWCy279cmhEN4hwDgo+SSm1U1Bx0JYRe9XoGmfcGcCz/eDRcTbbr7/tHhAPy1+4wuuwFWhhTNZSB9v4RLcvMgMbiHdn/TFOvDRS3NoL/BgJJrosw8/ChQ3AFQsqreZWB7k1yrxgxGGDwFarWCzHMw6y7ISbPZ7ttF1iA8wJP03AJWH6xauSlFcxlI3y/hqo6G9ERVDJCwCi4eBS4NBATIyddXFw3JNVFmikK2KRAA9zz7zD1blUmuVXMevnD3bG1A7dk9MGckWGzTiGIwKPRqHgbAsn1nbbJPvZCiWYgqLNs9BDWyi/bN7nkAGA0KiqI9VGDRV9EsRHnkG70BcCvIcnIkQriggDpwz4/agliH/oZlr9ps1z2vCwVgxf6zqHaY3s5ZpGgWoopTmw/U7uyea33MzaBVzWbp0yxcWIGbVjQbC7KdHIkQLio8BgbGaffXfFDsPFEZnRoE42kycDo1h72n022yTz2QolmIKs7S9DZQjHA63tpFwyhFs6gCCopams2ZpWwphLiqlndoi2KBNqNG0u5K79LTZOT6htpg3RUHqk4XDSmahajqfEIgspN2/7C2QqCbQUt9KZqFKzO7eQHgZpaWZiEq5abJ2jzO+VnaUtvZlR8ncH0jrWjecORCpfelF1I0C1EdNOyp3RYWzUUtzQVSNAsXZnbzAaRPsxCVZnSDId9AYCRcTIBfHoJKjnnpUD8IgC2JF6tMA40UzWUgU/MIl9fgRu32yCowF+i2T7PkmigPS2FLs0lamstNck2U4B0Ew364NDDwnw8rtbtmtf3x83AjI7eAvadtN6WdM0nRXAYyNY9weeFtwMMfclPh3L5LLc1mfRXNkmuiPMwmraXZZJGiubwk18QV1W4Nfd/W7i/7P0hcV+FdGQ0K7aJqALAxoWp00ZCiWYjqwGDUCmeAk1t029IsRHkUtTS7S9EshO3E3AethoJqhl8ehKyKF7yxUVoXjU1HpWgWQriSOu2025NbMBqL+jTLPM3Chbl5AGCw5Ds5ECGqEEWB296HoIaQdhLmPQYVnGu5TUQgALtOpdowQOeRolmI6uKyollmzxBVgWJw025VfS0HL4TL8/CDO6eD0QMOLIb1n1doNy3C/QE4fiGb1GzX/3ArRbMQ1UXtNtrtuX24K1qRIbNnCFemGE2AtDQLYRe1W0Pv17T7f78EZ/eVexeB3u7UCdS6Ue055fqDAaVoFqK68K8DJm+wFFAXbbJ5aWkWrkwxai3NBrXAyZEIUUXFPgSNeoE5F+Y+Aubyf0Atam3eXQW6aEjRLER1YTBAcEMAIi0nAWlpFi7OUNjSrErffCHsQlFgwCfgGaitKrv6/XLvokV4AAB7qsC0c1I0C1GdhDQBIFLVimazDAQULkwt7NPshrQ0C2E3/rWh77va/VVvw6n4cr28cZgvAAnnXX+5eymay0AmgRdVRnBjACIsJwAw66xmllwT5aEqRgCMyEDA8pJcE+XSagg0vx0sBbDgCTCX/YNq/RBtPnUpmqsJmQReVBlBDQAINSc5OZArk1wT5aEWds8wyuwZ5Sa5JspFUaDPO+ARAKe3w6YpZX5pVLBWNKdk5XMxM89eETqEFM1CVCf+tQEIsiQ7ORAhKq+opVm6ZwjhAH5h0Osl7f6y/4PUk2V6mZe7kfAATwCOnM+wU3COIUWzENWJXzggRbOoGtTC2TOkpVkIB2n7ANTtAHkZsPi5Mr+sfk2ttfnIOdfuoiFFsxDViV8tALzVbHzIRq3gKk9C6IGqFBbN0qdZCMcwGKD/h6AYYe8COLKyTC8r6qKRmJxlx+DsT4pmIaoTD1/w0ObMrKVccHIwQlSOzJ4hhBOEtYDYB7X7f70IltI/tIYXLnByOjXHnpHZnRTNQlQ3flq/5jDlopMDEaJyiopmaWkWwsG6T9QGBZ7ZCdt/LHXz2oV9ms+kZds7MruSolmI6sYnBIAauPaADCFURTuFKdLNSAjH8gmGbv/R7i97FfKu3Ve5VmHRLC3NQgjX4hkIQKCSgZQawpUpKM4OQYjqq8MjEBgJ6acxbJl6zU3DAwq7Z6TkuPRYGimahahuvAIBCMC1RzELIYRwIpMn9HgeAMP6zzCac6+6aVFLc3a+mbRs1x2DUO2K5t9//52mTZvSuHFjvv76a2eHI4TjFbY0+ytSNAshhKiEVkOhRhRK1nmizi+76maeJiM1vLXFiE6lum6/5mpVNBcUFDB+/HiWLVvGtm3beOedd0hOlvlqRTUjLc2iilGko5EQzmF0gxueAaDR2UWQf/Up5UL9tNbm5AzXXRWwWhXNGzdupEWLFtSpUwdfX1/69OnDX3/95eywhHCswpbmACUTF+5aJoS2tK8Qwrmi70INiMCzIBXD9llX3ayGj9bSfCFLimaHWLVqFf379yc8PBxFUZg3b16JbeLi4oiKisLT05OOHTuyceNG63OnTp2iTp061u/r1KnDyZNlWwZSiCrDMwAAf2lpFkIIUVlGE5aOjwNg2DwFLJYrblbD2x2Ai5muWzS7OTuA8sjMzCQ6OppRo0YxePDgEs/Pnj2b8ePH88UXX9CxY0c+/PBDevfuzf79+wkNDS338XJzc8nNvdSxPS0tDYD8/Hzy8/Mr/oPYQNHxJY5L9BKL3uNQDB64AZ5KPqnmAl28n/Wca0X08nctord4wPExWcza/MzqVY6p59+Rs0iulZ/e4gH9xZTf7A5Mf7+MKfkQBQeWoDa8scQ2AV5ayZmcnmP3uO21f5cqmvv06UOfPn2u+vz777/P6NGjGTlyJABffPEFCxcuZNq0aUycOJHw8PBiLcsnT56kQ4cOV93fG2+8wcsvv1zi8eXLl+Pt7V2Jn8R2lixZ4uwQAP3EAfqJRa9x1EzbSRfAi1y2bNlKwVGVrCznLm3qCrlWRC9/1yJ6iwccF1PO6X20BiwWC4sWLXJ6PGUhuVZ2evq7gf7iAX3F1DK4Gw3P/cX5Ra+xoWHJ+ZgvnDIABrbtPciinP12jcVeeaaoLjphnqIozJ07l4EDBwKQl5eHt7c3c+bMsT4GMGLECFJSUpg/fz4FBQU0a9aMFStWEBAQQLt27Vi7di3BwcFXPMaVPpFHRERw+vTpq77GUfLz81myZAm9evXCZDJV+zj0FIve41COrcXtuwEcttRm7+C/6d0ijLS0NEJCQkhNTcXf39/hseo514ro5e+q13jA8THt+WcB0StGcsRQj4jntzg9nrJITk6mdu3akmvXoLe/m97iAf3FlJ+fz9rfv+emPRNQUMl/fBPUqF9sm2/WJvL6H/vp16oWHw5tbdd47JVnLtXSfC3nz5/HbDYTFhZW7PGwsDD27dsHgJubG++99x49e/bEYrEwYcKEa/6T8PDwwMPDo8TjJpNJF29S0E8seokD9BOLbuPw9NNulDzc3Iy6iNMVcq2I3mLSWzzguJjcjMZix3R2PGXh7Dgk1ypOb/GAvmLK9AhDrd8dJWEFpt1zoOcLxZ4PKZw9IzW7wO4x22v/LjUQ0BYGDBjAgQMHOHToEA8//HCZXhMXF0fz5s2JjY21c3RCOIBJuwTrxdUnoncWyTUhHENyTdiDpfUw7c72WSUGBBYNBEzJdt2BgFWmaA4JCcFoNJKUlFTs8aSkJGrVqlWpfY8ZM4Y9e/awadOmSu1HCF0waZ/2PdHHAJLLSa4J4RiSa8Ie1Kb9wN0PUo7BsXXFnvPx0Do3ZOaanRGaTVSZotnd3Z127dqxdOlS62MWi4WlS5fSuXPnSu1bPpGLKqWwpdlbyUVvEzVLromKkNmay09yTdiFyRta3K7d/9eczT4eWneqjFxZRtshMjIyiI+PJz4+HoCEhATi4+M5duwYAOPHj2fKlCnMmDGDvXv38thjj5GZmWmdTaOi5BO5qFJMXta7illfXTQk10S5SLVcYZJrwm5a3and7l8E5ksFsq+1pdl1i2aXGgi4efNmevbsaf1+/PjxgDZDxvTp0xk2bBjnzp1j0qRJnDlzhjZt2rB48eISgwOFqNbcPK13jeaS0wIJIYQQFVbvem3l2axkOL4eoroCl7pnZOWZsVhUDAbX+9TrUkVzjx49KG2GvLFjxzJ27FibHjcuLo64uDjMZtfthyOEleFS2isWfb2nJdeEcAzJNWE3RhM07aN1z9i30Fo0F7U0A2TmFeDnqY9ZP8rDpbpnOItcxhJViqJgLkx9RdXXZTLJNSEcQ3JN2NV1t2m3e3+3jp3xcDNgLGxddtXBgFI0C1ENmSmc31Z1zX9cQgCo1k7N+hrQKkS11/BGMLpD6jFIPgxoi9L5uLv2YEApmoWohiyK9o9Lb90zhCgPRUYCCqFP7t4Q0VG7f2S59WFXHwwoRXMZyNQ8oqopamk26Kx7huSaEI4huSbsrkEP7fbICutD3kVFc56+zj1lJUVzGUjfL1HVFBXNis66Z0iuCeEYkmvC7hoUznaWsBoKr2p6uGllZ16B5Wqv0jUpmoWohoq6ZyDdM4QrU7TuGYr0aRZCf8LbaKsD5qbC2b0AuEvRLIRwNdaWZotrXiITQgihcwYj1InR7p/Qrmi4G7WyM1eK5qpL+n6JqsasFE05p6+WZsk1UR6KjAOsMMk14RB1C99fJzYD4GHSGmykpbkKk75foqqx6HQgoOSaEI4huSYcoqhoPqkVzUUtzXlmKZqFEC7i0jzN+iqahSgfaWoWQtfqtNduz+2D7BTrQMDcfH1d5SwrKZqFqIZknmZRlchAQCF0yrcmBERo98/uuTR7hrQ0CyFchV6nnBNCCFHFhDbTbs/uldkzqgMZMCGqGrOizz7NkmtCOIbkmnCYmtdpt+f2SdFcHciACVHVqEWpr7PuGZJrony0Ps3SOaP8JNeEw4Q2127P7r3Up1mKZiGEq1BlAJUQQghHCC1sab6se4YUzUIIIYQTyEdAIXQspCmgQNZ5AixpgAwEFEK4IlUubAvXJYubCOEC3L3Bvw4AQXknATCbXfPcI0WzENWQdM8QQgjhMDXqARCYewqAAosUzVWWjDIWwjEk10T5FH74c83zr1NJrgmHCtSK5oAcrWi2uOhVTimay0BGGYuqS1//uCTXREXI4iblJ7kmHKqwpdk/p7B7hrQ0CyGEEA4kvYyEcA2BkQD4FRXN0tIshHAZUmwIIYRwlMLuGb7ZpwGwSEuzEMLluOinfSEAmT5DCFfhVwsAr7zzgAwEFEK4EJk9QwghhMP4hgJgMmfjTY60NAshhBDO4ZonYCGqDXdfMHkDEKKkSp9mIYQrcs1/XEKAdM0XwmUoirW1uSYpMnuGEMKVSLkhqgJ5HwvhMnwKi2YlVeZprspkEnghHENyTQjHkFwTDlfY0hyipEpLc1Umk8CLKktnn/Yl14RwDMk14XC+YQDUVKR7hhDChbjmvyshrkw6aQjhAryDAAgkQ4pmIYTrUKTMEFWBvI2FcB2eAQD4K1mYXbNmlqJZiOrNRf9zCQFI1SyECykqmv+/vfuPibr+4wD+PH4cV8ovdwliCJUzy0xUBmHzV1FOnb+ma67NqDXrj/6okZlbbsysZEbqUmaNRCpbpM0fazgzmZaaxqbQEJwpQ9MMmmj8lFO51/eP5L6Chyefn+/zno/ts3Uf3p+7J8HzfN3d5w508HOaiSh4BOfdFRERBa2oGABAtKODp2cQERHZwcGHgUTqu+WZZg7NRBR0RLFPzyDqH56eQRQ0fOc0t/MvAhJREHFw2KB7AX+PiYLGzaE5Gld5TjMRERERkV/dQ7PjKsKky+Yw2oTk0Dxv3jzEx8djwYIFdkchshXPBSUiIkvcfCMgANyPDhuDaBeSQ/Obb76Jr776yu4YRLaRmy9rB+lpZUS98BeZSHkRTnjDnAAAl3TaHEabkByap0yZgujoaLtjEBGRDo6Q/BeMKHh5w6MAAE7x2JxEG+Xucn755RfMmjULSUlJcDgc2Llz521rCgsLkZqaCpfLhczMTFRUVFgflOiewGfoiIjIGt4IFwAgSq7ZnEQb5Ybm9vZ2jBkzBoWFhX6//t133yE3Nxd5eXk4fvw4xowZg2nTpuGff/7xrUlLS8MTTzxx23bx4kWrvg0ipQk/dYCIiCwm4f8NzU4E59AcYXeA3qZPn47p06f3+fU1a9Zg8eLFeOWVVwAAn332GcrKylBcXIxly5YBAKqqqgzJ4vF44PH8/yWE5uZmAMDly5cNuX49rl+/jo6ODjQ1NSEyMjLkc6iUJRhytHVeR8s1QWtrG5qamtDa2grAvs9tVrlr3VT5uaqaB7A+U3NLC1o8gjbcQFNTk+157kb37zS71jfVfm6q5QHUy3S3ea5fj8D9HoF4m/121iim9UwUBkB27Njhu+zxeCQ8PLzHPhGRl156SWbPnt2v696/f7/Mnz//jmvy8vIE/71+zY1bSGznz5/vV4+Mwq5xC7Wtrq6OXePGzeTN6J45RNR9/7zD4cCOHTswd+5cAMDFixcxdOhQ/Prrr8jKyvKtW7p0KX7++Wf89ttvd3W92dnZ+P3339He3o5BgwZh27ZtPa6vW+9H5P/++y9SUlLw559/IjY2Vt83p1NLSwuSk5Nx/vx5xMTEBD7gHs+hUpZgzCEiaG1tRVJSEsLCrD9rS+WudVPl56pqHkC9TKrlAf57ZnfYsGG4cuUK4uLiLL99di348wDqZVItj1k9U+70DCvs27fvrtZFRUUhKirqtv2xsbFK/FIAQExMjBJZVMkBqJMl2HLY+Q9mMHStmyo/126q5QHUy6RaHgC2PDgF2DU9VMsDqJdJtTxG90y5NwLeidvtRnh4OBobG3vsb2xsRGJiok2piIiIiOheF1RDs9PpxPjx41FeXu7b5/V6UV5e7vf0CiIiIiIiIyh3ekZbWxvOnDnju1xfX4+qqioMGjQIw4YNQ25uLnJycpCeno6MjAysW7cO7e3tvk/TMFNUVBTy8vL8vrRlNVWyqJJDpSzMoZ+K2VXLpFoeQL1MquUB1MukWh5AvUyq5QHUyxQqeZR7I+CBAwcwderU2/bn5OSgpKQEALBhwwZ8/PHHaGhoQFpaGj799FNkZmZanJSIiIiIQoVyQzMRERERkWqC6pxmIiIiIiI7cGgmIiIiIgqAQzMRERERUQAcmnspLCxEamoqXC4XMjMzUVFRccf127Ztw8iRI+FyuTB69Gjs3r3bliw1NTWYP38+UlNT4XA4sG7dOltyFBUVYeLEiYiPj0d8fDyys7MD/j80K8v27duRnp6OuLg4DBgwAGlpafj6668tz3Gr0tJSOBwO31+5tDJHSUkJHA5Hj83lchmSQwuVuqYlk5md05LH7O5pyWRmB7XkuZXRXdSayYpeqtY11XrW30xWdE21nvU3063M6potPTP0j3IHudLSUnE6nVJcXCw1NTWyePFiiYuLk8bGRr/rDx8+LOHh4bJ69Wqpra2V5cuXS2RkpFRXV1uepaKiQpYsWSLffvutJCYmytq1a3Vn0JLjxRdflMLCQqmsrJSTJ0/Kyy+/LLGxsXLhwgXLs+zfv1+2b98utbW1cubMGVm3bp2Eh4fLnj17LM3Rrb6+XoYOHSoTJ06UOXPm6MqgJcfmzZslJiZG/v77b9/W0NCgO4cWKnVNayazOqc1j5nd05rJrA5qzdPN6C7qyWR2L1Xrmmo905LJ7K6p1jMtmbqZ1TW7esah+RYZGRnyxhtv+C53dXVJUlKSrFq1yu/6F154QWbOnNljX2Zmprz++uuWZ7lVSkqKYXcsenKIiNy4cUOio6Plyy+/tD2LiMjYsWNl+fLllue4ceOGTJgwQb744gvJyckx5M6jvzk2b94ssbGxum/XCCp1TWumWxnZOSPyiBjbPaMyiRjTQT15zOiinkxm91K1rqnWM72ZRIzvmmo905rJzK7Z1TOennHTtWvXcOzYMWRnZ/v2hYWFITs7G0eOHPF7zJEjR3qsB4Bp06b1ud7MLGYwIkdHRweuX7+OQYMG2ZpFRFBeXo5Tp05h0qRJlud4//33MXjwYLz66quab9uIHG1tbUhJSUFycjLmzJmDmpoaQ/L0h0pd05PJTCp1z6hMRnVQbx6ju2hEJrN6qVrXVOuZUZmM7JpqPdOTyayu2dkzDs03Xbp0CV1dXUhISOixPyEhAQ0NDX6PaWho6Nd6M7OYwYgc7777LpKSkm67E7YqS3NzMwYOHAin04mZM2di/fr1eO655yzNcejQIWzatAlFRUWab9eIHI8++iiKi4uxa9cubNmyBV6vFxMmTMCFCxcMy3U3VOqankxmUql7ejMZ3UE9eczoot5MZvZSta6p1jOjMhnZNdV6pjWTmV2zs2fK/Rltunfk5+ejtLQUBw4csO0NZ9HR0aiqqkJbWxvKy8uRm5uLhx9+GFOmTLHk9ltbW7Fo0SIUFRXB7XZbcpt9ycrKQlZWlu/yhAkT8Nhjj+Hzzz/HypUrbUxGRlOhe93s7mA3lbp4K/YyuKnSNVV6BqjZNaN6xqH5JrfbjfDwcDQ2NvbY39jYiMTERL/HJCYm9mu9mVnMoCdHQUEB8vPzsW/fPjz55JO2ZQkLC8Pw4cMBAGlpaTh58iRWrVql+Y6kvznq6upw9uxZzJo1y7fP6/UCACIiInDq1Ck88sgjpufwJzIyEmPHjsWZM2f6fft6qNQ1PZnMpFL39GYyuoNa85jVRT2Z/DGyl6p1TbWe6c1kRtdU65mWTGZ3zc6e8fSMm5xOJ8aPH4/y8nLfPq/Xi/Ly8h6PTm6VlZXVYz0A/PTTT32uNzOLGbTmWL16NVauXIk9e/YgPT3d1iy9eb1eeDwey3KMHDkS1dXVqKqq8m2zZ8/G1KlTUVVVheTkZEty+NPV1YXq6moMGTJEUwatVOqankxmUql7ejP1preDWvOY1UU9mfwxspeqdU21nunJZFbXVOuZlkxmd83Wnul+K+E9pLS0VKKioqSkpERqa2vltddek7i4ON/HkixatEiWLVvmW3/48GGJiIiQgoICOXnypOTl5Rn6kXP9yeLxeKSyslIqKytlyJAhsmTJEqmsrJTTp09bmiM/P1+cTqd8//33PT7apbW1VVcOLVk++ugj2bt3r9TV1Ultba0UFBRIRESEFBUVWZqjN6PeRdzfHCtWrJAff/xR6urq5NixY7Jw4UJxuVxSU1OjO4vZ2c3smtZMZnVOax4zu6c1k1kd1JqnNzM+PUO1XqrWNdV6piWT2V1TrWdaMvVmdNfs6hmH5l7Wr18vw4YNE6fTKRkZGXL06FHf1yZPniw5OTk91m/dulVGjBghTqdTRo0aJWVlZbZkqa+vFwC3bZMnT7Y0R0pKit8ceXl5unP0N8t7770nw4cPF5fLJfHx8ZKVlSWlpaWW5+jNyDuP/uR46623fGsTEhJkxowZcvz4cUNyaKFS17RkMrNzWvKY3T0tmczsoJY8vZkxNPc3kxW9VK1rqvWsv5ms6JpqPetvpt7M6JodPXOIiPTvuWkiIiIiotDCc5qJiIiIiALg0ExEREREFACHZiIiIiKiADg0ExEREREFwKGZiIiIiCgADs1ERERERAFwaCYiIiIiCoBDMxERERFRAByaiYiIiIgC4NBMRERERBQAh2aylIhgzZo1eOihh3D//fdj7ty5aG5u7nN9U1MTBg8ejLNnz97xehcuXIhPPvnE4LREwYtdI7IGuxY6HCIidoeg0LFkyRLs2rULmzZtwoABAzB37lwsWLAAa9eu9bs+NzcXra2tKCoquuP1njhxApMmTUJ9fT1iY2PNiE4UVNg1ImuwayFEiCxy9OhRcTgccuzYMd++FStWyIgRI/yub29vl5iYGDly5MhdXX96erps2LDBkKxEwYxdI7IGuxZaeHoGWaagoADPPvssxo0b59uXkJCAS5cu+V2/e/duREVF4amnnrqr6581axZKS0sNyUoUzNg1Imuwa6GFQzNZwuPxoKysDPPmzeuxv7Ozs8+XnQ4ePIjx48ff9W1kZGSgoqICHo9HV1aiYMauEVmDXQs9HJrJEsePH8fVq1fx9ttvY+DAgb5t6dKlGDFihN9jzp07h6SkpLu+jaSkJFy7dg0NDQ1GxSYKOuwakTXYtdDDoZks8ccff2DAgAGorq5GVVWVb0tNTcXTTz/t95irV6/C5XL12PfNN9/0uHM6ePCg72v33XcfAKCjo8O8b4RIcewakTXYtdATYXcACg0tLS1wu90YPny4b9+5c+dw+vRpzJ8/3+8xbrcbV65c6bFv9uzZyMzM9F0eOnSo778vX74MAHjggQeMjE4UVNg1Imuwa6GHQzNZwu12o7m5GSICh8MBAPjwww8xY8YMPP74436PGTt2LLZs2dJjX3R0NKKjo/2uP3HiBB588EG43W5jwxMFEXaNyBrsWujh6RlkiWeeeQadnZ3Iz89HfX09PvjgA/zwww/YuHFjn8dMmzYNNTU1tz0q78vBgwfx/PPPGxWZKCixa0TWYNdCD4dmskRCQgJKSkqwceNGjBo1CkePHsWhQ4eQnJzc5zGjR4/GuHHjsHXr1oDX39nZiZ07d2Lx4sVGxiYKOuwakTXYtdDDvwhISisrK8M777yDEydOICys78d4GzduxI4dO7B3714L0xHdO9g1Imuwa8GL5zST0mbOnInTp0/jr7/+uuOj98jISKxfv97CZET3FnaNyBrsWvDiM81ERERERAHwnGYiIiIiogA4NBMRERERBcChmYiIiIgoAA7NREREREQBcGgmIiIiIgqAQzMRERERUQAcmomIiIiIAuDQTEREREQUAIdmIiIiIqIA/gd+o/JnJkRCGwAAAABJRU5ErkJggg==", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "loam_p = {\n", " \"k_s\": 1.915,\n", " \"theta_r\": 0.21,\n", " \"theta_s\": 0.52,\n", " \"n\": 1.704,\n", " \"alpha\": 0.015,\n", " \"xi\": 5.33,\n", "}\n", "genk_loam = pe.Kool(xi=loam_p.pop(\"xi\"), **loam_p)\n", "gen_loam = pe.Genuchten(**loam_p)\n", "\n", "siltysand = {\n", " \"k_s\": 25.92,\n", " \"theta_r\": 0.07,\n", " \"theta_s\": 0.38,\n", " \"n\": 2.02,\n", " \"alpha\": 0.02,\n", " \"xi\": 1.5,\n", "}\n", "genk_siltysand = pe.Kool(xi=siltysand.pop(\"xi\"), **siltysand)\n", "gen_siltysand = pe.Genuchten(**siltysand)\n", "\n", "siltyloam = {\n", " \"k_s\": 137.3,\n", " \"theta_r\": 0.03,\n", " \"theta_s\": 0.43,\n", " \"n\": 1.228,\n", " \"alpha\": 0.091,\n", " \"xi\": 5.53,\n", "}\n", "genk_siltyloam = pe.Kool(xi=siltyloam.pop(\"xi\"), **siltyloam)\n", "gen_siltyloam = pe.Genuchten(**siltyloam)\n", "\n", "f, axd = plt.subplot_mosaic(\n", " [[\"Loam\", \"SiltySand\", \"SiltyLoam\"]], figsize=(8, 4), sharey=True\n", ")\n", "axd[\"Loam\"].set_yscale(\"log\")\n", "pe.plot.swrc(gen_loam, ax=axd[\"Loam\"], label=rf\"$\\alpha$={gen_loam.alpha:.3f}\")\n", "pe.plot.swrc(genk_loam, ax=axd[\"Loam\"], label=rf\"$\\alpha_w$={genk_loam.alpha_w:.3f}\")\n", "pe.plot.swrc(\n", " gen_siltysand, ax=axd[\"SiltySand\"], label=rf\"$\\alpha$={gen_siltysand.alpha:.3f}\"\n", ")\n", "pe.plot.swrc(\n", " genk_siltysand,\n", " ax=axd[\"SiltySand\"],\n", " label=rf\"$\\alpha_w$={genk_siltysand.alpha_w:.3f}\",\n", ")\n", "pe.plot.swrc(\n", " gen_siltyloam, ax=axd[\"SiltyLoam\"], label=rf\"$\\alpha$={gen_siltyloam.alpha:.3f}\"\n", ")\n", "pe.plot.swrc(\n", " genk_siltyloam,\n", " ax=axd[\"SiltyLoam\"],\n", " label=rf\"$\\alpha_w$={genk_siltyloam.alpha_w:.3f}\",\n", ")\n", "for k in axd:\n", " axd[k].set_title(k, loc=\"left\")\n", " axd[k].set_xlim(0.0)\n", " axd[k].set_xticks([0.0, 0.1, 0.2, 0.3, 0.4, 0.5])\n", " axd[k].set_xlabel(r\"$\\theta$ (-)\")\n", " axd[k].legend(\n", " bbox_to_anchor=(1.05, 0.975), loc=\"lower right\", fontsize=\"small\", frameon=False\n", " )\n", "axd[\"Loam\"].set_ylabel(r\"$\\psi$ (cm)\")\n", "_ = axd[\"Loam\"].set_ylim(1e-1, 1e5)" ] }, { "cell_type": "markdown", "id": "c9588695", "metadata": {}, "source": [ "## References\n", "- Kool, J. B., and Parker, J. C. (1987). Development and evaluation of closed-form expressions for hysteretic soil hydraulic properties. Water Resources Research, 23(1), 105–114. https://doi.org/10.1029/WR023i001p00105\n", "- Bashir, R., Sharma, J., & Stefaniak, H. (2016). Effect of hysteresis of soil-water characteristic curves on infiltration under different climatic conditions. Canadian Geotechnical Journal, 53(2), 273–284. https://doi.org/10.1139/cgj-2015-0004" ] } ], "metadata": { "kernelspec": { "display_name": "pedon (3.13.5)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.13.5" } }, "nbformat": 4, "nbformat_minor": 5 }