{ "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_{ij}$ is non-dimensional and $n_{ij}$ has units of temperature because $A_{ij}$ has units of temperature." ] }, { "cell_type": "code", "execution_count": 1, "id": "de8cc3c3", "metadata": { "execution": { "iopub.execute_input": "2025-10-15T23:12:15.435640Z", "iopub.status.busy": "2025-10-15T23:12:15.435449Z", "iopub.status.idle": "2025-10-15T23:12:15.891631Z", "shell.execute_reply": "2025-10-15T23:12:15.891161Z" } }, "outputs": [ { "data": { "text/plain": [ "'0.23.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": "2025-10-15T23:12:15.893194Z", "iopub.status.busy": "2025-10-15T23:12:15.892975Z", "iopub.status.idle": "2025-10-15T23:12:15.898862Z", "shell.execute_reply": "2025-10-15T23:12:15.898452Z" } }, "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": "2025-10-15T23:12:15.900060Z", "iopub.status.busy": "2025-10-15T23:12:15.899942Z", "iopub.status.idle": "2025-10-15T23:12:16.290743Z", "shell.execute_reply": "2025-10-15T23:12:16.290263Z" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjUAAAG5CAYAAAB2j8WmAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjcsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvTLEjVAAAAAlwSFlzAAAPYQAAD2EBqD+naQAApBVJREFUeJzs3XdYU/f3wPF3wlYRRUDELSpuxL33VnDUWa1aZ9W2jrZabatVv87aatX+XNW66t6jagUHap04caDgVlBEBUR27u8PSmokCCgQCOf1PHmQm3vDSQvJyWeco1IURUEIIYQQIptTGzoAIYQQQoj0IEmNEEIIIYyCJDVCCCGEMAqS1AghhBDCKEhSI4QQQgijIEmNEEIIIYyCJDVCCCGEMAqS1AghhBDCKEhSI4QQQgijIEmNEEIIIYxClklqZsyYQc2aNbG2tsbBwYFOnTrh5+enc06TJk1QqVQ6t88++8xAEQshhBAiK8kySc3Ro0cZMWIEp06d4uDBg8TGxtKqVSsiIiJ0zhs8eDCBgYHa2+zZsw0UsRBCCCGyElNDB5Bo//79Ot+vXLkSBwcHfHx8aNSokfZ4rly5cHR0zOzwhBBCCJHFZZmk5m2hoaEA2Nra6hz/888/Wbt2LY6Ojri7u/PDDz+QK1cuvY8RHR1NdHS09nuNRsPz588pUKAAKpUq44IXQgghRLpRFIXw8HCcnJxQq5OfZFIpiqJkYlypotFo8PDw4OXLlxw/flx7fOnSpRQvXhwnJycuX77MuHHjqFWrFtu2bdP7OD/++COTJ0/OrLCFEEIIkYEePHhAkSJFkr0/SyY1w4YNY9++fRw/fvydwR86dIjmzZvj7++Ps7NzkvvfHqkJDQ2lWLFiPHjwgLx582ZI7EIIIYRIX2FhYRQtWpSXL19iY2OT7HlZbvrp888/Z8+ePXh7e78zoQGoXbs2QLJJjYWFBRYWFkmO582bV5IaIYQQIptJaelIlklqFEXhiy++YPv27Rw5coSSJUumeM3FixcBKFSoUAZHJ4QQQoisLsskNSNGjGDdunXs3LkTa2trgoKCALCxscHKyoqAgADWrVtHu3btKFCgAJcvX2b06NE0atSIKlWqGDh6IYQQQhhalllTk9yQ0h9//EH//v158OABffr0wdfXl4iICIoWLUrnzp35/vvvUz2VFBYWho2NDaGhoTL9JIQQQmQTqX3/zjIjNSnlVkWLFuXo0aOZFI0QQgghspssU1FYCCGEEOJDSFIjhBBCCKMgSY0QQgghjIIkNUIIIYQwCpLUCCGEEMIoSFIjhBBCCKMgSY0QQgghjIIkNUIIIYQwCpLUCCGEEMIoSFIjhBBCCKMgSY0QQgghjIIkNUIIIYQwCpLUCCGEEMIoSFIjhBBCCKMgSY0QQgghjIIkNUIIIYQwCpLUCCGEEMIoSFIjhBBCCKMgSY0QQgghjIIkNUIIIYQwCpLUCCGEEMIoSFIjhBBCCKMgSY0QQgghjIIkNUIIIYQwCpLUCCGEEMIoSFIjhBBCCKMgSY0QQgghjIIkNUIIIYQwCpLUCCGEEMIoSFIjhBBCCKMgSY0QQgghjIIkNUIIIYQwCpLUCCGEEMIoSFIjhBBCCKMgSY0QQgghjIIkNUIIIYQwCpLUCCGEEMIoSFIjhBBCCKMgSY0QQgghjIIkNUIIIYQwCpLUCCGEEMIoSFIjhBBCCKMgSY0QQgghjIIkNUIIIYQwCpLUCCGEEMIoSFIjhBBCCKMgSY0QQgghjIIkNUIIIYQwCpLUCCGEEMIoSFIjhBBCCKMgSY0QQgghjIIkNUIIIYQwCpLUCCGEEMIoSFIjhBBCCKMgSY0QQgghjIIkNUIIIYQwCpLUCCGEEMIoSFIjhBBCCKMgSY0QQgghjIIkNUIIIYQwCpLUCCGEEMIoSFIjhBBCCKMgSY0QQgghjEKWSWpmzJhBzZo1sba2xsHBgU6dOuHn56dzTlRUFCNGjKBAgQLkyZOHjz76iCdPnhgoYiGEEEJkJVkmqTl69CgjRozg1KlTHDx4kNjYWFq1akVERIT2nNGjR7N79242b97M0aNHefz4MV26dDFg1EIIIYTIKlSKoiiGDkKf4OBgHBwcOHr0KI0aNSI0NBR7e3vWrVtH165dAbhx4wbly5fn5MmT1KlTJ8XHDAsLw8bGhtDQUPLmzZvRT0EIIYQQ6SC1799ZZqTmbaGhoQDY2toC4OPjQ2xsLC1atNCeU65cOYoVK8bJkyf1PkZ0dDRhYWE6NyGEEEIYpyyZ1Gg0GkaNGkX9+vWpVKkSAEFBQZibm5MvXz6dcwsWLEhQUJDex5kxYwY2NjbaW9GiRTM6dCGEEEIYSJZMakaMGIGvry8bNmz4oMcZP348oaGh2tuDBw/SKUIhhBBCZDWmhg7gbZ9//jl79uzB29ubIkWKaI87OjoSExPDy5cvdUZrnjx5gqOjo97HsrCwwMLCIqNDFkIIIUQWkGVGahRF4fPPP2f79u0cOnSIkiVL6txfvXp1zMzM8PLy0h7z8/Pj/v371K1bN7PDFUIIIUQWk2VGakaMGMG6devYuXMn1tbW2nUyNjY2WFlZYWNjw8CBAxkzZgy2trbkzZuXL774grp166Zq55MQ2UF8fDzPnj0jKChIewsODiY8PJzw8HBevXql/RodHU18fDyHDh0CEhL/fPnyYWlpiaWlJVZWVtp/58uXjwIFCmBnZ4ednZ3234UKFcLS0tLAz1oIIdJHltnSrVKp9B7/448/6N+/P5BQfO+rr75i/fr1REdH07p1a/7v//4v2emnt8mWbpEVvHjxguvXr+Pv709AQID2671793j69CkajSZT4ylYsCDFihWjaNGiFCtWjDJlylCuXDnKlStHoUKFkv3bFEKIzJLa9+8sk9RkBklqRGbSaDQEBARw8eJFLl26pL2ltGBdpVLh4OCAo6Mjjo6O2NvbkzdvXvLkyUOePHmwtrYmd+7cWFpaYmpqyqFDh1CpVDRv3pz4+HiioqKIjIwkKipK+++XL1/y7Nkznj17RkhIiPbfUVFR74zF2toaFxcXypUrh5ubG9WrV8fNzU3+foQQmUqSGj0kqREZKTo6mrNnz3Ls2DGOHTvGiRMnkq2NVLRoUcqUKYOzszOlS5fG2dmZkiVL4uTkhJ2dHaamGT8zrCgKz58/5/79+zx48ID79+9z9+5dbt68yY0bN7h9+zbx8fFJrlOpVJQpU4bq1atTvXp16tatS40aNTA3N8/wmIUQOZMkNXpIUiPSU2xsLMePH8fT05Njx45x5swZoqOjdc6xtLSkcuXKuLq6am9VqlTBxsbGQFGnXnR0NAEBAfj5+eHr64uPjw8+Pj48fPgwyblWVlbUr1+fJk2a0KRJE2rWrClJjhAi3UhSo4ckNeJDPXv2jP3797N7924OHDigrXydyMHBgYYNG2pvVapUyZRRl8z09OlTbYJz7tw5Tpw4wbNnz3TOeTPJadOmDdWqVZO1OUKI9yZJjR6S1Ij3cf36dXbu3MmePXs4efKkzkJee3t72rRpQ+PGjWnYsCFlypTJcW/eGo2Ga9euceTIEY4ePcqRI0eSJDmFCxfGw8MDDw8PmjZtKvWjhBBpIkmNHpLUiNS6d+8eGzZsYN26dVy+fFnnPldXVzp06ECHDh2oVasWanWWKfeUJSiKok1yvLy8+Pvvv4mIiNDenydPHtq0aUPHjh1p166dtr+bEEIkR5IaPSSpEe/y6tUrNm3axIoVKzhx4oT2uJmZGS1atMDd3Z327dtTrFgxA0aZ/URFRXH48GF27tzJrl27CAwM1N5nYmJCo0aN6NWrFx999JEkOEIIvSSp0UOSGvE2RVE4deoUy5cvZ+PGjbx69QpI2OHTpEkTebNNZxqNBh8fH3bt2sXOnTu5cuWK9j4zMzPatm1Lr1698PDwIFeuXAaMVAiRlUhSo4ckNSJRZGQk69evZ8GCBVy8eFF7vEyZMgwcOJA+ffpQuHBhwwWYQ9y+fZvNmzcnmebLmzcvvXr1YsCAAdSsWTPHrVMSQuiSpEYPSWrEo0ePmD9/Pr///jvPnz8HErZd9+jRg4EDB9KgQQN5AzWQq1evsn79ev7880/u3r2rPV6xYkUGDhxIv379ZMRMiBxKkho9JKnJua5fv85PP/3E2rVriY2NBaBEiRIMHz6cgQMHyptlFqLRaDhy5AgrVqxg69at2qrHlpaW9OrVi+HDh1OjRg0DRymEyEyS1OghSU3Oc+bMGWbMmMGOHTu0xxo1asSYMWPo0KEDJiYmhgtOpOjly5ds2LCBxYsXc+nSJe3xmjVrMnz4cHr16iXbw4XIASSp0UOSmpzj/PnzTJw4kb1792qPdezYkXHjxlG3bl0DRibeR+KC7t9++43NmzcTExMDgKOjI59//jmfffYZBQoUMHCUQoiMIkmNHpLUGL/Lly8zadIk7ciMWq3mk08+Ydy4cZQvX96wwYl08fTpU5YvX85vv/3Go0ePgIQKxp9++imjR4+mdOnSBo5QCJHeUvv+LVXDhFF48OAB/fr1o2rVquzYsQOVSkXv3r25fv06K1eulITGiDg4ODB+/Hhu377NmjVrqFq1KpGRkfzf//0fZcuWpWfPnkkKJgohcgZJakS2Fh4eznfffUfZsmVZvXo1iqLQrVs3fH19Wbt2LWXLljV0iCKDmJub06dPH86fP8+hQ4do3749iqKwceNGXF1d6dixI2fOnDF0mEKITCRJjciWNBoNv//+O6VLl2b69OlERUXRsGFDzpw5w6ZNm6hQoYKhQxSZRKVS0bRpU/bs2cOlS5fo0aMHKpWKXbt2Ubt2bVq3bs2pU6cMHaYQIhNIUiOynfPnz1OvXj0GDx7M06dPKVOmDNu3b+fo0aPUrFnT0OEJA6pSpQobNmzg+vXr9OvXDxMTE/7++2/q1q2Lu7u7TqFFIYTxkaRGZBsvX77kiy++oGbNmpw+fRpra2t++eUXrl69SqdOnaRontBycXFh5cqV3Lp1iwEDBmBiYsKePXtwc3OjR48e3Lhxw9AhCiEygCQ1IstTFIUNGzbg4uLCwoUL0Wg09OzZkxs3bjB69GjMzMwMHaLIokqWLMny5cu5du0aPXv2BGDTpk1UrFiRTz/9VLt7SghhHCSpEVlaYGAgnTt3plevXjx9+hQXFxc8PT1Zv349Tk5Ohg5PZBNly5Zl/fr1XLp0CQ8PDzQaDStXrqRMmTJMnDiR8PBwQ4cohEgHktSILElRFFavXk3FihXZuXMnpqam/Pjjj1y+fJnmzZsbOjyRTVWpUoWdO3dy8uRJ6tevT2RkJFOnTqVMmTIsW7aMuLg4Q4cohPgAktSILCcwMBB3d3f69evHixcvqFatGj4+PkyaNAlzc3NDhyeMQJ06dTh27Bhbt26ldOnSPHnyhCFDhuDm5sbBgwcNHZ4Q4j1JUiOylF27dlG5cmX27t2Lubk506dP5/Tp01SpUsXQoQkjo1Kp6NKlC1evXmXevHnY2tri6+tLq1at+Oijj7h3756hQxRCpJEkNSJLeP36NcOHD6djx46EhIRQtWpVLly4wPjx4zE1NTV0eMKImZubM3LkSPz9/Rk5ciQmJiZs27aN8uXLM3XqVG2XcCFE1idJjTC4S5cuUaNGDRYtWgTAV199xalTp6SAnshU+fPnZ968eVy4cIHGjRsTGRnJxIkTqVixIrt37zZ0eEKIVJCkRhiMoij8/vvv1K5dm+vXr1OoUCH+/vtv5syZg4WFhaHDEzlU5cqVOXz4sHaH3e3bt/Hw8KBLly6yBVyILE6SGmEQkZGRDBw4kMGDBxMdHU2HDh24fPkyLVu2NHRoQqBSqejZsyd+fn6MHTsWU1NTtm/fTvny5fntt9+Ij483dIhCCD0kqRGZLiAggLp16/LHH3+gVquZMWMGO3fuxM7OztChCaEjT548zJo1i/Pnz1OnTh3Cw8P5/PPPqV+/vnQCFyILkqRGZKq9e/dSvXp1Ll26hIODA56ennz77beo1fKrKLKuypUrc/z4cX777Tfy5s3L6dOnqV69Ot999x3R0dGGDk8I8S95JxGZQlEUfvnlF9zd3QkNDaVevXqcP3+epk2bGjo0IVLFxMSE4cOHc+3aNbp06UJcXBzTp0+nevXqnDt3ztDhCSGQpEZkgpiYGAYNGsRXX32FoigMGTKEw4cPU7hwYUOHJkSaFS5cmK1bt7Jt2zYcHBy4evUqderU4fvvv5dRGyEMTJIakaGePXtGy5YtWbFiBWq1ml9//ZXFixdLZWCR7XXu3JmrV6/Ss2dP4uPjmTZtGjVq1OD8+fOGDk2IHEuSGpFhrl+/Tu3atfH29iZv3rzs3buXL7/8EpVKZejQhEgXdnZ2rF+/ni1btmBvb4+vry+1atVi0qRJxMbGGjo8IXIcSWpEhkhsGHj79m1KlSrFyZMnadOmjaHDEiJDfPTRR1y9epXu3bsTHx/PlClTqF+/Prdu3TJ0aELkKJLUiHS3Z88emjdvzosXL6hTpw6nT5+W6sDC6Nnb27Nx40Y2btxI/vz5OXv2LG5ubixfvhxFUQwdnhA5giQ1Il398ccfdOrUicjISNq3b4+np6fUnxE5Svfu3bl8+TJNmzYlIiKCQYMG0bVrV0JCQgwdmhBGT5IakS4URWHmzJkMGDCA+Ph4+vXrx/bt28mdO7ehQxMi0xUpUgRPT09mz56NmZkZ27Zto0qVKnh6eho6NCGMmiQ14oMpisK3337L+PHjARg3bhx//PEHZmZmBo5MCMNRq9V88803nDp1ChcXFx4/fkzLli0ZP348cXFxhg5PCKMkSY34IIqiMGbMGGbPng3Azz//zMyZM2WHkxD/qlatGufPn2fo0KEAzJw5k6ZNm/Lw4UMDRyaE8ZGkRrw3jUbDF198wbx58wBYtGgRY8aMMWxQQmRBuXLlYvHixWzatAlra2uOHz9O1apV2bdvX4rXxmsUTgaEsPPiI04GhBCvkUXHQiRHpeSgZflhYWHY2NgQGhpK3rx5DR1OtqbRaBg2bBhLly5FpVKxbNkyBjYpBfvGQdtZ4CztD4TQx9/fn+7du3PhwgUAvv32W6ZOnYqpqWmSc/f7BjJ59zUCQ6O0xwrZWDLJvQJtKhXKtJiFMLTUvn/LSI1IM41Gw+DBg7UJzR9//MHAAQPAazI880v4mnNyZSHSpHTp0vzzzz+MGDEC+G866tGjRzrn7fcNZNja8zoJDUBQaBTD1p5nv29gpsUsRHYhSY1IE0VRGDFihLbtwdq1a+nXrx8EeMHjhE+ePL6Q8L0QQi9LS0sWLlyoMx1VrVo1jh49CiRMOU3efQ19Hw0Sj03efU2mooR4iyQ1ItUURWHs2LEsXrwYlUrF2rVr+fjjjxNGZQ79D1QmCSeqTBK+l9EaId6pW7dunD9/HldXV54+fUrz5s2ZN28eZ+6EJBmheZMCBIZGcebO88wLVohsQJIakWrTpk1jzpw5ACxbtoxevXol3JE4SqPEJ3yvxMtojRCplDgd9fHHHxMfH8/o0aOZNPOXVF37NDz5xEeInEiSGpEqv/76Kz/88AMAc+fOZeDAgQl3vD1Kk0hGa4RItVy5crF27Vp+/fVXTE1N8T6wJ1XXOVhbZnBkQmQvktSIFP3xxx+MGjUKgMmTJ2v/DSQdpUkkozVCpIlKpeLLL7/Ey8uLfLHBxIUFoyga/eeSsAuqVknbzA1SiCxOkhrxTn/99ReDBw8G4KuvvtKO1gD/jdIk+2ukltEaIdKoUaNG+Jw7h92DI4AqSWKTWNZyknsFTNRS5FKIN0lSI5Ll4+ND9+7dtb2cfvrpJ91KwfExEPoI0P9pEjQQ9ijhPCFEqhUuXJjTW5bg+uoM8eG6jTAdbSxZ1KdainVqpGifyImk+J7Q6+7du9StW5eKVsEs71oAp4FrMXNpmfTE0IcQ8Sz5B8ptDzaFMy5QIYyYoijM+3U+381bgSpXPiqULMKu5b9QyLHgO6+Ton3C2KT2/VuSGpHEixcvqF+/PtevX+fKl3ZUyh8DTm4w+DBITychMt2BAwfo0aMHoaGhFC9enF27dlGlShW95yYW7Xv7hT3xLzc1ozxCZDVSUVi8l+joaDp16sT169f5uFbBhIQGZNGvEAbUunVrTp06RenSpbl37x716tVj586dSc6Ton0ip5OkRmglVgv29vYmb968LPu4uBTUEyKLKFeuHKdPn6Z58+ZERETQuXNnfv31V51zztx5LkX7RI4mSY3Q+u2331i+fDlqtRqvpRPI9eKGFNQTIguxtbVl3759DB06FEVRGDVqFKNGjSI+PuHvNLXF+KRonzBWktQIAI4cOaKtPzNr5kxqhO6TgnpCZEFmZmYsWrSIWbNmAQmFMbt27crr169TXYxPivYJYyVJjeDu3bt069aN+Ph4evfuzVed3KSgnhBZmEqlYuzYsWzcuBELCwt27NhB06ZNKZE7jkI2liS3nF+K9gljJ0lNDhcREUGnTp149uwZ1apVY9nSpagOS0E9IbKD7t274+npia2tLWfOnKF+vboMrGYDkCSxkaJ9IieQpCYHUxSFgQMHcunSJRwcHNi+fTtW5iZSUE+IbKRBgwacPHmSUqVKcefOHcb2asUX1axwtNGdYkpt0T4hsjOpU5ODLViwgC+//BJTU1MOHTpEw4YNE+6QgnpCZDvBwcF4eHhw6tQpLCwsWL9hI46VG/A0PAoH64QpJxmhEdmVFN/TQ5Ka//j4+FCvXj1iYmKYN28eI0eONHRIQogPFBkZSc+ePdm1axdqtZply5YxYMAAQ4clxAeT4nsiWWFhYfTo0YOYmBg6derEl19+aeiQhBDpwMrKiq1btzJgwAA0Gg0DBw5k5syZ5KDPriKHk6Qmh1EUhcGDBxMQEEDx4sVZsWKFbpNKIUS2Zmpqyu+//863334LwPjx4/nqq6/QaJJbJyeE8ZCkJodZunQpmzZtwtTUlA0bNpA/f35DhySESGcqlYoZM2bwyy+/ADB37lz69etHbGysgSMTImNlmaTG29sbd3d3nJycUKlU7NixQ+f+/v37o1KpdG5t2rQxTLDZ1KVLl7RrZ2bMmEGdOnUMHJEQIiONHj2aNWvWYGpqytq1a+nYsSMRERGGDkuIDJNlkpqIiAhcXV357bffkj2nTZs2BAYGam/r16/PxAizt6ioKD7++GOio6Np3749Y8aMMXRIQqTKs2fPeP36taHDyLb69OnDrl27sLKyYt++fbRu3ZrQ0FBDhyVEhjA1dACJ2rZtS9u2bd95joWFBY6OjpkUkXH57rvvuHbtGo6OjqxcuRK1Osvks0IkKywsDHt7e1QqlawJ+QBt27bFy8uLdu3aceLECVq0aMH+/fspUKCAoUMTIl1lq3e2I0eO4ODggIuLC8OGDSMkJOSd50dHRxMWFqZzy4mOHj3K3LlzAfj999+xs7MzcERCpM7Zs2eBhAXuBw4cMHA02VvdunU5fPgwdnZ2nDt3jiZNmhAUFGTosIRIV9kmqWnTpg2rV6/Gy8uLWbNmcfToUdq2bavtTqvPjBkzsLGx0d6KFi2aiRFnDeHh4fTv319bPbh9+/aGDkmIVCtdurT237KG7sNVrVoVb29vChUqhK+vL40aNeLBgweGDkuIdJMli++pVCq2b99Op06dkj3n9u3bODs74+npSfPmzfWeEx0dTXR0tPb7sLAwihYtmqOK7w0ZMoRly5ZRokQJLl26lGOet8h6oqOjefDgAffv3+fevXsEBgby/PlzQkJCtF9fvHhBdHQ0MTExxMTEEBYWRmRkpPYxHBwcMDc3195sbGywtbUlf/782q8FChSgSJEiFClShKJFi1KwYEFMTEzeEVnOExAQQPPmzbl37x7FixfHy8sLZ2dnQ4clRLJSW3wvy6ypSatSpUphZ2eHv79/skmNhYUFFhYWmRxZ1rF3716WLVuGSqVi5cqVktCITBEWFsbVq1e5cuWK9nbz5k2CgoI+uAjc06dP03yNqakpTk5OODs7U65cOVxcXChXrhzlypWjaNGiOXJ9mbOzM8eOHaN58+bcunWLhg0b4unpSYUKFQwdmhAfJNsmNQ8fPiQkJIRChaQ5mz4vXrxg0KBBAIwaNYrGjRsbOCJhjBRF4dq1a3h7e+Pt7c3Jkye5d+9esudbWVlRrFgxihcvTuHChSlQoAAFChTA1taWAgUKkD9/fqysrHRGY8qVK6e9/sqVK9pRnKioKF6+fMmLFy94/vy59mtwcDCPHj3i4cOHPH78mLi4OO7fv8/9+/c5fPhwknhcXFyoXr06NWvWpEaNGlSuXBlzc/MM+2+WVRQtWhRvb29atmyJr68vjRs3xsvLiypVqqTq+niNwpk7z6W3lMhSssz006tXr/D39wfAzc2NX375haZNm2Jra4utrS2TJ0/mo48+wtHRkYCAAMaOHUt4eDhXrlxJ9WiM0fd+CjgM+8ZB21kMnb2BpUuX4uLiwoULF7CysjJ0dMIIKIrCpUuXOHLkCN7e3hw7doxnz5I2P3VycqJSpUpUrlyZypUrU758eUqUKKHdyZQWDg4OBAcHv9cOqLi4OJ48ecL9+/e5desWfn5+3Lhxgxs3bnDr1i29xejMzc1xdXWlZs2a1KxZk1q1alG+fHmjrbwdEhJCmzZtOHfuHAUKFMDLywtXV9d3XrPfN5DJu68RGBqlPVbIxpJJ7hWkC7jIENmuoeWRI0do2rRpkuP9+vVj0aJFdOrUiQsXLvDy5UucnJxo1aoVU6dOpWDBgqn+GUad1CgKLGsKjy8Qbl2GvF/7AAlFDbXdt4V4D1FRURw+fJhdu3axe/duHj16pHO/lZUVdevWpVGjRjRs2BBXV9d03SocExPDhg0b6Nq1K7ly5Uq3x42Li+Pu3btcuXKFc+fOcfbsWc6dO8eLFy+SnGtnZ0ejRo1o0qQJjRs3plKlSkY1bfXy5Utat27NmTNnUkxs9vsGMmzted5+40hM+Rb1qSaJjUh32S6pyQxGndT4e8Laj7Tftl4bQbGmn7Js2TIDBiWyq7i4OA4fPsz69evZunWrTjmEXLly0bhxYxo1akSjRo2oUaOG0UzXKIrC7du3dZKcM2fO6CxWhoTRo5YtW9K6dWtatWqVpg9XWVVoaCitWrXizJkz2Nra4uXlRdWqVXXOidcoNJh1SGeE5k0qwNHGkuPjmslUlEhXktToYbRJTeIoTeBlUOKJ0yhceaam+NQAbKW4lkiDgIAAli1bxsqVK3ny5In2uJOTEx4eHnh4eNC0aVMsLS0NGGXmiomJ4dy5cxw9epQjR45w4sSJJK0GqlatioeHBx07dsTNzS3bTlWFhobSunVrTp8+ja2tLZ6enri5uWnvPxkQQq9lp1J8nPWD61DXWV57RPqRpEYPo01q3hql0eqzFUq3yPx4RLYSFxfHrl27WLx4MQcPHtQet7W1pXv37nz88cfUr1/fqKZbPkRMTAwnT57kwIEDHDhwgPPnz+vcX7RoUW2C06RJE8zMzAwU6fsJDQ2lTZs2nDp1ivz58+Pp6Um1atUA2HnxESM3XEzxMX7tWZWOVQtncKQiJ5GkRg+jTGr+HaVRAi+jUv4rRKioTFAVqgKDD0M2/dQoMtbr169ZsWIFc+bM0dmx1Lp1a4YOHUr79u2NZlopIz19+pR9+/axY8cO/v77b50+VQUKFOCjjz6iZ8+eNGrUKNvUywkLC6N169ZJEhsZqRGGIkmNHkaZ1CQ3SpNIRmvEW168eMGCBQtYsGCBdueSnZ0dgwYNYvDgwZQqVcrAEWZfkZGReHp6snPnTnbt2kVwcLD2PkdHR7p160avXr2oU6fOB09RZfSW6rCwMNq0acPJkyfJnz8/R44coWKlyjSYdYig0KgkC4VB1tSIjCNJjR5Gl9QkjtI8vohK70uMGpxcZbRGAAkjM/Pnz2fWrFm8fPkSgJIlSzLvy050MPFG3e4ncE66A1G8n7i4OI4cOcKGDRvYtm2bzq4qFxcX+vfvzyeffELhwmmfptnvG8iPu64RFPbfgl3HvJb86JG+W6rfHLGxt7fn6NGj3IvPx7C1CVNub77qyO4nkZFS+/4tk+TZWXwMhD5KJqEB0EDYo4TzRI4VGxvL4sWLKV26NOPHj+fly5dUqlSJ9evXc9PPDw+To6hDbsFfXyckyiJdmJqa0qJFC37//XeCgoLYs2cPvXv3JleuXPj5+TF+/HiKFStG27Zt2bRpk05Ll3fZ7xvIZ2vP6yQ0AEFhUXy29jz7fQPT7TnkzZuXffv2Ua1aNYKDg2nevDmlLSNY1Kcajja6i8UdbSwloREG90EjNdeuXeP+/fvExOi+aXp4eHxwYBnB6EZqgMfXz9C1XROioqKZOXMmrVq21D0htz3YyIK9nOrYsWMMHz4cX19fAEqUKMHUqVPp1atXwvqOW57w5xvTl723QhmZrsxI4eHhbN68mT/++IPjx49rj9vb2zNo0CCGDBlCiRIl9F4br1Go/r+DvHydtGhgovy5zDj3fct0nf4JCQmhSZMm+Pr6UqxYMby9vSlStJhUFBaZJkOnn27fvk3nzp25cuUKKpVK288lcY74XZ2zDckYk5revXuzbt06GjZsyNGjR7PtVlKRvp48ecLYsWNZvXo1kLBg9ccff2TIkCH/Lf5VFFhYHUIC/ruwgDN87iPTlZnE39+flStXsnLlSm1RQ5VKRfv27Rk+fDitW7fW2XV2wv8ZvX8/neLj/jmoNvVL26VrrE+ePKFx48b4+flRqlQpvL2932vqTIj3kaHTTyNHjqRkyZI8ffqUXLlycfXqVby9valRowZHjhx535hFGp08eZJ169ahUqmYN2+eJDQCRVFYsmQJ5cqVY/Xq1ahUKoYMGYKfnx+ff/657m4mfy/dhAYSvvf3ytygc7DSpUvzv//9j7t377Jt2zZatGiBoijs2bOHdu3aUbp0aebMmUNoaCiQUCcmNVJ7XloULFgQLy8vSpUqxe3bt2nevLlOLSMhsoL3SmpOnjzJlClTsLOzQ61Wo1aradCgATNmzODLL79M7xiFHhqNhlGjRgHw6aefautIiJwrMDCQ9u3b89lnn/Hy5UuqVavGqVOnWLJkSdK2BYoC+8fqf6D9Y2VtTSYzNTWlc+fOHDx4ED8/P0aPHk2+fPm4c+cO33zzDUWLFmXMmDGEhr5M5SNmzP+/woUL4+XlRdGiRfHz86NFixaEhKR/AiXE+3qvpCY+Ph5ra2sgYSvo48ePAShevDh+fn7pF51I1pYtWzhz5gx58uRh2rRphg5HGNjWrVupXLky+/btw9LSkrlz53LmzBlq1aql/wJ9ozSJZLTGoMqWLcsvv/zCo0eP+P3336lYsSLh4eHMnTuXX779LFWPUbdU+k49valEiRJ4eXlRqFAhfH19adeuHeHh4Rn284RIi/dKaipVqsSlS5cAqF27NrNnz+bEiRNMmTJFalxkgri4OH744QcAvvnmGxwdHQ0ckTCU169f079/f7p27UpISAhubm74+PgwatSo5Au9vWuUJpGM1hhcrly5GDhwIFeuXGHfvn00b96c13cvE/86lHcthcyXy4w6GVz4rkyZMnh6emJra8uZM2fo3LlzqndvCZGR3iup+f7779FoNABMmTKFO3fu0LBhQ/766y/mz5+frgGKpNasWcPNmzcpUKCAdgpK5Dy3b9+mbt26rFq1CrVazXfffcepU6eoUKHCuy+Mi4YX9959zot7CecJg1OpVLRp0wZPT08unPehUuQVgGQTm5ldKmfKLqQKFSqwb98+cufOjZeXF717986ym0REzpFuxfeeP39O/vz5s/RiVWPY/RQdHY2Liwv37t3jp59+4uuvvzZ0SMIA/v77b3r27MmLFy9wcHBg06ZNNG7cOHUXx0XDz+Uh8h1rIXIVgDHXwdQifQIW6WrtEV+m7vMjWvVfrRhVVCifVsnDD5+6Z+rrsKenJ+3btycmJoaBAweybNmyLP0+ILKnDNnSrdFo+Omnn9i1axcxMTE0b96cSZMmYWVllS5BZzRjSGp+++03Pv/8c5ycnPD39882/+1F+lAUhdmzZzN+/HgURaF27dps2bKFIkWKpO2BQh9CxLPk75f6RllevEZhv48/qzfv4MCOTYT6nwdFg6urK5MmTaJTp06Zllxs27aNbt26odFoGDt2LLNmzcqUnytyjlS/fytpMGXKFEWtViutWrVSOnbsqFhaWiqffvppWh7CoEJDQxVACQ0NNXQo7yUiIkJxdHRUAOX//u//DB2OyGRxcXHK8OHDFRK2tiiDBw9WoqKiDB2WyAKePn2qfPvtt0qePHm0vx+1atVSDh06lGkx/P7779qfPXPmzEz7uSJnSO37d5qSmtKlSyuLFy/Wfn/w4EHF3NxciY+Pf78oM1l2T2pmzZqlAEqJEiWU6OhoQ4cjMlFUVJTSvXt3BVBUKpUyf/58Q4cksqCQkBDlu+++U3LlyqVNMFq1aqWcO3cuU37+7NmztT936dKlmfIzRc6QIUmNubm5cv/+fZ1jFhYWyoMHD9IeoQFk56QmNDRUsbW1VQBl1apVhg5HZKKwsDClRYsWCqCYmZkpGzZsMHRIIosLCgpSPv/8c8XMzEybZHTv3l3x8/P74MeOi9co//g/U3ZceKj84/9MiYvX6Nw/btw4bfK9ZcuWD/55QihK6t+/07T7KS4uDktL3SZmZmZmxMYm34dEpI9Fixbx/PlzypUrR+/evQ0djkhPAYdhYa2Er28JCQmhWbNmeHp6kjt3bvbu3UuPHj0MEKTITgoWLMiCBQu4ceMGffr0QaVSsWnTJipUqMDQoUPfuxLwft9AGsw6RK9lpxi54SK9lp2iwaxDOk00Z8yYweDBg1EUhd69e+Pt7Z1eT0uIFKVpobBaraZt27ZYWPy3I2L37t00a9aM3Llza49t27YtfaNMJ9l1oXBUVBQlS5YkKCiIlStX0q9fP0OHJNKLosCypvD4Aji5weDD2r5LYWFhNG/enHPnzlGgQAH27dtHzZo1DRywyI4uX77Md999x549ewCwtrbm+++/Z+TIkTqv5++y3zeQYWvPJ6lVnLgU+c0O3fHx8XTt2pUdO3aQL18+jh8/TsWKFdPp2YicKEN6P/Xr1w8HBwdsbGy0tz59+uDk5KRzTKSv1atXExQURNGiRenVq5ehwxHpKcArIaGBhK8BCZV8IyMjcXd31yY0R48elYRGvLcqVaqwe/dujh07Ro0aNQgPD2fcuHFUrFiRnTt3vrOYHyTstJq8+5re5guJxybvvka8JuE7ExMT1q1bR7169Xj58iVt2rTh4cOH6fukhNAj3erUZAfZcaQmPj6ecuXK4e/vz9y5c6XYnjFJHKUJvAxKPKhMoFAVYvodoHOXLvz1119YW1tz+PBhqlevbuhohZHQaDSsWbOGb7/9lqCgIACaN2/OvHnzqFSpkt5rTgaE0GvZqRQfe/3gOtR9o5pxSEgI9evXx8/Pj0qVKnHs2DHy5cuXLs9D5CwZ2qVbZJ5t27bh7++Pra0tgwYNMnQ4Ij0ljtIo/1ZhVeLh8QXmDGvLX3/9hZWVFXv37pWERqQrtVpNv379uHnzJuPHj8fCwgIvLy9cXV0ZPny43gaVT8OjUvXYb59XoEAB9u/fj6OjI76+vtJOQWS4NI3UDBgwIFXnrVix4r0DykjZbaRGURRq1qyJj48PEydOZPLkyYYOSaSXt0dp/qVRVJx7HEuDVbHs3r2b1q1bGzBIkRMkdgLfunUrkNCkeM6cOfTt21dbvO99R2oSXbx4kUaNGhEeHk737t1Zv349arX+z9TxGoUzd57zNDwKB2tLapW0zZS2DyJry5CKwmq1muLFi+Pm5vbOOdjt27enLdpMkmWTmoDDsG8ctJ0Fzk21h728vGjRogVWVlbcv38fO7uM67wrMpm/J6z9KNm7jxT5kiaDpmZiQCKnO3LkCJ9//jlXr14FoHHjxixatIjy5csTr1FoMOsQQaFRetfVqABHG0uOj2uWbALi5eVF27ZtiY2NZfTo0fzyyy9JztnvG8jk3dcIDP1vxKeQjSWT3CtoFyGLnClDpp+GDRtGaGgod+7coWnTpixfvpzt27cnuYk0UBTwmgzP/BK+vpEszpw5E4BBgwZJQmNMFAUO/Y/k/vw0CjTRHJMu2SJTNWnShAsXLjBr1iysrKw4evQorq6ufP/998RERzHJPaFR6tspS+L3k9wrvHNEpXnz5vzxxx8AzJ07N0lSk7i76s2EBiAoNIpha8/rbBsXIjlpSmp+++03AgMDGTt2LLt376Zo0aJ0796dAwcOpLh6XiQjmd0v165dw9PTE7VazZgxYwwYoEh38TEQ+gjQ6L1brQLCHiWcJ0QmMjMzY+zYsVy7do327dsTGxvLtGnTqFSpEqpHl1nUpxqONrq1yhxtLHW2c79L7969mT17NgBff/219kNwWndXCZGcD9r9dO/ePVauXMnq1auJi4vj6tWr5MmTJz3jS1dZbvopmd0vDD7MF19+ycKFC+ncuXOWrfsjPsC/DSXj4uP5bOhQzl+4QIkSJVi1ahXWefJIQ0lhcIqisH37dr788ksePXoEQM+ePZn363xuh6vfe82LoiiMGDGCRYsWaUeE4mxLfdCaHWH8MmX3k1qtRqVSoSgK8fHxKV8gdCWz+yXSdy+rVq0CYPjw4QYMUGQYmyLgVJU5fx5k+T4f/CNyM2PlX1iXbQBOVSWhEQanUqno0qUL169fZ/To0ajVajZs2EDlShV5fPEIHasWpq5zgTQv4lWpVMyfP582bdpo6zFdu5O6Gjap3YUlcq40JzXR0dGsX7+eli1bUrZsWa5cucLChQu5f/9+lh6lyXIS11WoTHSPq0wI3zWO8PBwypYtS7NmzQwTn8hwly9fZuLEiQAsWLAAFxcXA0ckRFLW1tb88ssvnD59mkqVKhEcHEzXrl3p0aMHwcHB7/WYpqambNy4kSpVqvDkyRPmTJ2YquscrC1TPknkaGlKaoYPH06hQoWYOXMmHTp04MGDB2zevJl27doluz1PJOPtUZpESjwOsQ9p5WzCsGHD5L+rkYqJiaFv377Exsbi4eFB3759DR2SEO9Uo0YNzp07x3fffYeJiQmbNm2iYsWKbNmy5b0eL2/evOzZs4dChQpx49geTGPCkyxCTqQiYRdUrZK27x2/yBnSvKW7WLFiuLm5aesX6JNV14BkmTU12n4/l9C3WDReo3DhiYLzjHvkt5U/YmP0/fffM23aNOzs7PD19aVgwYKGDkmIVPPx8aF///74+voC0K1bN3777Tfs7e3f67EaNWqEUsQVh04TEpY0vHG/vt5SIufJkDU1ffv2pWnTpuTLl0+n19PbN5GCFHa/mKhVlHGwIn/e3HrvF9nb6dOnmTFjBgCLFy+WhEZkO9WrV+fcuXN8//33mJiYsHnzZipUqMDmzZvf67HWr19P1K1TPN0xndzqWJ3707K7Sgjp/WQo/+5+edPzFy+0xalWbt5LlQZtDBScyChxcXG4ubnh6+tL7969Wbt2raFDEuKD+Pj48Omnn3LlyhUg4cPvwoULsba2TtPjzJs3j9GjR4NKzczlWyjnVksqCgst6f2U1f27++XN2x/7z3PmQTSmRatLQmOkVqxYga+vL7a2tsyfP9/Q4YgsKjIykqdPn/Lo0SPWrl3L9evXCQwMJCQkhFevXmWp3aaJozbfffcdarWa1atX4+bmxunTp9P0OCNHjmTEiBGgaJg8ojfFlKfvtbtK5GwyUpOFuLm5cfHiRRYtWsRnn31m6HBEOgsPD6d06dI8ffqUX3/9lS+//NLQIQkDUBSFe/fuceHCBQICArh9+za3b9/m8ePHPH/+nOfPnxMZGZni49ja2uLg4ICDgwP29vYULFiQUqVKUbp0aUqXLo2zszOWlu+/W+h9ejAdO3aMPn36cP/+fUxMTJg8eTLffvstJiYm77wuUVxcHB06dODAgQMUKVKEs2fP4ujo+N7PQRiPDOn9lN1l5aTG19eXypUrY2ZmRmBgIAUKSIEpY5O4OLhMmTL4+vpibm5u6JBEJoiMjMTb25t//vmHs2fPcvbsWZ49e5byhSRsfY6Li9NuzEjLy7VKpaJIkSJUrVqV6tWrU61aNapXr45TpJ/eXnNv+pAeTC9fvuSzzz5j48aNADRq1Ii1a9dStGjRVMX98uVL6tSpg5+fH3Xq1OHw4cMflJwJ4yBJjR5ZOakZN24cs2fPplOnTtI/ywg9ePCAsmXLEhUVxfbt2+nUqZOhQxIZKCQkhC1btrBjxw6OHDlCVJRu0TgzMzMqV66Mi4sLpUqVomTJkhQtWpQCBQpga2tLgQIFyJMnT5KSDnFxcURHRxMREUFwcDBPnz7Vfg0MDCQgIAB/f39u3bpFWFiY3tjOD7PBzUHhMY4ENFtGrdq1sbCw0N6f2IPp7TeGtOxCiovX8L8lG1jw+yoingVi9eohy5YupWvXru+8LnF06PKtu3z/1Rc8u36aT/r0ZtWqVe/ccSuMX4YkNRMnTqRjx45Ur149XYLMbFk1qdFoNBQrVoxHjx6xdetWunTpYuiQRFok02X9TZ988glr166lUaNGHDlyRF6gjVB8fDz79u1j8eLFHDhwgLi4OO19RYoU4fN2lRhc5BZPq39NiWb9M3T0QVEUQkJCuH79OhcuXMDHxwcfHx+KRt9kX28r7Xmt10bg/ciMOnXq0KpVK1q2as3nB14QFKa/cm9qunHrG+WJCwvmuddSejUoz/z588mdO+nOTr3XhT/juecSfhzYibFjx77nfw1hDDIkqRkwYAB79uzB3Nwcd3d3PDw8aN68ebYZRs+qSc2hQ4do3rw5+fLlIygoSOdTk8jitDWHLoCTGww+DG8lLNevX6dChYQOx2fPnqVGjRqGiFRkkNevX7NkyRIWLlzI7du3tcfd3Nzo0aMHHTp0oEL58qh+b/bO35MMpyjEL2mMOugKKjTEK3A5WEW1RaHaUyyKVsbx4xkpPlRyPZiSG+UBBUWB4B3TKWUWxpYtWyhXrlyqr3u2YwbrfxqHh4dHap+tMDIZsvtpxYoVBAUFsX79eqytrRk1ahR2dnZ89NFHrF69mufPn39w4DnRmjVrAOjevbskNNlNMl3W3zRv3jwAOnbsKAmNEYmKimL+/PmUKlWKMWPGcPv2bfLnz8/XX3/N9evXOX/+POPGjaNixYqobh9K8fckwwV4YRJ0CdW/9bFMVODmoHD/0B8sXLgQd3d3chVI3aJcfT2Y3tVpG1SoVCrsWw3j6rXr1KxZk02bNqXuOiB/88H07tNHu21ciOSkeUu3Wq2mYcOGzJ49Gz8/P06fPk3t2rVZsmQJTk5ONGrUiDlz5mi7uop3e/36tbbM+CeffGLgaESavN2/S2WS8P0bg5/Pnj1j9erVAIwZM8YQUYoMcPDgQSpVqsTIkSN58uQJJUuWZOnSpTx8+JCffvpJZxQiNb8nGe4dveaK3lrJiOHD2bVrFzvWr0rVw+WzSDrKdObOc52pI31UuW2p3aE3r169okePHowcOZLjfkHvvk6lwjSvPbH5SuDh4fHe/aZEzvDBdWrKly/P2LFjOXHiBA8ePKBfv34cO3aM9evXp0d8Rm/Pnj28evWKEiVKUK9ePUOHI9IimS7rb34KX7x4MVFRUVSvXp2GDRsaKFCRXoKDg+nduzetWrWihOYufl/asGf+V/j5+TF48GBy5cqV9KJU/J5kuHf0mnszlvplHSlkY5lsDyZF0RAXFox7nfL0798fLy8vNJqEkZ/UdtAeN2ka48ePB2DZ/rP0X3Y8Vdc5OZfj7t27dO3aldjY2JQvEDlSuhbfs7e3Z+DAgezcuZOvv/46PR/aaCXudOrevbs0r8xO3vHJN/FTeHR0NAsXLgQSRmlkcXD2duzYMVxdXVm3bh1qtZoVHxejbH6F9hbnMDM11X9RKn5PMlxiDMm+3Ku1sZioVUxyT1j/pe+3VaVSYXp5O+FhYaxatYoWLVpQpkwZZsyYgWns61SF45gvF9OnT2fqH7ux7zQBjZlVyhcB038Yh7W1Nd7e3nz11VepukbkPPIuakAxMTH89ddfALLFN7tJxSffDRs28OTJEwoXLky3bt0ME6f4YIqiMGfOHJo2bUpgYCAVKlTg+u6FFDMNSTjhXaMuqRwhyVAp9JoDDYQ9SjgPaFOpEIv6VMPRRnd3ViEbSxb3qc5t7+0cO3aMoUOHkjdvXm7fvs2ECRPoVL8SZrGvkg3jzU7b8RqFvYFWqFSqFJP9xOu6NammXX+4YMECVq5cmbrnL3IUqVNjQH///TetW7emYMGCPH78WEZqsosUuqyDGpxcqb4kjPPnLzBz5kzGjRuX2VGKdKDRaBg1ahQLFiwAoHfv3ixZvJjc6zpA4OWE5ERlAoWqJN3RlMrfk0zZCaWn15yO3PZgU1jnUGoqCkdERLB582aWLVvGP//8g1XZuth3mpAwyvPGc3q7xs3JgBB6LTuVqtBV/17XsoIjZ+485//+WMvm1csh2J9j3kepWbNmqh5HZG+pff9OZsxUZIYdO3YA4OHhIQlNdpKKT75xz+/je+kOpqamDBo0KDOjE+kh4DDKvrH8ctWBBYv/QqVS8euvv/L555+jenPHG+iOupRu8d/xtIyQmGbwrkebIgm3NDBRq/Ru235T7ty56d+/P/379+fq1assXbqUtQfmYlXvE0zz2mvPs89jxpROlbVF+1K7/iZfLjNmdqkMQINZh/5dUFwax49nEBcWTJeRUzi3/XfpdC+0PiipSdzhVLhw4RTOFG/TaDTs2rULkKmnbMfUAoYcfucn3yWrNhETP4NWrZpJy4vsRlFQvCajenaTRnHXUKvVrFq1ij59+uiukXlzSilxjYxz8/9GKFLxe0Ju+4xPaDJJxYoV+fXXX5kaFsbyFX8wf/06noZHE//qBY+C/Fh7oxtFx4+nYsWKOFinrvDgAJeEiQR9NWxMre1QGg6lw2ff8c+mRZiZmaXzMxLZ0XtNP504cULbtAzAzs6O/v37891332WJaZ3kZKXpp7Nnz1KrVi1y587Ns2fPpLeJkalRowY+Pj4sWbKEIUOGGDockRb+nrD2I+23x0uMoUH/SXrvS6LPVt3RmhwsPj6e3bt3M2/ePI4ePao93rlzZ76dMIGRB0MJCo3SX59GUYgLf0bQ759R7uuNvNLo//ytKBo0ka9ooL7B2l9+lI7eRixDiu8lGjp0KOXLl+fs2bP4+fnx008/4enpSbVq1aQ+TSrt3LkTgDZt2khCk90EHIaFtRK+6nH37l18fHxQq9UyCpfdKArPt3xFnCbhrVaDigYxhxNGaNKwi0iAiYkJnTp14siRI/j4+NC1a1dUKhXbt2+nds2amFzajkLSXVYqEnZZVYm/ialj2WQTGgCVSo1JrryctKyF28Td7PcNzMinJLKB90pqAgICmDdvHtWqVaN06dL07duXc+fO4ebmxqhRo9I5ROOUuJ5G3vSyGUUBr8nwzC/hq543sK1btwIJ3YkdHBwyO0LxAe56rcA26i6m/37iV6P8t14mjbuIxH+qVavG5s2b8fX1pU+fPqjVak5s/I2n26ehjg7XOdfRxpJFfaqxd/E0+g39PNU/IzRWxWdrz0tik8O91/RTtWrVmD9/Pg0aNNA5fv36dWrVqkV4eHgyVxpWVpl+evDgAcWKFUOtVhMcHIytra3BYhFp9Pb0g57phnr16nHy5EkWLFjA55+n/kVZGJYmPp4bXztR1jpSm9QAurubwh6leReRSCogIIBZs2axcuVKYuPisShSkZqNmjOs/8f0aFZDO42Ull1SkLD93jGvBf+MbyFTUUYmQ6ef+vfvzxdffMGDBw90jhs6WcguDh06BCSsu5CEJhtJRbn7hw8fcvLkSQDptp7N7P9tLBVsonQTGtDd3WRTBJyqJn+ThCZVnJ2dWbp0Kf7+/gwaOIDYR1c5/uc8+rSqTb++n2gbg9YqaUshm9RPz6tUKp6Ex3AqQFop5FTvldSMGjWKS5cuUaZMGT7++GNmz57NjBkzGDhwILNnz07vGI1OYlLTrFkzA0ci0iQV5e4PHDgAQN26dXFycjJElOI9PH3yBMdrvxOvSW7gWtbLZIRixYqxbNkyrl27Rrdu3VAUhT///BMXFxdGjBhByLNgJrlXSFhnk4bHXbxK2vTkVO+V1AQGBrJv3z6mTJkCwMqVK5k4cSK3bt1i9uzZ9OnTh9mzZ7N///50DdYYKIoiSU12lMpy9+fOnQNIMjUrsra5P8/GKY/mHVMWsl4mI7m4uLBp0yZ8fHxo3bo1cXFx/N///R9lypTh6oF1zO9ZJUmF43fZvHq5tlo7JBQSPBkQws6LjzgZEPKO5FVkd+lWUTgqKoorV65w8eJFLl26xMWLF/H19eXly5fp8fDpIiusqfH396dMmTKYmZnx8uVL/Q3wRNaTyq28NWvW5Ny5c2zcuJHu3btnXnwi9QIOw75x0HYWODclIiKCIkWKkEcTyp9Lf6VRco1HZb1Mpjl69Chjxozh/PnzAJQpU4Y5P/9CgXK1GbHuAi8jY9A3dqMCLDSR+M3pQf58Nvj4+OAXYcnk3dd0OoHnszLj0/ol+LxZGVl7k02k9v1b2iRklLdeOBMtXbqUoUOH0qhRI53aDSILS2W5++i++7HOm5fY2FgCAgIoVapUZkcqUqL9f3kBnNxg8GGWLlvG0KFDcXZ25ubNm1LdO4uIj49n1apVjB8/nqdPnwLQqlUruo76H9OPPkV5a0N44r/m96jCtGHdOH36NBXb9CbCtZf+Wjj8V7E4sdKxyLoydKGwSME7tv3K1FM2lMqtvFcvXyA2Npb8+fNTsmTJzIxQpNabLQ4eX0Dx92L+/PkAjBgxQhKaLMTExIQBAwZw69Ytxo4di7m5OX///TfD3OtSPeoCBa3fqsQc+ZIprYvh7laUzZs3Y2dvz8uSLXjX5/aXr2MZJtvAjYqM1GSEZLb9KopCwYIFCQ4Oxtvbm4bJDXOLrCcVDQGXbtzL0KFDadmyJX///XfmxSZSJ3GU5o1GlK/yOmM95hy5c+fm4cOH5MuXz9BRimT4+/vz9ddfawuXFnR05Isp87F2KMzsyd/x6KI3ToUc+euvv3B1dWXhpr+Zcz42VY9dyMaS4+OayVRUFiYjNYbyjm2/V69eJTg4GCsrK2rXrm3YOEXapGIrb+Ii4Ro1ahgqSvEuenav5Qm9SStnE1q2bCkJTRZXunRpduzYwcGDB3FxceFJUBDfD+nOvsVT2TB/KhXKl+Px48c0bNgQT09PipatmOrHDgyN4syd5xkYvcgsktSkt3ds+02sX1K3bl3Mzc0NGKTICJLUZGHJ7F6LV2BqUwsaN2pkoMBEWrVo0YJLly4xefJkzM3N2b9/P61ataJbt240bNiQ8PBw2rZty6VTx9L0uPt8A2VnlBGQpCY9pbDt18cn4U2vZs2aBghOZCRFUfD19QXAzc3NwNGIJN7+sPEvExXUKmyKe4XcBgpMvA8LCwsmTpzI5cuXadq0KZGRkUyePJkXL17QokUL4uLimDS8N3lM4lL9mKtP3qPXslM0mHVI1thkY5LUpKdkXjgTR2tM73kD8kneGIWGhhIbmzB/L0X3shhtI0r94hUodXedFNbLhlxcXPDy8mLVqlUUKFAAX19fvLy8cHV1BUXD7S2z0vz/NSg0ShYPZ2NZJqnx9vbG3d0dJycnVCqVtuFjIkVRmDhxIoUKFcLKyooWLVpw69YtwwSrTwodfBVU9Cv6EIDq1atnYmAiM4SEhACQO3duLCwsUjhbZKr4GHh+O9m7TVSgCpfCetmVSqWib9++3Lhxg08//RRFUbh06RK5cuUi8uZJnu6YjhKV+n6EiSnQ5N3XZCoqG8oySU1ERASurq789ttveu+fPXs28+fPZ/HixZw+fZrcuXPTunVroqKi9J6f6VLY9qtCoYg1ONrlp0SJEpkamsh4z54l7IwqUKCAgSMRSZiYJyz01r7cqcHOhT9zD6HakldMuNcABh8BU0lGszM7OztWrFjBgQMHKFKkCK9fvwYg8uZJ7s/vTeTZzZjEpe79QkEWD2dXpoYOIFHbtm1p27at3vsURWHevHl8//33dOzYEYDVq1dTsGBBduzYQc+ePfVeFx0dTXR0tPb7sLCw9A88kakFDDmc7Lbfrdu2MvKXaVSp2xKVSrYNZivJFFJ8U+JIjZ2dXWZGJlIjwAueXH3jgAae+RETV5YLQRoaWDtLpWAj0qpVK65cucKoUaNYtWpVwkFFw9NDqwg7uRm78rWJsC1L3uruKT7W0/As8qFZpFqWGal5lzt37hAUFESLFi20x2xsbKhdu7Z2R5E+M2bMwMbGRnsrWrRoxgb6jm2/+y8F8Shckamn7OYdhRTflJjUyEhNFvOOxftNNScASUSNUb58+Vi5ciU7d+6kYMGC2uNRka8JvOSNEpm6D7gO1qnvNyWyhmyR1AQFBQHo/HImfp94nz7jx48nNDRUe3vw4EGGxvkuPj4+gCwSznbeqkD7ZkfuN8n0Uxb1jsX7Jcyf08rZRJIaI+bh4YGvry89evTQHovXKOSu0vqdlYYBHPNaUKukbUaHKNJZtkhq3peFhQV58+bVuRlCYrNPkEXC2co7Cim+TUZqsqAUFu9r/q1REx+X+m2/Ivuxs7Njw4YNbNy4EUtLSyyKVMQ0r12KywB61SomFYazoWyR1Dg6OgLw5MkTneNPnjzR3peV+fv7ExcXh42NDcWKFTN0OCK13lFI8W0yUpMFpbB4X62ConnVPAl8mLlxCYPo3r07/v7+ODmXS9X5JeykdlF2lC2SmpIlS+Lo6IiX139vJmFhYZw+fZq6desaMLLUSdx6XqZMGVkknF2kUEjx7dEa+f+aBSUu3h9yVO9tpcUAai6L4P6j5KewhXEpXLgwa5YuTNW5sp4me8oyu59evXqFv7+/9vs7d+5w8eJFbG1tKVasGKNGjeJ///sfZcqUoWTJkvzwww84OTnRqVMnwwWdSonPq0yZMgaORKTam2tp3vTmaE3p/xauJ47QJE5DiSzCpsi/27mTMi3qy6NwhcePH2dyUMKQ6jjbU8jGksDQSCDphxEV4GhjKetpsqksM1Jz7tw53NzctCXmx4wZg5ubGxMnTgRg7NixfPHFFwwZMoSaNWvy6tUr9u/fj6Vl1s+m3xypEdlACmsxQJ1ktCZxsWniNJTI+hIrPz969MjAkYjMZKJWMcm9AipU/FdqL4GiaFCASe4VZD1NNpVlkpomTZqgKEqS28qVK4GE4f0pU6YQFBREVFQUnp6elC1b1rBBp1JiUlO6dGkDRyJSJYW1GKCBMN0KtJLUZD+JSY2M1OQ8bSoVYlGfahSysdI5Hh8eQvD26ViF3DRQZOJDZZnpJ2Mm00/ZTAqFFAHIba9TgTYxqZHpp+yjSJEiqFQqwsLCePjwIUWK6J+mEsapTaVCtKzgyJk7zznhc4UfvvmSqAdXQdHQpEkT5syZw5gxY2S9XDaTZUZqjNXr1695+DBhd4UkNdnIOwop4lQ1SQXaxDU1MlKTfeTJk4c6deoA8Ndffxk4GmEIJmoVdZ0L8HX3JniuW4Tqjemor7/+mq5du6bYiideo3AyIISdFx9xMiBE+kUZmIzUZLCAgAAgocKlra0sPDNWMv2UTbzV8qJDhw6cPHmSPXv2MGTIEENHJwyofv367N+/n9atW2uPbdu2DTs7OyZMmEClSpUoVqwY5cuX1zat3e8byOTd1wgM/S/xccxrSa9axShhlwsH64QFx7I+J/OolJTKKhqRsLAwbGxsCA0NzbRCfNu3b6dLly7UrFmTM2fOZMrPFJkvPDxc+zsVERFBrly5DByRSEJRYFnThN1rTm4w+DCXr1zB1dUVKysrQkJCsLKySvlxhNF5+vQpc+fOZfHixbx8+fKd55qamlKpUiVKNOjEhdzVqa/25UfT1fwY15cTmspJzi9kY8kk9wq0qVSIeI3CmTvPeRoeJQlPGqX2/VtGajJYYmsG6cxt3PLkyYO5uTkxMTE8e/ZMiixmRXpaXlSu3JyiRYvy4MEDDh06RPv27Q0bo8h0p06dwsPDg+DgYCBhVD1fvnzcvXsXALVaTdPiKua3tWTsYdh7PYKLly4TXHckJorCWNONlFE/YqzpRjrGVOLtbeJBoVEMW3ueIY1KsutSoM6ojm1uMzpXLUyLCo6S4KQTWVOTwRKnI+zt7Q0cichIKpVKpqCysmRaXqhAm8js2bPHcPEJg3jy5Anu7u4EBwdTuXJldu7cybNnz7hz5w5fd6nB1eG5aVpcxZy21lSwN2FifYUuXTrzvyXrMc1rT2OTK7iqbwPgqr5NI/XlJD9D+fe2xPuOTkID8DwiluUn7tJr2SkazDrEft/ATHjWxk2SmgyW+AYnTfOMX+Ji4cRPfCILeUfLiw4dOgCwc+dOoqOjDRikyGwLFizg2bNnVK5cmZMnT+Lh4YGJiQkoCjNaWFLB3oTf2llS1T6hvEOtwqa8urSLjbv2AwpfmW4mTkl4G41T1Hxlupm3a9+kVuKIjiQ2H0aSmgwmSU3OkViH6Pz58waOROhIoeVF82bNcHJyIjAwkGXLlhkmRmEQBw4cAGDcuHHkzv1Gr6cAL0yfJIy6uNiZaHc0xSswo2Vu7gY+o5H6Mq7q25iqEhIeU5Um2dGa1EhMhSbvvkZMnEZ2VL0nSWoymCQ1RiLgMCyslfA1GU2bNgXg8OHkzxEG8PYoTaJ/R2ssH53g+++/B2DatGm8fv3aAEEKQwgMTBgVKVOmjO7f+FtJcOJaFxMVVCsInVrV0RmlSfShozUKEBgaRZ0ZXvRadoqRGy7Sa9kpqvx4gGFrfTjh/0wSnBRIUpPBJKkxAooCXpPhmV/C12Q2DCYmNSdOnCAmJkbvOSKTpbLlxcABAyhRogRBQUEsXJi6hoci+3N0dAQgKDAQdgxL+BvfMkB/EvyvOEXNj6ardEZpEn3oaE2i5xG6rx8RMfHs8w2i9++nqf6/gzJF9Q6S1GQwSWqMgJ5dM/pUrFgRe3t7Xr9+Ldv3s4pUtrwwN4Eff/wRgFmzZhEaGppZEQoDqlChAgCPj6+G8H8Thcjn77zGVKWhlDqI5AZMNIrqg0ZrUvLydSyfydqbZElSkx6SmZpQFEWSmuwumV0z+kZrVCoVTZo0AWQKKstIbHkx5Gjyt8FHwNSCPn36UK5cOZ4/f87cuXMNHbnIBN27dwfAQ/FM03WKAsntvlarFAqpQjAnDkioUzO0UUk9/cA/zOTd12QqSg8pvveh9BT04t9eIYk/DxLaJUhhr2zI3xPWfpT0eJ+tULpFksOLFi1i+PDhNG3alEOHDmVCgCI9bd68me7du2Ntbc2dO3e0O9qEcYqNjaVfAyfWtU37dPFzJQ+fxnxD3Bvl3l4eXUWFAiZ88t2vvLAsolNgT1/14TcpigaVKm3jDOsH16Guc874HZXie5lF39TEv292L168AMDS0lISmuzozVGaN+fXE0drnJtrE9hEietq/vnnH6KiorC0tMzMiMUH+uijj6hatSoXL15k1qxZzJ4929AhiQxkZmrKsg5WCdOUqfBmIhOi5CUI3YQi7LkJ/5w8zbFb/fjrr78oVeq/HnFvNtA8eC2IHRcf66yd0bwOwyR3vjTF/zT83X2pciKZfvoQKUxNxMcnvBGamJgk9wgiK0th14y+tTUuLi44OjoSHR3NqVOnMilQkV7UajX/+9//AFi4cCGPHz82cEQiQ906SO74d6yfaj0LhhzlYtudtI+eRrvoGVxSynBVKZkkoQH49ddfKVKkCH5+ftSpUyfJ2rrEBpoT3Sty9rsWrB9ch197VmX94NrMrGeKEvEcRUlu/VdSDtbyoeltktR8iHcU9ALQaBJ+OdVq+c+c7aRy18zba2tUKpVs7c7m2rVrR7169YiMjKRv377aDyfCyCgK7Pry3eecmAeFXLlnkXwi8yYz6wKcPn0aNzc3goODad68eZJp6MSu3nsuJyTMHao4UdfZjt69evBrv4aoVCpSsyqkkE3C1JbQJe+27yuFgl4oivYXU6WSfh7ZTip3zegbtm7WrBkAXl76d0mJLCzgMKrfarPuf4PJnTs3Xl5eTJo0ydBRiYwQGwWvUthB9CoQYqNSPSIyb+YU7Ozs8Pb2pkWLFrx69Yp27dqxc+dOIKGrd4NZh3Rq0LzZHqFTjRIs7lOdvBYpj+5Pcq8gvaL0kIXC7yu5BaSJ+mzlpqYYLi4u2NjYpNj5VWRBoQ8h4h19nHLbg03hJIfv379PiRIlUBSFa9euUb58+QwMUqSbtxb9r88zhI979wYS+kJJs0sjExUGM4vx7q3XKvj2PvHm1jSYdYig0KhkzlaIDw/h4aIB9P2kDytXriQmJoZevXqxfft2TExM+PrXtWx8YJ3k+sS0ZFGfarSpVAhIGM1Z4HWT3w75EavoJjj5c5kxo0tl7bk5RWrfv2Wk5n2kcmpC8++wtUw/ZVM2RcCpavI3PQkNQLFixfDw8AAS5thFNvHWov9etRwYMWIEAH369OHOnTsGDE6kO8u8MOggdFmmvcV3WsK8B5X4eOtrPtkRhWfJ8WCZFxO1iknuCTVt3h4bURQNigJ9KlhiolaxevVq/ve//2FhYcGmTZvo168f8RqFtVejUPSkRG+2R0jcom2iVjGqpQs3prkzr2Mp7J6e4+WJDQStn0C+I7OokDc2A//DZG/ybvs+0jg1IdNPOc/o0aMBWL16NSEhIQaORqQomUX/P8+ZQ61atXj58iVdu3YlKkp2mxiVIjWhSnftzaRqT0YsPoKqSnfWXoqh9YAJrFq1CkjYvbSoTzUcbXSnokyiwwneMZ1zO35nwYIFAEycOJF169ZhamrKihUr+HjkD5jmtSdpSpQgsT3CmTu6hf9M1Co61S3P2eUTmf1JQ0ye+XPk8CGqVKnCn3/+maq1NzmNTD+9r1RMTVx7FErFihUpUKCAtgifyBkURaF69epcuHCBadOmMWHCBEOHJN7lHfWI7puXpVq1aoSEhDBkyBCWLFmS+fGJTBUfH89nn33G77//DsBvv/3G8OHDE+7TKJy585yn4QlrbfLFPqN6NTciIyNZsGAB9+7dY86cOZibm3PQ0xMzpwrsuxLI6lP3Uvy5v/asSseq+keAAfz9/fnkk0+0Oyu7d+/O4sWLyZ8/fzo866xNpp8yWiqmJhJ3P8lIjRFLppq0SqXSjtYsXLhQekFlZSks+i9WtCh//vknKpWKpUuXsnr1asPEKTKNiYkJS5cuZeTIkQCMGDGC4cOHEx0drd2W3bFqYeo6F6B8ORd++uknAMaOHcuAAQPo0qULJiWq02fjHXotO5WqhAZS3qJdunRpjh07xpQpUzAxMWHTpk24ublJW5Y3SFKTgWT3k5FLodFljx49KFSoEIGBgWzatMlAQYoUpaIeUevWrZk4cSIAQ4cOlWrROYBKpWLu3LlMnToVlUrFokWLaNy4MQ8fPkxy7vDhw2nVqhWRkZH069eP3uN+wqHTBFS5U7flWkXqt2ibmpryww8/cPLkSZydnbl37x4NGjTg119/lekoJKnJULly5QIgIiLCwJGIDJFCo0tzc3PtQtO5c+fKC05WlIZ6RD/88AMdOnQgKiqKDh064OmZtn5BIvtRqVR8//337Nmzh3z58nH69GmqVauWpAaVSqVixYoV5MuXj7PnfPhh+yVQqVL3gfbf14W0btGuWbMmPj4+fPTRR8TGxjJq1Ci6du2a43faSlKTgRKbWL5+/ZrXr18bOBqRrlLZ6HLo0KFYWVlx/vx5jh07ZoBAxTulYdG/iYkJmzdvpn379kRGRuLu7s6BAwcyM1phIO3atcPHx4eqVasSHBxMixYtmD59uk5hxsKFC7NgwQIsilQkQjFP9WPHhT+jzNOjtChnn+a4bGxs2Lx5M/Pnz8fMzIxt27Zp1/LlVLJQOAMpioKFhQWxsbHcu3ePYsWKZfjPFJkkDY0uP/vsM5YsWULHjh3ZsWNH5sQnUi+N9Yiio6Pp1q0bu3fvxsLCgu3bt9O2bdtMCFQYWmRkJMOGDdPuiKpTpw4rV67ExcUFSHjNb9L3a+4VbpbiY/WtW5y8obeZMLg7sTHR9OjRg7Vr12Jq+n4tGc+ePUv37t25e/cuVlZW/P7773z88cc657y9yDmx2WZ2kNr3b0lqMljhwoV5/Pgx586do3r16pnyM0UGSyzSFng5aaPLQlV0OrUD3Lhxg/Lly6NSqThx4gR169Y1QNAiPcXExNCjRw927NiBubk5W7dupUOHDoYOS2QCRVFYtWoVI0eOJCwsDEtLS9ZMHshHuU6jajubvc+LMGLrzRQfJ7HD9s6dO+nWrRuxsbEfnNi8ePGC3r17s2/fPgC++uorZs6ciampqd4u4YVsLJnkXiFbFPKT3U9ZROIUlGzpNiJpbHRZrlw5+vbti6IoDBo0iOjo6EwMVmQEc3NzNm3axEcffURMTAxdunTRlsIXxk2lUtG/f398fX1p3bo1UVFRFLv1B6pnN4naO4E21ZyxMdMk25hSUTTY5TLRLgru2LEjW7ZswczMjI0bN9K7d2/i4uKA//pE7bz4iJMBIdrifMnJnz8/u3fvZvz48QD8/PPPtGnThs0nbzFs7XmdhAYgKDSKYWvPa9s0GANJajKYvX3CPGlwcLCBIxHp4j0bXf7yyy84ODhw7do1pk+fnuFhioxndv84mxvfY2r/ZsTGxtK1a1e2bdtm6LBEJilatCj79u3jrwVfUatwwsiK5fNrrPmxH5M7VkaFKmlioyiAimd/L+J5yH8fdD08PNi6dStmZmZs2rSJ3r17s/fSo3f2iUqOiYkJ06dPZ9OmTeTKlQuvQ4f5Zv1pve0d9FUzzu4kqclgiSM1ktQYifdsdFmgQAEWLlwIwPTp07l8+XLGxiky1r/b+VUhN/muTjy9evUkLi6Obt26MWvWLNnplkOogLbm51BUCW+lcRqFCkHbGPlRI7oXDkOJeKFzfsG8Fpie+oN7J3bRvXt3YmP/a3fg7u6uTWx2X3zAiHUXPmhkpVu3bpw6dYriNZpBruSL8yVXzTi7er+JO5FqiSM1Mv1kJEwtYMjhlBeWmlokOdy1a1c6derEjh07GDhwICdPnnzvuXNhYG9s51cFXmDNjxPIk8eaZcuW8e2333L27Fn++OMPrK2tDRyoyFD//h4krqAzVauoVdiUcmaP+Wlkb6q4VsXveTxmeQuwYPY0PmlTF7/Ohald+2+OHDnCmDFjtK0VICGx2bxlKyMOPCdhTEdX4rHJu6/RsoJjiot8K1euzJSffmXivpT7lj0NN44WIDJSk8Fk+skIvWejS5VKxW+//YaNjQ3nzp2TZpfZlZ7t/CZHprNk8WKWLFmCmZkZW7dupXbt2vj5+Rk0VJGBkqlErahMWN67JObm5ly+dJGYh768unqEGaM/JS42hgoVKvDnn38CCdXGE1sxJHKoWA9Ta7tka9ykdWSlTJGCqTovpWrG2YUkNRksMakJCgoycCQiK3BycuLnn38G4IcffsDf39/AEYk0e3uh+L8LxFW3DzFkyBC8vb1xcnLi+vXr1KpVi127dhk2XpExktkwoFLiKaJ6yp2Dy3B3d9dORd68eZMWLVoQHh6Oh4cHU6dOBRKqEf/zzz/a61M7YpLa82qVtKWQjWUyrTT/q2ZcvXj+NC1KzqokqclgifULrl69auBIRFYxYMAAmjdvTmRkJIMHD9b2CBPZQAp9olAU6tSpw/nz52nYsCFhYWF07NiRiRMn6hRqE9lcKjYMON1Ywa6dOzl69ChlypQB4Pjx45QoUYIdv4xigs0OJvdrSmxsLF26dNG2X0jtiElqzzNRq5jkXgHQ0yP834TLw7UQjX86nOZFyVmRJDUZzNXVFYDbt28TFhZm4GhEVpDYGDFXrlwcOXKEZcuWGTokkVqp3M5fsGBBvLy8tA0Rp06diru7Oy9evHj7EUV2lIYNA40aNcLPz09bn+r58+c4XV+GOuQWQ0oFUqlSRZ48eUKnTp2IjIxM9chKavpEJWpTqRCL+lTD0UY3EYoLf4bZbW+WeN8xmu3eUnwvExQpUoRHjx5x/Phx6tevn2k/V2Rt8+bNY/To0VhaWuLl5UW9evUMHZJ4l8Sii48vof/NTA1OrkmKL65du5YhQ4YQGRlJ8eLFWbZsGS1btsy0sEUGSWMl6qdPn1K+fHlq5A/lQJ/c2uO99lqw+/prIiIi+HTAAAZPmI3ntSCWn7ib5CEVRYNKpWJxn+rvVTDvzYrCr4If8XX/zph0nIaptZ3O72wiFeBoY8nxcc0MXnlYiu9lIYmjNZcuXTJwJMKgAg7DwloJX4EvvvhC2yDR3d1dFpVmde+5nb9Pnz78888/lCpVinv37tGqVSsGDRqU4xsPZntp3DDg4ODAwgULmNrUgrh/16vEaWB01QgiIiKwKluXvy0a0mvZKW1C83YeER8ewovdszF/ev29QjZRq6jrXICOVQvTu2UtFm7Yh2lee70JDWTP7d6S1GSCKlWqAEhtkpzs37omPPNL+KoomJiYsGHDBmrWrMnz589p27YtT548MXSkIjmJ2/mHHE3+NviI3u38VatW5dKlS3zxxReoVCqWL19OxYoV2b17d+Y/D2EwPWs6UKuwKab/ZiumaqhV2JTeHepj32kCJtZ2OucnzqMMqF+CPwfWpE7wX4RdO0anTp3SpWmlqXWBVJ2XnbZ7S1KTCWSkRrxZ1+TNtRe5c+dmz549ODs7c+fOHdq3b8+rV68MGKh4p/fczg+QJ08e5s+fj7e3N2XLluXx48d4eHjQu3dvqWOVEygKqsP/0xbq0x5WmTChxmtUKiXJNu7EtSH7fIOo42zP2jWradKkCeHh4bRt25bbt29/UEjpvSg5K5CkJhMkJjVXrlyRnS45kZ66Jm+2UnBwcGD//v3Y2dnh4+NDjx49tL1fhPFp0KABFy9eZOzYsajVatatW0eFChXYvHmzVCI2ZomF+t5qnaBS4qmgBNBIfSXZSwNDo/hp5TZMTEzYsWMHrq6uPHnyhNatW/P06dM094hKlNpFyRqNkm22estC4fQWcBj2jYO2s8C5KQBxcXFYW1sTFRXFzZs3tdv7RA7h7wlrP0p6vM9WKN1C++3p06dp2rSpdqv3kiVLki3AJYzD2bNn+fTTT3F8fYP5bS3ZFF6DobPWU6hQ1u+aLNIghUXmGkXFFaUkHWOmomfjNQDBu2aTPyyAzz77jI4dO+Lh4YGz6gFLuhZkrvUY9kSU056blu7b+30DGbb2fEKYOiEnLErOl8ucl6//a+dgqM7eslDYEPSsmwAwNTWlYsWKgExB5TipqGuSqHbt2mzYsAG1Ws2yZcuYNm1aJgcrMlvNmjXxOXeOVZ+UooK9Ce0szlKmTBmmTJlCRESEocMT6SWFReZqlUIhVQjmJD9Cm9dM4dGjR/zwww/UqFGDKlUqM9vdHudcrxgcu5Y3U5K0bMdObru3JjIcRYGXr3UXvmf1rd6S1KSnZNZNANSqVQsALy8vfVcKY5XKuiaJPDw8tI0vf/jhB1auXJlJgQpDsXh4nMKqhAXitQqbUt8xikmTJlG2bFlWrFghRfuMwRuLzOMHH6Gf2U+0j56mc/OI/h8xmCW5VFE0OOQxI+DU36xZs4ZatWoRExND1LX9VMv/GgBX9W0aqf/biJLW7tttKhXi+LhmrB9ch197VmXNgJpYWSTGon+dT1bt7C1JTXpJYd2Eu7s7ALt27ZJ1NTlFKqqOvj1aAzBs2DC+/fZbAAYPHsyqVasyNk5hOG+9bigqE9YNqkDJkiV4/PgxAwcOpFq1ahw8eNDAgYoP9u8i8zNRxTgaXpirSkmdWxB6diIpCS0sY0+vx8zUhD59+nD69GlOnzrFgm5FiVP+7Q6uqPnKdDNvjtakdTv2m9u9TU1MiDHJlW79pzKTJDXpJZl+MImfxJs1a0aePHl4/PgxPj4+BgxUZJr3rGsCMG3aNPr3709cXBz9+/fnxx9/lEWkxuit1w2VEk+BqHvc2Pt/zJkzh3z58nH58mVatWpF27ZtuXIl+cWkIntIy/Zo+zxmvP57Hhf3rGTixIna47UKhFPW/BmmqoTXFlOVJslozfv8vLRekxW3ektSkx5SsW7CwsKCNm3aALBz504DBCky3QfUNVGr1Sxfvpzx48cDMHnyZPr160dMTNIESGRT73jdMD82i6/GjMHf359Ro0ZhZmbG/v37qVq1KoMGDSIwMGuuZxApS+326B/al+fUd61Y/P0wAGbOnMmBAwe0vzdvbw3XN1oD8OvMKezYsSNNrx3Zeau3JDXpIZXrJjp27AhIUpOjfEBdE7VazfTp01m6dCkmJiasWbOG1q1bS/8gY5GK140CBQowd+5crl+/TteuXdFoNCxfvpzSpUszZswYHj16ZJjYxXtLaRs1ikIhG0v61y+JiVrFRx99xPDhwwH45JNPCDm7We/W8CSjNYpCXFgwhzf9TufOnXF0dGTo0KF4e3unuAQixRhJe/+pzCJJzYdKw7qJdu3aYWJigq+vLwEBAZkZpcjGBg8ezN69e7G2tubIkSPUq1ePO3fuGDos8SHSuN7K2dmZzZs3c+LECerWrcvr16+ZO3cupUqVYujQofJ6ko28q2u2omhQgI5FY3V6Lf3888+4uroSHBxM8IaRKMmkGxpFxVemm1H+TZQntCnDuLHfUKhQIV68eMHSpUtp3LgxxYsXZ+zYsVy8eFHvtHaKMSoKnYsn/Iz3qY+TkSSp+VBpWDdha2tLo0aNABmtEWnTunVrjh07RuHChblx4wZ16tThzJkzhg5LvK/3XG9Vr149Tpw4wf79+2nYsCExMTEsXbqUsmXL0rt3b3x9fTM8dPHhkttGnUuJJnjHdJZOHK6zpd/S0pKNGzeSP29u8qlfoUJ/8qDdGv76GU93TOfA8tlMnz6dBw8e4OnpyYABA8ibNy8PHz7kp59+ws3NjUqVKjFt2rT/EuN/e9S1sbqhN0bzuNcE75jOjBnTqTllH72WnWLkhov0WnaKBrMOGXyrtxTfSw9p6NY6f/58Ro4cSaNGjTh69Gj6xSByhEePHtG+fXsuXbqElZUV69ato1OnToYOS7yPNHZ51ufYsWPMmDGDffv2aY917NiRCRMmaMtIiKzrza7ZDtaWlLczo0rlSjx48ICvvvqKOXPm6Jy/Zs0aJnzej4J5TFi0eDEWBUvz/HUM+awStl+/jIzF2rYQsZERNGncmOjoaCZMmKBT8yoqKoq//vqLdevWsWfPHqKjo7X3Va3qyt6OkTgRBE5uMPgw8Qo6Mbo65ab1gG+4X6w1gM4OqcR/LepTLd2L86X2/VuSmkx27949SpQogVqt5smTJ9jZ2aV8kRBveO37F8GrBzBwyzMO3dXw888/M2rUKKk+nINduHCB6dOns3XrVu10QosWLZgwYQJNmjSR341sZO/evXTo0AG1Ws3p06epUaOGzv39+/dn08mb2Lcahir3f2ta3q70++eff9KnTx8A1q1bR69evZL8rJcvX7J9+3bWrVvH4cOHaV4CDvTJrb1/tdIJt65fU6lSJe3vULxGof5ML4JCo/R291YBjjaWHB/XTGcK7UNJReEsqnjx4ri6uqLRaNi6dauhwxHZjaKQ65/ZFM8Vye89i6IoCmPGjKFz587S4TsHc3NzY/PmzVy7do1+/fphYmKCp6cnzZo1o2rVqixdulQapWYT7du3p1evXmg0GgYNGpSkD1zHzyfh0GkC5Mqvc/ztSr+9e/dm3LhxAAwYMIBz584l+Vn58uXj008/5eDBgwQFBrL203LE/zvMEadRKPd4C1WqVMHFxYUJEybg4+PDmTshBIVF601owPA1bCSpMYC+ffsCCYu/pFqoSJM3qlaXMH/O1lnDMTMzY+fOnVSsWJEtW7YYOEBhSOXKlWPlypX4+/szfPhwrKysuHz5MkOHDqVw4cKMHDkSPz8/Q4cpUjBv3jzy58/PpUuXWLBggfZ4vEZh9sHboFIl29H7zUq/06ZNo3379kRFRdGpU6d3lgKwC72EfewDTP59WFO1ilqFTWlfzpJbt24xY8YMatSoQa8Bn6XqORiqho0kNQYwZMgQ8ufPz61bt9i2bZuhwxHZhZ6q1V3yXeXc2bO4uroSEhJCt27d6N27N8+fZ71KnyLzlChRgt9++41Hjx7xyy+/ULp0acLCwpg/fz7lypWjZcuWbN++XbrBZ1EODg7MmjULSGiX8uDBAyBhbUtgaPLJwtujJCYmJqxbt47y5cvz6NEjOnfuTFSUnuvfUTNp55c1Wb9uHV27diVXrlwE3bmZuudgoBo2ktQYQJ48efjyyy8BmD59ulSKFamTTNXqKrmDOXPmDN999x1qtZp169ZRqVIl/vrrL8PGKwwuf/78jB49Gj8/P/bv34+HhwdqtRpPT0+6dOlCqVKlmDZtmkxdZkEDBw6kXr16REREMHLkSOD9Kv3mzZuXXbt2kT9/fk6fPs3gIUM5GfBMdxv2O2ommTy5RM+a9mzevJng4GDW/ToFi/jXKEoyO/cUBfvcpgarYSMLhQ0kJCSE4sWLExERwb59+7TVhoXQS1FgWVMIvKz7wqMygUJVYPBhUKk4ffo0/fr1004xDB48mJ9//hlra2sDBS6ymrt377JkyRJ+//13nj1L2H1lZmZGly5d+KpzNWoEb0bVdjY4NzVwpOLKlStUq1aNuLg4du3ahV2FevRadirF69YPrkNdZ91eUl5eXnT8fBL5mg7CNK+99nihvBbsyz2JfC+uob/EgBqcXLWvMQD7fQMZtvZ8ko3lCYmOCttrW7mw+4+0PdkUyELhLK5AgQIMGTIEgBkzZhg4GpHlpbJqde3atblw4QKjRo0CYNmyZVSpUkXKBwitEiVKMGPGDB48eMCaNWuoU6cOsbGxbNy4EcVzMqpnN7m/chCXL10ydKg5XuXKlRkzZgwAn3/+ORUdLN670m9swQrYdfwWE2vdHbfPw14R+/wBaamZlFhnp9BbNWyslGhe7f+Fn0d/kpqnlyFkpMaAHj58SKlSpYiNjeX48ePUr1/f0CGJrChxlObxJVL7SQrgyJEjfPrpp9y9exeAUaNGMX36dKysrDIlbJF9XLhwgVNrpzHM+r9u4K3XRvDUujJ9+/bl448/pmDBggaMMOeKiIigQoUK3L9/n7Fjx9L0k1EMW3se0O3ylDhKMqeTC93qltF5jHiNQoNZh5Jdj+NECGWso1jRvyYm+nY1JVMz6e06O7VK2qJCQa1O//ESGanJBooUKUK/fv0AGa0R7/Ce1WebNGnC5cuXGTx4MJCwo6JChQqsW7cuxd4vImdxq1qVYS7PUP5dKBqvwLRmVly8eJExY8ZQuHBhOnTowKZNm/QvNBUZJnfu3CxcuBCAX375hcLKM72VflWRoQTvmM6GOd8mWaeZ0gLjxxTgaHhhzkQVS1OPOhO1irrOBehYtTB1nQtgolZlSEKTFjJSY2C3bt2iXLlyaDQaLl26RJUqVQwdksiKPrD67L59+xg8eLC2AWK1atWYPXs2zZs3T+9IRXbk7wlrP0pyeE++Afxvw0lOnz6tPWZjY8NHH31Et27daNasGebm5pkZaY7VpUsXtm/fTr169Th27BgKKp1REtMXd6lXrx4mji58MfY7urRtQa2StpioVey8+IiRGy6m+DN+7VmVjlXfXcXaUKSisB5ZMakB6NGjB5s2baJHjx5s2LDB0OEIIxUREcG8efM4vX4WMxtr+HJfFGYuLZk1a5Yk0zlZKhah+928yZo1a1izZg3379/XnpIvXz46duxI165dadmyJRYWFgZ4AjnDgwcPqFChAq9evWLZsmUMGjRI5/79voF89ecpIpT/kszEKsM2VubvvcA4q5CkRo+smtRcvHgRNzc3AP766y/atm1r4IiE0VIUYhc1xOzpFc4+jqfWsghUKhV9+/ZlypQpFCtWzNARisyWzCiNVp+tULoFABqNBm9vbzZv3szWrVt1toLnzZsXd3d3unbtSuvWrWXtVgaYO3cuY8aMIX/+/Ny4cQMHBwcg+d1Iiatjfvu4GlP3XiMoNEpvK0xF0eCQ24xT37dO19YG6UnW1GQjVatW1dYhGDhwICEhIQaOSBitAC/Mnl4BoKaTCdM+bYqiKKxatYqyZcsybtw4Xrx4YeAgRaZJLLqW7FuBOuH+fz/7qtVqmjRpoi3s5+3tzZdffknhwoUJCwvjzz//pHPnztjb29OzZ0+2bNmi021afJgvvvgCV1dXXrx4wTfffAMkLNadvPua/mTl369T917jh/blAZLunFIUQMXLIys4eu2hbv2abCjbJDU//vgjqn9LQyfeypUrZ+iw0s2MGTMoV64cgYGBDB8+XAryifSnpyLxhDoaTp86ReN/O/rOnj0bZ2dnfv75Z1kQmhO85yJ0SKhW27BhQ3799Vfu37/PP//8w5gxYyhWrBgRERFs3LiRbt26UaBAAdq1a8f//d//ce/evYSLAw7DwloJX0WqmZqasmTJElQqFatXr+bw4cOprjKcP7eF3gXGBfNaoLl6gMhy7Riw9jIjN1yk17JTNJh1SNtHKjvJNtNPP/74I1u2bMHT01N7zNTUNE1drrPq9FOic+fOUbduXeLi4vjzzz/5+OOPDR2SMCbJTTP02Yri3Jy//vqLcePGcfXqVSCh+erkyZPp1auXLAY1Zh+4CP1tiqJw7tw5tmzZwpYtW7h9+7bO/ZUqVeSvzlEUVQejFHJDNeRwss0RhX7Dhw9n0aJFuLi4MG3tAb7a4pviNYmLgN/ehv0iIprh6y6gKIpOP6nEfy3qU03b+duQjG5NzY8//siOHTu4ePFiqq+Jjo4mOjpa+31YWBhFixbNskkNwJQpU5g0aRL58uXjypUrFClSxNAhCWOQyorE8fHxrFq1ih9++IHHjx8D4OTkxOeff87QoUOxtTVM6XORPSmKwrVr19i7dy979uzhxIkTtCip4kCf3NpzZj2uTfFmA2jdujX58+d/x6OJRC9fvqRcuXI8efKEYZPm8ldUmRSv0bcIOKX6NSrA0caS4+OaGXytjVGuqbl16xZOTk6UKlWK3r1766zC12fGjBnY2Nhob0WLFs2kSN/fhAkTqFWrFi9fvuTTTz+VeiIifaSyIrGJiQkDBgzg1q1bzJw5E0dHRx4/fsyECRMoUqQIw4YNky7PItVUKhUVK1Zk7NixeHt7E/z0KesGViD+34/ScRqFpsoJevXqhb29PY0aNWLq1KmcPHlSmm2+Q758+Zg7dy4AK2ZNwC63abJVhhVFQ3xYMLkiHie5L60NMrODbDNSs2/fPl69eoWLiwuBgYFMnjyZR48e4evrm2xfm+w4UgPg5+eHm5sbkZGRLFiwgM8//9zQIYns7D0rEkPC39DGjRuZO3euzihpu3btGD16NM2bN9cZshbinZKZAh163J6lXgE6x2xsbGjatCktWrSgZcuWlClTRn7X3qAoCq1atcLT05O63YYRWKp9wvE3zlH9e97THdOpmDeWkydPYmpqqr0/O9WvMbrpp7e9fPmS4sWL88svvzBw4MBUXZPV19S8aeHChXzxxRdYWVlx4cIFXFxcDB2SyK7iomFuJYh4mvw5eRxglC+Y6q8zoigKR48eZe7cuezevVu7kL1y5cqMHj2aXr16YWlpqfdaIYAUp0DvtvyDvw8e5ODBg3h5eSXZhVesWDFtgtO8eXPs7e3J6W7dukXlypWJjo7m+8Wb8XyeX2fkpZCNJV/UL8SIjvV5GRrG55Pn0sL9I21LgzN3nmeb+jVGn9QA1KxZkxYtWqS6xUB2Smo0Gg1t2rTh4MGD1KxZkxMnTmBmZmbosER2lY6LQW/dusX8+fP5448/tNt1HRwc+GlYez7OfxHTDj9Lh2eRVBrq4cTHx3P+/Hk8PT05ePAgJ06cICZGdweWq6srTZo0oXHjxjRs2DBNm0aMSeI6TEdHR3yvXuPmCw1BoZE8j4jBNo8Fjnkt2f7XQdZdj9Ltzm1jyQ/tK6RcvyaPOae+a5Vt1tRk26Tm1atXFCtWjB9//JEvv/wyVddkm6Qm4DDsG0dwzW8o22YIL1++pHv37vz55586Q4dCGNKLFy9YtmwZCxYs4OHDh5welItahU25E52fF1234latmkwXiAQfMAUKCdWwjx8/zsF/R3IuX76c5JyKFSvSuHFj7S2nNOCMjo6mSpUq3Lx5ky+//JK2g79l8u5rSdfKKIrOf9vEfw1pVJKl3ncSTnnrfAWwubIJn53LDf6h2uiSmq+//hp3d3eKFy/O48ePmTRpEhcvXuTatWupHobMFkmN9o//Aji5sa/IWDp26kRsbCy9e/dm1apVmJiYGDpKIbRiY2M5sWoKTR7O1x5rvTaChxYufPLJJ/Tu3TtbLNIXGSgdpkDf9OTJE44cOcLRo0c5evQo165dS3KOi4uLNsFp2LChUf8Oenp60rJlS3K71Meu0/hUX5e4u+mH9uWZuve6TiLkkMeMu9vmEHzhIBMmTGDatGkZEHnqGV1S07NnT7y9vQkJCcHe3p4GDRowbdo0nJ2dU/0Y2SKpeXuIts9Wdl6NoGvXrsTFxdG/f3+WL19u8E6oQmj9m4grgZdRKfHEK3AhUEPNZa+AhB0wjRs35pNPPqFr165Z929PZKx0rofzpuDgYI4dO6ZNci5fvpykgGnRokWpV6+e9ubq6mrw0Yf09FHXbpx2aI+ptV2a6/6sH1xHu8YmsX5NrZK2bN+2lW7duqFSqTh0+DAWRSrp3J+ZU1JGl9Skhyyf1Ly9kO6NGiJbtm6lZ8+exMfHM3jwYBYvXiyJjcgaklkrsd9uCDO3nOXo0aPaY5aWlnTs2JFPPvmEVq1aGdWbisg6Xrx4wbFjx/D29ubo0aNcuHCB+HjdcgZWVlbUrFlTm+TUrVs36bqcf5cC0HZWll8ntvPkNUbuvPNe175rd9OAAQPYcOIG9q2Go8r9Xx2hxGaZmVWYT5IaPbJ8UvOOiq+UbsH69evp06cPGo2G4cOHs3DhQlmzIAwrFUX97t2/z7p161izZg3Xr1/XnmJvb0+vXr345JNPqF69uvwuiwwTERHB2bNn+eeff7Q3fT3OypYtS+3atalVqxa1atakxsWxqAMvgpNbsut9sorUbs/W5127m7afvcPoLVdRwKAVhyWp0SNLJzWprPi6evVq+vfvj6IojBw5krlz58qbgTCcNOxoURSF8+fPs2bNGtavX8/Tp/+tryhXrhx9+vShT58+FC9ePKOjFjmcRqPh5s2bOknOmwk3QCtnE53Kx56OwyncuC8uLi5ZcpT8ZEBIqrZnv0lRNDhaW/DPhJZ6p5KyUsVhSWr0yNJJTRreHJYvX86gQYMA+Oabb5g1a5YkNiLzfcCOlri4OP7++2/WrFnDjh07dJpn1qtXD3d3d9zd3alQoYL8botM8fz5c06dOsWZM2c4e/YMU4v9QxV7DaZqFXEahfOB8dT+/TV58+alRo0a1KxZk2rVqlGtWjVKlSpl8EQnpQQkiX93N5V5cgTPlXP0npLaRCkz6thIUqNHlk1q3uPNYfHixQwbNgyA7777jqlTp8qLv8hc6bSjJSwsjG3btrFmzRoOHz6ss8CzZMmS2gSnUaNG0lhTZI5kPmR23BzPrmsRSY7b2Njg5uamTXKqVatG2bJlM32n6n7fQIatPZ/wN5TC+4GdlZobG6YRceMEO3bsoGPHjknOyUoVhyWp0SPLJjXv+eawYMECbY2eYcOG8csvv0hVV5G50nlHy6NHj9i9eze7du3i0KFDOm1O8ubNS+vWrXF3d6ddu3YUKGDYCqfCSL1jKYDiWIUrdX7lzNmz+Pj4cP78eS5duqTze5ood+7cVK1aVXtzdXWlUqVKWFlZZWj4+30D+WH7ZYIj/uudlVhoL39uc53dS99NGM+sWbNwcnLi2rVr2NjY6DyWjNRkcVk2qYH3fnOYO3cuY8aMAaBq1aps3LiRsmXLZlSUQmSaiIgIDh48yO7du9mzZ4/OGhy1Wq0zTVWuXDkZqRTpIw1LASChTtP169c5f/4858+fx8fHh4sXL/L69eskl6rValxcXHB1ddXeqlatiqOjY7r+/sZrFIZNnMOf2/ZQIJcpV4/swjpP7iTnRUZGUrlyZQICAhg2bBj/93//l+RxGsw6lGzFYUhImGRNjYFk6aTmA+zbt4++ffvy7Nkz8uTJw+LFi+ndu7ehwxIi3Wg0Gs6ePcvu3bvZvXt3koqyzs7O2gSnYcOGslVcvJ8PrHycKD4+nps3b+Lj48OlS5e4ePEily5dIjg4WO/59vb2VKxYMcntQ0YjX79+TYUKFbh37x7ff/89U6dO1Xve4cOHadasGQDHjh2jQYMGOvcnTmmBbsVhRdGgQsWCXq64V834woaS1OhhrEkNwOPHj/n444+1NUE+/fRTFixYQO7cSbNzIbK7e/fusWfPHnbt2sWRI0d0+gLlzZuXRo0a0aRJE5o2bYqrq6tU4Rapk86Vj9+kKApBQUHaBCfx5ufnh0ajL4ECR0dHvclOvnz5UvUzt2/fTpcuXTA3N+fq1auULl1a73mDBg1i+fLllCtXjosXL2Jhofvc9vsGJmm9oHkVwrODi/m6RwumTJmSuv8IH0CSGj2MOamBhE8HU6dOZcqUKSiKQvny5dm0aROVKlUydGhCZJjw8HD+/vtvdu/ezd69e3n2THcaN1++fDRu3Fib5FSuXNngO1VEFpaBlY/1ef36NdevX+fq1av4+vpy9epVrl69yr1795K9pnDhwjpJTqVKlahQoQLW1tY65ymKQps2bfj7779p3749e/bs0ft4L168oHz58jx58oT//e9/fPfdd0nOidcoOhWH7/scomeP7piYmHD69GmqV6/+Yf8hUiBJjR7GntQkOnz4ML179yYwMBBLS0vmz5/PoEGDZM2BMHrx8fFcvHiRw4cPc/jwYY4dO0Z4eLjOOba2tjRu3JimTZvStGlTKlSoIEmOyHLCw8O5fv26TqJz9epVHj58mOw1xYoV00l2EmvqNGzYkNjYWHbv3k2HDh30Xrtu3Tp69+6NpaUlV69epVSpUinG2KNHDzZt2kTFihXx8fFJMsKTniSp0SOnJDUAT58+pW/fvhw4cABI6J21ZMkSo3/eQrwpLi6O8+fPc/jwYY4cOcKxY8eIiNDdkmtnZ6cdxWnSpAnly5eXDwAiy3r58iXXrl3TSXR8fX0JCgpK9horKysiIyOxtrZm3LhxVKpUCRcXF5ydnbXrzxRFoWXLlnh5edG2bVv27t2b4t/Bs2fPqFixIk+fPmXcuHHMnDkzXZ/rmySp0SMnJTWQsLhyzpw5TJgwgfj4eJydnVm5cmWShWBC5BSxsbH4+PhoR3JOnDiRZJdKwYIFadKkCfXq1aNWrVq4ubll6CdQIdLD8+fPdZKcGzdu4Ofn986RHRMTE0qVKkW5cuVwcXHBxsaGyZMnExcXx+bNm+natWuKP3fHjh107twZtVrNiRMnqFOnTno+LS1JavTIaUlNopMnT9KzZ0/u378PQOfOnZkxYwYuLi4GjkwIw4qJieHs2bPakZwTJ07oVDcGMDMzo2rVqtSuXVt7K126tIzmiGzh1atX3Lx5k5UrV7JgwQJMTEwoX748d+7cSTJq+SaVSkXVqlUpW7YspUuXpkyZMpQuXZrSpUvj4OCg8/v/ySefsHbtWsqWLcvFixczpBaPJDV65NSkBhKy+G+//Zbly5ej0WgwMTFh8ODBTJo0CUdHR0OHJ0SWEB0dzZkzZzh69CinT5/m9OnTerfh2traJjQ9rFVL2wAxSYdnIbIQRVFo0qQJ3t7e9OrViz///JNHjx7h5+enc7tx48Y7FykD5MmTR5vglClThkKFCvHjjz/y/PlzRo0axdy5c9M9fklq9MjJSU2ia9euMX78eHbt2gUkVL385ptv+Oqrr8iTJ4+BoxMia1EUhbt372oTnNOnT3P+/Hm9FWSdnZ21IzkybSWyogsXLlC9enUUReH48ePUr19f73m7du2iY8eOqNVqRo8eTUREBP7+/vj7+3Pv3j1SShvGjBnDzz//nK6xS1KjhyQ1//H29uabb77hzJkzQMI6gh9//JGBAwdK4TIh3iEmJobLly9z+vRpzpw5w+nTp/Hz80tyXuK0VbVq1bTVYytXrpxk260QmWnIkCEsW7aM6tWrc+bMmWR3/nXr1o0tW7ZQp04dTpw4oT0vOjqau3fvcuvWLW2ik3i7ffs2iqLQq1cv1q1bl65xS1KjhyQ1uhRFYevWrYwfPx5/f38AXFxcmDlzJh07dpQ1A0Kk0osXLzh79qzOiM7b9XISOTs765TJd3V1pXjx4vL3JjLF06dPKVOmDGFhYaxYsYJPP/1U73n3HzzEtcVHxJhY8vWIwUwc1jvFVgjBwcF4enrStm3bVBcITC1JavSQpEa/mJgYli5dypQpU7TrB+rVq8dPP/1EvXr1DBydENlP4rTVmTNndCrIPn78WO/5NjY2VKlSRSfRSbb5YcBh2DcO2s4C56YZ/EyEMfr555/5+uuvKViwIDdv3kzyfqivgnBBa3Mmd6xEm0qFMjtcQJIavSSpebewsDB++uknfv75ZyIjIwFo3bo1w4cPp127dpiamho4QiGyt2fPnumUyL906RLXrl0jNjY2yblqtZqyZcvqTF2VL1eOUl6DUAVeACe3FHsQCaFPTEwMlSpV4tatW0nqyyT2ekqSGCgKKpWKRX2qGSSxkaRGD0lqUufx48dMmjSJFStWaHuSFClShMGDBzNw4EAKF06/EuFC5HQxMTHcuHEjSbKjb9dVK2cTDvT5r5/bWrpgXqEt5cqVo2zZslhaWmZm6CIb27NnD+7u7jp9oRK7cr85QvO2zOrK/TZJavSQpCZtAgICWLp0KStWrNCuDzAxMcHDw4PPPvuMFi1aSHl5ITKAoigEBgbqJDlXr/qyos5dqjqqMFWriNMonA+Mp/bvCcUD1Wo1JUuWpHz58tpbuXLlKF++fLqvbxDZn6IotG3blgMHDtCxY0d27NjByYAQei07leK16wfXoa7z+3cQfx+S1OghSc37iY6OZuvWrSxevBjzh/8wv60lX+6L4q66BEOHDqV///7Y29sbOkwhjJu/J6z9KMnhUeeKsfL4A0JDQ5O91NHRUSfZSSykVrx4cZlWzsGuXbtGlSpViI+P5+DBg/wTWYgVJ+6meN2vPavSsWrmjthLUqOHJDUfSFF4Pb8uuV5cxycIaiwJA8Dc3JyuXbvy2Wef0aBBA9nFIUR6UxRY1hQCL4MS/99xlQkUqoIy6BBPnj7l+vXr2tuNGze4fv06jx49SvZhTU1NKVGihLaQmrOzs/bfJUuWlDo7OcDIkSOZP38+Ls27E1Wjb6quWdHHlWaVimRwZLokqdFDkpoP9NYnxQP2Q/l+5WHOnTunPVahQgU+++wzPvnkExnyFiK9JDNKo9VnK5RuofeusLAwbty4oU1yrl+/TkBAAP7+/klaQrxJpVJRrFixJMlO6dKlKVWqFLlz5072WpF9PH/+nDJly2LVfQ6mee2Ad3woVRTiwp/RMm8QnXp+goO1JbVK2mbK+hpJavSQpOYDvP1J8d9PiAw+zDkfH5YsWcK6deu0zQHNzc1p1qwZHh4eeHh4yOJiId5X4t/e40uARs8JanByTfNOKI1GQ2BgoE7xtMRkx9/fn/Dw8Hde7+TkRKlSpShRooT2Vrx4cUqUKEHRokXTPsojW9UN5ttflrPhaWra5SjEvw7DJJeN9kghG0smuVfI8B1RktToIUnNB0juk+IbnxBDQ0NZu3YtixcvxtfXV+e06tWr07FjRzw8PKhSpYpMUQmRWnHRMLcSRDxN/pw8DjDKF0zTZ7pIURSCg4P1Jjv+/v48f/78nderVCqcnJy0Sc6bCU+JEiUoVqyY7k4tbeImW9UNYbvPA0Zvvpy6kxVF5/9N4r8yequ3JDV6SFLznlKYz3/7BUhRFG7cuMHOnTvZtWsXp06d0ukVUrx4cTw8POjYsSONGjWStgxCpCT0IUTor1AMQG57sMm80dAXL17g7+/PnTt3uHv3rvZ279497t69qx2xfRdHR0dtktOipIqB5nu1973ssAKb6l3kw08mSe2up+SoAMcM3uotSY0ektS8pw+Yzwd48uQJe/bsYefOnRw8eFBnHt/GxoZ27drh4eFB27ZtsbGxSfZxhBBZn6IoPHv2LEmi8+YtIiJC55rTg3JRrZCJzlb1hqvjKFKkCEWLFqVIkSI6/078am9vL4lPOtDWp3kZ+UEjZBm51VuSGj0kqXkP6Tyf//r1azw9Pdm5cye7d+/WKTBmampKkyZN6NixI+7u7hQvXjz9nocQIktQFIXnz59rExzF35OukUmbH7ZeG8HfAfF6HuE/5ubmSRIeJycnnJycKFSokParFCVM2X7fQD5bex5F0aBS/Vd/TAVJqwsnIyO3ektSo4ckNe8hA+fz4+PjOX36NLt27WLXrl1cv35d5/6yZcvSoEED6tevT4MGDShTpox8KhPCmCQzta2oTIixLccZ15k8fPSIhw8f8uDBA52vQUFBqf4x+fPnT5LoSPKT1H7fQEauOkG0yX89xwrZWNKzZjHmet5M8XoZqclkktS8p0yaz7916xa7du1i586dnDhxQtuiIZG9vT3169fXJjnVqlXD3Nz8g3+uEMJAPmBqOyYmhsePH2uTnMSEJzAwkMePH2u/RkdHpzqc/Pnz6yQ8jo6OFCxYMMnNzs7OaKupP3ocSKWmnYgxsWT8qOGMH9QdgAazDhEUGqV/1EZRKJTPStbUZDZJarKPFy9e8M8//3DixAmOHz/OmTNnkrw4WVpaUqtWLW2SU7duXfLnz2+giIUQaZJBW9V1f4TCy5cvefz4sU6io+/ru2r2JIlMrcbe3j7ZpOfNm729PSYmJu8Vv6FMmzaN77//nmLFinHjxg2srKy0jS5BdzpKUTSAih9bFObTlm4ZFpMkNXpIUpN9RUdHc/78eY4fP65NdEJCQpKcV6lSJW2SU79+fUqUKCFTVkJkRQbYqp6cN5OfxETn8ePHPHnyRHsLCgriyZMnel933kWlUmFnZ4eDgwP29vbY29tr/63vmK2trcFHgV6/fo2LiwsPHz5kxowZfPvtt0DC9NTk3dd0Gl6axoTzeO982lR0ZNu2bRkWkyQ1ekhSYzwUReHmzZs6Sc6tW7eSnOfk5ET9+vWpWbMmlStXpkqVKhQqVEgSHSGygiy2VT01YmNjCQ4O1kl4krsFBweT1rdYtVqNnZ2d9lagQAGdr/qO2djYpHsitGbNGvr27Yu1tTX+/v44ODgACTul/r+9e4+Kss7/AP4eLjODyB3lrgaoaYh3CMWk1KXVzHbb9VKRkbseTU1ly7tRmuJhXTNTcyNTt2OxR8XWvKAbRa2XFgMxBNQAXUyRi9xBhmHm+/ujH3McGJAZmBmY3q9znjPPPPN95vl8ZtB5n2eeeZ60G+UoqWlAXwc57OuLMGrkCKhUKnz99dd48knjnDiRoUYHhhrLVlxcjPPnz2uCTnp6OpqamlqNc3V1RXBwsCbkDBs2DEFBQTztOxF1KZVKhbKyMk3AKSkpQWlpaZvzFRUVBm3HysoKbm5u7QaflrfOzs7tfi2mVqsREhKC9PR0LFiwAB9++GGbYxcvXoxdu3Zh+PDhSE9PN8rXbQw1OjDU/LrU19fj4sWLOHfuHC5fvowff/wR169fb3UAMvDLLmJ/f38EBwdrBR5/f/8e9304EfVMSqUSZWVlKC0tRVlZGe7du4eysjKt+ZbLHnY5i7ZIJBI4OTnB1dUVLi4ucHV11Zp3cXFBaWkp4uPjIZFIkJSUhLFjx8LFxQV2dnZae7vLysoQGBiIqqoqfPzxx5g3b15XvSQaDDU6MNTQ/fv3kZubi6ysLPz444+a2+LiYp3j7ezsEBQUpLVXJzg4GO7u7iaunIioNYVCgfLy8g6HoLKyMlRXV3dqmzKZrFUQunv3Li5evAh7e3ucOHECEydO7KIOf8FQowNDDbWlpKQEWVlZWmHnypUrbf4iwtPTE0FBQZorGDdfxZhXLyai7q6xsRGVlZUoLy9HeXk5KioqWs033965cweZmZkAfvmaS9ee7paefvppnDp1qktrZqjRgaGG9KFSqZCXl9dqr05BQUG763l6emqCTsvJ3d2dBykTUY+ydOlS7NixA0FBQfjPf/6DqqoqnUEoLS0Nly5dQnx8PP7whz90aQ0MNTow1FBXqKmpQXZ2NnJycpCfn681PexAP0dHxzYDj6+vL4/fIaJu5969ewgMDERlZaXRjpl5GIYaHRhqyNgqKiq0Qk5eXp5m/vbt2+2uK5VKMWDAAE3I8ff3h6+vL3x8fODr6wsvLy9e0ZyIzGLbtm34y1/+Ak9PT/z000/o3bu3SbfPUKMDQw2Z0/3793Hjxg2doefmzZtQKpXtri+RSODh4QEfHx9N0Gk57+vra/L/bIjI8ikUCgwdOhQFBQWIjY3F22+/bdLtM9TowFBD3ZVKpcKtW7e0As/Nmzdx+/8v5nfnzp2Hhp5mjo6OrYJOy/BjydeuISLjOHToEGbOnAl7e3vk5eXB09PTZNtmqNGBoYZ6KrVajdLSUty+fVsTdFrO//zzzx0+Z4VUKoW3t7cm7DRfo+bB07c3zzs7O/PgZiKCEAKPP/440tLSsHDhQuzevdtk22ao0YGhhixdTU1Nu6Hn9u3bKCkp0evU7ba2tppr1+gKPS2X9e7dmyGIyEKlpqbiySefhLW1NXJycjBo0CCTbJehRgeGGqJfzlFRVFSkCT23b9/WOlV7SUmJZt6Qk3TJZLJ2w0+fPn3g5uYGFxcXzSSVSo3QKREZw7Rp03Dy5Ek8//zzOHz4sEm2yVCjA0MNkX4aGhpQVlamM/Douq2rqzNoO/b29lohp/kspS3nW953dnaGjY1NF3dNRO3JysrC8OHDIYTAhQsX8Pjjjxt9mww1OjDUEBlXfX19m3t9HrytqKhARUUFKisr9b6KcUuOjo5thp4H7zs7O8PR0VFr6t27Nw+YJjJAdHQ09u/fjwkTJuDbb781+lfODDU6MNQQdS9qtVrr7KTNYaet+QfvG3ohv5YcHBzg5OTUKvC0N7Uc7+DgwBMn0q/KrVu3MGjQIDQ0NODYsWOYPn26UbfHUKMDQw2R5VAqlaisrHxoEGqer6qqQnV1Naqrq1FVVQWVStWl9djb2+sMQL1799aa7O3tdc63vC+Xy3nANXVrK1euRHx8PIYOHYrLly8b9atghhodGGqICPjlp6kNDQ2akPPg9GD46cikUCiMUqOVlVW7oae9+7169Wp3YmCirlBRUYGAgABUVFQY/fIJDDU6MNQQUVdTKBSoqalpMwzV1taitrYWdXV1Oudb3r9//75J6n5Y8OnIZGdnBzs7O8jlcsjlcp3zMpmMxy1ZsL/97W9444034OPjg+vXr6NXr15G2Q5DjQ4MNUTU3alUKtTX17cZetoLRA8Go/r6+lZTY2OjWXqSyWSasNNW+Gk53xXjZDIZ90gZWUNDAwYPHozCwkLExcVh1apVRtkOQ40ODDVE9GvW1NTUZuAxdGpoaNBM9+/f19yq1WpztwsAOvccPThJpdKHznd0XEfWkUqlFhe0Pv30U7z88stwcnJCfn4+3NzcunwbDDU6MNQQEZmGUqnUGXjam++qcd39Y83W1tag8GRra6vztr3HOjKmvcc68tWhWq3GqFGjcPnyZSxfvhzbtm3r8teso5/fPGsVERF1OVtbW9ja2sLBwcGk2xVCaAUqXeGnsbERCoVCMz14vyseazmu5S/tlEollEolamtrTfraGMLKyqpDwan5grvbt29HUFAQXn31VbPUy1BDREQWQyKRaD6Eu8seeZVKpQk6hoYmpVKJxsZGze2D8/o81t4YXcdcqdVqTSjsCCEEPvvsM4YaIiIiS2Rtba35pVh3JoSASqUyODhdv34dOTk5WLlypdl6YKghIiIiSCQS2NjYwMbGptsHsLbw5AFERERkERhqiIiIyCIw1BAREZFFYKghIiIii9DjQs2uXbswYMAAyOVyhIaGIi0tzdwlERERUTfQo0LNP//5T8TExCA2NhYZGRkYPnw4IiMjUVJSYu7SiIiIyMx61GUSQkNDMXbsWOzcuRPALycF8vPzw5IlS3ReRKv5pEXNqqqq0K9fP9y6davbnJSJiIiI2lddXQ0/Pz9UVlbCycmpzXE95jw1jY2NSE9Px+rVqzXLrKysMHnyZFy4cEHnOnFxcXjnnXdaLffz8zNanURERGQcNTU1lhFqysrKoFKp4OHhobXcw8MDV69e1bnO6tWrERMTo7mvVqtRXl4ONzc3k14ltTlhWuIeIkvuDbDs/iy5N8Cy+2NvPZcl92fM3oQQqKmpgbe3d7vjekyoMUTz1U0f5OzsbJ5iADg6OlrcH3EzS+4NsOz+LLk3wLL7Y289lyX3Z6ze2ttD06zHHCjs7u4Oa2trFBcXay0vLi6Gp6enmaoiIiKi7qLHhBqpVIrRo0cjJSVFs0ytViMlJQVhYWFmrIyIiIi6gx719VNMTAzmzp2LMWPGICQkBNu3b0ddXR2io6PNXVq7ZDIZYmNjW30VZgksuTfAsvuz5N4Ay+6PvfVcltxfd+itR/2kGwB27tyJv/71r7h79y5GjBiBHTt2IDQ01NxlERERkZn1uFBDREREpEuPOaaGiIiIqD0MNURERGQRGGqIiIjIIjDUEBERkUVgqOkiu3btwoABAyCXyxEaGoq0tLQ2xyYkJGDChAlwcXGBi4sLJk+e3O54c9Ont6SkJIwZMwbOzs6wt7fHiBEj8Omnn5qwWv3p09+DEhMTIZFI8Nxzzxm3wE7Qp7f9+/dDIpFoTXK53ITV6kff962yshKLFi2Cl5cXZDIZBg0ahJMnT5qoWv3p019ERESr904ikWDatGkmrLjj9H3vtm/fjsGDB8POzg5+fn5Yvnw5GhoaTFSt/vTpT6lUYsOGDQgICIBcLsfw4cORnJxswmo77rvvvsP06dPh7e0NiUSCL7744qHrpKamYtSoUZDJZAgMDMT+/fuNW6SgTktMTBRSqVR88sknIjs7W/z5z38Wzs7Oori4WOf4F154QezatUtcunRJ5ObmildeeUU4OTmJn3/+2cSVP5y+vX3zzTciKSlJ5OTkiLy8PLF9+3ZhbW0tkpOTTVx5x+jbX7MbN24IHx8fMWHCBDFjxgzTFKsnfXvbt2+fcHR0FEVFRZrp7t27Jq66Y/TtTaFQiDFjxoipU6eKs2fPihs3bojU1FSRmZlp4so7Rt/+7t27p/W+XblyRVhbW4t9+/aZtvAO0Le3gwcPCplMJg4ePChu3LghTp8+Lby8vMTy5ctNXHnH6NvfihUrhLe3tzhx4oTIz88Xu3fvFnK5XGRkZJi48oc7efKkWLt2rUhKShIAxNGjR9sdX1BQIHr16iViYmJETk6O+OCDD4z+ecBQ0wVCQkLEokWLNPdVKpXw9vYWcXFxHVq/qalJODg4iAMHDhirRIN1tjchhBg5cqRYt26dMcrrNEP6a2pqEuPGjRMff/yxmDt3brcNNfr2tm/fPuHk5GSi6jpH394+/PBD4e/vLxobG01VYqd09t/de++9JxwcHERtba2xSjSYvr0tWrRIPPXUU1rLYmJixPjx441ap6H07c/Ly0vs3LlTa9nvf/978eKLLxq1zs7qSKhZsWKFeOyxx7SWzZo1S0RGRhqtLn791EmNjY1IT0/H5MmTNcusrKwwefJkXLhwoUPPUV9fD6VSCVdXV2OVaZDO9iaEQEpKCq5du4YnnnjCmKUaxND+NmzYgL59+2LevHmmKNMghvZWW1uL/v37w8/PDzNmzEB2drYpytWLIb0dO3YMYWFhWLRoETw8PBAUFITNmzdDpVKZquwO64r/U/bu3YvZs2fD3t7eWGUaxJDexo0bh/T0dM1XOAUFBTh58iSmTp1qkpr1YUh/CoWi1de8dnZ2OHv2rFFrNYULFy5ovRYAEBkZ2eG/Y0P0qMskdEdlZWVQqVTw8PDQWu7h4YGrV6926DlWrlwJb2/vVm++uRnaW1VVFXx8fKBQKGBtbY3du3djypQpxi5Xb4b0d/bsWezduxeZmZkmqNBwhvQ2ePBgfPLJJwgODkZVVRW2bt2KcePGITs7G76+vqYou0MM6a2goABff/01XnzxRZw8eRJ5eXl47bXXoFQqERsba4qyO6yz/6ekpaXhypUr2Lt3r7FKNJghvb3wwgsoKytDeHg4hBBoamrCggULsGbNGlOUrBdD+ouMjMS2bdvwxBNPICAgACkpKUhKSuqWgVtfd+/e1flaVFdX4/79+7Czs+vybXJPjZlt2bIFiYmJOHr0aLc+KFMfDg4OyMzMxMWLF7Fp0ybExMQgNTXV3GV1Wk1NDaKiopCQkAB3d3dzl9PlwsLC8PLLL2PEiBGYOHEikpKS0KdPH/z97383d2mdplar0bdvX3z00UcYPXo0Zs2ahbVr12LPnj3mLq3L7d27F8OGDUNISIi5S+kSqamp2Lx5M3bv3o2MjAwkJSXhxIkT2Lhxo7lL6xLvv/8+Bg4ciEcffRRSqRSLFy9GdHQ0rKz48WwI7qnpJHd3d1hbW6O4uFhreXFxMTw9Pdtdd+vWrdiyZQu++uorBAcHG7NMgxjam5WVFQIDAwEAI0aMQG5uLuLi4hAREWHMcvWmb3/5+fm4efMmpk+frlmmVqsBADY2Nrh27RoCAgKMW3QHdebvspmtrS1GjhyJvLw8Y5RoMEN68/Lygq2tLaytrTXLhgwZgrt376KxsRFSqdSoNeujM+9dXV0dEhMTsWHDBmOWaDBDelu/fj2ioqLwpz/9CQAwbNgw1NXVYf78+Vi7dm23+vA3pL8+ffrgiy++QENDA+7duwdvb2+sWrUK/v7+pijZqDw9PXW+Fo6OjkbZSwNwT02nSaVSjB49GikpKZplarUaKSkpCAsLa3O9+Ph4bNy4EcnJyRgzZowpStWbob21pFaroVAojFFip+jb36OPPoqsrCxkZmZqpmeffRZPPvkkMjMz4efnZ8ry29UV751KpUJWVha8vLyMVaZBDOlt/PjxyMvL04RQALh+/Tq8vLy6VaABOvfeHTp0CAqFAi+99JKxyzSIIb3V19e3Ci7N4VR0s0sXdua9k8vl8PHxQVNTE44cOYIZM2YYu1yjCwsL03otAODf//63Xp8fejPaIci/IomJiUImk4n9+/eLnJwcMX/+fOHs7Kz5OWxUVJRYtWqVZvyWLVuEVCoVhw8f1voZZk1NjblaaJO+vW3evFmcOXNG5Ofni5ycHLF161ZhY2MjEhISzNVCu/Ttr6Xu/OsnfXt75513xOnTp0V+fr5IT08Xs2fPFnK5XGRnZ5urhTbp21thYaFwcHAQixcvFteuXRPHjx8Xffv2Fe+++665WmiXoX+X4eHhYtasWaYuVy/69hYbGyscHBzE559/LgoKCsSZM2dEQECAmDlzprlaaJe+/X3//ffiyJEjIj8/X3z33XfiqaeeEo888oioqKgwUwdtq6mpEZcuXRKXLl0SAMS2bdvEpUuXxP/+9z8hhBCrVq0SUVFRmvHNP+l+8803RW5urti1axd/0t1TfPDBB6Jfv35CKpWKkJAQ8f3332semzhxopg7d67mfv/+/QWAVlNsbKzpC+8AfXpbu3atCAwMFHK5XLi4uIiwsDCRmJhohqo7Tp/+WurOoUYI/XpbtmyZZqyHh4eYOnVqtzxXRjN937fz58+L0NBQIZPJhL+/v9i0aZNoamoycdUdp29/V69eFQDEmTNnTFyp/vTpTalUirffflsEBAQIuVwu/Pz8xGuvvdYtP/Sb6dNfamqqGDJkiJDJZMLNzU1ERUWJ27dvm6Hqh/vmm290fnY19zN37lwxceLEVuuMGDFCSKVS4e/vb/RzJ0mE6Gb774iIiIgMwGNqiIiIyCIw1BAREZFFYKghIiIii8BQQ0RERBaBoYaIiIgsAkMNERERWQSGGiIiIrIIDDVERERkERhqiIiIyCIw1BCRxYmIiMCyZcs69RxCCMyfPx+urq6QSCTIzMzsktqIyHgYaojIqKKjo7Fu3Tpzl6G35ORk7N+/H8ePH0dRURGCgoLMXRIRPYSNuQsgIsulUqlw/PhxnDhxwtyl6C0/Px9eXl4YN25cm2MaGxshlUpNWBURtYd7aohI4/PPP4ednR2Kioo0y6KjoxEcHIyqqiq9n+/8+fOwtbXF2LFjdT4eERGBJUuWYNmyZXBxcYGHhwcSEhJQV1eH6OhoODg4IDAwEKdOndJaT6FQ4PXXX0ffvn0hl8sRHh6OixcvtlmHWq1GXFwcHnnkEdjZ2WH48OE4fPhwm+NfeeUVLFmyBIWFhZBIJBgwYICm3sWLF2PZsmVwd3dHZGQkgF/26oSHh8PZ2Rlubm545plnkJ+fr7X9+Ph4BAYGQiaToV+/fti0aVNHX0Yi6iCGGiLSmD17NgYNGoTNmzcDAGJjY/HVV1/h1KlTcHJy0vv5jh07hunTp0MikbQ55sCBA3B3d0daWhqWLFmChQsX4o9//CPGjRuHjIwM/OY3v0FUVBTq6+s166xYsQJHjhzBgQMHkJGRgcDAQERGRqK8vFznNuLi4vCPf/wDe/bsQXZ2NpYvX46XXnoJ3377rc7x77//PjZs2ABfX18UFRVpBaYDBw5AKpXi3Llz2LNnDwCgrq4OMTEx+OGHH5CSkgIrKyv87ne/g1qtBgCsXr0aW7Zswfr165GTk4PPPvsMHh4eer+eRPQQgojoAV9++aWQyWTi3XffFS4uLuLKlSuax5577jnh7Owsnn/++Q4918CBA8Xx48fbfHzixIkiPDxcc7+pqUnY29uLqKgozbKioiIBQFy4cEEIIURtba2wtbUVBw8e1IxpbGwU3t7eIj4+XvO8S5cuFUII0dDQIHr16iXOnz+vte158+aJOXPmtFnbe++9J/r379+q3pEjR7bftBCitLRUABBZWVmiurpayGQykZCQ8ND1iKhzeEwNEWl55plnMHToUGzYsAFnzpzBY489pnls6dKlePXVV3HgwIGHPk9ubi7u3LmDSZMmtTsuODhYM29tbQ03NzcMGzZMs6x5j0ZJSQmAX451USqVGD9+vGaMra0tQkJCkJub2+r58/LyUF9fjylTpmgtb2xsxMiRIx/aR0ujR49uteynn37CW2+9hf/+978oKyvT7KEpLCxEfX09FArFQ18HIuo8hhoi0pKcnIyrV69CpVK1+ookIiICqampHXqeY8eOYcqUKZDL5e2Os7W11bovkUi0ljV/ddUcFPRVW1sLADhx4gR8fHy0HpPJZHo/n729fatl06dPR//+/ZGQkABvb2+o1WoEBQWhsbERdnZ2BtVNRPrjMTVEpJGRkYGZM2di7969mDRpEtavX2/wc/3rX//CjBkzurC6XwQEBGiOaWmmVCpx8eJFDB06tNX4oUOHQiaTobCwEIGBgVqTn59fp+u5d+8erl27hnXr1mHSpEkYMmQIKioqNI8PHDgQdnZ2SElJ6fS2iKh93FNDRACAmzdvYtq0aVizZg3mzJkDf39/hIWFISMjA6NGjdLruUpKSvDDDz/g2LFjXV6nvb09Fi5ciDfffBOurq7o168f4uPjUV9fj3nz5rUa7+DggDfeeAPLly+HWq1GeHg4qqqqcO7cOTg6OmLu3LmdqsfFxQVubm746KOP4OXlhcLCQqxatUrzuFwux8qVK7FixQpIpVKMHz8epaWlyM7O1tS7c+dOHD16lMGHqJMYaogI5eXlePrppzFjxgzNB3JoaCh++9vfYs2aNUhOTtbr+b788kuEhITA3d3dGOViy5YtUKvViIqKQk1NDcaMGYPTp0/DxcVF5/iNGzeiT58+iIuLQ0FBAZydnTFq1CisWbOm07VYWVkhMTERr7/+OoKCgjB48GDs2LEDERERmjHr16+HjY0N3nrrLdy5cwdeXl5YsGCB5vGysjKtn4ATkWEkQghh7iKIqOdITU3Fzp072z3Py7PPPovw8HCsWLHChJUR0a8d99QQUYdNnjwZly9fRl1dHXx9fXHo0CGEhYW1GhceHo45c+aYoUIi+jXjnhoiIiKyCPz1ExEREVkEhhoiIiKyCAw1REREZBEYaoiIiMgiMNQQERGRRWCoISIiIovAUENEREQWgaGGiIiILAJDDREREVkEhhoiIiKyCAw1REREZBH+D2bJ8d+aQKuGAAAAAElFTkSuQmCC", "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.14" } }, "nbformat": 4, "nbformat_minor": 5 }