{ "cells": [ { "cell_type": "markdown", "id": "92b3c17a", "metadata": {}, "source": [ "# Advanced cubic mixing rules\n", "\n", "\n", "In the advanced cubic mixing rules, the conventional cubic EOS is taken as the basis for the method (usually Peng-Robinson), but different rules are used for the attractive term $a_m$. The formulation reads:\n", "\n", "$$\n", "\\frac{a_m}{b_m} = \\sum_i z_i\\frac{a_i}{b_i} + \\frac{a^{E,\\gamma}_{\\rm res}}{CEoS}\n", "$$\n", "\n", "where $CEoS$ is a scaling parameter that is in principle linked with the EOS coefficients, but can also be allowed to be an adjustable parameter. The $a_i$ and $b_i$ are the pure fluid values of component $i$. The $z_i$ are mole fractions. The mixture covolume is given by \n", "$$\n", "b_m = \\sum_i\\sum_j z_iz_jb_{ij}\n", "$$\n", "with \n", "$$\n", "b_{ij} = \\left(\\frac{b_i^{1/s} + b_j^{1/s}}{2}\\right)^s\n", "$$\n", "\n", "The heart of the method is the definition of $a^{E,\\gamma}_{\\rm res}$, the residual contribution (not in the conventional thermodynamic sense) to the excess Helmholtz energy. There are many possible models here, but one that seems to work well is that of Wilson, for which the expression reads:\n", "\n", "$$\n", "\\frac{a^{E,\\gamma}_{\\rm res}}{RT} = -\\sum_i z_i\\ln\\left(\\sum_j z_j \\Omega_{ji}(T)\\right) - \\sum_iz_i\\ln\\left(\\frac{\\phi_i}{z_i}\\right)\n", "$$\n", "with \n", "$$\n", "\\Omega_{ji}=\\frac{v_j}{v_i}\\exp(-A_{ij}/T)\n", "$$\n", "and \n", "$$\n", "\\frac{\\phi_i}{z_i} = \\frac{v_i}{\\sum_kz_kv_k}\n", "$$\n", "with the $v_i=b_i$. The parameter $A_{ij}\\neq A_{ji}$ in general, and is also given temperature dependence, which is also not supposed to be present according to the derivation. Thus, the models for $A_{ij}$ read something like this here:\n", "$$\n", "A_{ij} = m_{ij}T + n_{ij}\n", "$$\n", "so $m$ is non-dimensional and $n$ has units of temperature." ] }, { "cell_type": "code", "execution_count": 1, "id": "de8cc3c3", "metadata": { "execution": { "iopub.execute_input": "2024-03-15T22:41:25.443825Z", "iopub.status.busy": "2024-03-15T22:41:25.443664Z", "iopub.status.idle": "2024-03-15T22:41:25.901290Z", "shell.execute_reply": "2024-03-15T22:41:25.900689Z" } }, "outputs": [ { "data": { "text/plain": [ "'0.19.1'" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import numpy, matplotlib.pyplot as plt, numpy as np, pandas\n", "import teqp\n", "teqp.__version__" ] }, { "cell_type": "code", "execution_count": 2, "id": "7e7e2f96", "metadata": { "execution": { "iopub.execute_input": "2024-03-15T22:41:25.903639Z", "iopub.status.busy": "2024-03-15T22:41:25.903261Z", "iopub.status.idle": "2024-03-15T22:41:25.909075Z", "shell.execute_reply": "2024-03-15T22:41:25.908670Z" } }, "outputs": [], "source": [ "# Four isotherms of experimental data from doi: 10.1016/j.fluid.2016.05.015\n", "import io, pandas\n", "dat = pandas.read_csv(io.StringIO(\"\"\"PointID y1 uy1 x1 ux1 p/bar up T/K\n", "1 0.0274 0.0007 0.0068 0.0002 59.830 0.053 293.10 \n", "2 0.0664 0.0014 0.0183 0.0004 64.864 0.080 293.10 \n", "3 0.0978 0.0020 0.0298 0.0007 69.772 0.080 293.10 \n", "4 0.1199 0.0024 0.0424 0.0009 74.737 0.080 293.10 \n", "5 0.1219 0.0028 0.1132 0.0023 89.869 0.080 293.10 \n", "6 0.1339 0.0024 0.0995 0.0022 89.198 0.080 293.10 \n", "7 0.1399 0.0026 0.0943 0.0020 88.853 0.080 293.10 \n", "8 0.1461 0.0027 0.0823 0.0019 86.962 0.080 293.10 \n", "9 0.1466 0.0028 0.0778 0.0017 85.942 0.080 293.10 \n", "10 0.1466 0.0028 0.0772 0.0016 85.868 0.080 293.10 \n", "1 0.1378 0.0027 0.0159 0.0004 42.667 0.051 273.08 \n", "2 0.2143 0.0038 0.0297 0.0007 49.547 0.051 273.08 \n", "3 0.2612 0.0043 0.0411 0.0009 55.238 0.051 273.08 \n", "4 0.3209 0.0049 0.0609 0.0013 65.069 0.088 273.08 \n", "5 0.3554 0.0051 0.0786 0.0016 73.395 0.088 273.08 \n", "6 0.3758 0.0052 0.0978 0.0019 81.061 0.088 273.08 \n", "7 0.3903 0.0053 0.1190 0.0023 90.706 0.088 273.08 \n", "8 0.3914 0.0053 0.1477 0.0028 100.966 0.088 273.08 \n", "9 0.3879 0.0053 0.1614 0.0030 104.806 0.088 273.08 \n", "10 0.3724 0.0052 0.1875 0.0033 110.846 0.088 273.08\n", "11 0.3550 0.0051 0.2068 0.0036 114.105 0.088 273.08\n", "12 0.2727 0.0044 0.2531 0.0041 118.020 0.088 273.08\n", "13 0.3343 0.0049 0.2268 0.0038 116.295 0.088 273.08\n", "1 0.2048 0.0038 0.0106 0.0003 25.754 0.050 253.05 \n", "2 0.3019 0.0049 0.0217 0.0005 30.479 0.050 253.05 \n", "3 0.4638 0.0056 0.0436 0.0010 45.352 0.050 253.05 \n", "4 0.5319 0.0056 0.0647 0.0014 58.188 0.050 253.05 \n", "5 0.5854 0.0054 0.1077 0.0021 78.315 0.084 253.05 \n", "6 0.5979 0.0054 0.1497 0.0028 98.276 0.084 253.05\n", "7 0.5898 0.0054 0.1801 0.0032 109.241 0.084 253.05\n", "8 0.5042 0.0057 0.0570 0.0012 51.343 0.084 253.05\n", "9 0.5644 0.0055 0.0861 0.0017 67.594 0.084 253.05\n", "10 0.5949 0.0054 0.1267 0.0024 86.883 0.084 253.05\n", "11 0.5826 0.0054 0.2015 0.0035 116.614 0.084 253.05\n", "12 0.5537 0.0055 0.2431 0.0040 129.873 0.084 253.05\n", "13 0.4973 0.0055 0.2971 0.0046 139.161 0.084 253.05\n", "14 0.4971 0.0055 0.2972 0.0046 139.261 0.084 253.05\n", "1 0.7076 0.0050 0.0257 0.0006 27.983 0.056 223.10\n", "2 0.7774 0.0041 0.0522 0.0011 44.918 0.056 223.10\n", "3 0.8077 0.0036 0.0930 0.0019 64.906 0.081 223.10\n", "4 0.8131 0.0035 0.1261 0.0024 84.799 0.081 223.10\n", "5 0.8057 0.0035 0.1584 0.0029 104.410 0.081 223.10\n", "6 0.7843 0.0038 0.1982 0.0035 125.782 0.081 223.10\n", "7 0.7533 0.0041 0.2380 0.0040 144.287 0.081 223.10\n", "8 0.7150 0.0045 0.2813 0.0044 159.015 0.081 223.10\n", "9 0.6942 0.0047 0.3064 0.0047 165.347 0.081 223.10\n", "\"\"\"), sep='\\s+', engine='python')" ] }, { "cell_type": "code", "execution_count": 3, "id": "49c01826", "metadata": { "execution": { "iopub.execute_input": "2024-03-15T22:41:25.910898Z", "iopub.status.busy": "2024-03-15T22:41:25.910738Z", "iopub.status.idle": "2024-03-15T22:41:26.343925Z", "shell.execute_reply": "2024-03-15T22:41:26.343490Z" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjUAAAG5CAYAAAB2j8WmAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/H5lhTAAAACXBIWXMAAA9hAAAPYQGoP6dpAACkFUlEQVR4nOzdd1hT9/fA8XfCVhFFQMQtKm7EvfdWcNRZrVpn1baOtlptq1W/ztpq1f5c1brq3qNqBQdqnThxoOBWUEQFRHbu7w9KaiQIKBAI5/U8eZCbe8NJC8nJZ5yjUhRFQQghhBAim1MbOgAhhBBCiPQgSY0QQgghjIIkNUIIIYQwCpLUCCGEEMIoSFIjhBBCCKMgSY0QQgghjIIkNUIIIYQwCpLUCCGEEMIoSFIjhBBCCKMgSY0QQgghjEKWSWpmzJhBzZo1sba2xsHBgU6dOuHn56dzTpMmTVCpVDq3zz77zEARCyGEECIryTJJzdGjRxkxYgSnTp3i4MGDxMbG0qpVKyIiInTOGzx4MIGBgdrb7NmzDRSxEEIIIbISU0MHkGj//v06369cuRIHBwd8fHxo1KiR9niuXLlwdHTM7PCEEEIIkcVlmaTmbaGhoQDY2trqHP/zzz9Zu3Ytjo6OuLu788MPP5ArVy69jxEdHU10dLT2e41Gw/PnzylQoAAqlSrjghdCCCFEulEUhfDwcJycnFCrk59kUimKomRiXKmi0Wjw8PDg5cuXHD9+XHt86dKlFC9eHCcnJy5fvsy4ceOoVasW27Zt0/s4P/74I5MnT86ssIUQQgiRgR48eECRIkWSvT9LJjXDhg1j3759HD9+/J3BHzp0iObNm+Pv74+zs3OS+98eqQkNDaVYsWI8ePCAvHnzZkjsQgghhEhfYWFhFC1alJcvX2JjY5PseVlu+unzzz9nz549eHt7vzOhAahduzZAskmNhYUFFhYWSY7nzZtXkhohhBAim0lp6UiWSWoUReGLL75g+/btHDlyhJIlS6Z4zcWLFwEoVKhQBkcnhBBCiKwuyyQ1I0aMYN26dezcuRNra2uCgoIAsLGxwcrKioCAANatW0e7du0oUKAAly9fZvTo0TRq1IgqVaoYOHohhBBCGFqWWVOT3JDSH3/8Qf/+/Xnw4AF9+vTB19eXiIgIihYtSufOnfn+++9TPZUUFhaGjY0NoaGhMv0khBBCZBOpff/OMiM1KeVWRYsW5ejRo5kUjRBCCCGymyxTUVgIIYQQ4kNIUiOEEEIIoyBJjRBCCCGMgiQ1QgghhDAKktQIIYQQwihIUiOEEEIIoyBJjRBCCCGMgiQ1QgghhDAKktQIIYQQwihIUiOEEEIIoyBJjRBCCCGMgiQ1QgghhDAKktQIIYQQwihIUiOEEEIIoyBJjRBCCCGMgiQ1QgghhDAKktQIIYQQwihIUiOEEEIIoyBJjRBCCCGMgiQ1QgghhDAKktQIIYQQwihIUiOEEEIIoyBJjRBCCCGMgiQ1QgghhDAKktQIIYQQwihIUiOEEEIIoyBJjRBCCCGMgiQ1QgghhDAKktQIIYQQwihIUiOEEEIIoyBJjRBCCCGMgiQ1QgghhDAKktQIIYQQwihIUiOEEEIIoyBJjRBCCCGMgiQ1QgghhDAKktQIIYQQwihIUiOEEEIIoyBJjRBCCCGMgiQ1QgghhDAKktQIIYQQwihIUiOEEEIIoyBJjRBCCCGMgiQ1QgghhDAKktQIIYQQwihIUiOEEEIIoyBJjRBCCCGMgiQ1QgghhDAKktQIIYQQwihIUiOEEEIIoyBJjRBCCCGMgiQ1QgghhDAKktQIIYQQwihIUiOEEEIIoyBJjRBCCCGMgiQ1QgghhDAKktQIIYQQwihIUiOEEEIIoyBJjRBCCCGMgiQ1QgghhDAKktQIIYQQwihIUiOEEEIIoyBJjRBCCCGMQpZJambMmEHNmjWxtrbGwcGBTp064efnp3NOVFQUI0aMoECBAuTJk4ePPvqIJ0+eGChiIYQQQmQlWSapOXr0KCNGjODUqVMcPHiQ2NhYWrVqRUREhPac0aNHs3v3bjZv3szRo0d5/PgxXbp0MWDUQgghhMgqVIqiKIYOQp/g4GAcHBw4evQojRo1IjQ0FHt7e9atW0fXrl0BuHHjBuXLl+fkyZPUqVMnxccMCwvDxsaG0NBQ8ubNm9FPQQghhBDpILXv31lmpOZtoaGhANja2gLg4+NDbGwsLVq00J5Trlw5ihUrxsmTJ/U+RnR0NGFhYTo3IYQQQhinLJnUaDQaRo0aRf369alUqRIAQUFBmJubky9fPp1zCxYsSFBQkN7HmTFjBjY2Ntpb0aJFMzp0IYQQQhhIlkxqRowYga+vLxs2bPigxxk/fjyhoaHa24MHD9IpQiGEEEJkNaaGDuBtn3/+OXv27MHb25siRYpojzs6OhITE8PLly91RmuePHmCo6Oj3seysLDAwsIio0MWQgghRBaQZUZqFEXh888/Z/v27Rw6dIiSJUvq3F+9enXMzMzw8vLSHvPz8+P+/fvUrVs3s8MVQgghRBaTZUZqRowYwbp169i5cyfW1tbadTI2NjZYWVlhY2PDwIEDGTNmDLa2tuTNm5cvvviCunXrpmrnkxDZQXx8PM+ePSMoKEh7Cw4OJjw8nPDwcF69eqX9Gh0dTXx8PIcOHQISEv98+fJhaWmJpaUlVlZW2n/ny5ePAgUKYGdnh52dnfbfhQoVwtLS0sDPWggh0keW2dKtUqn0Hv/jjz/o378/kFB876uvvmL9+vVER0fTunVr/u///i/Z6ae3yZZukRW8ePGC69ev4+/vT0BAgPbrvXv3ePr0KRqNJlPjKViwIMWKFaNo0aIUK1aMMmXKUK5cOcqVK0ehQoWS/dsUQojMktr37yyT1GQGSWpEZtJoNAQEBHDx4kUuXbqkvaW0YF2lUuHg4ICjoyOOjo7Y29uTN29e8uTJQ548ebC2tiZ37txYWlpiamrKoUOHUKlUNG/enPj4eKKiooiMjCQqKkr775cvX/Ls2TOePXtGSEiI9t9RUVHvjMXa2hoXFxfKlSuHm5sb1atXx83NTf5+hBCZSpIaPSSpERkpOjqas2fPcuzYMY4dO8aJEyeSrY1UtGhRypQpg7OzM6VLl8bZ2ZmSJUvi5OSEnZ0dpqYZPzOsKArPnz/n/v37PHjwgPv373P37l1u3rzJjRs3uH37NvHx8UmuU6lUlClThurVq1O9enXq1q1LjRo1MDc3z/CYhRA5kyQ1ekhSI9JTbGwsx48fx9PTk2PHjnHmzBmio6N1zrG0tKRy5cq4urpqb1WqVMHGxsZAUadedHQ0AQEB+Pn54evri4+PDz4+Pjx8+DDJuVZWVtSvX58mTZrQpEkTatasKUmOECLdSFKjhyQ14kM9e/aM/fv3s3v3bg4cOKCtfJ3IwcGBhg0bam9VqlTJlFGXzPT06VNtgnPu3DlOnDjBs2fPdM55M8lp06YN1apVk7U5Qoj3JkmNHpLUiPdx/fp1du7cyZ49ezh58qTOQl57e3vatGlD48aNadiwIWXKlMlxb94ajYZr165x5MgRjh49ypEjR5IkOYULF8bDwwMPDw+aNm0q9aOEEGkiSY0ektSI1Lp37x4bNmxg3bp1XL58Wec+V1dXOnToQIcOHahVqxZqdZYp95QlKIqiTXK8vLz4+++/iYiI0N6fJ08e2rRpQ8eOHWnXrp22v5sQQiRHkho9JKkR7/Lq1Ss2bdrEihUrOHHihPa4mZkZLVq0wN3dnfbt21OsWDEDRpn9REVFcfjwYXbu3MmuXbsIDAzU3mdiYkKjRo3o1asXH330kSQ4Qgi9JKnRQ5Ia8TZFUTh16hTLly9n48aNvHr1CkjY4dOkSRN5s01nGo0GHx8fdu3axc6dO7ly5Yr2PjMzM9q2bUuvXr3w8PAgV65cBoxUCJGVSFKjhyQ1IlFkZCTr169nwYIFXLx4UXu8TJkyDBw4kD59+lC4cGHDBZhD3L59m82bNyeZ5subNy+9evViwIAB1KxZM8etUxJC6JKkRg9JasSjR4+YP38+v//+O8+fPwcStl336NGDgQMH0qBBA3kDNZCrV6+yfv16/vzzT+7evas9XrFiRQYOHEi/fv1kxEyIHEqSGj0kqcm5rl+/zk8//cTatWuJjY0FoESJEgwfPpyBAwfKm2UWotFoOHLkCCtWrGDr1q3aqseWlpb06tWL4cOHU6NGDQNHKYTITJLU6CFJTc5z5swZZsyYwY4dO7THGjVqxJgxY+jQoQMmJiaGC06k6OXLl2zYsIHFixdz6dIl7fGaNWsyfPhwevXqJdvDhcgBJKnRQ5KanOP8+fNMnDiRvXv3ao917NiRcePGUbduXQNGJt5H4oLu3377jc2bNxMTEwOAo6Mjn3/+OZ999hkFChQwcJRCiIwiSY0ektQYv8uXLzNp0iTtyIxareaTTz5h3LhxlC9f3rDBiXTx9OlTli9fzm+//cajR4+AhArGn376KaNHj6Z06dIGjlAIkd5S+/4tVcOEUXjw4AH9+vWjatWq7NixA5VKRe/evbl+/TorV66UhMaIODg4MH78eG7fvs2aNWuoWrUqkZGR/N///R9ly5alZ8+eSQomCiFyBklqRLYWHh7Od999R9myZVm9ejWKotCtWzd8fX1Zu3YtZcuWNXSIIoOYm5vTp08fzp8/z6FDh2jfvj2KorBx40ZcXV3p2LEjZ86cMXSYQohMJEmNyJY0Gg2///47pUuXZvr06URFRdGwYUPOnDnDpk2bqFChgqFDFJlEpVLRtGlT9uzZw6VLl+jRowcqlYpdu3ZRu3ZtWrduzalTpwwdphAiE0hSI7Kd8+fPU69ePQYPHszTp08pU6YM27dv5+jRo9SsWdPQ4QkDqlKlChs2bOD69ev069cPExMT/v77b+rWrYu7u7tOoUUhhPGRpEZkGy9fvuSLL76gZs2anD59Gmtra3755ReuXr1Kp06dpGie0HJxcWHlypXcunWLAQMGYGJiwp49e3Bzc6NHjx7cuHHD0CEKITKAJDUiy1MUhQ0bNuDi4sLChQvRaDT07NmTGzduMHr0aMzMzAwdosiiSpYsyfLly7l27Ro9e/YEYNOmTVSsWJFPP/1Uu3tKCGEcJKkRWVpgYCCdO3emV69ePH36FBcXFzw9PVm/fj1OTk6GDk9kE2XLlmX9+vVcunQJDw8PNBoNK1eupEyZMkycOJHw8HBDhyiESAeS1IgsSVEUVq9eTcWKFdm5cyempqb8+OOPXL58mebNmxs6PJFNValShZ07d3Ly5Enq169PZGQkU6dOpUyZMixbtoy4uDhDhyiE+ACS1IgsJzAwEHd3d/r168eLFy+oVq0aPj4+TJo0CXNzc0OHJ4xAnTp1OHbsGFu3bqV06dI8efKEIUOG4ObmxsGDBw0dnhDiPUlSI7KUXbt2UblyZfbu3Yu5uTnTp0/n9OnTVKlSxdChCSOjUqno0qULV69eZd68edja2uLr60urVq346KOPuHfvnqFDFEKkkSQ1Ikt4/fo1w4cPp2PHjoSEhFC1alUuXLjA+PHjMTU1NXR4woiZm5szcuRI/P39GTlyJCYmJmzbto3y5cszdepUbZdwIUTWJ0mNMLhLly5Ro0YNFi1aBMBXX33FqVOnpICeyFT58+dn3rx5XLhwgcaNGxMZGcnEiROpWLEiu3fvNnR4QohUkKRGGIyiKPz+++/Url2b69evU6hQIf7++2/mzJmDhYWFocMTOVTlypU5fPiwdofd7du38fDwoEuXLrIFXIgsTpIaYRCRkZEMHDiQwYMHEx0dTYcOHbh8+TItW7Y0dGhCoFKp6NmzJ35+fowdOxZTU1O2b99O+fLl+e2334iPjzd0iEIIPSSpEZkuICCAunXr8scff6BWq5kxYwY7d+7Ezs7O0KEJoSNPnjzMmjWL8+fPU6dOHcLDw/n888+pX7++dAIXIguSpEZkqr1791K9enUuXbqEg4MDnp6efPvtt6jV8qsosq7KlStz/PhxfvvtN/Lmzcvp06epXr063333HdHR0YYOTwjxL3knEZlCURR++eUX3N3dCQ0NpV69epw/f56mTZsaOjQhUsXExIThw4dz7do1unTpQlxcHNOnT6d69eqcO3fO0OEJIZCkRmSCmJgYBg0axFdffYWiKAwZMoTDhw9TuHBhQ4cmRJoVLlyYrVu3sm3bNhwcHLh69Sp16tTh+++/l1EbIQxMkhqRoZ49e0bLli1ZsWIFarWaX3/9lcWLF0tlYJHtde7cmatXr9KzZ0/i4+OZNm0aNWrU4Pz584YOTYgcS5IakWGuX79O7dq18fb2Jm/evOzdu5cvv/wSlUpl6NCESBd2dnasX7+eLVu2YG9vj6+vL7Vq1WLSpEnExsYaOjwhchxJakSGSGwYePv2bUqVKsXJkydp06aNocMSIkN89NFHXL16le7duxMfH8+UKVOoX78+t27dMnRoQuQoktSIdLdnzx6aN2/OixcvqFOnDqdPn5bqwMLo2dvbs3HjRjZu3Ej+/Pk5e/Ysbm5uLF++HEVRDB2eEDmCJDUiXf3xxx906tSJyMhI2rdvj6enp9SfETlK9+7duXz5Mk2bNiUiIoJBgwbRtWtXQkJCDB2aEEZPkhqRLhRFYebMmQwYMID4+Hj69evH9u3byZ07t6FDEyLTFSlSBE9PT2bPno2ZmRnbtm2jSpUqeHp6Gjo0IYyaJDXigymKwrfffsv48eMBGDduHH/88QdmZmYGjkwIw1Gr1XzzzTecOnUKFxcXHj9+TMuWLRk/fjxxcXGGDk8IoyRJjfggiqIwZswYZs+eDcDPP//MzJkzZYeTEP+qVq0a58+fZ+jQoQDMnDmTpk2b8vDhQwNHJoTxkaRGvDeNRsMXX3zBvHnzAFi0aBFjxowxbFBCZEG5cuVi8eLFbNq0CWtra44fP07VqlXZt29fitfGaxROBoSw8+IjTgaEEK+RRcdCJEel5KBl+WFhYdjY2BAaGkrevHkNHU62ptFoGDZsGEuXLkWlUrFs2TIGNikF+8ZB21ngLO0PhNDH39+f7t27c+HCBQC+/fZbpk6diqmpaZJz9/sGMnn3NQJDo7THCtlYMsm9Am0qFcq0mIUwtNS+f8tIjUgzjUbD4MGDtQnNH3/8wcABA8BrMjzzS/iac3JlIdKkdOnS/PPPP4wYMQL4bzrq0aNHOuft9w1k2NrzOgkNQFBoFMPWnme/b2CmxSxEdiFJjUgTRVEYMWKEtu3B2rVr6devHwR4weOET548vpDwvRBCL0tLSxYuXKgzHVWtWjWOHj0KJEw5Td59DX0fDRKPTd59TaaihHiLJDUi1RRFYezYsSxevBiVSsXatWv5+OOPE0ZlDv0PVCYJJ6pMEr6X0Roh3qlbt26cP38eV1dXnj59SvPmzZk3bx5n7oQkGaF5kwIEhkZx5s7zzAtWiGxAkhqRatOmTWPOnDkALFu2jF69eiXckThKo8QnfK/Ey2iNEKmUOB318ccfEx8fz+jRo5k085dUXfs0PPnER4icSJIakSq//vorP/zwAwBz585l4MCBCXe8PUqTSEZrhEi1XLlysXbtWn799VdMTU3xPrAnVdc5WFtmcGRCZC+S1IgU/fHHH4waNQqAyZMna/8NJB2lSSSjNUKkiUql4ssvv8TLy4t8scHEhQWjKBr955KwC6pWSdvMDVKILE6SGvFOf/31F4MHDwbgq6++0o7WAP+N0iT7a6SW0Roh0qhRo0b4nDuH3YMjgCpJYpNY1nKSewVM1FLkUog3SVIjkuXj40P37t21vZx++ukn3UrB8TEQ+gjQ/2kSNBD2KOE8IUSqFS5cmNNbluD66gzx4bqNMB1tLFnUp1qKdWqkaJ/IiaT4ntDr7t271K1bl4pWwSzvWgCngWsxc2mZ9MTQhxDxLPkHym0PNoUzLlAhjJiiKMz7dT7fzVuBKlc+KpQswq7lv1DIseA7r5OifcLYpPb9W5IakcSLFy+oX78+169f58qXdlTKHwNObjD4MEhPJyEy3YEDB+jRowehoaEUL16cXbt2UaVKFb3nJhbte/uFPfEvNzWjPEJkNVJRWLyX6OhoOnXqxPXr1/m4VsGEhAZk0a8QBtS6dWtOnTpF6dKluXfvHvXq1WPnzp1JzpOifSKnk6RGaCVWC/b29iZv3rws+7i4FNQTIosoV64cp0+fpnnz5kRERNC5c2d+/fVXnXPO3HkuRftEjiZJjdD67bffWL58OWq1Gq+lE8j14oYU1BMiC7G1tWXfvn0MHToURVEYNWoUo0aNIj4+4e80tcX4pGifMFaS1AgAjhw5oq0/M2vmTGqE7pOCekJkQWZmZixatIhZs2YBCYUxu3btyuvXr1NdjE+K9gljJUmN4O7du3Tr1o34+Hh69+7NV53cpKCeEFmYSqVi7NixbNy4EQsLC3bs2EHTpk0pkTuOQjaWJLecX4r2CWMnSU0OFxERQadOnXj27BnVqlVj2dKlqA5LQT0hsoPu3bvj6emJra0tZ86coX69ugysZgOQJLGRon0iJ5CkJgdTFIWBAwdy6dIlHBwc2L59O1bmJlJQT4hspEGDBpw8eZJSpUpx584dxvZqxRfVrHC00Z1iSm3RPiGyM6lTk4MtWLCAL7/8ElNTUw4dOkTDhg0T7pCCekJkO8HBwXh4eHDq1CksLCxYv2EjjpUb8DQ8CgfrhCknGaER2ZUU39NDkpr/+Pj4UK9ePWJiYpg3bx4jR440dEhCiA8UGRlJz5492bVrF2q1mmXLljFgwABDhyXEB5PieyJZYWFh9OjRg5iYGDp16sSXX35p6JCEEOnAysqKrVu3MmDAADQaDQMHDmTmzJnkoM+uIoeTpCaHURSFwYMHExAQQPHixVmxYoVuk0ohRLZmamrK77//zrfffgvA+PHj+eqrr9BoklsnJ4TxkKQmh1m6dCmbNm3C1NSUDRs2kD9/fkOHJIRIZyqVihkzZvDLL78AMHfuXPr160dsbKyBIxMiY2WZpMbb2xt3d3ecnJxQqVTs2LFD5/7+/fujUql0bm3atDFMsNnUpUuXtGtnZsyYQZ06dQwckRAiI40ePZo1a9ZgamrK2rVr6dixIxEREYYOS4gMk2WSmoiICFxdXfntt9+SPadNmzYEBgZqb+vXr8/ECLO3qKgoPv74Y6Kjo2nfvj1jxowxdEhCpMqzZ894/fq1ocPItvr06cOuXbuwsrJi3759tG7dmtDQUEOHJUSGMDV0AInatm1L27Zt33mOhYUFjo6OmRSRcfnuu++4du0ajo6OrFy5ErU6y+SzQiQrLCwMe3t7VCqVrAn5AG3btsXLy4t27dpx4sQJWrRowf79+ylQoIChQxMiXWWrd7YjR47g4OCAi4sLw4YNIyQk5J3nR0dHExYWpnPLiY4ePcrcuXMB+P3337GzszNwREKkztmzZ4GEBe4HDhwwcDTZW926dTl8+DB2dnacO3eOJk2aEBQUZOiwhEhX2SapadOmDatXr8bLy4tZs2Zx9OhR2rZtq+1Oq8+MGTOwsbHR3ooWLZqJEWcN4eHh9O/fX1s9uH379oYOSYhUK126tPbfsobuw1WtWhVvb28KFSqEr68vjRo14sGDB4YOS4h0kyWL76lUKrZv306nTp2SPef27ds4Ozvj6elJ8+bN9Z4THR1NdHS09vuwsDCKFi2ao4rvDRkyhGXLllGiRAkuXbqUY563yHqio6N58OAB9+/f5969ewQGBvL8+XNCQkK0X1+8eEF0dDQxMTHExMQQFhZGZGSk9jEcHBwwNzfX3mxsbLC1tSV//vzarwUKFKBIkSIUKVKEokWLUrBgQUxMTN4RWc4TEBBA8+bNuXfvHsWLF8fLywtnZ2dDhyVEslJbfC/LrKlJq1KlSmFnZ4e/v3+ySY2FhQUWFhaZHFnWsXfvXpYtW4ZKpWLlypWS0IhMERYWxtWrV7ly5Yr2dvPmTYKCgj64CNzTp0/TfI2pqSlOTk44OztTrlw5XFxcKFeuHOXKlaNo0aI5cn2Zs7Mzx44do3nz5ty6dYuGDRvi6elJhQoVDB2aEB8k2yY1Dx8+JCQkhEKFpDmbPi9evGDQoEEAjBo1isaNGxs4ImGMFEXh2rVreHt74+3tzcmTJ7l3716y51tZWVGsWDGKFy9O4cKFKVCgAAUKFMDW1pYCBQqQP39+rKysdEZjypUrp73+ypUr2lGcqKgoXr58yYsXL3j+/Ln2a3BwMI8ePeLhw4c8fvyYuLg47t+/z/379zl8+HCSeFxcXKhevTo1a9akRo0aVK5cGXNz8wz7b5ZVFC1aFG9vb1q2bImvry+NGzfGy8uLKlWqpOr6eI3CmTvPpbeUyFKyzPTTq1ev8Pf3B8DNzY1ffvmFpk2bYmtri62tLZMnT+ajjz7C0dGRgIAAxo4dS3h4OFeuXEn1aIzR934KOAz7xkHbWQydvYGlS5fi4uLChQsXsLKyMnR0wggoisKlS5c4cuQI3t7eHDt2jGfPkjY/dXJyolKlSlSuXJnKlStTvnx5SpQood3JlBYODg4EBwe/1w6ouLg4njx5wv3797l16xZ+fn7cuHGDGzducOvWLb3F6MzNzXF1daVmzZrUrFmTWrVqUb58eaOtvB0SEkKbNm04d+4cBQoUwMvLC1dX13des983kMm7rxEYGqU9VsjGkknuFaQLuMgQ2a6h5ZEjR2jatGmS4/369WPRokV06tSJCxcu8PLlS5ycnGjVqhVTp06lYMGCqf4ZRp3UKAosawqPLxBuXYa8X/sACUUNtd23hXgPUVFRHD58mF27drF7924ePXqkc7+VlRV169alUaNGNGzYEFdX13TdKhwTE8OGDRvo2rUruXLlSrfHjYuL4+7du1y5coVz585x9uxZzp07x4sXL5Kca2dnR6NGjWjSpAmNGzemUqVKRjVt9fLlS1q3bs2ZM2dSTGz2+wYybO153n7jSEz5FvWpJomNSHfZLqnJDEad1Ph7wtqPtN+2XhtBsaafsmzZMgMGJbKruLg4Dh8+zPr169m6datOOYRcuXLRuHFjGjVqRKNGjahRo4bRTNcoisLt27d1kpwzZ87oLFaGhNGjli1b0rp1a1q1apWmD1dZVWhoKK1ateLMmTPY2tri5eVF1apVdc6J1yg0mHVIZ4TmTSrA0caS4+OayVSUSFeS1OhhtElN4ihN4GVQ4onTKFx5pqb41ABspbiWSIOAgACWLVvGypUrefLkifa4k5MTHh4eeHh40LRpUywtLQ0YZeaKiYnh3LlzHD16lCNHjnDixIkkrQaqVq2Kh4cHHTt2xM3NLdtOVYWGhtK6dWtOnz6Nra0tnp6euLm5ae8/GRBCr2WnUnyc9YPrUNdZXntE+pGkRg+jTWreGqXR6rMVSrfI/HhEthIXF8euXbtYvHgxBw8e1B63tbWle/fufPzxx9SvX9+opls+RExMDCdPnuTAgQMcOHCA8+fP69xftGhRbYLTpEkTzMzMDBTp+wkNDaVNmzacOnWK/Pnz4+npSbVq1QDYefERIzdcTPExfu1ZlY5VC2dwpCInkaRGD6NMav4dpVECL6NS/itEqKhMUBWqAoMPQzb91Cgy1uvXr1mxYgVz5szR2bHUunVrhg4dSvv27Y1mWikjPX36lH379rFjxw7+/vtvnT5VBQoU4KOPPqJnz540atQo29TLCQsLo3Xr1kkSGxmpEYYiSY0eRpnUJDdKk0hGa8RbXrx4wYIFC1iwYIF255KdnR2DBg1i8ODBlCpVysARZl+RkZF4enqyc+dOdu3aRXBwsPY+R0dHunXrRq9evahTp84HT1Fl9JbqsLAw2rRpw8mTJ8mfPz9HjhyhYqXKNJh1iKDQqCQLhUHW1IiMI0mNHkaX1CSO0jy+iErvS4wanFxltEYACSMz8+fPZ9asWbx8+RKAkiVLMu/LTnQw8Ubd7idwTroDUbyfuLg4jhw5woYNG9i2bZvOrioXFxf69+/PJ598QuHCaZ+m2e8byI+7rhEU9t+CXce8lvzokb5bqt8csbG3t+fo0aPci8/HsLUJU25vvurI7ieRkVL7/i2T5NlZfAyEPkomoQHQQNijhPNEjhUbG8vixYspXbo048eP5+XLl1SqVIn169dz088PD5OjqENuwV9fJyTKIl2YmprSokULfv/9d4KCgtizZw+9e/cmV65c+Pn5MX78eIoVK0bbtm3ZtGmTTkuXd9nvG8hna8/rJDQAQWFRfLb2PPt9A9PtOeTNm5d9+/ZRrVo1goODad68OaUtI1jUpxqONrqLxR1tLCWhEQb3QSM1165d4/79+8TE6L5penh4fHBgGcHoRmqAx9fP0LVdE6Kiopk5cyatWrbUPSG3PdjIgr2c6tixYwwfPhxfX18ASpQowdSpU+nVq1fC+o5bnvDnG9OXvbdCGZmuzEjh4eFs3ryZP/74g+PHj2uP29vbM2jQIIYMGUKJEiX0XhuvUaj+v4O8fJ20aGCi/LnMOPd9y3Sd/gkJCaFJkyb4+vpSrFgxvL29KVK0mFQUFpkmQ6efbt++TefOnbly5QoqlUrbzyVxjvhdnbMNyRiTmt69e7Nu3ToaNmzI0aNHs+1WUpG+njx5wtixY1m9ejWQsGD1xx9/ZMiQIf8t/lUUWFgdQgL+u7CAM3zuI9OVmcTf35+VK1eycuVKbVFDlUpF+/btGT58OK1bt9bZdXbC/xm9fz+d4uP+Oag29UvbpWusT548oXHjxvj5+VGqVCm8vb3fa+pMiPeRodNPI0eOpGTJkjx9+pRcuXJx9epVvL29qVGjBkeOHHnfmEUanTx5knXr1qFSqZg3b54kNAJFUViyZAnlypVj9erVqFQqhgwZgp+fH59//rnubiZ/L92EBhK+9/fK3KBzsNKlS/O///2Pu3fvsm3bNlq0aIGiKOzZs4d27dpRunRp5syZQ2hoKJBQJyY1UnteWhQsWBAvLy9KlSrF7du3ad68uU4tIyGygvdKak6ePMmUKVOws7NDrVajVqtp0KABM2bM4Msvv0zvGIUeGo2GUaNGAfDpp59q60iInCswMJD27dvz2Wef8fLlS6pVq8apU6dYsmRJ0rYFigL7x+p/oP1jZW1NJjM1NaVz584cPHgQPz8/Ro8eTb58+bhz5w7ffPMNRYsWZcyYMYSGvkzlI2bM/7/ChQvj5eVF0aJF8fPzo0WLFoSEpH8CJcT7eq+kJj4+HmtrayBhK+jjx48BKF68OH5+fukXnUjWli1bOHPmDHny5GHatGmGDkcY2NatW6lcuTL79u3D0tKSuXPncubMGWrVqqX/An2jNIlktMagypYtyy+//MKjR4/4/fffqVixIuHh4cydO5dfvv0sVY9Rt1T6Tj29qUSJEnh5eVGoUCF8fX1p164d4eHhGfbzhEiL90pqKlWqxKVLlwCoXbs2s2fP5sSJE0yZMkVqXGSCuLg4fvjhBwC++eYbHB0dDRyRMJTXr1/Tv39/unbtSkhICG5ubvj4+DBq1KjkC729a5QmkYzWGFyuXLkYOHAgV65cYd++fTRv3pzXdy8T/zqUdy2FzJfLjDoZXPiuTJkyeHp6Ymtry5kzZ+jcuXOqd28JkZHeK6n5/vvv0Wg0AEyZMoU7d+7QsGFD/vrrL+bPn5+uAYqk1qxZw82bNylQoIB2CkrkPLdv36Zu3bqsWrUKtVrNd999x6lTp6hQocK7L4yLhhf33n3Oi3sJ5wmDU6lUtGnTBk9PTy6c96FS5BWAZBObmV0qZ8oupAoVKrBv3z5y586Nl5cXvXv3zrKbRETOkW7F954/f07+/Pmz9GJVY9j9FB0djYuLC/fu3eOnn37i66+/NnRIwgD+/vtvevbsyYsXL3BwcGDTpk00btw4dRfHRcPP5SHyHWshchWAMdfB1CJ9Ahbpau0RX6bu8yNa9V+tGFVUKJ9WycMPn7pn6uuwp6cn7du3JyYmhoEDB7Js2bIs/T4gsqcM2dKt0Wj46aef2LVrFzExMTRv3pxJkyZhZWWVLkFnNGNIan777Tc+//xznJyc8Pf3zzb/7UX6UBSF2bNnM378eBRFoXbt2mzZsoUiRYqk7YFCH0LEs+Tvl/pGWV68RmG/jz+rN+/gwI5NhPqfB0WDq6srkyZNolOnTpmWXGzbto1u3bqh0WgYO3Yss2bNypSfK3KOVL9/K2kwZcoURa1WK61atVI6duyoWFpaKp9++mlaHsKgQkNDFUAJDQ01dCjvJSIiQnF0dFQA5f/+7/8MHY7IZHFxccrw4cMVEra2KIMHD1aioqIMHZbIAp4+fap8++23Sp48ebS/H7Vq1VIOHTqUaTH8/vvv2p89c+bMTPu5ImdI7ft3mpKa0qVLK4sXL9Z+f/DgQcXc3FyJj49/vygzWXZPambNmqUASokSJZTo6GhDhyMyUVRUlNK9e3cFUFQqlTJ//nxDhySyoJCQEOW7775TcuXKpU0wWrVqpZw7dy5Tfv7s2bO1P3fp0qWZ8jNFzpAhSY25ubly//59nWMWFhbKgwcP0h6hAWTnpCY0NFSxtbVVAGXVqlWGDkdkorCwMKVFixYKoJiZmSkbNmwwdEgiiwsKClI+//xzxczMTJtkdO/eXfHz8/vgx46L1yj/+D9Tdlx4qPzj/0yJi9fo3D9u3Dht8r1ly5YP/nlCKErq37/TtPspLi4OS0vdJmZmZmbExibfh0Skj0WLFvH8+XPKlStH7969DR2OSE8Bh2FhrYSvbwkJCaFZs2Z4enqSO3du9u7dS48ePQwQpMhOChYsyIIFC7hx4wZ9+vRBpVKxadMmKlSowNChQ9+7EvB+30AazDpEr2WnGLnhIr2WnaLBrEM6TTRnzJjB4MGDURSF3r174+3tnV5PS4gUpWmhsFqtpm3btlhY/LcjYvfu3TRr1ozcuXNrj23bti19o0wn2XWhcFRUFCVLliQoKIiVK1fSr18/Q4ck0ouiwLKm8PgCOLnB4MPavkthYWE0b96cc+fOUaBAAfbt20fNmjUNHLDIji5fvsx3333Hnj17ALC2tub7779n5MiROq/n77LfN5Bha88nqVWcuBT5zQ7d8fHxdO3alR07dpAvXz6OHz9OxYoV0+nZiJwoQ3o/9evXDwcHB2xsbLS3Pn364OTkpHNMpK/Vq1cTFBRE0aJF6dWrl6HDEekpwCshoYGErwEJlXwjIyNxd3fXJjRHjx6VhEa8typVqrB7926OHTtGjRo1CA8PZ9y4cVSsWJGdO3e+s5gfJOy0mrz7mt7mC4nHJu++Rrwm4TsTExPWrVtHvXr1ePnyJW3atOHhw4fp+6SE0CPd6tRkB9lxpCY+Pp5y5crh7+/P3LlzpdieMUkcpQm8DEo8qEygUBVi+h2gc5cu/PXXX1hbW3P48GGqV69u6GiFkdBoNKxZs4Zvv/2WoKAgAJo3b868efOoVKmS3mtOBoTQa9mpFB97/eA61H2jmnFISAj169fHz8+PSpUqcezYMfLly5cuz0PkLBnapVtknm3btuHv74+trS2DBg0ydDgiPSWO0ij/VmFV4uHxBeYMa8tff/2FlZUVe/fulYRGpCu1Wk2/fv24efMm48ePx8LCAi8vL1xdXRk+fLjeBpVPw6NS9dhvn1egQAH279+Po6Mjvr6+0k5BZLg0jdQMGDAgVeetWLHivQPKSNltpEZRFGrWrImPjw8TJ05k8uTJhg5JpJe3R2n+pVFUnHscS4NVsezevZvWrVsbMEiREyR2At+6dSuQ0KR4zpw59O3bV1u8731HahJdvHiRRo0aER4eTvfu3Vm/fj1qtf7P1PEahTN3nvM0PAoHa0tqlbTNlLYPImvLkIrCarWa4sWL4+bm9s452O3bt6ct2kySZZOagMOwbxy0nQXOTbWHvby8aNGiBVZWVty/fx87u4zrvCsymb8nrP0o2buPFPmSJoOmZmJAIqc7cuQIn3/+OVevXgWgcePGLFq0iPLlyxOvUWgw6xBBoVF619WoAEcbS46Pa5ZsAuLl5UXbtm2JjY1l9OjR/PLLL0nO2e8byOTd1wgM/W/Ep5CNJZPcK2gXIYucKUOmn4YNG0ZoaCh37tyhadOmLF++nO3btye5iTRQFPCaDM/8Er6+kSzOnDkTgEGDBklCY0wUBQ79j+T+/DQKNNEcky7ZIlM1adKECxcuMGvWLKysrDh69Ciurq58//33xERHMck9oVHq2ylL4veT3Cu8c0SlefPm/PHHHwDMnTs3SVKTuLvqzYQGICg0imFrz+tsGxciOWlKan777TcCAwMZO3Ysu3fvpmjRonTv3p0DBw6kuHpeJCOZ3S/Xrl3D09MTtVrNmDFjDBigSHfxMRD6CNDovVutAsIeJZwnRCYyMzNj7NixXLt2jfbt2xMbG8u0adOoVKkSqkeXWdSnGo42urXKHG0sdbZzv0vv3r2ZPXs2AF9//bX2Q3Bad1cJkZwP2v107949Vq5cyerVq4mLi+Pq1avkyZMnPeNLV1lu+imZ3S8MPswXX37JwoUL6dy5c5at+yM+wL8NJePi4/ls6FDOX7hAiRIlWLVqFdZ58khDSWFwiqKwfft2vvzySx49egRAz549mffrfG6Hq997zYuiKIwYMYJFixZpR4TibEt90JodYfwyZfeTWq1GpVKhKArx8fEpXyB0JbP7JdJ3L6tWrQJg+PDhBgxQZBibIuBUlTl/HmT5Ph/8I3IzY+VfWJdtAE5VJaERBqdSqejSpQvXr19n9OjRqNVqNmzYQOVKFXl88QgdqxamrnOBNC/iValUzJ8/nzZt2mjrMV27k7oaNqndhSVyrjQnNdHR0axfv56WLVtStmxZrly5wsKFC7l//36WHqXJchLXVahMdI+rTAjfNY7w8HDKli1Ls2bNDBOfyHCXL19m4sSJACxYsAAXFxcDRyREUtbW1vzyyy+cPn2aSpUqERwcTNeuXenRowfBwcHv9ZimpqZs3LiRKlWq8OTJE+ZMnZiq6xysLVM+SeRoaUpqhg8fTqFChZg5cyYdOnTgwYMHbN68mXbt2iW7PU8k4+1RmkRKPA6xD2nlbMKwYcPkv6uRiomJoW/fvsTGxuLh4UHfvn0NHZIQ71SjRg3OnTvHd999h4mJCZs2baJixYps2bLlvR4vb9687Nmzh0KFCnHj2B5MY8KTLEJOpCJhF1StkrbvHb/IGdK8pbtYsWK4ublp6xfok1XXgGSZNTXafj+X0LdYNF6jcOGJgvOMe+S3lT9iY/T9998zbdo07Ozs8PX1pWDBgoYOSYhU8/HxoX///vj6+gLQrVs3fvvtN+zt7d/rsRo1aoRSxBWHThMSljS8cb++3lIi58mQNTV9+/aladOm5MuXT6fX09s3kYIUdr+YqFWUcbAif97ceu8X2dvp06eZMWMGAIsXL5aERmQ71atX59y5c3z//feYmJiwefNmKlSowObNm9/rsdavX0/UrVM83TGd3OpYnfvTsrtKCOn9ZCj/7n550/MXL7TFqVZu3kuVBm0MFJzIKHFxcbi5ueHr60vv3r1Zu3atoUMS4oP4+Pjw6aefcuXKFSDhw+/ChQuxtrZO0+PMmzeP0aNHg0rNzOVbKOdWSyoKCy3p/ZTV/bv75c3bH/vPc+ZBNKZFq0tCY6RWrFiBr68vtra2zJ8/39DhiCwqMjKSp0+f8ujRI9auXcv169cJDAwkJCSEV69eZandpomjNt999x1qtZrVq1fj5ubG6dOn0/Q4I0eOZMSIEaBomDyiN8WUp++1u0rkbDJSk4W4ublx8eJFFi1axGeffWbocEQ6Cw8Pp3Tp0jx9+pRff/2VL7/80tAhCQNQFIV79+5x4cIFAgICuH37Nrdv3+bx48c8f/6c58+fExkZmeLj2Nra4uDggIODA/b29hQsWJBSpUpRunRpSpcujbOzM5aW779b6H16MB07dow+ffpw//59TExMmDx5Mt9++y0mJibvvC5RXFwcHTp04MCBAxQpUoSzZ8/i6Oj43s9BGI8M6f2U3WXlpMbX15fKlStjZmZGYGAgBQpIgSljk7g4uEyZMvj6+mJubm7okEQmiIyMxNvbm3/++YezZ89y9uxZnj17lvKFJGx9jouL027MSMvLtUqlokiRIlStWpXq1atTrVo1qlevjlOkn95ec2/6kB5ML1++5LPPPmPjxo0ANGrUiLVr11K0aNFUxf3y5Uvq1KmDn58fderU4fDhwx+UnAnjIEmNHlk5qRk3bhyzZ8+mU6dO0j/LCD148ICyZcsSFRXF9u3b6dSpk6FDEhkoJCSELVu2sGPHDo4cOUJUlG7RODMzMypXroyLiwulSpWiZMmSFC1alAIFCmBra0uBAgXIkydPkpIOcXFxREdHExERQXBwME+fPtV+DQwMJCAgAH9/f27dukVYWJje2M4Ps8HNQeExjgQ0W0at2rWxsLDQ3p/Yg+ntN4a07EKKi9fwvyUbWPD7KiKeBWL16iHLli6la9eu77wucXTo8q27fP/VFzy7fppP+vRm1apV79xxK4xfhiQ1EydOpGPHjlSvXj1dgsxsWTWp0Wg0FCtWjEePHrF161a6dOli6JBEWiTTZf1Nn3zyCWvXrqVRo0YcOXJEXqCNUHx8PPv27WPx4sUcOHCAuLg47X1FihTh83aVGFzkFk+rf02JZv0zdPRBURRCQkK4fv06Fy5cwMfHBx8fH4pG32Rfbyvtea3XRuD9yIw6derQqlUrWrZqzecHXhAUpr9yb2q6cesb5YkLC+a511J6NSjP/PnzyZ076c5OvdeFP+O55xJ+HNiJsWPHvud/DWEMMiSpGTBgAHv27MHc3Bx3d3c8PDxo3rx5thlGz6pJzaFDh2jevDn58uUjKChI51OTyOK0NYcugJMbDD4MbyUs169fp0KFhA7HZ8+epUaNGoaIVGSQ169fs2TJEhYuXMjt27e1x93c3OjRowcdOnSgQvnyqH5v9s7fkwynKMQvaYw66AoqNMQrcDlYRbVFodpTLIpWxvHjGSk+VHI9mJIb5QEFRYHgHdMpZRbGli1bKFeuXKqve7ZjBut/GoeHh0dqn60wMhmy+2nFihUEBQWxfv16rK2tGTVqFHZ2dnz00UesXr2a58+ff3DgOdGaNWsA6N69uyQ02U0yXdbfNG/ePAA6duwoCY0RiYqKYv78+ZQqVYoxY8Zw+/Zt8ufPz9dff83169c5f/4848aNo2LFiqhuH0rx9yTDBXhhEnQJ1b/1sUxU4OagcP/QHyxcuBB3d3dyFUjdolx9PZje1WkbVKhUKuxbDePqtevUrFmTTZs2pe46IH/zwfTu00e7bVyI5KR5S7daraZhw4bMnj0bPz8/Tp8+Te3atVmyZAlOTk40atSIOXPmaLu6ind7/fq1tsz4J598YuBoRJq83b9LZZLw/RuDn8+ePWP16tUAjBkzxhBRigxw8OBBKlWqxMiRI3ny5AklS5Zk6dKlPHz4kJ9++klnFCI1vycZ7h295oreWsmI4cPZtWsXO9avStXD5bNIOsp05s5znakjfVS5bandoTevXr2iR48ejBw5kuN+Qe++TqXCNK89sflK4OHh8d79pkTO8MF1asqXL8/YsWM5ceIEDx48oF+/fhw7doz169enR3xGb8+ePbx69YoSJUpQr149Q4cj0iKZLutvfgpfvHgxUVFRVK9enYYNGxooUJFegoOD6d27N61ataKE5i5+X9qwZ/5X+Pn5MXjwYHLlypX0olT8nmS4d/SaezOW+mUdKWRjmWwPJkXREBcWjHud8vTv3x8vLy80moSRn9R20B43aRrjx48HYNn+s/RfdjxV1zk5l+Pu3bt07dqV2NjYlC8QOVK6Ft+zt7dn4MCB7Ny5k6+//jo9H9poJe506t69uzSvzE7e8ck38VN4dHQ0CxcuBBJGaWRxcPZ27NgxXF1dWbduHWq1mhUfF6NsfoX2FucwMzXVf1Eqfk8yXGIMyb7cq7WxmKhVTHJPWP+l77dVpVJhenk74WFhrFq1ihYtWlCmTBlmzJiBaezrVIXjmC8X06dPZ+ofu7HvNAGNmVXKFwHTfxiHtbU13t7efPXVV6m6RuQ88i5qQDExMfz1118AssU3u0nFJ98NGzbw5MkTChcuTLdu3QwTp/hgiqIwZ84cmjZtSmBgIBUqVOD67oUUMw1JOOFdoy6pHCHJUCn0mgMNhD1KOA9oU6kQi/pUw9FGd3dWIRtLFvepzm3v7Rw7doyhQ4eSN29ebt++zYQJE+hUvxJmsa+SDePNTtvxGoW9gVaoVKoUk/3E67o1qaZdf7hgwQJWrlyZuucvchSpU2NAf//9N61bt6ZgwYI8fvxYRmqyixS6rIManFypviSM8+cvMHPmTMaNG5fZUYp0oNFoGDVqFAsWLACgd+/eLFm8mNzrOkDg5YTkRGUChaok3dGUyt+TTNkJpafXnI7c9mBTWOdQaioKR0REsHnzZpYtW8Y///yDVdm62HeakDDK88ZzervGzcmAEHotO5Wq0FX/XteygiNn7jzn//5Yy+bVyyHYn2PeR6lZs2aqHkdkb6l9/05mzFRkhh07dgDg4eEhCU12kopPvnHP7+N76Q6mpqYMGjQoM6MT6SHgMMq+sfxy1YEFi/9CpVLx66+/8vnnn6N6c8cb6I66lG7x3/G0jJCYZvCuR5siCbc0MFGr9G7bflPu3Lnp378//fv35+rVqyxdupS1B+ZiVe8TTPPaa8+zz2PGlE6VtUX7Urv+Jl8uM2Z2qQxAg1mH/l1QXBrHj2cQFxZMl5FTOLf9d+l0L7Q+KKlJ3OFUuHDhFM4Ub9NoNOzatQuQqadsx9QChhx+5yffJas2ERM/g1atmknLi+xGUVC8JqN6dpNGcddQq9WsWrWKPn366K6ReXNKKXGNjHPz/0YoUvF7Qm77jE9oMknFihX59ddfmRoWxvIVfzB//TqehkcT/+oFj4L8WHujG0XHj6dixYo4WKeu8OAAl4SJBH01bEyt7VAaDqXDZ9/xz6ZFmJmZpfMzEtnRe00/nThxQtu0DMDOzo7+/fvz3XffZYlpneRkpemns2fPUqtWLXLnzs2zZ8+kt4mRqVGjBj4+PixZsoQhQ4YYOhyRFv6esPYj7bfHS4yhQf9Jeu9Los9W3dGaHCw+Pp7du3czb948jh49qj3euXNnvp0wgZEHQwkKjdJfn0ZRiAt/RtDvn1Hu64280uj//K0oGjSRr2igvsHaX36Ujt5GLEOK7yUaOnQo5cuX5+zZs/j5+fHTTz/h6elJtWrVpD5NKu3cuROANm3aSEKT3QQchoW1Er7qcffuXXx8fFCr1TIKl90oCs+3fEWcJuGtVoOKBjGHE0Zo0rCLSICJiQmdOnXiyJEj+Pj40LVrV1QqFdu3b6d2zZqYXNqOQtJdVioSdllVib+JqWPZZBMaAJVKjUmuvJy0rIXbxN3s9w3MyKcksoH3SmoCAgKYN28e1apVo3Tp0vTt25dz587h5ubGqFGj0jlE45S4nkbe9LIZRQGvyfDML+GrnjewrVu3AgndiR0cHDI7QvEB7nqtwDbqLqb/fuJXo/y3XiaNu4jEf6pVq8bmzZvx9fWlT58+qNVqTmz8jafbp6GODtc519HGkkV9qrF38TT6Df081T8jNFbFZ2vPS2KTw73X9FO1atWYP38+DRo00Dl+/fp1atWqRXh4eDJXGlZWmX568OABxYoVQ61WExwcjK2trcFiEWn09vSDnumGevXqcfLkSRYsWMDnn6f+RVkYliY+nhtfO1HWOlKb1AC6u5vCHqV5F5FIKiAggFmzZrFy5Upi4+KxKFKRmo2aM6z/x/RoVkM7jZSWXVKQsP3eMa8F/4xvIVNRRiZDp5/69+/PF198wYMHD3SOGzpZyC4OHToEJKy7kIQmG0lFufuHDx9y8uRJAOm2ns3s/20sFWyidBMa0N3dZFMEnKomf5OEJlWcnZ1ZunQp/v7+DBo4gNhHVzn+5zz6tKpNv76faBuD1ippSyGb1E/Pq1QqnoTHcCpAWinkVO+V1IwaNYpLly5RpkwZPv74Y2bPns2MGTMYOHAgs2fPTu8YjU5iUtOsWTMDRyLSJBXl7g8cOABA3bp1cXJyMkSU4j08ffIEx2u/E69JbuBa1stkhGLFirFs2TKuXbtGt27dUBSFP//8ExcXF0aMGEHIs2AmuVdIWGeThsddvEra9ORU75XUBAYGsm/fPqZMmQLAypUrmThxIrdu3WL27Nn06dOH2bNns3///nQN1hgoiiJJTXaUynL3586dA0gyNSuytrk/z8Ypj+YdUxayXiYjubi4sGnTJnx8fGjdujVxcXH83//9H2XKlOHqgXXM71klSYXjd9m8erm2WjskFBI8GRDCzouPOBkQ8o7kVWR36VZROCoqiitXrnDx4kUuXbrExYsX8fX15eXLl+nx8OkiK6yp8ff3p0yZMpiZmfHy5Uv9DfBE1pPKrbw1a9bk3LlzbNy4ke7du2defCL1Ag7DvnHQdhY4NyUiIoIiRYqQRxPKn0t/pVFyjUdlvUymOXr0KGPGjOH8+fMAlClThjk//0KBcrUZse4CLyNj0Dd2owIsNJH4zelB/nw2+Pj44BdhyeTd13Q6geezMuPT+iX4vFkZWXuTTaT2/VvaJGSUt144Ey1dupShQ4fSqFEjndoNIgtLZbn76L77sc6bl9jYWAICAihVqlRmRypSov1/eQGc3GDwYZYuW8bQoUNxdnbm5s2bUt07i4iPj2fVqlWMHz+ep0+fAtCqVSu6jvof048+RXlrQ3jiv+b3qMK0Yd04ffo0Fdv0JsK1l/5aOPxXsTix0rHIujJ0obBIwTu2/crUUzaUyq28Vy9fIDY2lvz581OyZMnMjFCk1pstDh5fQPH3Yv78+QCMGDFCEposxMTEhAEDBnDr1i3Gjh2Lubk5f//9N8Pc61I96gIFrd+qxBz5kimti+HuVpTNmzdjZ2/Py5IteNfn9pevYxkm28CNiozUZIRktv0qikLBggUJDg7G29ubhskNc4usJxUNAZdu3MvQoUNp2bIlf//9d+bFJlIncZTmjUaUr/I6Yz3mHLlz5+bhw4fky5fP0FGKZPj7+/P1119rC5cWdHTkiynzsXYozOzJ3/HoojdOhRz566+/cHV1ZeGmv5lzPjZVj13IxpLj45rJVFQWJiM1hvKObb9Xr14lODgYKysrateubdg4RdqkYitv4iLhGjVqGCpK8S56dq/lCb1JK2cTWrZsKQlNFle6dGl27NjBwYMHcXFx4UlQEN8P6c6+xVPZMH8qFcqX4/HjxzRs2BBPT0+Klq2Y6scODI3izJ3nGRi9yCyS1KS3d2z7TaxfUrduXczNzQ0YpMgIktRkYcnsXotXYGpTCxo3amSgwERatWjRgkuXLjF58mTMzc3Zv38/rVq1olu3bjRs2JDw8HDatm3LpVPH0vS4+3wDZWeUEZCkJj2lsO3XxyfhTa9mzZoGCE5kJEVR8PX1BcDNzc3A0Ygk3v6w8S8TFdQqbIp7hdwGCky8DwsLCyZOnMjly5dp2rQpkZGRTJ48mRcvXtCiRQvi4uKYNLw3eUziUv2Yq0/eo9eyUzSYdUjW2GRjktSkp2ReOBNHa0zveQPySd4YhYaGEhubMH8vRfeyGG0jSv3iFSh1d50U1suGXFxc8PLyYtWqVRQoUABfX1+8vLxwdXUFRcPtLbPS/P81KDRKFg9nY1kmqfH29sbd3R0nJydUKpW24WMiRVGYOHEihQoVwsrKihYtWnDr1i3DBKtPCh18FVT0K/oQgOrVq2diYCIzhISEAJA7d24sLCxSOFtkqvgYeH472btNVKAKl8J62ZVKpaJv377cuHGDTz/9FEVRuHTpErly5SLy5kme7piOEpX6foSJKdDk3ddkKiobyjJJTUREBK6urvz222967589ezbz589n8eLFnD59mty5c9O6dWuioqL0np/pUtj2q0KhiDU42uWnRIkSmRqayHjPniXsjCpQoICBIxFJmJgnLPTWvtypwc6FP3MPodqSV0y41wAGHwFTSUazMzs7O1asWMGBAwcoUqQIr1+/BiDy5knuz+9N5NnNmMSl7v1CQRYPZ1emhg4gUdu2bWnbtq3e+xRFYd68eXz//fd07NgRgNWrV1OwYEF27NhBz5499V4XHR1NdHS09vuwsLD0DzyRqQUMOZzstt+t27Yy8pdpVKnbEpVKtg1mK8kUUnxT4kiNnZ1dZkYmUiPAC55cfeOABp75ERNXlgtBGhpYO0ulYCPSqlUrrly5wqhRo1i1alXCQUXD00OrCDu5GbvytYmwLUve6u4pPtbT8CzyoVmkWpYZqXmXO3fuEBQURIsWLbTHbGxsqF27tnZHkT4zZszAxsZGeytatGjGBvqObb/7LwXxKFyRqafs5h2FFN+UmNTISE0W847F+001JwBJRI1Rvnz5WLlyJTt37qRgwYLa41GRrwm85I0SmboPuA7Wqe83JbKGbJHUBAUFAej8ciZ+n3ifPuPHjyc0NFR7e/DgQYbG+S4+Pj6ALBLOdt6qQPtmR+43yfRTFvWOxfslzJ/TytlEkhoj5uHhga+vLz169NAei9co5K7S+p2VhgEc81pQq6RtRoco0lm2SGrel4WFBXnz5tW5GUJis0+QRcLZyjsKKb5NRmqyoBQW72v+rVETH5f6bb8i+7Gzs2PDhg1s3LgRS0tLLIpUxDSvXYrLAHrVKiYVhrOhbJHUODo6AvDkyROd40+ePNHel5X5+/sTFxeHjY0NxYoVM3Q4IrXeUUjxbTJSkwWlsHhfrYKiedU8CXyYuXEJg+jevTv+/v44OZdL1fkl7KR2UXaULZKakiVL4ujoiJfXf28mYWFhnD59mrp16xowstRJ3HpepkwZWSScXaRQSPHt0Rr5/5oFJS7eH3JU722lxQBqLovg/qPkp7CFcSlcuDBrli5M1bmyniZ7yjK7n169eoW/v7/2+zt37nDx4kVsbW0pVqwYo0aN4n//+x9lypShZMmS/PDDDzg5OdGpUyfDBZ1Kic+rTJkyBo5EpNqba2ne9OZoTen/Fq4njtAkTkOJLMKmyL/buZMyLerLo3CFx48fZ3JQwpDqONtTyMaSwNBIIOmHERXgaGMp62myqSwzUnPu3Dnc3Ny0JebHjBmDm5sbEydOBGDs2LF88cUXDBkyhJo1a/Lq1Sv279+PpWXWz6bfHKkR2UAKazFAnWS0JnGxaeI0lMj6Eis/P3r0yMCRiMxkolYxyb0CKlT8V2ovgaJoUIBJ7hVkPU02lWWSmiZNmqAoSpLbypUrgYTh/SlTphAUFERUVBSenp6ULVvWsEGnUmJSU7p0aQNHIlIlhbUYoIEw3Qq0ktRkP4lJjYzU5DxtKhViUZ9qFLKx0jkeHx5C8PbpWIXcNFBk4kNlmeknYybTT9lMCoUUAchtr1OBNjGpkemn7KNIkSKoVCrCwsJ4+PAhRYron6YSxqlNpUK0rODImTvPOeFzhR+++ZKoB1dB0dCkSRPmzJnDmDFjZL1cNpNlRmqM1evXr3n4MGF3hSQ12cg7CiniVDVJBdrENTUyUpN95MmThzp16gDw119/GTgaYQgmahV1nQvwdfcmeK5bhOqN6aivv/6arl27ptiKJ16jcDIghJ0XH3EyIET6RRmYjNRksICAACChwqWtrSw8M1Yy/ZRNvNXyokOHDpw8eZI9e/YwZMgQQ0cnDKh+/frs37+f1q1ba49t27YNOzs7JkyYQKVKlShWrBjly5fXNq3d7xvI5N3XCAz9L/FxzGtJr1rFKGGXCwfrhAXHsj4n86iUlMoqGpGwsDBsbGwIDQ3NtEJ827dvp0uXLtSsWZMzZ85kys8UmS88PFz7OxUREUGuXLkMHJFIQlFgWdOE3WtObjD4MJevXMHV1RUrKytCQkKwsrJK+XGE0Xn69Clz585l8eLFvHz58p3nmpqaUqlSJUo06MSF3NWpr/blR9PV/BjXlxOayknOL2RjyST3CrSpVIh4jcKZO895Gh4lCU8apfb9W0ZqMlhiawbpzG3c8uTJg7m5OTExMTx79kyKLGZFelpeVK7cnKJFi/LgwQMOHTpE+/btDRujyHSnTp3Cw8OD4OBgIGFUPV++fNy9excAtVpN0+Iq5re1ZOxh2Hs9gouXLhNcdyQmisJY042UUT9irOlGOsZU4u1t4kGhUQxbe54hjUqy61KgzqiObW4zOlctTIsKjpLgpBNZU5PBEqcj7O3tDRyJyEgqlUqmoLKyZFpeqECbyOzZs8dw8QmDePLkCe7u7gQHB1O5cmV27tzJs2fPuHPnDl93qcHV4blpWlzFnLbWVLA3YWJ9hS5dOvO/JesxzWtPY5MruKpvA+Cqvk0j9eUkP0P597bE+45OQgPwPCKW5Sfu0mvZKRrMOsR+38BMeNbGTZKaDJb4BidN84xf4mLhxE98Igt5R8uLDh06ALBz506io6MNGKTIbAsWLODZs2dUrlyZkydP4uHhgYmJCSgKM1pYUsHehN/aWVLVPqG8Q63Cpry6tIuNu/YDCl+ZbiZOSXgbjVPUfGW6mbdr36RW4oiOJDYfRpKaDCZJTc6RWIfo/PnzBo5E6Eih5UXzZs1wcnIiMDCQZcuWGSZGYRAHDhwAYNy4ceTO/UavpwAvTJ8kjLq42JlodzTFKzCjZW7uBj6jkfoyrurbmKoSEh5TlSbZ0ZrUSEyFJu++RkycRnZUvSdJajKYJDVGIuAwLKyV8DUZTZs2BeDw4eTPEQbw9ihNon9HaywfneD7778HYNq0abx+/doAQQpDCAxMGBUpU6aM7t/4W0lw4loXExVUKwidWtXRGaVJ9KGjNQoQGBpFnRle9Fp2ipEbLtJr2Smq/HiAYWt9OOH/TBKcFEhSk8EkqTECigJek+GZX8LXZDYMJiY1J06cICYmRu85IpOlsuXFwAEDKFGiBEFBQSxcmLqGhyL7c3R0BCAoMBB2DEv4G98yQH8S/K84Rc2Ppqt0RmkSfehoTaLnEbqvHxEx8ezzDaL376ep/r+DMkX1DpLUZDBJaoyAnl0z+lSsWBF7e3tev34t2/ezilS2vDA3gR9//BGAWbNmERoamlkRCgOqUKECAI+Pr4bwfxOFyOfvvMZUpaGUOojkBkw0iuqDRmtS8vJ1LJ/J2ptkSVKTHpKZmlAURZKa7C6ZXTP6RmtUKhVNmjQBZAoqy0hseTHkaPK3wUfA1II+ffpQrlw5nj9/zty5cw0ducgE3bt3B8BD8UzTdYoCye2+VqsUCqlCMCcOSKhTM7RRST39wD/M5N3XZCpKDym+96H0FPTi314hiT8PEtolSGGvbMjfE9Z+lPR4n61QukWSw4sWLWL48OE0bdqUQ4cOZUKAIj1t3ryZ7t27Y21tzZ07d7Q72oRxio2NpV8DJ9a1Tft08XMlD5/GfEPcG+XeXh5dRYUCJnzy3a+8sCyiU2BPX/XhNymKBpUqbeMM6wfXoa5zzvgdleJ7mUXf1MS/b3YvXrwAwNLSUhKa7OjNUZo359cTR2ucm2sT2ESJ62r++ecfoqKisLS0zMyIxQf66KOPqFq1KhcvXmTWrFnMnj3b0CGJDGRmasqyDlYJ05Sp8GYiE6LkJQjdhCLsuQn/nDzNsVv9+OuvvyhV6r8ecW820Dx4LYgdFx/rrJ3RvA7DJHe+NMX/NPzdfalyIpl++hApTE3Exye8EZqYmCT3CCIrS2HXjL61NS4uLjg6OhIdHc2pU6cyKVCRXtRqNf/73/8AWLhwIY8fPzZwRCJD3TpI7vh3rJ9qPQuGHOVi2520j55Gu+gZXFLKcFUpmSShAfj1118pUqQIfn5+1KlTJ8nausQGmhPdK3L2uxasH1yHX3tWZf3g2sysZ4oS8RxFSW79V1IO1vKh6W2S1HyIdxT0AtBoEn451Wr5z5ztpHLXzNtra1QqlWztzubatWtHvXr1iIyMpG/fvtoPJ8LIKArs+vLd55yYB4VcuWeRfCLzJjPrApw+fRo3NzeCg4Np3rx5kmnoxK7eey4nJMwdqjhR19mO3r168Gu/hqhUKlKzKqSQTcLUltAl77bvK4WCXiiK9hdTpZJ+HtlOKnfN6Bu2btasGQBeXvp3SYksLOAwqt9qs+5/g8mdOzdeXl5MmjTJ0FGJjBAbBa9S2EH0KhBio1I9IjJv5hTs7Ozw9vamRYsWvHr1inbt2rFz504goat3g1mHdGrQvNkeoVONEizuU528FimP7k9yryC9ovSQhcLvK7kFpIn6bOWmphguLi7Y2Nik2PlVZEGhDyHiHX2cctuDTeEkh+/fv0+JEiVQFIVr165Rvnz5DAxSpJu3Fv2vzzOEj3v3BhL6QkmzSyMTFQYzi/Hurdcq+PY+8ebWNJh1iKDQqGTOVogPD+HhogH0/aQPK1euJCYmhl69erF9+3ZMTEz4+te1bHxgneT6xLRkUZ9qtKlUCEgYzVngdZPfDvkRq+gmOPlzmTGjS2XtuTlFat+/ZaTmfaRyakLz77C1TD9lUzZFwKlq8jc9CQ1AsWLF8PDwABLm2EU28dai/161HBgxYgQAffr04c6dOwYMTqQ7y7ww6CB0Waa9xXdawrwHlfh462s+2RGFZ8nxYJkXE7WKSe4JNW3eHhtRFA2KAn0qWGKiVrF69Wr+97//YWFhwaZNm+jXrx/xGoW1V6NQ9KREb7ZHSNyibaJWMaqlCzemuTOvYynsnp7j5YkNBK2fQL4js6iQNzYD/8Nkb/Ju+z7SODUh0085z+jRowFYvXo1ISEhBo5GpCiZRf8/z5lDrVq1ePnyJV27diUqSnabGJUiNaFKd+3NpGpPRiw+gqpKd9ZeiqH1gAmsWrUKSNi9tKhPNRxtdKeiTKLDCd4xnXM7fmfBggUATJw4kXXr1mFqasqKFSv4eOQPmOa1J2lKlCCxPcKZO7qF/0zUKjrVLc/Z5ROZ/UlDTJ75c+TwIapUqcKff/6ZqrU3OY1MP72vVExNXHsUSsWKFSlQoIC2CJ/IGRRFoXr16ly4cIFp06YxYcIEQ4ck3uUd9Yjum5elWrVqhISEMGTIEJYsWZL58YlMFR8fz2effcbvv/8OwG+//cbw4cMT7tMonLnznKfhCWtt8sU+o3o1NyIjI1mwYAH37t1jzpw5mJubc9DTEzOnCuy7EsjqU/dS/Lm/9qxKx6r6R4AB/P39+eSTT7Q7K7t3787ixYvJnz9/OjzrrE2mnzJaKqYmEnc/yUiNEUummrRKpdKO1ixcuFB6QWVlKSz6L1a0KH/++ScqlYqlS5eyevVqw8QpMo2JiQlLly5l5MiRAIwYMYLhw4cTHR2t3ZbdsWph6joXoHw5F3766ScAxo4dy4ABA+jSpQsmJarTZ+Mdei07laqEBlLeol26dGmOHTvGlClTMDExYdOmTbi5uUlbljdIUpOBZPeTkUuh0WWPHj0oVKgQgYGBbNq0yUBBihSloh5R69atmThxIgBDhw6VatE5gEqlYu7cuUydOhWVSsWiRYto3LgxDx8+THLu8OHDadWqFZGRkfTr14/e437CodMEVLlTt+VaReq3aJuamvLDDz9w8uRJnJ2duXfvHg0aNODXX3+V6SgkqclQuXLlAiAiIsLAkYgMkUKjS3Nzc+1C07lz58oLTlaUhnpEP/zwAx06dCAqKooOHTrg6Zm2fkEi+1GpVHz//ffs2bOHfPnycfr0aapVq5akBpVKpWLFihXky5ePs+d8+GH7JVCpUveB9t/XhbRu0a5ZsyY+Pj589NFHxMbGMmrUKLp27Zrjd9pKUpOBEptYvn79mtevXxs4GpGuUtnocujQoVhZWXH+/HmOHTtmgEDFO6Vh0b+JiQmbN2+mffv2REZG4u7uzoEDBzIzWmEg7dq1w8fHh6pVqxIcHEyLFi2YPn26TmHGwoULs2DBAiyKVCRCMU/1Y8eFP6PM06O0KGef5rhsbGzYvHkz8+fPx8zMjG3btmnX8uVUslA4AymKgoWFBbGxsdy7d49ixYpl+M8UmSQNjS4/++wzlixZQseOHdmxY0fmxCdSL431iKKjo+nWrRu7d+/GwsKC7du307Zt20wIVBhaZGQkw4YN0+6IqlOnDitXrsTFxQVIeM1v0vdr7hVuluJj9a1bnLyht5kwuDuxMdH06NGDtWvXYmr6fi0Zz549S/fu3bl79y5WVlb8/vvvfPzxxzrnvL3IObHZZnaQ2vdvSWoyWOHChXn8+DHnzp2jevXqmfIzRQZLLNIWeDlpo8tCVXQ6tQPcuHGD8uXLo1KpOHHiBHXr1jVA0CI9xcTE0KNHD3bs2IG5uTlbt26lQ4cOhg5LZAJFUVi1ahUjR44kLCwMS0tL1kweyEe5TqNqO5u9z4swYuvNFB8nscP2zp076datG7GxsR+c2Lx48YLevXuzb98+AL766itmzpyJqamp3i7hhWwsmeReIVsU8pPdT1lE4hSUbOk2ImlsdFmuXDn69u2LoigMGjSI6OjoTAxWZARzc3M2bdrERx99RExMDF26dNGWwhfGTaVS0b9/f3x9fWndujVRUVEUu/UHqmc3ido7gTbVnLEx0yTbmFJRNNjlMtEuCu7YsSNbtmzBzMyMjRs30rt3b+Li4oD/+kTtvPiIkwEh2uJ8ycmfPz+7d+9m/PjxAPz888+0adOGzSdvMWzteZ2EBiAoNIpha89r2zQYA0lqMpi9fcI8aXBwsIEjEeniPRtd/vLLLzg4OHDt2jWmT5+e4WGKjGd2/zibG99jav9mxMbG0rVrV7Zt22bosEQmKVq0KPv27eOvBV9Rq3DCyIrl82us+bEfkztWRoUqaWKjKICKZ38v4nnIfx90PTw82Lp1K2ZmZmzatInevXuz99Kjd/aJSo6JiQnTp09n06ZN5MqVC69Dh/lm/Wm97R30VTPO7iSpyWCJIzWS1BiJ92x0WaBAARYuXAjA9OnTuXz5csbGKTLWv9v5VSE3+a5OPL169SQuLo5u3boxa9Ys2emWQ6iAtubnUFQJb6VxGoUKQdsY+VEjuhcOQ4l4oXN+wbwWmJ76g3sndtG9e3diY/9rd+Du7q5NbHZffMCIdRc+aGSlW7dunDp1iuI1mkGu5IvzJVfNOLt6v4k7kWqJIzUy/WQkTC1gyOGUF5aaWiQ53LVrVzp16sSOHTsYOHAgJ0+efO+5c2Fgb2znVwVeYM2PE8iTx5ply5bx7bffcvbsWf744w+sra0NHKjIUP/+HiSuoDNVq6hV2JRyZo/5aWRvqrhWxe95PGZ5C7Bg9jQ+aVMXv86FqV37b44cOcKYMWO0rRUgIbHZvGUrIw48J2FMR1fiscm7r9GygmOKi3wrV67MlJ9+ZeK+lPuWPQ03jhYgMlKTwWT6yQi9Z6NLlUrFb7/9ho2NDefOnZNml9mVnu38Jkems2TxYpYsWYKZmRlbt26ldu3a+Pn5GTRUkYGSqUStqExY3rsk5ubmXL50kZiHvry6eoQZoz8lLjaGChUq8OeffwIJ1cYTWzEkcqhYD1Nru2Rr3KR1ZKVMkYKpOi+lasbZhSQ1GSwxqQkKCjJwJCIrcHJy4ueffwbghx9+wN/f38ARiTR7e6H4vwvEVbcPMWTIELy9vXFycuL69evUqlWLXbt2GTZekTGS2TCgUuIponrKnYPLcHd3105F3rx5kxYtWhAeHo6HhwdTp04FEqoR//PPP9rrUztiktrzapW0pZCNZTKtNP+rZly9eP40LUrOqiSpyWCJ9QuuXr1q4EhEVjFgwACaN29OZGQkgwcP1vYIE9lACn2iUBTq1KnD+fPnadiwIWFhYXTs2JGJEyfqFGoT2VwqNgw43VjBrp07OXr0KGXKlAHg+PHjlChRgh2/jGKCzQ4m92tKbGwsXbp00bZfSO2ISWrPM1GrmOReAdDTI/zfhMvDtRCNfzqc5kXJWZEkNRnM1dUVgNu3bxMWFmbgaERWkNgYMVeuXBw5coRly5YZOiSRWqnczl+wYEG8vLy0DRGnTp2Ku7s7L168ePsRRXaUhg0DjRo1ws/PT1uf6vnz5zhdX4Y65BZDSgVSqVJFnjx5QqdOnYiMjEz1yEpq+kQlalOpEIv6VMPRRjcRigt/htltb5Z43zGa7d5SfC8TFClShEePHnH8+HHq16+faT9XZG3z5s1j9OjRWFpa4uXlRb169QwdkniXxKKLjy+h/81MDU6uSYovrl27liFDhhAZGUnx4sVZtmwZLVu2zLSwRQZJYyXqp0+fUr58eWrkD+VAn9za4732WrD7+msiIiL4dMAABk+Yjee1IJafuJvkIRVFg0qlYnGf6u9VMO/NisKvgh/xdf/OmHSchqm1nc7vbCIV4GhjyfFxzQxeeViK72UhiaM1ly5dMnAkwqACDsPCWglfgS+++ELbINHd3V0WlWZ177mdv0+fPvzzzz+UKlWKe/fu0apVKwYNGpTjGw9me2ncMODg4MDCBQuY2tSCuH/Xq8RpYHTVCCIiIrAqW5e/LRrSa9kpbULzdh4RHx7Ci92zMX96/b1CNlGrqOtcgI5VC9O7ZS0WbtiHaV57vQkNZM/t3pLUZIIqVaoASG2SnOzfuiY880v4qiiYmJiwYcMGatasyfPnz2nbti1PnjwxdKQiOYnb+YccTf42+Ije7fxVq1bl0qVLfPHFF6hUKpYvX07FihXZvXt35j8PYTA9azpQq7Appv9mK6ZqqFXYlN4d6mPfaQIm1nY65yfOowyoX4I/B9akTvBfhF07RqdOndKlaaWpdYFUnZedtntLUpMJZKRGvFnX5M21F7lz52bPnj04Oztz584d2rdvz6tXrwwYqHin99zOD5AnTx7mz5+Pt7c3ZcuW5fHjx3h4eNC7d2+pY5UTKAqqw//TFurTHlaZMKHGa1QqJck27sS1Ift8g6jjbM/aNatp0qQJ4eHhtG3bltu3b39QSOm9KDkrkKQmEyQmNVeuXJGdLjmRnromb7ZScHBwYP/+/djZ2eHj40OPHj20vV+E8WnQoAEXL15k7NixqNVq1q1bR4UKFdi8ebNUIjZmiYX63mqdoFLiqaAE0Eh9JdlLA0Oj+GnlNkxMTNixYweurq48efKE1q1b8/Tp0zT3iEqU2kXJGo2SbbZ6y0Lh9BZwGPaNg7azwLkpAHFxcVhbWxMVFcXNmze12/tEDuHvCWs/Snq8z1Yo3UL77enTp2natKl2q/eSJUuSLcAljMPZs2f59NNPcXx9g/ltLdkUXoOhs9ZTqFDW75os0iCFReYaRcUVpSQdY6aiZ+M1AMG7ZpM/LIDPPvuMjh074uHhgbPqAUu6FmSu9Rj2RJTTnpuW7tv7fQMZtvZ8Qpg6IScsSs6Xy5yXr/9r52Cozt6yUNgQ9KybADA1NaVixYqATEHlOKmoa5Kodu3abNiwAbVazbJly5g2bVomBysyW82aNfE5d45Vn5Sigr0J7SzOUqZMGaZMmUJERIShwxPpJYVF5mqVQiFVCOYkP0Kb10zh0aNH/PDDD9SoUYMqVSoz290e51yvGBy7ljdTkrRsx05uu7cmMhxFgZevdRe+Z/Wt3pLUpKdk1k0A1KpVCwAvLy99Vwpjlcq6Jok8PDy0jS9/+OEHVq5cmUmBCkOxeHicwqqEBeK1CptS3zGKSZMmUbZsWVasWCFF+4zBG4vM4wcfoZ/ZT7SPnqZz84j+HzGYJblUUTQ45DEj4NTfrFmzhlq1ahETE0PUtf1Uy/8aAFf1bRqp/9uIktbu220qFeL4uGasH1yHX3tWZc2AmlhZJMaif51PVu3sLUlNeklh3YS7uzsAu3btknU1OUUqqo6+PVoDMGzYML799lsABg8ezKpVqzI2TmE4b71uKCoT1g2qQMmSJXj8+DEDBw6kWrVqHDx40MCBig/27yLzM1HFOBpemKtKSZ1bEHp2IikJLSxjT6/HzNSEPn36cPr0aU6fOsWCbkWJU/7tDq6o+cp0M2+O1qR1O/ab271NTUyIMcmVbv2nMpMkNeklmX4wiZ/EmzVrRp48eXj8+DE+Pj4GDFRkmvesawIwbdo0+vfvT1xcHP379+fHH3+URaTG6K3XDZUST4Goe9zY+3/MmTOHfPnycfnyZVq1akXbtm25ciX5xaQie0jL9mj7PGa8/nseF/esZOLEidrjtQqEU9b8GaaqhNcWU5UmyWjN+/y8tF6TFbd6S1KTHlKxbsLCwoI2bdoAsHPnTgMEKTLdB9Q1UavVLF++nPHjxwMwefJk+vXrR0xM0gRIZFPveN0wPzaLr8aMwd/fn1GjRmFmZsb+/fupWrUqgwYNIjAwa65nEClL7fboH9qX59R3rVj8/TAAZs6cyYEDB7S/N29vDdc3WgPw68wp7NixI02vHdl5q7ckNekhlesmOnbsCEhSk6N8QF0TtVrN9OnTWbp0KSYmJqxZs4bWrVtL/yBjkYrXjQIFCjB37lyuX79O165d0Wg0LF++nNKlSzNmzBgePXpkmNjFe0tpGzWKQiEbS/rXL4mJWsVHH33E8OHDAfjkk08IObtZ79bwJKM1ikJcWDCHN/1O586dcXR0ZOjQoXh7e6e4BCLFGEl7/6nMIknNh0rDuol27dphYmKCr68vAQEBmRmlyMYGDx7M3r17sba25siRI9SrV487d+4YOizxIdK43srZ2ZnNmzdz4sQJ6taty+vXr5k7dy6lSpVi6NCh8nqSjbyra7aiaFCAjkVjdXot/fzzz7i6uhIcHEzwhpEoyaQbGkXFV6abUf5NlCe0KcO4sd9QqFAhXrx4wdKlS2ncuDHFixdn7NixXLx4Ue+0dooxKgqdiyf8jPepj5ORJKn5UGlYN2Fra0ujRo0AGa0RadO6dWuOHTtG4cKFuXHjBnXq1OHMmTOGDku8r/dcb1WvXj1OnDjB/v37adiwITExMSxdupSyZcvSu3dvfH19Mzx08eGS20adS4kmeMd0lk4crrOl39LSko0bN5I/b27yqV+hQn/yoN0a/voZT3dM58Dy2UyfPp0HDx7g6enJgAEDyJs3Lw8fPuSnn37Czc2NSpUqMW3atP8S43971LWxuqE3RvO41wTvmM6MGdOpOWUfvZadYuSGi/RadooGsw4ZfKu3FN9LD2no1jp//nxGjhxJo0aNOHr0aPrFIHKER48e0b59ey5duoSVlRXr1q2jU6dOhg5LvI80dnnW59ixY8yYMYN9+/Zpj3Xs2JEJEyZoy0iIrOvNrtkO1paUtzOjSuVKPHjwgK+++oo5c+bonL9mzRomfN6PgnlMWLR4MRYFS/P8dQz5rBK2X7+MjMXathCxkRE0adyY6OhoJkyYoFPzKioqir/++ot169axZ88eoqOjtfdVrerK3o6ROBEETm4w+DDxCjoxujrlpvWAb7hfrDWAzg6pxH8t6lMt3Yvzpfb9W5KaTHbv3j1KlCiBWq3myZMn2NnZpXyREG947fsXwasHMHDLMw7d1fDzzz8zatQoqT6cg124cIHp06ezdetW7XRCixYtmDBhAk2aNJHfjWxk7969dOjQAbVazenTp6lRo4bO/f3792fTyZvYtxqGKvd/a1rervT7559/0qdPHwDWrVtHr169kvysly9fsn37dtatW8fhw4dpXgIO9MmtvX+10gm3rl9TqVIl7e9QvEah/kwvgkKj9Hb3VgGONpYcH9dMZwrtQ0lF4SyqePHiuLq6otFo2Lp1q6HDEdmNopDrn9kUzxXJ7z2LoigKY8aMoXPnztLhOwdzc3Nj8+bNXLt2jX79+mFiYoKnpyfNmjWjatWqLF26VBqlZhPt27enV69eaDQaBg0alKQPXMfPJ+HQaQLkyq9z/O1Kv71792bcuHEADBgwgHPnziX5Wfny5ePTTz/l4MGDBAUGsvbTcsT/O8wRp1Eo93gLVapUwcXFhQkTJuDj48OZOyEEhUXrTWjA8DVsJKkxgL59+wIJi7+kWqhIkzeqVpcwf87WWcMxMzNj586dVKxYkS1bthg4QGFI5cqVY+XKlfj7+zN8+HCsrKy4fPkyQ4cOpXDhwowcORI/Pz9DhylSMG/ePPLnz8+lS5dYsGCB9ni8RmH2wdugUiXb0fvNSr/Tpk2jffv2REVF0alTp3eWArALvYR97ANM/n1YU7WKWoVNaV/Oklu3bjFjxgxq1KhBrwGfpeo5GKqGjSQ1BjBkyBDy58/PrVu32LZtm6HDEdmFnqrVXfJd5dzZs7i6uhISEkK3bt3o3bs3z59nvUqfIvOUKFGC3377jUePHvHLL79QunRpwsLCmD9/PuXKlaNly5Zs375dusFnUQ4ODsyaNQtIaJfy4MEDIGFtS2Bo8snC26MkJiYmrFu3jvLly/Po0SM6d+5MVJSe699RM2nnlzVZv24dXbt2JVeuXATduZm652CgGjaS1BhAnjx5+PLLLwGYPn26VIoVqZNM1eoquYM5c+YM3333HWq1mnXr1lGpUiX++usvw8YrDC5//vyMHj0aPz8/9u/fj4eHB2q1Gk9PT7p06UKpUqWYNm2aTF1mQQMHDqRevXpEREQwcuRI4P0q/ebNm5ddu3aRP39+Tp8+zeAhQzkZ8Ex3G/Y7aiaZPLlEz5r2bN68meDgYNb9OgWL+NcoSjI79xQF+9ymBqthIwuFDSQkJITixYsTERHBvn37tNWGhdBLUWBZUwi8rPvCozKBQlVg8GFQqTh9+jT9+vXTTjEMHjyYn3/+GWtrawMFLrKau3fvsmTJEn7//XeePUvYfWVmZkaXLl34qnM1agRvRtV2Njg3NXCk4sqVK1SrVo24uDh27dqFXYV69Fp2KsXr1g+uQ11n3V5SXl5edPx8EvmaDsI0r732eKG8FuzLPYl8L66hv8SAGpxcta8xAPt9Axm29nySjeUJiY4K22tbubD7j7Q92RTIQuEsrkCBAgwZMgSAGTNmGDgakeWlsmp17dq1uXDhAqNGjQJg2bJlVKlSRcoHCK0SJUowY8YMHjx4wJo1a6hTpw6xsbFs3LgRxXMyqmc3ub9yEJcvXTJ0qDle5cqVGTNmDACff/45FR0s3rvSb2zBCth1/BYTa90dt8/DXhH7/AFpqZmUWGen0Fs1bKyUaF7t/4WfR3+SmqeXIWSkxoAePnxIqVKliI2N5fjx49SvX9/QIYmsKHGU5vElUvtJCuDIkSN8+umn3L17F4BRo0Yxffp0rKysMiVskX1cuHCBU2unMcz6v27grddG8NS6Mn379uXjjz+mYMGCBoww54qIiKBChQrcv3+fsWPH0vSTUQxbex7Q7fKUOEoyp5ML3eqW0XmMeI1Cg1mHkl2P40QIZayjWNG/Jib6djUlUzPp7To7tUraokJBrU7/8RIZqckGihQpQr9+/QAZrRHv8J7VZ5s0acLly5cZPHgwkLCjokKFCqxbty7F3i8iZ3GrWpVhLs9Q/l0oGq/AtGZWXLx4kTFjxlC4cGE6dOjApk2b9C80FRkmd+7cLFy4EIBffvmFwsozvZV+VZGhBO+YzoY53yZZp5nSAuPHFOBoeGHORBVLU486E7WKus4F6Fi1MHWdC2CiVmVIQpMWMlJjYLdu3aJcuXJoNBouXbpElSpVDB2SyIo+sPrsvn37GDx4sLYBYrVq1Zg9ezbNmzdP70hFduTvCWs/SnJ4T74B/G/DSU6fPq09ZmNjw0cffUS3bt1o1qwZ5ubmmRlpjtWlSxe2b99OvXr1OHbsGAoqnVES0xd3qVevHiaOLnwx9ju6tG1BrZK2mKhV7Lz4iJEbLqb4M37tWZWOVd9dxdpQpKKwHlkxqQHo0aMHmzZtokePHmzYsMHQ4QgjFRERwbx58zi9fhYzG2v4cl8UZi4tmTVrliTTOVkqFqH73bzJmjVrWLNmDffv39eeki9fPjp27EjXrl1p2bIlFhYWBngCOcODBw+oUKECr169YtmyZQwaNEjn/v2+gXz15ykilP+SzMQqwzZW5u+9wDirkKRGj6ya1Fy8eBE3NzcA/vrrL9q2bWvgiITRUhRiFzXE7OkVzj6Op9ayCFQqFX379mXKlCkUK1bM0BGKzJbMKI1Wn61QugUAGo0Gb29vNm/ezNatW3W2gufNmxd3d3e6du1K69atZe1WBpg7dy5jxowhf/783LhxAwcHByD53UiJq2N++7gaU/deIyg0Sm8rTEXR4JDbjFPft07X1gbpSdbUZCNVq1bV1iEYOHAgISEhBo5IGK0AL8yeXgGgppMJ0z5tiqIorFq1irJlyzJu3DhevHhh4CBFpkksupbsW4E64f5/P/uq1WqaNGmiLezn7e3Nl19+SeHChQkLC+PPP/+kc+fO2Nvb07NnT7Zs2aLTbVp8mC+++AJXV1devHjBN998AyQs1p28+5r+ZOXfr1P3XuOH9uUBku6cUhRAxcsjKzh67aFu/ZpsKNskNT/++COqf0tDJ97KlStn6LDSzYwZMyhXrhyBgYEMHz5cCvKJ9KenIvGEOhpOnzpF4387+s6ePRtnZ2d+/vlnWRCaE7znInRIqFbbsGFDfv31V+7fv88///zDmDFjKFasGBEREWzcuJFu3bpRoEAB2rVrx//93/9x7969hIsDDsPCWglfRaqZmpqyZMkSVCoVq1ev5vDhw6muMpw/t4XeBcYF81qguXqAyHLtGLD2MiM3XKTXslM0mHVI20cqO8k2008//vgjW7ZswdPTU3vM1NQ0TV2us+r0U6Jz585Rt25d4uLi+PPPP/n4448NHZIwJslNM/TZiuLcnL/++otx48Zx9epVIKH56uTJk+nVq5csBjVmH7gI/W2KonDu3Dm2bNnCli1buH37ts79lSpV5K/OURRVB6MUckM15HCyzRGFfsOHD2fRokW4uLgwbe0Bvtrim+I1iYuA396G/SIimuHrLqAoik4/qcR/LepTTdv525CMbk3Njz/+yI4dO7h48WKqr4mOjiY6Olr7fVhYGEWLFs2ySQ3AlClTmDRpEvny5ePKlSsUKVLE0CEJY5DKisTx8fGsWrWKH374gcePHwPg5OTE559/ztChQ7G1NUzpc5E9KYrCtWvX2Lt3L3v27OHEiRO0KKniQJ/c2nNmPa5N8WYDaN26Nfnz53/Ho4lEL1++pFy5cjx58oRhk+byV1SZFK/Rtwg4pfo1KsDRxpLj45oZfK2NUa6puXXrFk5OTpQqVYrevXvrrMLXZ8aMGdjY2GhvRYsWzaRI39+ECROoVasWL1++5NNPP5V6IiJ9pLIisYmJCQMGDODWrVvMnDkTR0dHHj9+zIQJEyhSpAjDhg2TLs8i1VQqFRUrVmTs2LF4e3sT/PQp6wZWIP7fj9JxGoWmygl69eqFvb09jRo1YurUqZw8eVKabb5Dvnz5mDt3LgArZk3ALrdpslWGFUVDfFgwuSIeJ7kvrQ0ys4NsM1Kzb98+Xr16hYuLC4GBgUyePJlHjx7h6+ubbF+b7DhSA+Dn54ebmxuRkZEsWLCAzz//3NAhiezsPSsSQ8Lf0MaNG5k7d67OKGm7du0YPXo0zZs31xmyFuKdkpkCHXrcnqVeATrHbGxsaNq0KS1atKBly5aUKVNGftfeoCgKrVq1wtPTk7rdhhFYqn3C8TfOUf173tMd06mYN5aTJ09iamqqvT871a8xuumnt718+ZLixYvzyy+/MHDgwFRdk9XX1Lxp4cKFfPHFF1hZWXHhwgVcXFwMHZLIruKiYW4liHia/Dl5HGCUL5jqrzOiKApHjx5l7ty57N69W7uQvXLlyowePZpevXphaWmp91ohgBSnQO+2/IO/Dx7k4MGDeHl5JdmFV6xYMW2C07x5c+zt7cnpbt26ReXKlYmOjub7xZvxfJ5fZ+SlkI0lX9QvxIiO9XkZGsbnk+fSwv0jbUuDM3eeZ5v6NUaf1ADUrFmTFi1apLrFQHZKajQaDW3atOHgwYPUrFmTEydOYGZmZuiwRHaVjotBb926xfz58/njjz+023UdHBz4aVh7Ps5/EdMOP0uHZ5FUGurhxMfHc/78eTw9PTl48CAnTpwgJkZ3B5arqytNmjShcePGNGzYME2bRoxJ4jpMR0dHfK9e4+YLDUGhkTyPiME2jwWOeS3Z/tdB1l2P0u3ObWPJD+0rpFy/Jo85p75rlW3W1GTbpObVq1cUK1aMH3/8kS+//DJV12SbpCbgMOwbR3DNbyjbZggvX76ke/fu/PnnnzpDh0IY0osXL1i2bBkLFizg4cOHnB6Ui1qFTbkTnZ8XXbfiVq2aTBeIBB8wBQoJ1bCPHz/OwX9Hci5fvpzknIoVK9K4cWPtLac04IyOjqZKlSrcvHmTL7/8kraDv2Xy7mtJ18oois5/28R/DWlUkqXedxJOeet8BbC5sgmfncsN/qHa6JKar7/+Gnd3d4oXL87jx4+ZNGkSFy9e5Nq1a6kehswWSY32j/8COLmxr8hYOnbqRGxsLL1792bVqlWYmJgYOkohtGJjYzmxagpNHs7XHmu9NoKHFi588skn9O7dO1ss0hcZKB2mQN/05MkTjhw5wtGjRzl69CjXrl1Lco6Li4s2wWnYsKFR/w56enrSsmVLcrvUx67T+FRfl7i76Yf25Zm697pOIuSQx4y72+YQfOEgEyZMYNq0aRkQeeoZXVLTs2dPvL29CQkJwd7engYNGjBt2jScnZ1T/RjZIql5e4i2z1Z2Xo2ga9euxMXF0b9/f5YvX27wTqhCaP2biCuBl1Ep8cQrcCFQQ81lr4CEHTCNGzfmk08+oWvXrln3b09krHSuh/Om4OBgjh07pk1yLl++nKSAadGiRalXr5725urqavDRh/T0UddunHZoj6m1XZrr/qwfXEe7xiaxfk2tkrZs37aVbt26oVKpOHT4MBZFKuncn5lTUkaX1KSHLJ/UvL2Q7o0aIlu2bqVnz57Ex8czePBgFi9eLImNyBqSWSux324IM7ec5ejRo9pjlpaWdOzYkU8++YRWrVoZ1ZuKyDpevHjBsWPH8Pb25ujRo1y4cIH4eN1yBlZWVtSsWVOb5NStWzfpupx/lwLQdlaWXye28+Q1Ru68817Xvmt304ABA9hw4gb2rYajyv1fHaHEZpmZVZhPkho9snxS846Kr5Ruwfr16+nTpw8ajYbhw4ezcOFCWbMgDCsVRf3u3b/PunXrWLNmDdevX9eeYm9vT69evfjkk0+oXr26/C6LDBMREcHZs2f5559/tDd9Pc7Kli1L7dq1qVWrFrVq1qTGxbGoAy+Ck1uy632yitRuz9bnXbubtp+9w+gtV1HAoBWHJanRI0snNams+Lp69Wr69++PoiiMHDmSuXPnypuBMJw07GhRFIXz58+zZs0a1q9fz9On/62vKFeuHH369KFPnz4UL148o6MWOZxGo+HmzZs6Sc6bCTdAK2cTncrHno7DKdy4Ly4uLllylPxkQEiqtme/SVE0OFpb8M+ElnqnkrJSxWFJavTI0klNGt4cli9fzqBBgwD45ptvmDVrliQ2IvN9wI6WuLg4/v77b9asWcOOHTt0mmfWq1cPd3d33N3dqVChgvxui0zx/PlzTp06xZkzZzh79gxTi/1DFXsNpmoVcRqF84Hx1P79NXnz5qVGjRrUrFmTatWqUa1aNUqVKmXwRCelBCSJf3c3lXlyBM+Vc/SektpEKTPq2EhSo0eWTWre481h8eLFDBs2DIDvvvuOqVOnyou/yFzptKMlLCyMbdu2sWbNGg4fPqyzwLNkyZLaBKdRo0bSWFNkjmQ+ZHbcHM+uaxFJjtvY2ODm5qZNcqpVq0bZsmUzfafqft9Ahq09n/A3lML7gZ2VmhsbphFx4wQ7duygY8eOSc7JShWHJanRI8smNe/55rBgwQJtjZ5hw4bxyy+/SFVXkbnSeUfLo0eP2L17N7t27eLQoUM6bU7y5s1L69atcXd3p127dhQoYNgKp8JIvWMpgOJYhSt1fuXM2bP4+Phw/vx5Ll26pPN7mih37txUrVpVe3N1daVSpUpYWVllaPj7fQP5YftlgiP+652VWGgvf25znd1L300Yz6xZs3BycuLatWvY2NjoPJaM1GRxWTapgfd+c5g7dy5jxowBoGrVqmzcuJGyZctmVJRCZJqIiAgOHjzI7t272bNnj84aHLVarTNNVa5cORmpFOkjDUsBIKFO0/Xr1zl//jznz5/Hx8eHixcv8vr16ySXqtVqXFxccHV11d6qVq2Ko6Njuv7+xmsUhk2cw5/b9lAglylXj+zCOk/uJOdFRkZSuXJlAgICGDZsGP/3f/+X5HEazDqUbMVhSEiYZE2NgWTppOYD7Nu3j759+/Ls2TPy5MnD4sWL6d27t6HDEiLdaDQazp49y+7du9m9e3eSirLOzs7aBKdhw4ayVVy8nw+sfJwoPj6emzdv4uPjw6VLl7h48SKXLl0iODhY7/n29vZUrFgxye1DRiNfv35NhQoVuHfvHt9//z1Tp07Ve97hw4dp1qwZAMeOHaNBgwY69ydOaYFuxWFF0aBCxYJerrhXzfjChpLU6GGsSQ3A48eP+fjjj7U1QT799FMWLFhA7txJs3Mhsrt79+6xZ88edu3axZEjR3T6AuXNm5dGjRrRpEkTmjZtiqurq1ThFqmTzpWP36QoCkFBQdoEJ/Hm5+eHRqMvgQJHR0e9yU6+fPlS9TO3b99Oly5dMDc35+rVq5QuXVrveYMGDWL58uWUK1eOixcvYmGh+9z2+wYmab2geRXCs4OL+bpHC6ZMmZK6/wgfQJIaPYw5qYGETwdTp05lypQpKIpC+fLl2bRpE5UqVTJ0aEJkmPDwcP7++292797N3r17efZMdxo3X758NG7cWJvkVK5c2eA7VUQWloGVj/V5/fo1169f5+rVq/j6+nL16lWuXr3KvXv3kr2mcOHCOklOpUqVqFChAtbW1jrnKYpCmzZt+Pvvv2nfvj179uzR+3gvXrygfPnyPHnyhP/973989913Sc6J1yg6FYfv+xyiZ4/umJiYcPr0aapXr/5h/yFSIEmNHsae1CQ6fPgwvXv3JjAwEEtLS+bPn8+gQYNkzYEwevHx8Vy8eJHDhw9z+PBhjh07Rnh4uM45tra2NG7cmKZNm9K0aVMqVKggSY7IcsLDw7l+/bpOonP16lUePnyY7DXFihXTSXYSa+o0bNiQ2NhYdu/eTYcOHfReu27dOnr37o2lpSVXr16lVKlSKcbYo0cPNm3aRMWKFfHx8UkywpOeJKnRI6ckNQBPnz6lb9++HDhwAEjonbVkyRKjf95CvCkuLo7z589z+PBhjhw5wrFjx4iI0N2Sa2dnpx3FadKkCeXLl5cPACLLevnyJdeuXdNJdHx9fQkKCkr2GisrKyIjI7G2tmbcuHFUqlQJFxcXnJ2dtevPFEWhZcuWeHl50bZtW/bu3Zvi38GzZ8+oWLEiT58+Zdy4ccycOTNdn+ubJKnRIyclNZCwuHLOnDlMmDCB+Ph4nJ2dWblyZZKFYELkFLGxsfj4+GhHck6cOJFkl0rBggVp0qQJ9erVo1atWri5uWXoJ1Ah0sPz5891kpwbN27g5+f3zpEdExMTSpUqRbly5XBxccHGxobJkycTFxfH5s2b6dq1a4o/d8eOHXTu3Bm1Ws2JEyeoU6dOej4tLUlq9MhpSU2ikydP0rNnT+7fvw9A586dmTFjBi4uLgaOTAjDiomJ4ezZs9qRnBMnTuhUNwYwMzOjatWq1K5dW3srXbq0jOaIbOHVq1fcvHmTlStXsmDBAkxMTChfvjx37txJMmr5JpVKRdWqVSlbtiylS5emTJkylC5dmtKlS+Pg4KDz+//JJ5+wdu1aypYty8WLFzOkFo8kNXrk1KQGErL4b7/9luXLl6PRaDAxMWHw4MFMmjQJR0dHQ4cnRJYQHR3NmTNnOHr0KKdPn+b06dN6t+Ha2tomND2sVUvbADFJh2chshBFUWjSpAne3t706tWLP//8k0ePHuHn56dzu3HjxjsXKQPkyZNHm+CUKVOGQoUK8eOPP/L8+XNGjRrF3Llz0z1+SWr0yMlJTaJr164xfvx4du3aBSRUvfzmm2/46quvyJMnj4GjEyJrURSFu3fvahOc06dPc/78eb0VZJ2dnbUjOTJtJbKiCxcuUL16dRRF4fjx49SvX1/vebt27aJjx46o1WpGjx5NREQE/v7++Pv7c+/ePVJKG8aMGcPPP/+crrFLUqOHJDX/8fb25ptvvuHMmTNAwjqCH3/8kYEDB0rhMiHeISYmhsuXL3P69GnOnDnD6dOn8fPzS3Je4rRVtWrVtNVjK1eunGTbrRCZaciQISxbtozq1atz5syZZHf+devWjS1btlCnTh1OnDihPS86Opq7d+9y69YtbaKTeLt9+zaKotCrVy/WrVuXrnFLUqOHJDW6FEVh69atjB8/Hn9/fwBcXFyYOXMmHTt2lDUDQqTSixcvOHv2rM6Iztv1chI5OzvrlMl3dXWlePHi8vcmMsXTp08pU6YMYWFhrFixgk8//VTvefcfPMS1xUfEmFjy9YjBTBzWO8VWCMHBwXh6etK2bdtUFwhMLUlq9JCkRr+YmBiWLl3KlClTtOsH6tWrx08//US9evUMHJ0Q2U/itNWZM2d0Ksg+fvxY7/k2NjZUqVJFJ9FJtvlhwGHYNw7azgLnphn8TIQx+vnnn/n6668pWLAgN2/eTPJ+qK+CcEFrcyZ3rESbSoUyO1xAkhq9JKl5t7CwMH766Sd+/vlnIiMjAWjdujXDhw+nXbt2mJqaGjhCIbK3Z8+e6ZTIv3TpEteuXSM2NjbJuWq1mrJly+pMXZUvV45SXoNQBV4AJ7cUexAJoU9MTAyVKlXi1q1bSerLJPZ6SpIYKAoqlYpFfaoZJLGRpEYPSWpS5/Hjx0yaNIkVK1Zoe5IUKVKEwYMHM3DgQAoXTr8S4ULkdDExMdy4cSNJsqNv11UrZxMO9Pmvn9taumBeoS3lypWjbNmyWFpaZmboIhvbs2cP7u7uOn2hErtyvzlC87bM6sr9Nklq9JCkJm0CAgJYunQpK1as0K4PMDExwcPDg88++4wWLVpIeXkhMoCiKAQGBuokOVev+rKizl2qOqowVauI0yicD4yn9u8JxQPVajUlS5akfPny2lu5cuUoX758uq9vENmfoii0bduWAwcO0LFjR3bs2MHJgBB6LTuV4rXrB9ehrvP7dxB/H5LU6CFJzfuJjo5m69atLF68GPOH/zC/rSVf7ovirroEQ4cOpX///tjb2xs6TCGMm78nrP0oyeFR54qx8vgDQkNDk73U0dFRJ9lJLKRWvHhxmVbOwa5du0aVKlWIj4/n4MGD/BNZiBUn7qZ43a89q9KxauaO2EtSo4ckNR9IUXg9vy65XlzHJwhqLAkDwNzcnK5du/LZZ5/RoEED2cUhRHpTFFjWFAIvgxL/33GVCRSqgjLoEE+ePuX69eva240bN7h+/TqPHj1K9mFNTU0pUaKEtpCas7Oz9t8lS5aUOjs5wMiRI5k/fz4uzbsTVaNvqq5Z0ceVZpWKZHBkuiSp0UOSmg/01ifFA/ZD+X7lYc6dO6c9VqFCBT777DM++eQTGfIWIr0kM0qj1WcrlG6h966wsDBu3LihTXKuX79OQEAA/v7+SVpCvEmlUlGsWLEkyU7p0qUpVaoUuXPnTvZakX08f/6cMmXLYtV9DqZ57YB3fChVFOLCn9EybxCden6Cg7UltUraZsr6Gklq9JCk5gO8/Unx30+IDD7MOR8flixZwrp167TNAc3NzWnWrBkeHh54eHjI4mIh3lfi397jS4BGzwlqcHJN804ojUZDYGCgTvG0xGTH39+f8PDwd17v5OREqVKlKFGihPZWvHhxSpQoQdGiRdM+yiNb1Q3m21+Ws+FpatrlKMS/DsMkl432SCEbSya5V8jwHVGS1OghSc0HSO6T4hufEENDQ1m7di2LFy/G19dX57Tq1avTsWNHPDw8qFKlikxRCZFacdEwtxJEPE3+nDwOMMoXTNNnukhRFIKDg/UmO/7+/jx//vyd16tUKpycnLRJzpsJT4kSJShWrJjuTi1t4iZb1Q1hu88DRm++nLqTFUXn/03ivzJ6q7ckNXpIUvOeUpjPf/sFSFEUbty4wc6dO9m1axenTp3S6RVSvHhxPDw86NixI40aNZK2DEKkJPQhROivUAxAbnuwybzR0BcvXuDv78+dO3e4e/eu9nbv3j3u3r2rHbF9F0dHR22S06KkioHme7X3veywApvqXeTDTyZJ7a6n5KgAxwze6i1JjR6S1LynD5jPB3jy5Al79uxh586dHDx4UGce38bGhnbt2uHh4UHbtm2xsbFJ9nGEEFmfoig8e/YsSaLz5i0iIkLnmtODclGtkInOVvWGq+MoUqQIRYsWpUiRIjr/Tvxqb28viU860NaneRn5QSNkGbnVW5IaPSSpeQ/pPJ//+vVrPD092blzJ7t379YpMGZqakqTJk3o2LEj7u7uFC9ePP2ehxAiS1AUhefPn2sTHMXfk66RSZsftl4bwd8B8Xoe4T/m5uZJEh4nJyecnJwoVKiQ9qsUJUzZft9APlt7HkXRoFL9V39MBUmrCycjI7d6S1KjhyQ17yED5/Pj4+M5ffo0u3btYteuXVy/fl3n/rJly9KgQQPq169PgwYNKFOmjHwqE8KYJDO1rahMiLEtxxnXmTx89IiHDx/y4MEDna9BQUGp/jH58+dPkuhI8pPUft9ARq46QbTJfz3HCtlY0rNmMeZ63kzxehmpyWSS1LynTJrPv3XrFrt27WLnzp2cOHFC26Ihkb29PfXr19cmOdWqVcPc3PyDf64QwkA+YGo7JiaGx48fa5OcxIQnMDCQx48fa79GR0enOpz8+fPrJDyOjo4ULFgwyc3Ozs5oq6k/ehxIpaadiDGxZPyo4Ywf1B2ABrMOERQapX/URlEolM9K1tRkNklqso8XL17wzz//cOLECY4fP86ZM2eSvDhZWlpSq1YtbZJTt25d8ufPb6CIhRBpkkFb1XV/hMLLly95/PixTqKj7+u7avYkiUytxt7ePtmk582bvb09JiYm7xW/oUybNo3vv/+eYsWKcePGDaysrLSNLkF3OkpRNICKH1sU5tOWbhkWkyQ1ekhSk31FR0dz/vx5jh8/rk10QkJCkpxXqVIlbZJTv359SpQoIVNWQmRFBtiqnpw3k5/EROfx48c8efJEewsKCuLJkyd6X3feRaVSYWdnh4ODA/b29tjb22v/re+Yra2twUeBXr9+jYuLCw8fPmTGjBl8++23QML01OTd13QaXprGhPN473zaVHRk27ZtGRaTJDV6SFJjPBRF4ebNmzpJzq1bt5Kc5+TkRP369alZsyaVK1emSpUqFCpUSBIdIbKCLLZVPTViY2MJDg7WSXiSuwUHB5PWt1i1Wo2dnZ32VqBAAZ2v+o7Z2NikeyK0Zs0a+vbti7W1Nf7+/jg4OAAJO6X+v717j4qyzv8A/h4uM4PIHeWuBqhpiHcIxaTUpdXMdtv1UpGRux5NTWXLu1Ga4mFdM1NzI1O3Y7FHxda8oBtFrZcWAzEE1ABdTJGL3EGGYeb7+6MfcxwYkBmYGZjer3OeM888833m+Xxm0HmfZ555nrQb5SipaUBfBzns64swauQIqFQqfP3113jySeOcOJGhRgeGGstWXFyM8+fPa4JOeno6mpqaWo1zdXVFcHCwJuQMGzYMQUFBPO07EXUplUqFsrIyTcApKSlBaWlpm/MVFRUGbcfKygpubm7tBp+Wt87Ozu1+LaZWqxESEoL09HQsWLAAH374YZtjFy9ejF27dmH48OFIT083ytdtDDU6MNT8utTX1+PixYs4d+4cLl++jB9//BHXr19vdQAy8MsuYn9/fwQHB2sFHn9//x73fTgR9UxKpRJlZWUoLS1FWVkZ7t27h7KyMq35lssedjmLtkgkEjg5OcHV1RUuLi5wdXXVmndxcUFpaSni4+MhkUiQlJSEsWPHwsXFBXZ2dlp7u8vKyhAYGIiqqip8/PHHmDdvXle9JBoMNTow1ND9+/eRm5uLrKws/Pjjj5rb4uJinePt7OwQFBSktVcnODgY7u7uJq6ciKg1hUKB8vLyDoegsrIyVFdXd2qbMpmsVRC6e/cuLl68CHt7e5w4cQITJ07sog5/wVCjA0MNtaWkpARZWVlaYefKlStt/iLC09MTQUFBmisYN1/FmFcvJqLurrGxEZWVlSgvL0d5eTkqKipazTff3rlzB5mZmQB++ZpL157ulp5++mmcOnWqS2tmqNGBoYb0oVKpkJeX12qvTkFBQbvreXp6aoJOy8nd3Z0HKRNRj7J06VLs2LEDQUFB+M9//oOqqiqdQSgtLQ2XLl1CfHw8/vCHP3RpDQw1OjDUUFeoqalBdnY2cnJykJ+frzU97EA/R0fHNgOPr68vj98hom7n3r17CAwMRGVlpdGOmXkYhhodGGrI2CoqKrRCTl5enmb+9u3b7a4rlUoxYMAATcjx9/eHr68vfHx84OvrCy8vL17RnIjMYtu2bfjLX/4CT09P/PTTT+jdu7dJt89QowNDDZnT/fv3cePGDZ2h5+bNm1Aqle2uL5FI4OHhAR8fH03QaTnv6+tr8v9siMjyKRQKDB06FAUFBYiNjcXbb79t0u0z1OjAUEPdlUqlwq1bt7QCz82bN3H7/y/md+fOnYeGnmaOjo6tgk7L8GPJ164hIuM4dOgQZs6cCXt7e+Tl5cHT09Nk22ao0YGhhnoqtVqN0tJS3L59WxN0Ws7//PPPHT5nhVQqhbe3tybsNF+j5sHTtzfPOzs78+BmIoIQAo8//jjS0tKwcOFC7N6922TbZqjRgaGGLF1NTU27oef27dsoKSnR69Tttra2mmvX6Ao9LZf17t2bIYjIQqWmpuLJJ5+EtbU1cnJyMGjQIJNsl6FGB4Yaol/OUVFUVKQJPbdv39Y6VXtJSYlm3pCTdMlksnbDT58+feDm5gYXFxfNJJVKjdApERnDtGnTcPLkSTz//PM4fPiwSbbJUKMDQw2RfhoaGlBWVqYz8Oi6raurM2g79vb2WiGn+SylLedb3nd2doaNjU0Xd01E7cnKysLw4cMhhMCFCxfw+OOPG32bDDU6MNQQGVd9fX2be30evK2oqEBFRQUqKyv1vopxS46Ojm2GngfvOzs7w9HRUWvq3bs3D5gmMkB0dDT279+PCRMm4NtvvzX6V84MNTow1BB1L2q1WuvspM1hp635B+8beiG/lhwcHODk5NQq8LQ3tRzv4ODAEyfSr8qtW7cwaNAgNDQ04NixY5g+fbpRt8dQowNDDZHlUCqVqKysfGgQap6vqqpCdXU1qqurUVVVBZVK1aX12Nvb6wxAvXv31prs7e11zre8L5fLecA1dWsrV65EfHw8hg4disuXLxv1q2CGGh0YaogI+OWnqQ0NDZqQ8+D0YPjpyKRQKIxSo5WVVbuhp737vXr1andiYKKuUFFRgYCAAFRUVBj98gkMNTow1BBRV1MoFKipqWkzDNXW1qK2thZ1dXU651vev3//vknqfljw6chkZ2cHOzs7yOVyyOVynfMymYzHLVmwv/3tb3jjjTfg4+OD69evo1evXkbZDkONDgw1RNTdqVQq1NfXtxl62gtEDwaj+vr6VlNjY6NZepLJZJqw01b4aTnfFeNkMhn3SBlZQ0MDBg8ejMLCQsTFxWHVqlVG2Q5DjQ4MNUT0a9bU1NRm4DF0amho0Ez379/X3KrVanO3CwA69xw9OEml0ofOd3RcR9aRSqUWF7Q+/fRTvPzyy3ByckJ+fj7c3Ny6fBsMNTow1BARmYZSqdQZeNqb76px3f1jzdbW1qDwZGtrq/O2vcc6Mqa9xzry1aFarcaoUaNw+fJlLF++HNu2bevy16yjn988axUREXU5W1tb2NrawsHBwaTbFUJoBSpd4aexsREKhUIzPXi/Kx5rOa7lL+2USiWUSiVqa2tN+toYwsrKqkPBqfmCu9u3b0dQUBBeffVVs9TLUENERBZDIpFoPoS7yx55lUqlCTqGhialUonGxkbN7YPz+jzW3hhdx1yp1WpNKOwIIQQ+++wzhhoiIiJLZG1trfmlWHcmhIBKpTI4OF2/fh05OTlYuXKl2XpgqCEiIiJIJBLY2NjAxsam2wewtvDkAURERGQRGGqIiIjIIjDUEBERkUVgqCEiIiKL0ONCza5duzBgwADI5XKEhoYiLS3N3CURERFRN9CjQs0///lPxMTEIDY2FhkZGRg+fDgiIyNRUlJi7tKIiIjIzHrUZRJCQ0MxduxY7Ny5E8AvJwXy8/PDkiVLdF5Eq/mkRc2qqqrQr18/3Lp1q9uclImIiIjaV11dDT8/P1RWVsLJyanNcT3mPDWNjY1IT0/H6tWrNcusrKwwefJkXLhwQec6cXFxeOedd1ot9/PzM1qdREREZBw1NTWWEWrKysqgUqng4eGhtdzDwwNXr17Vuc7q1asRExOjua9Wq1FeXg43NzeTXiW1OWFa4h4iS+4NsOz+LLk3wLL7Y289lyX3Z8zehBCoqamBt7d3u+N6TKgxRPPVTR/k7OxsnmIAODo6WtwfcTNL7g2w7P4suTfAsvtjbz2XJfdnrN7a20PTrMccKOzu7g5ra2sUFxdrLS8uLoanp6eZqiIiIqLuoseEGqlUitGjRyMlJUWzTK1WIyUlBWFhYWasjIiIiLqDHvX1U0xMDObOnYsxY8YgJCQE27dvR11dHaKjo81dWrtkMhliY2NbfRVmCSy5N8Cy+7Pk3gDL7o+99VyW3F936K1H/aQbAHbu3Im//vWvuHv3LkaMGIEdO3YgNDTU3GURERGRmfW4UENERESkS485poaIiIioPQw1REREZBEYaoiIiMgiMNQQERGRRWCo6SK7du3CgAEDIJfLERoairS0tDbHJiQkYMKECXBxcYGLiwsmT57c7nhz06e3pKQkjBkzBs7OzrC3t8eIESPw6aefmrBa/enT34MSExMhkUjw3HPPGbfATtCnt/3790MikWhNcrnchNXqR9/3rbKyEosWLYKXlxdkMhkGDRqEkydPmqha/enTX0RERKv3TiKRYNq0aSasuOP0fe+2b9+OwYMHw87ODn5+fli+fDkaGhpMVK3+9OlPqVRiw4YNCAgIgFwux/Dhw5GcnGzCajvuu+++w/Tp0+Ht7Q2JRIIvvvjioeukpqZi1KhRkMlkCAwMxP79+41bpKBOS0xMFFKpVHzyySciOztb/PnPfxbOzs6iuLhY5/gXXnhB7Nq1S1y6dEnk5uaKV155RTg5OYmff/7ZxJU/nL69ffPNNyIpKUnk5OSIvLw8sX37dmFtbS2Sk5NNXHnH6Ntfsxs3bggfHx8xYcIEMWPGDNMUqyd9e9u3b59wdHQURUVFmunu3bsmrrpj9O1NoVCIMWPGiKlTp4qzZ8+KGzduiNTUVJGZmWniyjtG3/7u3bun9b5duXJFWFtbi3379pm28A7Qt7eDBw8KmUwmDh48KG7cuCFOnz4tvLy8xPLly01cecfo29+KFSuEt7e3OHHihMjPzxe7d+8WcrlcZGRkmLjyhzt58qRYu3atSEpKEgDE0aNH2x1fUFAgevXqJWJiYkROTo744IMPjP55wFDTBUJCQsSiRYs091UqlfD29hZxcXEdWr+pqUk4ODiIAwcOGKtEg3W2NyGEGDlypFi3bp0xyus0Q/pramoS48aNEx9//LGYO3dutw01+va2b98+4eTkZKLqOkff3j788EPh7+8vGhsbTVVip3T23917770nHBwcRG1trbFKNJi+vS1atEg89dRTWstiYmLE+PHjjVqnofTtz8vLS+zcuVNr2e9//3vx4osvGrXOzupIqFmxYoV47LHHtJbNmjVLREZGGq0ufv3USY2NjUhPT8fkyZM1y6ysrDB58mRcuHChQ89RX18PpVIJV1dXY5VpkM72JoRASkoKrl27hieeeMKYpRrE0P42bNiAvn37Yt68eaYo0yCG9lZbW4v+/fvDz88PM2bMQHZ2tinK1YshvR07dgxhYWFYtGgRPDw8EBQUhM2bN0OlUpmq7A7riv9T9u7di9mzZ8Pe3t5YZRrEkN7GjRuH9PR0zVc4BQUFOHnyJKZOnWqSmvVhSH8KhaLV17x2dnY4e/asUWs1hQsXLmi9FgAQGRnZ4b9jQ/SoyyR0R2VlZVCpVPDw8NBa7uHhgatXr3boOVauXAlvb+9Wb765GdpbVVUVfHx8oFAoYG1tjd27d2PKlCnGLldvhvR39uxZ7N27F5mZmSao0HCG9DZ48GB88sknCA4ORlVVFbZu3Ypx48YhOzsbvr6+pii7QwzpraCgAF9//TVefPFFnDx5Enl5eXjttdegVCoRGxtrirI7rLP/p6SlpeHKlSvYu3evsUo0mCG9vfDCCygrK0N4eDiEEGhqasKCBQuwZs0aU5SsF0P6i4yMxLZt2/DEE08gICAAKSkpSEpK6paBW193797V+VpUV1fj/v37sLOz6/Jtck+NmW3ZsgWJiYk4evRotz4oUx8ODg7IzMzExYsXsWnTJsTExCA1NdXcZXVaTU0NoqKikJCQAHd3d3OX0+XCwsLw8ssvY8SIEZg4cSKSkpLQp08f/P3vfzd3aZ2mVqvRt29ffPTRRxg9ejRmzZqFtWvXYs+ePeYurcvt3bsXw4YNQ0hIiLlL6RKpqanYvHkzdu/ejYyMDCQlJeHEiRPYuHGjuUvrEu+//z4GDhyIRx99FFKpFIsXL0Z0dDSsrPjxbAjuqekkd3d3WFtbo7i4WGt5cXExPD09211369at2LJlC7766isEBwcbs0yDGNqblZUVAgMDAQAjRoxAbm4u4uLiEBERYcxy9aZvf/n5+bh58yamT5+uWaZWqwEANjY2uHbtGgICAoxbdAd15u+yma2tLUaOHIm8vDxjlGgwQ3rz8vKCra0trK2tNcuGDBmCu3fvorGxEVKp1Kg166Mz711dXR0SExOxYcMGY5ZoMEN6W79+PaKiovCnP/0JADBs2DDU1dVh/vz5WLt2bbf68Dekvz59+uCLL75AQ0MD7t27B29vb6xatQr+/v6mKNmoPD09db4Wjo6ORtlLA3BPTadJpVKMHj0aKSkpmmVqtRopKSkICwtrc734+Hhs3LgRycnJGDNmjClK1ZuhvbWkVquhUCiMUWKn6Nvfo48+iqysLGRmZmqmZ599Fk8++SQyMzPh5+dnyvLb1RXvnUqlQlZWFry8vIxVpkEM6W38+PHIy8vThFAAuH79Ory8vLpVoAE6994dOnQICoUCL730krHLNIghvdXX17cKLs3hVHSzSxd25r2Ty+Xw8fFBU1MTjhw5ghkzZhi7XKMLCwvTei0A4N///rdenx96M9ohyL8iiYmJQiaTif3794ucnBwxf/584ezsrPk5bFRUlFi1apVm/JYtW4RUKhWHDx/W+hlmTU2NuVpok769bd68WZw5c0bk5+eLnJwcsXXrVmFjYyMSEhLM1UK79O2vpe786yd9e3vnnXfE6dOnRX5+vkhPTxezZ88WcrlcZGdnm6uFNunbW2FhoXBwcBCLFy8W165dE8ePHxd9+/YV7777rrlaaJehf5fh4eFi1qxZpi5XL/r2FhsbKxwcHMTnn38uCgoKxJkzZ0RAQICYOXOmuVpol779ff/99+LIkSMiPz9ffPfdd+Kpp54SjzzyiKioqDBTB22rqakRly5dEpcuXRIAxLZt28SlS5fE//73PyGEEKtWrRJRUVGa8c0/6X7zzTdFbm6u2LVrF3/S3VN88MEHol+/fkIqlYqQkBDx/fffax6bOHGimDt3ruZ+//79BYBWU2xsrOkL7wB9elu7dq0IDAwUcrlcuLi4iLCwMJGYmGiGqjtOn/5a6s6hRgj9elu2bJlmrIeHh5g6dWq3PFdGM33ft/Pnz4vQ0FAhk8mEv7+/2LRpk2hqajJx1R2nb39Xr14VAMSZM2dMXKn+9OlNqVSKt99+WwQEBAi5XC78/PzEa6+91i0/9Jvp019qaqoYMmSIkMlkws3NTURFRYnbt2+boeqH++abb3R+djX3M3fuXDFx4sRW64wYMUJIpVLh7+9v9HMnSYToZvvviIiIiAzAY2qIiIjIIjDUEBERkUVgqCEiIiKLwFBDREREFoGhhoiIiCwCQw0RERFZBIYaIiIisggMNURERGQRGGqIiIjIIjDUEJHFiYiIwLJlyzr1HEIIzJ8/H66urpBIJMjMzOyS2ojIeBhqiMiooqOjsW7dOnOXobfk5GTs378fx48fR1FREYKCgsxdEhE9hI25CyAiy6VSqXD8+HGcOHHC3KXoLT8/H15eXhg3blybYxobGyGVSk1YFRG1h3tqiEjj888/h52dHYqKijTLoqOjERwcjKqqKr2f7/z587C1tcXYsWN1Ph4REYElS5Zg2bJlcHFxgYeHBxISElBXV4fo6Gg4ODggMDAQp06d0lpPoVDg9ddfR9++fSGXyxEeHo6LFy+2WYdarUZcXBweeeQR2NnZYfjw4Th8+HCb41955RUsWbIEhYWFkEgkGDBggKbexYsXY9myZXB3d0dkZCSAX/bqhIeHw9nZGW5ubnjmmWeQn5+vtf34+HgEBgZCJpOhX79+2LRpU0dfRiLqIIYaItKYPXs2Bg0ahM2bNwMAYmNj8dVXX+HUqVNwcnLS+/mOHTuG6dOnQyKRtDnmwIEDcHd3R1paGpYsWYKFCxfij3/8I8aNG4eMjAz85je/QVRUFOrr6zXrrFixAkeOHMGBAweQkZGBwMBAREZGory8XOc24uLi8I9//AN79uxBdnY2li9fjpdeegnffvutzvHvv/8+NmzYAF9fXxQVFWkFpgMHDkAqleLcuXPYs2cPAKCurg4xMTH44YcfkJKSAisrK/zud7+DWq0GAKxevRpbtmzB+vXrkZOTg88++wweHh56v55E9BCCiOgBX375pZDJZOLdd98VLi4u4sqVK5rHnnvuOeHs7Cyef/75Dj3XwIEDxfHjx9t8fOLEiSI8PFxzv6mpSdjb24uoqCjNsqKiIgFAXLhwQQghRG1trbC1tRUHDx7UjGlsbBTe3t4iPj5e87xLly4VQgjR0NAgevXqJc6fP6+17Xnz5ok5c+a0Wdt7770n+vfv36rekSNHtt+0EKK0tFQAEFlZWaK6ulrIZDKRkJDw0PWIqHN4TA0RaXnmmWcwdOhQbNiwAWfOnMFjjz2meWzp0qV49dVXceDAgYc+T25uLu7cuYNJkya1Oy44OFgzb21tDTc3NwwbNkyzrHmPRklJCYBfjnVRKpUYP368ZoytrS1CQkKQm5vb6vnz8vJQX1+PKVOmaC1vbGzEyJEjH9pHS6NHj2617KeffsJbb72F//73vygrK9PsoSksLER9fT0UCsVDXwci6jyGGiLSkpycjKtXr0KlUrX6iiQiIgKpqakdep5jx45hypQpkMvl7Y6ztbXVui+RSLSWNX911RwU9FVbWwsAOHHiBHx8fLQek8lkej+fvb19q2XTp09H//79kZCQAG9vb6jVagQFBaGxsRF2dnYG1U1E+uMxNUSkkZGRgZkzZ2Lv3r2YNGkS1q9fb/Bz/etf/8KMGTO6sLpfBAQEaI5paaZUKnHx4kUMHTq01fihQ4dCJpOhsLAQgYGBWpOfn1+n67l37x6uXbuGdevWYdKkSRgyZAgqKio0jw8cOBB2dnZISUnp9LaIqH3cU0NEAICbN29i2rRpWLNmDebMmQN/f3+EhYUhIyMDo0aN0uu5SkpK8MMPP+DYsWNdXqe9vT0WLlyIN998E66urujXrx/i4+NRX1+PefPmtRrv4OCAN954A8uXL4darUZ4eDiqqqpw7tw5ODo6Yu7cuZ2qx8XFBW5ubvjoo4/g5eWFwsJCrFq1SvO4XC7HypUrsWLFCkilUowfPx6lpaXIzs7W1Ltz504cPXqUwYeokxhqiAjl5eV4+umnMWPGDM0HcmhoKH77299izZo1SE5O1uv5vvzyS4SEhMDd3d0Y5WLLli1Qq9WIiopCTU0NxowZg9OnT8PFxUXn+I0bN6JPnz6Ii4tDQUEBnJ2dMWrUKKxZs6bTtVhZWSExMRGvv/46goKCMHjwYOzYsQMRERGaMevXr4eNjQ3eeust3LlzB15eXliwYIHm8bKyMq2fgBORYSRCCGHuIoio50hNTcXOnTvbPc/Ls88+i/DwcKxYscKElRHRrx331BBRh02ePBmXL19GXV0dfH19cejQIYSFhbUaFx4ejjlz5pihQiL6NeOeGiIiIrII/PUTERERWQSGGiIiIrIIDDVERERkERhqiIiIyCIw1BAREZFFYKghIiIii8BQQ0RERBaBoYaIiIgsAkMNERERWQSGGiIiIrIIDDVERERkEf4PZsnx35pAq4YAAAAASUVORK5CYII=", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# Model from Lasala, FPE, 2016: https://doi.org/10.1016/j.fluid.2016.05.015\n", "j = {\n", " \"kind\": \"advancedPRaEres\",\n", " \"model\": {\n", " \"Tcrit / K\": [304.21, 126.19],\n", " \"pcrit / Pa\": [7.383e6, 3395800.0],\n", " \"alphas\": [{\"type\": \"PR78\", \"acentric\": 0.22394}, {\"type\": \"PR78\", \"acentric\": 0.0372}],\n", " \"aresmodel\": {\"type\": \"Wilson\", \"m\": [[0.0, -3.4768], [3.5332, 0.0]], \"n\": [[0.0, 825], [-585, 0.0]]},\n", " \"options\": {\"s\": 2.0, \"brule\": \"Quadratic\", \"CEoS\": -0.52398}\n", " }\n", "}\n", "\n", "model = teqp.make_model(j)\n", "for T in [223.15, 253.05, 273.08, 293.1]:\n", " ipure = 0\n", "\n", " [rhoL0, rhoV0] = model.superanc_rhoLV(T, ipure)\n", "\n", " rhovecL0 = np.array([0.0, 0.0]); rhovecL0[ipure] = rhoL0\n", " rhovecV0 = np.array([0.0, 0.0]); rhovecV0[ipure] = rhoV0\n", "\n", " J = model.trace_VLE_isotherm_binary(T, rhovecL0, rhovecV0)\n", " df = pandas.DataFrame(J)\n", " plt.plot(df['xL_0 / mole frac.'], df['pL / Pa']/1e6,'k')\n", " plt.plot(df['xV_0 / mole frac.'], df['pV / Pa']/1e6,'k')\n", "\n", "plt.plot(1-dat['x1'], dat['p/bar']/10, 'o')\n", "plt.plot(1-dat['y1'], dat['p/bar']/10, '^')\n", "\n", "plt.gca().set(xlabel='$x_1$ / mole frac.', ylabel='$p$ / MPa', ylim=(0, 25))\n", "plt.show()" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.11.0" } }, "nbformat": 4, "nbformat_minor": 5 }