{ "cells": [ { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "from feos.dft import ExternalPotential, HelmholtzEnergyFunctional, Adsorption1D, Geometry, Pore1D, State\n", "from feos.pcsaft import PcSaftParameters\n", "from feos.si import ANGSTROM, KELVIN, BAR, NAV, KILO, METER, MOL, SIArray1\n", "\n", "import pandas as pd\n", "import numpy as np\n", "import os\n", "import json\n", "from scipy.spatial import Voronoi, distance\n", "from itertools import permutations" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "class SolidStructure:\n", "\n", " N_atom = 0\n", " dimensions = {\n", " \"Lx\": 0.0,\n", " \"Ly\": 0.0,\n", " \"Lz\": 0.0\n", " }\n", " forcefield = ''\n", " total_mass = 0.0\n", "\n", " def __init__(self,structure):\n", " \n", " self.name = structure\n", "\n", " database = os.path.join(os.getcwd(),'structure_parameters','solid_database.json')\n", " \n", " with open(database) as f:\n", " data = json.load(f)\n", "\n", " for i in range(len(data)):\n", " if data[i]['Name']==self.name:\n", " self.N_atom = data[i]['N_atom']\n", " self.dimensions = data[i]['Dimensions']\n", " self.forcefield = data[i]['Forcefield']\n", " \n", " \n", " def read_structure(self):\n", " struct_param = os.path.join(os.getcwd(),'structure_parameters')\n", " structure_df = pd.read_csv(os.path.join(struct_param,'{}.dat'.format(self.name)),names=['x','y','z','Type'], delim_whitespace=True)\n", " \n", " filename = os.path.join(struct_param,'{}.dat'.format(self.forcefield))\n", " forcefield_df = pd.read_csv(filename,names=['Type','sigma','epsilon','mass'], delim_whitespace=True)\n", "\n", " coordinates = np.array([structure_df[\"x\"], structure_df[\"y\"], structure_df[\"z\"]])\n", " \n", " sigma_ss = np.zeros(len(structure_df[\"x\"]))\n", " epsilon_k_ss = np.zeros_like(sigma_ss)\n", " \n", " for i in range(len(sigma_ss)):\n", " sigma_ss[i] = forcefield_df.sigma[forcefield_df.Type==structure_df[\"Type\"][i]]\n", " epsilon_k_ss[i] = forcefield_df.epsilon[forcefield_df.Type==structure_df[\"Type\"][i]]\n", " \n", " self.total_mass = np.sum(np.array([forcefield_df.mass[forcefield_df.Type==t]for t in structure_df.Type])) * 1.66054e-27\n", "\n", " return coordinates, sigma_ss, epsilon_k_ss " ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "structure = SolidStructure('LTA')\n", "coordinates, sigma_ss, epsilon_ss = structure.read_structure()\n", "system_size = [structure.dimensions[\"Lx\"] * ANGSTROM, structure.dimensions[\"Lx\"] * ANGSTROM, structure.dimensions[\"Lx\"] * ANGSTROM]" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "6.84534814234455 3220\n", "[5.9595, 5.9594999999999985, 5.9595]\n", "CPU times: user 328 ms, sys: 3.58 ms, total: 331 ms\n", "Wall time: 331 ms\n" ] } ], "source": [ "%%time\n", "vor = Voronoi(coordinates.T)\n", "\n", "minimum_distance = 0.0\n", "dist = np.zeros(int(vor.vertices.size / 2))\n", "for i,point in enumerate(vor.vertices):\n", " if(point[0] >= 0.0 and point[0] <= structure.dimensions[\"Lx\"] and point[1] >= 0.0 and point[1] <= structure.dimensions[\"Ly\"]):\n", " dist[i] = np.min(distance.cdist(np.array([point]),coordinates.T))\n", "pore_radius = np.max(dist)\n", "print(pore_radius,np.argmax(dist))\n", "pore_center = [vor.vertices[np.argmax(dist)][0], vor.vertices[np.argmax(dist)][1], vor.vertices[np.argmax(dist)][2]]\n", "print(pore_center)" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "params = PcSaftParameters.from_json(['argon'],'structure_parameters/noble_gases.json')\n", "func = HelmholtzEnergyFunctional.pcsaft(params)" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "potential = ExternalPotential.FreeEnergyAveraged(coordinates * ANGSTROM, sigma_ss, epsilon_ss, pore_center, system_size, [51, 51])\n", "pore = Pore1D(Geometry.Spherical, 5.9595 * ANGSTROM, potential, 128)" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [], "source": [ "isotherm = Adsorption1D.adsorption_isotherm(func, 298.0 * KELVIN, SIArray1.linspace(1.0e-3 * BAR, 5.0 * BAR, 51), pore)" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAA30AAAFGCAYAAADeqPb+AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAACpWElEQVR4nOzdd3jUVdbA8e+d9N57o4QSOklogoAgguiigCjWxbIurGvfVfddV9d1d9W1rw2xwOpaYVEsFEUBKdJ7CiQBQhLSe09m5r5/BCIlZRKSTMr5PA+PzMz93d8ZVDJn7r3nKK01QgghhBBCCCG6J4O1AxBCCCGEEEII0X4k6RNCCCGEEEKIbkySPiGEEEIIIYToxiTpE0IIIYQQQohuTJI+IYQQQgghhOjGJOkTQgghhBBCiG6s3ZM+pdT7SqkcpdThRl5XSql/K6WSlVIHlVLR7R2TEEII0ZUopR5USsUppQ4rpT5RSjlaOyYhhBBdR0es9C0DZjTx+pVAv9O/7gbe6oCYhBBCiC5BKRUC3AfEaq2HADbAfOtGJYQQoitp96RPa/0TUNDEkGuAD3Sd7YCnUiqoveMSQgghuhBbwEkpZQs4A6esHI8QQogupDOc6QsB0s56nH76OSGEEKLH01pnAC8AJ4FMoFhr/Z11oxJCCNGV2Fo7AEA18JxucKBSd1O3BRQXF5eYgQMHtmdcQgghOok9e/bkaa39rB2HNSilvKjbFdMbKAKWK6Vu0Vr/97xx8jNSCCF6GEt/PnaGpC8dCDvrcSiNbFvRWi8BlgDExsbq3bt3t390QgghrE4plWrtGKzocuC41joXQCm1ErgEOCfpk5+RQgjR81j687EzbO/8CrjtdBXPsdRtW8m0dlBCCCFEJ3ESGKuUclZKKWAqkGDlmIQQQnQh7b7Sp5T6BJgM+Cql0oEnATsArfViYDUwE0gGKoDb2zsmIYQQoqvQWu9QSq0A9gJGYB+nV/SEEEIIS7R70qe1vrGZ1zVwT3vHIYQQQnRVWusnqfvSVAghhGixznCmTwghhBBCCNHN1NbWkp6eTlVVlbVD6fIcHR0JDQ3Fzs6uVddL0ieEEEIIIYRoc+np6bi5udGrVy/qjiSL1tBak5+fT3p6Or17927VHJ2hkIsQQgghhBCim6mqqsLHx0cSvouklMLHx+eiVkwl6RNCCCGEEEK0C0n42sbF/jlK0ieEEEIIIYToltauXcuAAQOIjIzk2WefbXCM1pr77ruPyMhIhg0bxt69e5u9fvny5QwePBiDwUBjfVE3btzI1Vdf3bZvqJUk6RNCCCGEEEJ0OyaTiXvuuYc1a9YQHx/PJ598Qnx8/AXj1qxZQ1JSEklJSSxZsoRFixY1e/2QIUNYuXIlEydObNf424okfUIIIYQQQohuZ+fOnURGRtKnTx/s7e2ZP38+q1atumDcqlWruO2221BKMXbsWIqKisjMzGzy+qioKAYMGNBsDCUlJcyePZtBgwaxcOFCzGYzAIsWLSI2NpbBgwfz5JO/dOTp1asXf/vb35gwYQLLly9voz8Jqd4phBBCCCGE6IYyMjIICwurfxwaGsqOHTssGpeRkWHx9U3ZuXMn8fHxREREMGPGDFauXMl1113HP/7xD7y9vTGZTEydOpWDBw8ybNgwoK49w5YtW1r6dpskSZ8QQgghhBDtIDf1OD//7xNiZl5LyMBB1g7Hqp76Oo74UyVtOuegYHee/NXgRl/XWl/wXEMFURobZ+n1TRk9ejR9+vQB4MYbb2TLli1cd911fP755yxZsgSj0UhmZibx8fH1Sd8NN9zQontYQpI+IYQQQggh2ljyru2sfu0FaqurSNm9k6l3LGTY5TOsHVaPEhoaSlpaWv3j9PR0goODLR5XU1Nj0fVNOT9JVEpx/PhxXnjhBXbt2oWXlxcLFiw4px2Di4tLi+5hCUn6hBBCCCGEaCNaa3Z+uZwtn35AYN9+TF/0AJv++z7fv/M6OSdSmHLHQgwGG2uH2eGaWpFrL6NGjSIpKYnjx48TEhLCp59+yscff3zBuFmzZvH6668zf/58duzYgYeHB0FBQfj5+Vl0fVN27tzJ8ePHiYiI4LPPPuPuu++mpKQEFxcXPDw8yM7OZs2aNUyePLmN3nXDJOkTQgghhBCijaQe2MuWTz9g4PhJXLHwPuzsHZj96BNs/vg/7P56JWGDhzFg3KXWDrNHsLW15fXXX2f69OmYTCbuuOMOBg+uSz4XL14MwMKFC5k5cyarV68mMjISZ2dnli5d2uz1X3zxBffeey+5ublcddVVjBgxgnXr1l0Qw7hx43jsscc4dOgQEydOZPbs2RgMBkaOHMngwYPp06cP48ePb/c/C9XQXtWuIDY2VjfWE0MIIUT3opTao7WOtXYcXYX8jBTCej7762MU5WRx17/fwcbWrv55bTbzzu/vxC+iF7MffbKJGbqPhIQEoqKirB1Gt9HQn6elPx+lZYMQQgghhBBtIONIAukJh4m9avY5CR+AMhiImjCJ4/v3UF5UaKUIRU8lSZ8QQgghhBBtYOeXn+Po5s6wqdMbfH3QxClos5kj237q4MhETydJnxBCCCGEEBcpN/U4x/buIvrKX2Hn6NjgGJ/QcAL6RBL3048dHJ3o6STpE0IIIYQQ4iLtXLUCO0cnRk7/VZPjBk2cQs7xFPLSUjsoMiEk6RNCCCGEEOKilObncWTbZoZPuxJHV9cmxw68ZCLKYCB+84YOik4ISfqEEEIIIYS4KEk7t6G1maFTGj7LdzZnD096j4ghYfMGzGZTB0QnhCR9QgghhBBCXJSj27fiG94L7+AQi8ZHTZhMWUE+WclH2zkycccdd+Dv78+QIUMaHaO15r777iMyMpJhw4axd+/eBse5NrOK25lJ0ieEEEIIIUQrlRcVknEknn6jL7H4mtCougREkr72t2DBAtauXdvkmDVr1pCUlERSUhJLlixh0aJFbR6H0Whs8zlbQpI+IYQQQgghWilp58+gNf3Hjrf4GldvH1y9vMlKSWrHyATAxIkT8fb2bnLMqlWruO2221BKMXbsWIqKisjMzGxw7MMPP0x0dDRTp04lNzcXgHfeeYdRo0YxfPhw5s6dS0VFBVCXcD700ENcdtllPProo237xlpIkj4hhBBCCCFaKWnHFryDQ/EJDW/RdQF9+0vS10lkZGQQFhZW/zg0NJSMjIwLxpWXlxMdHc3evXuZNGkSTz31FABz5sxh165dHDhwgKioKN577736a44ePcr69et58cUX2/+NNMHWqncXQgghhBCii6ooKSYt7jCjr52HUqpF1wb27UfK7u1UlZfh6NJ1z4pZbM1jkHWobecMHApXPnvR02itL3iuoX+fBoOBG264AYBbbrmFOXPmAHD48GEef/xxioqKKCsrY/r0Xwr6zJs3Dxsbm4uO8WLJSp8QQgghhBCtkLxrO1qbW7S184zAvv0AyD6W3NZhiRYKDQ0lLS2t/nF6ejrBwcHNXncmMVywYAGvv/46hw4d4sknn6Sqqqp+jIuLS9sH3Aqy0ieEEEJ0YkqpAcBnZz3VB3hCa/2KdSISQpyRtGMrngFB+EX0bvG1AaeTvqyUJCKGjmjjyDqhNliRay+zZs3i9ddfZ/78+ezYsQMPDw+CgoIuGGc2m1mxYgXz58/n448/ZsKECQCUlpYSFBREbW0tH330ESEhllVx7UiS9AkhhBCdmNb6CDACQCllA2QAX1gzJiEEVFeUc/LwAWKuurbFWzsBnFzd8AwIIlvO9bWrG2+8kY0bN5KXl0doaChPPfUUd955J4sXLwZg4cKFzJw5k9WrVxMZGYmzszNLly5tcC4XFxfi4uKIiYnBw8ODzz6r+z7u6aefZsyYMURERDB06FBKS0s77P1ZSpI+IYQQouuYCqRorVOtHYgQPd3Jwwcwm0z0GTmq1XME9O3HqSMJbRiVON8nn3zS4PMLFy6s/71SijfeeKPZucrKyoC6JO9sixYtarDNw7Jly1oQafuSM31CCCFE1zEfaPgTjBCiQ6Ue3I+doxNB/Qe0eo7Avv0ozc+lvKiwDSMT4kKS9AkhhBBdgFLKHpgFLG/k9buVUruVUrvP9I4SQrSf1EP7CBs0BBtbu1bPEXjWuT4h2pMkfUIIIUTXcCWwV2ud3dCLWuslWutYrXWsn59fB4cmRM9SnJNFUVYmEcNGXtQ8/r37opRBkj7R7iTpE0IIIbqGG5GtnUJ0CqkH9wNcdNJn7+iEd0go2SlH2yAqIRonSZ8QQgjRySmlnIFpwEprxyKEgNSD+3D18cU7OPSi5wrs25+slKQGG4QL0VYk6RNCCCE6Oa11hdbaR2tdbO1YhOjpzGYTJw8fIGLoiFa1ajhfYN9+VJaWUJonZ3FF+5GkTwghhBBCCAtlH0umqryMXhe5tfMMn7BwAApOpbfJfOIXaWlpXHbZZURFRTF48GBeffXVBsdprbnvvvuIjIxk2LBh7N27t8Fxrq6u7Rluu5KkTwghhBBCCAudOc8XPnREm8znFRgMQGFmRpvMJ35ha2vLiy++SEJCAtu3b+eNN94gPj7+gnFr1qwhKSmJpKQklixZ0mDPvYtlNBrbfM6WkKRPCCGEEEIIC6Ue3Id/7744u3u0yXwuXt7YOThSmHWqTeYTvwgKCiI6OhoANzc3oqKiyMi4MLletWoVt912G0opxo4dS1FREZmZmQ3O+fDDDxMdHc3UqVM50x7nnXfeYdSoUQwfPpy5c+dSUVEBwIIFC3jooYe47LLLePTRR9vpXVpGkj4hhBBCCCEsUFtVxamjiUS00SofgFIKz6BgijIl6WtPJ06cYN++fYwZM+aC1zIyMggLC6t/HBoa2mByWF5eTnR0NHv37mXSpEk89dRTAMyZM4ddu3Zx4MABoqKieO+99+qvOXr0KOvXr+fFF19sh3dlOVur3l0IIYQQQoguIjP5CGaTkbBBQ9t0Xq/AYHJOpLTpnJ3NczufI7EgsU3nHOg9kEdHN7+CVlZWxty5c3nllVdwd3e/4PWGKqc2VKTHYDBwww03AHDLLbcwZ84cAA4fPszjjz9OUVERZWVlTJ8+vf6aefPmYWNjY/F7ai+y0ieEEEIIIYQF0hPiQCmCB0S16bxeQcEU52RjsvK5r+6otraWuXPncvPNN9cnaecLDQ0lLS2t/nF6ejrBwcHNzn0mMVywYAGvv/46hw4d4sknn6Sqqqp+jIuLy0W+g7YhK31CCCGEEEJYICMxDr+I3jg4t+0Hec/AYLTZTHFONt7BIW06d2dhyYpcW9Nac+eddxIVFcVDDz3U6LhZs2bx+uuvM3/+fHbs2IGHhwdBQUEXjDObzaxYsYL58+fz8ccfM2HCBABKS0sJCgqitraWjz76iJCQzvfvUJI+IYQQQgghmmEyGjmVlMjQy65o87m9guqShKKsU9026bOGrVu38uGHHzJ06FBGjBgBwD//+U9mzpzJ4sWLAVi4cCEzZ85k9erVREZG4uzszNKlSxucz8XFhbi4OGJiYvDw8OCzzz4D4Omnn2bMmDFEREQwdOhQSktLO+T9tYQkfUIIIYQQQjQj50QKxupqQgYObvO5vYLOtG2QYi5tacKECQ2e14O6ZO8MpRRvvPFGs/OVlZUBdUne2RYtWtRgm4dly5a1INr2JWf6hBBCCCGEaEZGQhwAoVFtn/Q5ubnj4OwibRtEu5GkTwghhBBCiGakJ8bjGRiEi6dXm8+tlMIzMFgatIt2I0mfEEIIIYQQTdBmMxlH4ttla+cZXkHBFMlKn2gnHZL0KaVmKKWOKKWSlVKPNfC6h1Lqa6XUAaVUnFLq9o6ISwghhBBCiObkZ6RRVVpCaDsnfSV5uRhratrtHqLnavekTyllA7wBXAkMAm5USg06b9g9QLzWejgwGXhRKWXf3rEJIYQQQgjRnIzEuvN8Ie1wnu8Mr8Bg0JrinKx2u4fouTpipW80kKy1Pqa1rgE+Ba45b4wG3FRdh0NXoACQ7pRCCCGEEMLq0hPicPH0wjPgwt5tbcXzdAXPAjnXJ9pBRyR9IUDaWY/TTz93tteBKOAUcAi4X2ttPn8ipdTdSqndSqndubm57RWvEEKITqK0qpaXvz9q7TCEED1cRmLdeb669Yn24RV4uleftG1oM1VVVYwePZrhw4czePBgnnzyyQbHaa257777iIyMZNiwYezdu7fBca6uru0ZbrvqiD59Df3fcX7DjOnAfmAK0Bf4Xim1WWtdcs5FWi8BlgDExsY23HRDCCFEl1dZY+I/P59g8aYUiipqrR2OEKIHK8nNoTQ/l9hfzWnX+zi6uuLk5i5tG9qQg4MDP/74I66urtTW1jJhwgSuvPJKxo4de864NWvWkJSURFJSEjt27GDRokXs2LGjTWMxGo3Y2lqvRXpHrPSlA2FnPQ6lbkXvbLcDK3WdZOA4MLADYhNCCNGJVBtN/GfbCSY+v4Fn1yQyPNSTr38/wdphCSF6sPTE9uvPdz7PoGBZ6WtDSqn61bna2lpqa2sbXK1dtWoVt912G0opxo4dS1FREZmZmQ3O+fDDDxMdHc3UqVM5s/PwnXfeYdSoUQwfPpy5c+dSUVEBwIIFC3jooYe47LLLePTRR9vpXVqmI5K+XUA/pVTv08VZ5gNfnTfmJDAVQCkVAAwAjnVAbEIIIToBo8nM57vSmPLCJp78Ko7evi4sXziO/9wxmqGhHtYOTwjRg2UkxGHv5IxveES738tLevW1OZPJxIgRI/D392fatGmMGTPmgjEZGRmEhf2yRhUaGkpGxoX/HsrLy4mOjmbv3r1MmjSJp556CoA5c+awa9cuDhw4QFRUFO+99179NUePHmX9+vW8+OKL7fDuLNfua4xaa6NS6vfAOsAGeF9rHaeUWnj69cXA08AypdQh6raDPqq1zmvv2IQQQliX2az59lAmL39/lGN55QwL9eCZOUO5tJ9vu56dEUIIS6UnxhEyIAqDwabd7+UVFEL8Tz9SW1WFnaNju9+vI2X9859UJyS26ZwOUQMJ/L//a3KMjY0N+/fvp6ioiNmzZ3P48GGGDBlyzhitLzw11tDPIIPBwA033ADALbfcwpw5dVt+Dx8+zOOPP05RURFlZWVMnz69/pp58+ZhY9P+/+00p0M2lmqtVwOrz3tu8Vm/PwVc0RGxCCGEsD6tNT8m5vDCd0dJyCyhf4Ari2+JYfrgAEn2hBCdRkVJMQUZaQy69LIOuZ9HQCAAJXk5+ISGd8g9ewpPT08mT57M2rVrL0j6QkNDSUv7pe5keno6wcHBzc555ufVggUL+PLLLxk+fDjLli1j48aN9WNcXFza5g1cJOudJhRCCNEj/ZySz/PrEtl7sohwb2deuWEEvxoejI1Bkj0hROeScSQeaN/+fGdz9/EDoCQvt9slfc2tyLWH3Nxc7Ozs8PT0pLKykvXr1zd4tm7WrFm8/vrrzJ8/nx07duDh4UFQ0IXtOcxmMytWrGD+/Pl8/PHHTJhQd+a8tLSUoKAgamtr+eijjwgJOb9RgfVJ0ieEEKJDHEgr4oXvjrA5KY9Ad0f+OXso82JDsbPpiOPlQgjRchkJcdjY2RHYt3+rrk+Ny2ffdycZMjGEyBj/Zse7+dYlfaV50pqsLWRmZvLrX/8ak8mE2Wzm+uuv5+qrrwZg8eK6TYcLFy5k5syZrF69msjISJydnVm6dGmD87m4uBAXF0dMTAweHh589tlnADz99NOMGTOGiIgIhg4dSmlpace8wRaQpE8IIUS7OppdyovfHWFdXDbeLvY8flUUt4yNwNHO+mcchBCiKRmJcQT27Y+tnV2Lrqsqr2XriiQSf87C1s5AxpFCju3zZ+L8ATi6Nj6Xq5c3ymCgRJK+NjFs2DD27dvX4GsLFy6s/71SijfeeKPZ+crKyoC6JO9sixYtYtGiRReMX7ZsWQuibV+S9AkhhGgXJ/MreHn9Ub7cn4GrvS0PTevPHRN64+ogP3qEEJ1fTVUl2cdTGH3NdS26rqywmhXP7aaipIaYGRHEXNmLAz+ksevb46QfLeK6R2Jw93Vq8FqDjQ2uXj6U5kvSJ9qW/OQVQgjRprJLqnjtxyQ+3ZmGrY3i7ol9WDixL14u9tYOrctSSnkC7wJDAA3cobX+2apBCdHNZR49gjabCRnYsvN821YmU1VWy9xHYgjo5Q5A7MxeRAz1YcVzu9m/Po2J8xvfLurm6yfbO0Wbk6RPCCFEmygsr2HxphSWbTuByay5cXQ4v58SSYB79yo7biWvAmu11ted7nnrbO2AhOju0hPjUMpAcP8oi685lVxE0q5sYmf2qk/4zvALc6P/qAAStp1i9K964+jS8DZPd18/MpOPXFTsQpxPkj4hhBAXpazayHubj/PO5mOU1xiZPTKEB6b2J9xH8pK2oJRyByYCCwC01jVAjTVjEqInyEiMwy+iNw7Olv1dZjZrNn92FFcvB6JnNNzIffjUcBJ/ziJucwYxM3o1OMbN14+j27eizWaUQQpdibYhSZ8QQohWqao18d/tqby5MYWC8hqmDw7g4SsG0D/AzdqhdTd9gFxgqVJqOLAHuF9rXW7dsITovkzGWjKTjjB0quVtpOO3nCIvrYzpvxmCnX3Dhap8Q10JHejFoQ3pjLg8HBvbC5M6Nx9fzCYjFSXFuHh6tfo9CHE2+fpACCFEixhNZj7bdZIpL2zk798mMCjInVX3jOftW2Ml4WsftkA08JbWeiRQDjx2/iCl1N1Kqd1Kqd25uXIeSIiLkX0sBWNNNaEWnuerqTSyfVUKIf096Rvt1+TYEdPCKS+uIXl3doOvu/ue6dWX07KghWiCJH1CCCEsYjZrvjl4iite/olH/3cIP3dHPr5rDP+9awzDwzytHV53lg6ka613nH68grok8Bxa6yVa61itdayfX9MfOoUQTctIjAOwuIjL0Z1ZVJcbGTu7L0qpJseGD/LGK8iF/T+kobW+4HU3H+nV19ZMJhMjR46s79F3Pq019913H5GRkQwbNoy9e/c2OM7V1bU9w2xXsr1TCCFEk7TWbDqay/PrjhB3qoT+Aa4suTWGaYMCmv1wIy6e1jpLKZWmlBqgtT4CTAXirR2XEN1ZemIcXkHBFm2v1FpzePMpfMNcLyje0hClFMOnhLLxoyPkpJZecI27b10Td+nV13ZeffVVoqKiKCkpafD1NWvWkJSURFJSEjt27GDRokXs2LGjwbGtZTQasbW1XuolK31CCCEatftEATe8vZ0FS3dRXFnLS9cPZ839E7licKAkfB3rXuAjpdRBYATwT+uGI0T3pc1mTiXGW7zKl32ihPz0MgZfGmLx34u9h9et5qUnFlzwmoOLC3YOjpTm51ketGhUeno63377LXfddVejY1atWsVtt92GUoqxY8dSVFREZmZmg2MffvhhoqOjmTp1Kme20r/zzjuMGjWK4cOHM3fuXCoqKgBYsGABDz30EJdddhmPPvpo27+5FpCkTwghxAUSMku4c9kurlv8M8fyyvnbNYP58eHJzIkOxcYgyV5H01rvP711c5jW+lqtdaG1YxKiu8pPP0lVeZnFSV/cTxnYOdjQf3SAxfdwdrfHO9iFjKNFF7ymlJJefW3ogQce4F//+heGJiqhZmRkEBYWVv84NDSUjIyMC8aVl5cTHR3N3r17mTRpEk899RQAc+bMYdeuXRw4cICoqCjee++9+muOHj3K+vXrefHFF9vwXbWcbO8UQghRLzW/nJe+P8pXB07h6mDLH6cP4PbxvXC2lx8XQoieIT2xbve0JUVcqitqSd6dw4Cxgdg7tuzvyZD+XiRsO4XJaL6giqe7r1+32965+fOj5KWVtemcvmGuXHp9443uv/nmG/z9/YmJiWHjxo2NjmvobGVDq7YGg4EbbrgBgFtuuYU5c+YAcPjwYR5//HGKioooKytj+vTp9dfMmzcPG5uGq7l2JPkpLoQQguySKv79QxKf7UrD1kbx24l9WTipD57O9tYOTQghOlRGYhwuXt54BAQ2O/bIjiyMtWYGXxrS4vuEDPDk0MZ0ck6UEBTpec5rbr5+5Jw41uI5xbm2bt3KV199xerVq6mqqqKkpIRbbrmF//73v+eMCw0NJS0trf5xeno6wcHBzc5/JjFcsGABX375JcOHD2fZsmXnJJguLi5t82YukiR9QgjRgxVX1LL4pxSWbj2O0aSZPzqM+6b0w9/d0dqhCSFEh9Nak54YR8jAwc2ez9Nac/inU/hHuOEX3vJ2NSH9vEBBxtHCC5M+H18qiosw1tRga989vnxrakWuvTzzzDM888wzAGzcuJEXXnjhgoQPYNasWbz++uvMnz+fHTt24OHhQVBQ0AXjzGYzK1asYP78+Xz88cdMmDABgNLSUoKCgqitreWjjz4iJKTlXwK0N0n6hBCiB6qsMbF023EWb0yhtNrINcODeXBafyJ8Osc3kkIIYQ0luTmU5ecROnBQs2NzT5ZSmFnO5JsHtOpejq52+IS4kn6kiNiZ5752poJnaUEeXoHNrziJllu8eDEACxcuZObMmaxevZrIyEicnZ1ZunRpg9e4uLgQFxdHTEwMHh4efPbZZwA8/fTTjBkzhoiICIYOHUppaWmHvQ9LSdInhBA9SK3JzKe70vj3D0nkllYzZaA/f7hiAIOCmy8zLoQQ3V1L+vOl7M3BYFD0jfZv9f1C+3txeHMGplozNna/nOs7u1efJH1tY/LkyUyePLn+8cKFC+t/r5TijTfeaHaOsrK6M4lPP/30Oc8vWrSIRYsWXTB+2bJlrQu2HUjSJ4QQPYDZrPn64Cle+v4oqfkVjOrlxZs3RzOql7e1QxNCiE4jPTEOB2cXfMMjmhyntSZ5Tw6hUV44uti1+n4hAzw58GMa2SeKCe73S09Ad9+6pK+7FXMR1iNJnxBCdGNaazYezeVfa4+QkFnCwEA3li4YxeQBftJnTwghzpOREEfwgCgMhqarLeaeLKUkr4qYK3td1P2C+3miFKQfKTon6XP18QWQtg2izUjSJ4QQ3dSe1EKeW5vIzuMFhHs78+r8EfxqWDAG6bMnhBAXqCgppuBUOoMmTW12bPKeuq2dfUb4XdQ9HZzt8A1zI+NIIVzdu/55Wzs7nD08Kc2XpE+0DUn6hBCimzmaXcrz647wfXw2vq4OPH3NYG4YFY69beONaYUQoqc7c56vuf58WmtS9l781s4zQgZ4cXBDGsYaE7b2v6wwdsdefcJ6JOkTQohuIr2wgpe/T2LlvnRc7aWxuhBCtERGYhw2dnYE9O3X5LgzWztjZ/Zqk/sGR3qw//uT5KWXEdjHo/55N18/8tNOtsk9hJBPAkII0cXll1XzxoYU/rs9FRT85tI+LJrUFy+X7tHbSQghOkJ6QjxBkQOwtWt69e7M1s7ewy9ua+cZPiGuAORnnJv0ufv6cXz/HrTWcgZbXDRJ+oQQoosqrzby7ubjvLP5GBU1RubFhPHAtH4EeThZOzQhhOhSaqoqyTmRwuhr5jU57petnd5tsrUTwM3HETtHG/LTy8573g9jdTVV5WU4uba8+buo06tXL9zc3LCxscHW1pbdu3dfMEZrzf3338/q1atxdnZm2bJlREdHXzDO1dW1vm1DVyNJnxBCdDE1RjOf7DzJaz8mkVdWw4zBgfxhen8i/eVDgRBCtMapo4los7nZpuz5GeV1VTtn9Gqzeyul8A1xJf9U+TnPu3r7AFBWkC9J30XasGEDvr6+jb6+Zs0akpKSSEpKYseOHSxatIgdO3a0aQxGoxFbW+ulXpL0CSFEF2E2a746cIoXvz9CWkElY/t4885tAxkZ7tX8xUIIIRqVkRiHUgaC+kc1OS71cB4AEUN92vT+3iGuJO/OPmcrp4tXXR/V8oJ8/MJ7ten9xLlWrVrFbbfdhlKKsWPHUlRURGZmJkFBQReMffjhh9mwYQNeXl58+umn+Pn58c4777BkyRJqamqIjIzkww8/xNnZmQULFuDt7c2+ffuIjo7mxRdftMK7qyOl3IQQopPTWrPhSA5XvbaFBz7bj5uDHf+5YzSf/GasJHxCCNEGMhLi8OvVGwdn5ybHpR7Kxy/cDRcPhza9v2+IC9UVRsoKq+ufczu90ldamN+m9+pplFJcccUVxMTEsGTJkgbHZGRkEBYWVv84NDSUjIyMC8aVl5cTHR3N3r17mTRpEk899RQAc+bMYdeuXRw4cICoqCjee++9+muOHj3K+vXrrZrwgaz0CSFEp7b3ZCHPrUlkh/TaE0KIdmEy1pKZdIRhl89oclxVWS1Zx4ovuiF7Q84u5uLm7QiAi9cv2zu7gw3LlpCTeqxN5/SP6MNlC+5ucszWrVsJDg4mJyeHadOmMXDgQCZOnHjOGK31Bdc1VDzHYDBwww03AHDLLbcwZ84cAA4fPszjjz9OUVERZWVlTJ8+vf6aefPmYWNjc8FcHU2SPiGE6ISSc8p4fl0i6+Ky8XW152/XDGa+9NoTQog2l30sGWNtDSFRTffnOxmfj9Ztv7UT6rZ3Ql3S12to3dkzWzs7HN3cKS8saPP79STBwcEA+Pv7M3v2bHbu3HlB0hcaGkpaWlr94/T09PrrmnImMVywYAFffvklw4cPZ9myZWzcuLF+jIuLSxu8i4snSZ8QQnQimcWVvPJ9Esv3pOFkZ8ODl/fnrkt74+Igf10LIUR7SE+oa8oeMqDpIi4nDuXj5GZHQIR7m8fg4GSLm7cj+RnnFnNx8/KmtJus9DW3ItceysvLMZvNuLm5UV5eznfffccTTzxxwbhZs2bx+uuvM3/+fHbs2IGHh0eD5/nMZjMrVqxg/vz5fPzxx0yYMAGA0tJSgoKCqK2t5aOPPiIkJKTd31tLyacIIYToBIoranlzUzLLtp7ArDW/vqQXv78sEh/Xtj03IoQQ4lwZiXF4BYXg4tn4GWmzWXMyPp/eQ31R7bS93ifUlfyMc9sBuHj7yErfRcjOzmb27NlAXfXMm266iRkz6rbxLl68GICFCxcyc+ZMVq9eTWRkJM7OzixdurTB+VxcXIiLiyMmJgYPDw8+++wzAJ5++mnGjBlDREQEQ4cOpbS0tAPeXctI0ieEEFZUVWviP9tO8MaGZEqrjVw7IoSHpvUnzLvpYgJCCCEuntlsIuNIPP3HjG9yXPaxYqrLjUQMbbzs/8XyCXbh5OF8TLVmbOzqtvK7evmQe6Jtz8H1JH369OHAgQMNvrZw4cL63yuleOONN5qd70yPvqeffvqc5xctWsSiRYsuGL9s2bIWRNu+JOkTQggrMJrMrNybwcvrj5JZXMXkAX48Mn0gg4LbftuQEEKIhuWmnqC6vJywQUObHHfiUD4GgyJskHe7xeIT6orZrCnMLsc3tK4vn6u3D+XFRZiMRmys2ONNdH3yX48QQnQgrTXfx2fzr3VHSM4pY3iYJy9dP4Jxfdu+MIAQQoimpccfBiC0maQv9XAeQZEeODi130fn+gqe6WW/JH1e3qA1FcVFuPm03yqj6P4k6RNCiA6y+0QBz6xJZE9qIX18XXjr5mhmDAlssCy0EEKI9pcWfwjPgKAmE6rSgiryM8q5ZE5ku8bi6e+Eja3hnGIurt6/tG2QpE9cDEn6hBCinSVll/Lc2iOsT8jGz82Bf8wewvWxYdjZSPsFIYSwFm02k5FwmMjRlzQ5LvVwXfXM9mjVcDaDjQGvIOdzirmcnfR1VVpr+XKzDTTUS7AlJOkTQoh2cnb7BRd7W/44fQC3j++Fs7381SuEENaWe/IEVeVlhA1uZmvnoTzcfR3xCmz/Alu+Ia6cTPilWqerV90ZwrLCrpn0OTo6kp+fj4+PjyR+F0FrTX5+Po6Ojq2eQz55CCFEGyuurOWtjSks3XocreH28b2557JIvF3srR2a6KKUUieAUsAEGLXWsdaNSIiuLz3+EAChUUMaHWOsMZGeWEjU+OAOSVq8Q1xJ3J5FVVktjq52OLt7YLCx7bIrfaGhoaSnp5Obm2vtULo8R0dHQkNDW329JH1CCNFGqmpNfPhzKq9vSKakqpbZI0J4UNoviLZzmdY6z9pBCNFdpMUfwiMgEHdfv0bHZBwtwlhrbvetnWecWU0syqkg0NUDZTDg4uVFWRft1WdnZ0fv3r2tHYZAkj4hhLhoJrPmy30ZvPT9UTKKKpnU349HZgxgcLCHtUMTQgjRAG02kx5/mMjR45ocl3ooD1t7AyH9PTskLg8/JwCKcyoI7FP3M8TVy7vLrvSJzkOSPiGEaCWtNZuO5vLsmkQSs0oZGuLB89cN45JIqbAm2pwGvlNKaeBtrfUSawckRFdWf56viVYNWmtOHMondKA3tnY2HRKXu68TSkFRTmX9c67ePuSnneyQ+4vuS5I+IYRohYPpRTyzOpGfj+UT7u3MazeO5KqhQRgMclBdtIvxWutTSil/4HulVKLW+qezByil7gbuBggPD7dGjEJ0GekJZ/rzNX6eryCznNKCKmKujOiosLCxNeDm40hx7llJn5cPqQf3d1gMonuSpE8IIVogNb+c59cd4ZuDmXi72PPXXw3ipjER2NtK+wXRfrTWp07/M0cp9QUwGvjpvDFLgCUAsbGxF1fbW4huLi3uEB7+Abj7+jc6JvXQ6VYNQzrmPN8ZHv7OFOdU1D929fahprKCmqpK7B2dOjQW0X1I0ieEEBbIL6vmtR+T+WhHKrYGA/dOieTuiX1wc7Szdmiim1NKuQAGrXXp6d9fAfzNymEJ0WWZzSbS4w/Rd9TYJselHs7HJ9QVV6/Wl8lvDQ8/J44eL6nvb3d2rz7v4NZXbxQ9W4ckfUqpGcCrgA3wrtb62QbGTAZeAeyAPK31pI6ITQghmlJRY+S9zcd5+6djVNaauD42jAcv74e/e8d+CBBdj1JqMFALDANOaK13t3KqAOCL0+XibYGPtdZr2yZKIXqe3BPHqSovI2LoiEbHVJXXkplSTPQVHb9V2tPfmZpKI1XltTi52v/Sq6+gQJI+0WrtnvQppWyAN4BpQDqwSyn1ldY6/qwxnsCbwAyt9cnTZxaEEMJqjCYzy/ek8/L3R8kpreaKQQE8MmMgkf6u1g5NdB1PA2HAP0///srWTKK1PgYMb8O4hOjRUg/tByB8SOP/W6UlFKDNmoihHV+Y65cKnpV1Sd+Zlb4u2qBddA4dsdI3Gkg+/UMLpdSnwDVA/FljbgJWaq1PQt2ZhQ6ISwghLqC15vv4bP617gjJOWXERHjx5s3RxPbytnZooutJBTZorb9QSk21djBCiDonDx/AJzQcF0+vRsekHsrH0cWOgN7uHRhZHQ//c9s2/LLSJ0mfaL2OSPpCgLSzHqcDY84b0x+wU0ptBNyAV7XWH3RAbEIIUW/vyUKeWZ3ArhOF9PFz4e1bY7hiUACnt9UJ0VIfaq33nv79BqtGIoQAwFhbS0ZiPEOnXtHoGLNZkxqXT/hgb6tUZD6/bYO9kzP2Tk6y0icuSkckfQ3933J+VTFbIAaYCjgBPyultmutj54zkZSjFkK0g+N55Ty/LpHVh7LwdXXgH7OHcENsGLY2UpFTtN5ZCR9a6/9ZMxYhRJ3MowkYa6oJHzKi0TE5J0qoKqullxW2dkLjbRtkpU9cjI5I+tKpO9NwRihwqoExeVrrcqBcKfUTdecXzkn6pBy1EKIt5ZVV8+8fkvh4x0nsbQ08eHl/7rq0Ny4OUthYtB2lVCzwZyCCup+7CtBa62FWDUyIHujk4QMoZSCsif58Jw7loQyKsEFts61f19ai7FpW6bmhtg1lhQVtEo/omTrik80uoJ9SqjeQAcyn7gzf2VYBryulbAF76rZ/vtwBsQkheqAzFTkXb0qhymjmxtFh3D+1P35uDtYOTXRPHwF/BA4BZivHIkSPlnr4AIF9++Hg7NLEmHwC+7jj6HJxLXmq4uPJe+stSr9fj+ukSfg9cD+OUVEWXXtB2wYvb9JON5QXojXaPenTWhuVUr8H1lHXsuF9rXWcUmrh6dcXa60TlFJrgYPU/UB8V2st/2ULIdqU0WRmxZ50XjpdkXP64LqKnH39pCKnaFe5WuuvrB2EED1ddUUFWclHGX3NdY2OKSusJi+tjHGz+7b6PqayMk49+hhlP/yAwc0Nj7lzKF3/A8dnz8F95pUE/u1pbFwbTzqhgbYN3j6UFxagzWaUQY4eiJbrkD1MWuvVwOrznlt83uPngec7Ih4hRM+itebHxByeXZNIUk4Z0eGeUpFTdKQnlVLvAj8A1Wee1FqvtF5IQvQ86QmH0WZzk60aUg/nARAxxKdV99Bak/l/f6Zs40Z877sX71tvxcbNDdOjj5L//vvkv70E24BAAh59pMl5zm/b4OLlg9lkorK0BGcPz1bFJno2ObgihOjWDqQV8c/VCew4XkBvXxcW3xLN9MGBUpFTdKTbgYGAHb9s79SAJH1CdKCThw9ga2dPcP/Gt1imHs7H1dsB7+CmV+IaU7DsP5R+9x3+jzyCzx231z9v4+6O/wMPYMzLo+DDD/GcNw+HPr0bnef8tg1up3v1lRbkS9InWkWSPiFEt3Qyv4J/rUvkm4OZ+LjY8/Q1g5k/Ohw7qcgpOt5wrfVQawchRE+XenAfwQMHYWtv3+DrplozaYmFDBzTui8GK3bvJueFF3CbNg3v2xc0OMb/wQcpXbuO7GefIXzJkkbnqm/bcLqC55kG7eWFBdC79VtPRc8lSZ8QolspLK/htR+T+XD7CWwMinunRHL3xD64OV7cgXwhLsJ2pdQgrXW8tQMRoqcqycslP/0kgydf3uiYjKRCjNUmIoa2fGunqbiY9AcfxD40lKBn/tlo0mjr44PvPfeQ89xzlG7ciNvkyQ2Oq2/bcLpXn4s0aBcXSZI+IUS3UFVr4j/bTvD6hmTKq43MiwnjoSv6E+DuaO3QhJgA/FopdZy6M33SskGIDpZ6cB8AvYdHNz7mUD42dgZCBni1eP78d9/DlJtH2P8WY+PadHEw75tvoujzz8l55llcL7kE1cjK49ltG1w8vUApSiXpE60kSZ8QokszmzWrDmTwwrqjZBRVctkAPx67MooBgW7WDk30cEqpccB2YIa1YxGipzuxfw+u3j74hEU0+LrWmhOH8ggd6IWdvU2L5q7NzqHgww9xv/pqnAYPbna8srfH/49/IP1391C6cSPuV1zR4DgPPyeSTtS1bbCxtcXFw5PyQkn6ROtI0ieE6LK2JefxzzUJHM4oYUiIO89fN4xLIn2tHZYQZ/waeAM4CqwF1mqts6wbkhA9j9lkIvXwfvqNvqTRbZdF2RWU5FUxclp4i+fPe+tNtNGI3/33WXyN68SJ2Hh5UbruuyaTvuoKI9XlRhxd7XDx8pbtnaLVJOkTQnQ5R7NLeWZ1AhuO5BLi6cQrN4xg1vBgDAapyCk6D631QgCl1EDgSmCZUsoD2EBdErhVa22yYohC9AiZyUepLi+n1/CYRsecOFSXTEUMbdkXhzUnTlC0fAVeN9yAfViYxdcpW1vcLr+ckm+/xVxVhcHxwqMI7r51FTxL8itxdLXD1duH0tycFsUnxBlSxk4I0WXklFTx2P8OMuOVn9idWsifrhzIDw9P4tqRIZLwiU5La52otX5Zaz0DmAJsAeYBO6wbmRA9w4kDe1DKQMTQEY2OST2ch3ewC27eLTsHnvvv11D29vguWtjiuNxmTMdcUUH5li0Nvu7uWxdLSV4VAK5e3pQVFrT4PkKArPQJIbqA8mojS346xpKfjmE0m1lwSW/unRKJl0vDh9+F6Ky01pXA6tO/hBAd4MSBvQT2649jIwVWqiuNZCYVM2Ka5St1ANXHjlGyejU+v/0ttn5+LY7LZcwYbDw9KVm7DrfLL6wq6ubzy0of1LVtqCwtwVhbi62dVKQWLSNJnxCi0zKazCzfk85L3x8lt7Saq4YG8ciMAUT4tK5prhAdSSlVSl0Tdqir2Hk2rbV27+CQhOhxKkqKyUpJYtzcGxsdkxZfgNmsiRjSsq2dBR98gLK3x/u2W1sVm7K1xW3a5ZSsXoO5uhqDg8M5rzs42eLgbEvpmZW+s3r1efgHtOqeoueSpE8I0elordl4JJd/rk4gKaeM2Agv3r41hujwlpfRFsJatNZSQlYIKzt5aD9oTe8RjZ/nSz2ch4OzLYF9LP8exlRURPGXq3D/1dXY+rS8r98ZbtNnULR8BeVbtuA2deoFr7v7OlGSf2Z7Z919ygryJekTLSZJnxCiUzmcUcw/VyewLSWf3r4uLL4lmumDAxutuCZEV6CUGg5cevrhT1rrg9aMR4ie4sSBvTi6uhHQN7LB181mTerhfMIHeWOwsbzUReHny9FVVXjf9uuLis9lzGhsPDzqtng2lPT5OFKQWQ78stJXJm0bRCtI0ieE6BROFVXywndH+GJfBp5Odjw1azA3jQnHrgU/hIXojJRS9wO/AVaefuojpdQSrfVrVgxLiG7PbDZxbN9ueg2PxmBouPde9vESKktr6T3c8jN5uraWwo8+wnncWBwH9L+oGJWdHa7TLqd0zdoGt3i6+Tpx4nA+Wutfkr4CKeYiWk6SPiGEVZVW1fLWxhTe23IcDfx2Yl9+d1lf3B3lkLroNu4ExmitywGUUs8BPwOS9AnRjrKSj1JZUkyfmNGNjjlxMBeDQRE+2NvieUvWfYcxO5vAp/7aBlGC29SpFK/4H5X7D+Ay5txY3X0cMdWaqSipwdndFRs7O1npE60iSZ8QwipqTWY+3XmSV9YnkV9ew+yRITx8RX9CvZytHZoQbU0BZ/fjM3FhYRchRBtL2bMTZTDQu4n+fMcP5BHc3xMHZ8u/aCz44APse/XCdeLEtggT5+hoACr37rkg6XPzqWvbUJpfhYuHA67ePtKgXbSKJH1CiA6ltWZ9Qg7PrEngWG45Y3p7s/SqKIaFelo7NCHay1Jgh1Lqi9OPrwXes144QvQMx/buInTg4EZbNRRlV1CYVcGQSSEWz1l5OI6qgwcJePxxlKFtjh/YeHjg0K8fFXv2XvBafYP2vEoC+3jg6uUjK32iVSTpE0J0mEPpxfxjdTzbjxXQx8+Fd2+LZWqUvxRpEd2a1volpdQmYDx1K3y3a633WTksIbq1ktwc8k6eYNItdzQ65vjBPAB6DbW8VUPR55+jHB3xuGbWRcd4NqeYaEq++RZtMqFsfjl/eGal7+wG7dnHk9v03qJnkKRPCNHuMooqeWFdXZEWbxd7nr5mMPNHS5EW0XNorfcAe6wdhxA9RcqeHQD0iRnT6JgTB/PwCXGtX01rjrm8nJJvvsH9yiuxcWvbjizO0dEUffoZ1UlJOA4cWP+8nb0NTu725zRoT9m7E621fGEqWkSSPiFEuymtquXN00VaABZN7suiyVKkRfQsSqlY4M9ABHU/dxV1zdmHtXAeG2A3kKG1vrrNAxWiG0nZsxOvoBC8gxveullVVktmchExV/ayeM7i1asxV1Tgef28NoryF07RdecOK/bsOSfpg7piLqX5v6z0Gaurqa4ox9Gl4W2rQjSk2aRPKRWutT7ZEcEIIboHo8nMJ2cVabl2RDB/nDGQEE/Lvk0Vopv5CPgjcAgwX8Q89wMJgOUdpIXogWoqK0iPP8SIGb9qdMyJw3loDb2Ht2Rr53Ic+vXDacSINojyXHYhwdgGBFC5dx/cfPM5r7n7OpF9vBg4q1dfQb4kfaJFLFnpW6uU8qfuB80h4ODpfx7SWpe0Z3BCiK5Fa82PiTn8c3UCKbnljJYiLUIA5Gqtv7qYCZRSocBVwD+Ah9okKiG6qdSD+zEZjfRtqlXDgTxcPOzxC7Nsm2ZVQgJVhw4R8H//1y7bKpVSOEWPpGLvhcVc3HwcSdmTg9l8Vq++wgJ8wyLaPA7RfTWb9GmtByml7IHBwFBgGHWVx4Yppaq11r3bN0QhRFdwOKOYf65OYFtKPn18XVhyawzTBgXImQMh4Eml1LvAD0D1mSe11isbv+QCrwCPAG17kEiIbihlzw4cXFwIGTCowddra0ykxuUzYGwQymDZz6ii5ctRDg54zGp89fBiOUfHULpmLbWnTmEXHFz/vLuPI2azpqywClevX1b6hGgJi870aa1rgH1KqWSgEvAF+lG34ieE6MGyiqt4ft0RVu5Lx9PJjqdmDeamMVKkRYiz3A4MBOz4ZXunBixK+pRSVwM5Wus9SqnJTYy7G7gbIDw8/CLCFaLrMhmNpOzeQd/o0RjOqoJ5trT4Aow1ZvqO9LNoTnNlJcVffY3b9Cuw8fRsw2jP5RxT16+vYu8+PM5O+k4XminNr8K/V10TeUn6REtZcqZvAHVbSq4G/IDvqTufcPfpZFAI0QOVVxt5e1MKSzYfw2yGuy/tw+8ui8TDSYq0CHGe4VrroRdx/XhgllJqJuAIuCul/qu1vuXsQVrrJcASgNjYWH0R9xOiy0qLO0hVeRn9xoxvdEzKvhwcXGwJ7u9p0Zyl33+PuawMz7nXtVGUDXPo3x+DszOVe/fgcfVV9c+f3bYhpL8Xjq5ulBUWtGssovuxZKUvAdgHPAt8pbWubma8EKIbM5k1K/ak8cJ3R8ktreZXw4N5ZPoAwrydrR2aEJ3VdqXUIK11fGsu1lr/CfgTwOmVvj+cn/AJIeok7diGnYMjEcNHNvi6yWjmxMF8+ozwxcbCHSlFK7/ALjQU51GxbRnqBZStLU4jRlzQpN3N2xEUv7Rt8PKWlT7RYpYkfYuoO8v3e+ANpVQ+pwu5UFfM5cv2C08I0ZlsTsrlH98mkJhVSnS4J2/fGkN0uJe1wxKis5sA/FopdZy6M32tatkghGia2Wwiefd2ekePws7eocEx6UcKqak00nekv0Vz1qRnULF9O773/h5laP9jC04x0eS9/gam0tL6XoA2tgZcPR1+advg7UN5oSR9omUsKeTy9tmPT1cQG0ZdIjgX+LJdIhNCdBpJ2aX8c3UCG47kEubtxOs3jeSqoUFSpEUIy8xoq4m01huBjW01nxDdSUZiPBXFRfQfc0mjY47tzcHOwYbQKMu+sCz+8ktQCs9rr22bIJvhNHQoaE11YiLOo0bVP+/u60RJ3i8N2vNOnuiQeET30eLm7FrrdCAdWN324QghOpO8smpeWX+UT3am4Wxvw//NHMivL+mFg23Dh+OFEA3y1VrvOfsJpdSvgFQrxSNEt5S0Yxu2dvb0HtnwNkyzWXPsQB69hvpga9f8zzFtNlP8xRc4jx2DXUjDTd7bmsPpxuxVCecmfW4+jmQcKQTqtneWFxVhNpkaLVYjxPlanPQJIbq/qloTS7ee4M0NyVTUmrhlTDj3X94fbxd7a4cmRFf0jlLq11rrQwBKqRuBB4CvrRqVEN2INptJ2rmNiOHR2Ds6NTgmM6mIqrJa+li4tbNi5y5qMzLwe+D+tgy1SbZ+ftj4+FCVmHjO8+4+jhwpqsZUa8bV2wetzVQUF9X37ROiOZL0CSHqaa355mAmz65JJKOokqkD/fnTzCgi/V2tHZoQXdl1wAql1M3Une+7DbjCuiEJ0b1kJh+lrCCfS2/8daNjUvblYmNnIHywt0VzFn+xEoOrK26XX95WYTZLKYXjwIFUJSac87y7rxNoKC2swuWsXn2S9AlLSdInhABg78lC/v5NPHtPFhEV5M6/rhvG+Ehfa4clRJentT6mlJpP3Rn4NOAKrXWldaMSontJ2rkNg40tfWJGN/i62axJ2ZdD+CBv7B2b//hrKiujZN13eMyahcGp4ZXD9uIYNZCC/3yArqlB2dftsDnTtqE0rwq304leaWE+gR0amejKJOkToodLL6zgubVH+PrAKfzcHPjX3GHMjQnFxiBFWoS4GEqpQ9Q1YT/DG7ABdiilkOqdQrQNbTaTuO0neg0fiaNLwztTMpOKqCiuod+oAIvmLFmzBl1VhefcOW0ZqkUcBkaha2upPn4cxwEDgF8atJfkV9JraF3SV14gvfqE5STpE6KHKq2q5a2NKby75TgGBfdOiWThpL64OMhfC0K0kautHYAQPUF6Yhxl+XlMvPn2Rsck7c7G1t5Ar6GW7WApXvkF9n374jis47+bcYyqK+ZSnZhYn/S5eDpgsFGU5FXh5B6IMhgok7YNogXk050QPYzRZObz3em89P0R8spqmDMyhD9MH0CwZ8duXxGiu9NaS3VOITpA4pZN2Do4EBkzpsHXTSYzKXtz6T3MFzuH5qtdVh87TuW+ffj/8Q9WaU1kHxGBcnCgKiERj2uuAcBgULh6O1KaX4nBYIOLlzdlstInWkCSPiF6kM1Jufz9mwSOZJcyqpcX7y8YxbBQT2uHJYQQQrSKyVjL0e1biIwdi52jY4Nj0hMLqSqvJTLWsq2dxV98ATY2eMya1ZahWkzZ2uLQv3+DFTxLTjdod/PykZU+0SKS9AnRAyTnlPKPb+uaq4d7O/PWzdHMGBIozdWFEEJ0aScO7KWqvIyoCZMbHZO8Oxt7J1siBjdf6VKbTBSvWoXrpZdi6+fXhpG2jOPAgZR+9x1a6/qf1e4+jhw/mAeAi5c3hZkZVotPdD2S9AnRjRWU1/DK+qN8tOMkznbSXF0IIUT3krBlE45u7kQMG9ng66ZaM8f25dJnpB82doZm5yvfuhVjTg4ej/+5rUNtEYeogRQtX44xKwu7oCAA3HydqCytpbbahKu3D2nxB60ao+haJOkTohuqMZr54OcTvPpDEhU1Jm4aHc4Dl/fDx9XB2qEJIYQQbaKmqpKU3TsYPGkKNrYNf6RNjcunpspEPwu3dhat/AIbLy/cJk9uw0hbznFgFABVCYn1SZ/7mbYN+VW4evtQXV5ObXUVdg4Nb2sV4myS9AnRjWitWReXzTNrEkjNr2DyAD/+PDOKfgFu1g5NCCGEaFMpu7ZjrKlm4PhJjY5J2p2No4sdIQO9mp3PWFhI2Q8/4Hnj/Pr+eNbi0L8/KEVVYgJuUy4Dzm3b4OpV12C+rCAfr6AQq8Upug5J+oToJg5nFPP0N/HsOF5AP39X/nPHaCb1t955BCF6OqXUOGC71lo3O1gI0WJxP/2Im68fIQMGNfh6daWR4wfyGHRJEDY2zW/tLP5yFbq2Fq9589o61BazcXXBPjyc6sQj9c+dadBekleFh2/d+cSywgJJ+oRFJOkToovLKani+XVHWLE3HS9ne/5+7RDmjwrD1oIfcEKIdvVr4A2l1FFgLbBWa51l5ZiE6BZKcnNIPbSfcXPnowwN/7xL2ZODqdbMgLFBzc6ntaZo+XKcRozAoV+/tg63VRyioqiKj69/7Oxuj42dgdL8SkL6nU76CqSCp7CMJH1CdFFVtSbe3XyMNzemUGsy85tL+3DPZZF4ONlZOzQhBKC1XgiglBoIXAksU0p5ABuoSwK3aq1NVgxRiC4rbtMPoDWDJ13e6JjE7Zl4BTrj36v5Iw6V+/ZRc+wYQf/4e1uGeVEcBw6gdO1aTGXl2Li6oJSqb9vg6l23uldWKL36hGUk6ROii9Fa89WBU/xr7REyiiqZMTiQP80cSISPi7VDE0I0QGudCCQCLyulnIDLgHnAS0CsNWMToivSZjOHN64nfMhwPPwbLtBSnFtBZnIxY6/tY1F7oqLPl2NwccH9yivbOtxWO7PiWJOSjNPw4QC4+ThRkleJg7Mzdo5OstInLCZJnxBdyN6ThTz9TTz7ThYxONidF+YNZ1zf5vsOCSE6B611JbD69C8hRCukxR+iJDebCfNvbXTMke1ZoKD/6MBm5zOVlFCydi0e11yDwdm5LUO9KA6RkQBUJ/+S9Ln7OpJ9vBgAVy9vSfqExSTpE6ILOFVUyb/WJvLl/lP4uTnwr+uGMTc6FBuDNFcXQgjRsxze8D0Ozi5Ejh7X4Otaa47syCJ0gBdu3s23Myj59lt0VRWe113X1qFeFLvQUJSDA9VJyfXPufk4Ul1hpLqiFldvH9neKSzWIUmfUmoG8CpgA7yrtX62kXGjgO3ADVrrFR0RmxCdWUWNkcWbjrHkpxTMGu65rC+LJkfi6iDf1wghhOh5qsrLSNqxjcGTL8fOvuHes5kpxZTkVTH66t4WzVm4fDkOUVE4DhnclqFeNGVjg33fPlQn/5L0ufucadtQ16svIzG+scuFOEe7f3JUStkAbwDTgHRgl1LqK611fAPjngPWtXdMQnR2ZrPmy/0ZPLc2keySaq4eFsSjMwYS5t15tp0IIYQQHS1x608Ya2sYctm0xsf8nImtgw29RzTftqjycBzV8QkEPPEXi87+dTSHyEgqdu6qf+zue7pBe14Vrl7elBfmo7XulLGLzqUjlgtGA8la62MASqlPgWuA87+auBf4HzCqA2ISotPak1rA376O50B6McNCPXjjpmhie3lbOywhhBDCqrTWHPjuW/x69SGgT2SDY2oqjSTtzqFfrD/2js1/zC1avhzl6IjH1Ve3dbhtwiGyHyVffY2ptBQbN7dzG7R7+2AyGqksLcHZ3cPKkYrOriOSvhAg7azH6cCYswcopUKA2cAUJOkTPVRGUSXPrknk6wOnCHB34MV5w5k9MgSDnNsTostTShkAg9baaO1YhOiqMhLjyEtLZdrd9za6snV0VzbGahODJzTfsNxcUUHJN9/gPn06Nu7ubR0uNaYatmduZ0/2HqL9o7kk5BLsDC1rq3R2MRfnkSNxcLbFztGGkvwqAnv90qtPkj7RnI5I+hr6v1Kf9/gV4FGttamp5Wml1N3A3QDh4eFtFZ8QVlVebeTtTSm8/dMxAO6bEsnCyX1xtpdze0J0B0qp3wNPAjVKqTzgNa31u1YOS4guZ/+6b3FwcSFqwqRGx8RvOYVPqKtFvflK1qzFXF6O5/Xz2jJMiqqKeGH3C/xw8gfKassAeJ/38XLwYmafmdw78l5c7Cxrs+TQ79ykr65XnxOleZVEjqzbBVReWAC9+rTpexDdT0d8qkwHws56HAqcOm9MLPDp6YTPF5iplDJqrb88e5DWegmwBCA2Nvb8xFGILsVs1nyxL4N/ras7tzdreDCPXjmQEE8na4cmhGhbDwNDtdZZSqkg4J9KqVCt9V+tHJcQXUZZYQFJO7cxcsbV2Dk0XJEzJ7WE3JOlTJzf37LefMuXY9+3L07R0W0W58mSk/zuh9+RWZbJVX2uYlrENGICYtiVtYtvjn3Dp4mfklqSymtTXsPW0PzHcLuQEJSjIzXJKfXPufs6Upxbiat3EACl0rZBWKAjkr5dQD+lVG8gA5gP3HT2AK11fXklpdQy4JvzEz4hupOzz+0ND/XgzZujiYmQc3tCdFNlQA6A1jpTKXUnsB/4qxVjEqJLOfTDOswmE8OnzWx0TNyWU9jaGeg/pvnefNVJSVTu34//o4+2WRGU/Tn7ue/H+9Bo3p3+LiP9R9a/NilsEpPCJvH5kc95evvTPL/ref405k/NzqkMBhz6nFvB083HkbSEApw9PAGkV5+wSLsnfVpr4+mtLeuoa9nwvtY6Tim18PTri9s7BiE6i/PP7b10/XCuHSHn9oTo5t4CliulHtVaJwPhQIWlFyulHIGfAAfqfm6v0Fo/2S6RCtEJmYxGDq5fQ6/h0XgFNXxWr6bKSNLObCJHBeDgZEEBlxUrwM4Oj2tmtUmMcflx3PXdXfg7+/PW5W8R4R7R4LjrB1zPyZKT/Cf+P4S7h3Nz1M3Nzu3QL5Ly7TvqH7v7OGGsMVNbBc4ennXbO4VoRoccGtJarwZWn/dcg8me1npBR8QkREeqqDGyeKOc2xOiJ9Jav6mUygLeVUoNA9yBj5RS84D9WuukZqaoBqZorcuUUnbAFqXUGq319nYOXYhOIWX3dsoKC7j8N/c0OiZpVza11SYGTwhudj5zZSVFX3yJ+7TLsfW++F02hVWFPLjhQbwdvfnwyg/xcfJpcvyDMQ9ysvQk/9r1L6L9o4nyiWpyvH1kJMWrvsJUUoKNu3t924aSvCpcvXwoK5SVPtE8g7UDEKI7qzu3l86UFzbx7x+TmTYogB8ensRDVwyQhE+IHkRrvVJrPRnwB6KBH4FLgLctuFZrrctOP7Q7/UvOtYseQWvN7m++wCMgkN4jYxsdc2hjBj4hrgT0br4KZ/E332AuKcHr5uZX2ZpjMpt45KdHyK/M5+XJLzeb8AHYGGz4x4R/4GLnwtsHm/0r4KwKnnXn+s5t2+AtZ/qERSTpE6Kd7DtZyJy3tvHgZwfwc3NgxcJxvH5TNKFe0mBdiJ5Ka23UWh/UWv9Ha/2g1nqKJdcppWyUUvupOxv4vdZ6RzOXCNEtZByJJzPpCDFXXYvBYNPgmFNHi8jPKGPYlNBmz+dprSn870c4DBzYJgVcXtv3Gtszt/P42McZ7DvY4uvc7N24JeoWfjj5A0cLjzY51qFfPwCqk+s2Bbj5nG7Qnn96pS8/r5XRi55Ekj4h2lhWcRUPfraf2W9uI6OokuevG8aqe8ZLg3UhRKtprU1a6xHUVcAerZQacv4YpdTdSqndSqndubm5HR6jEO1h99df4OjmzpDJlzc65sCPaTi62NF/VECz81Xu2UP1kSN43XzTRRdw2ZW1i/cOv8d1/a9jdr/ZLb7+5qibcbFz4Z2D7zQ5zi44GOXkVF/Mxd7RFkcXO0ryKnH386eytITaqqpWvQfRc0jSJ0Qbqao18e8fkrjshY18eyiT303uy4Y/TGZebJgUahFCtAmtdRGwEZjRwGtLtNaxWutYPz+/jg5NiDZXcCqdlD07GHHFzEbbNJTkVXL8YB6DLw3G1r7hlcBz5vzoIwweHnhcffVFxVZprOTJbU8S7hbOI6MeadUcHg4ezB8wn3Un1nGs+Fij45TBgEPfvtScVcHT3deRkvwqPPzrEt3i3OxWxSB6Dkn6hLhIWmu+PZjJ1Bc38dL3R5k8wI8fHprEIzMG4uog5/aEEBdHKeWnlPI8/Xsn4HIg0apBCdEBdn/zBba2doyc3niCdnBjOkophkxquKrn2Wqzsyn97ns8587F4HRxPXHf3P8maaVp/PWSv+Jk2/q5bht8G462jrx78N0mxzlERlKddG7bhpK8Sjz869pTFOdI0ieaJkmfEBfhcEYxN7y9nXs+3ou7kx2f/GYsb90SQ5i3nNsTQrSZIGCDUuogdb1vv9daf2PlmIRoV+VFhcT/9CODJk2p70d3vpoqIwlbM+kb7YerV8MrgWcr+uwzMJvxuunGi4rtcN5hPoj/gHn95zEqcNRFzeXt6M28/vNYfXw1mWWZjY5ziOyLMTcXU3ExUNe2obSgCjefulV9SfpEcyTpE6IVckureex/B/nV61tIyS3jn7OH8s29ExjXt/mqXUII0RKnC7+M1FoP01oP0Vr/zdoxCdHe9qxehcloJOaqxs/KHdmeRU2lkeFTwpqdz1xVReEnn+I6eTL2oaGtjqvWXMtftv4FXydfHox5sNXznG3+wPmYtIm1J9Y2Osb+TAXPlDMVPB0xGzUaZ2wdHCjJzWqTWET3JXvPhGiBGqOZZduO8+8fkqmqNXHXhN7cO7Uf7o521g5NCCGE6BYqSorZv/YbBl4yEe/ghrdtmk1m9q8/SUBvd8vaNKz6ClNhId63L7io2D5N/JTkomReuewV3OzdLmquM8LcwhjqO5Q1x9dw+5DbGxzjEHm6gmdSMs7R0bidbttQWlCFh18AxTk5bRKL6L4k6RPCAlprfkzM4e/fJnA8r5wpA/3581VR9PVztXZoQgghRLey55svqK2pZuyc+Y2OSd6TQ0leFeOv69d8mwazmYKlS3EcPBjnUa3fjplfmc+b+99kfPB4poRZ1G3FYjN6zeD53c9zovgEvTx6XfC6XXAQytm5voKn+5m2DXmVePgHSCEX0SzZ3ilEM5KyS7nt/Z3c+Z/dGBQsu30U7y8YJQmfEEII0cYqSorZt/YbBoy7FJ/Qhrdtaq3Zuy4Vr0Bneg/zbXbOso2bqDlxAu/bb7+oNg3/3vdvqoxVPDr60Ytu93C+6b2mo1CNbvE8U8Hz/F59JflVuPv5UyJn+kQzJOkTohHFFbX89as4Zry6mf1pRTxx9SDWPjCRyQP8rR2aEEII0S3t+fZLamuqGTe38VW+1MP55GeUEz09AmVBS6SCpUuxDQ7CffoVrY7rcN5hvkj6glsG3UJvj96tnqcxAS4BRAdEs+b4GrTWDY5xiIysX+mztbPB2cO+rm2DXwDVFeVUlZW1eVyi+5CkT4jzGE1mPtyeyuQXNvDBzyeYPyqMjX+YzB0TemNnI//LCCGEEO2hsrSkbpVv7AR8QsMbHbd3XSquXg70G21BM/ZDh6nYtQvvW29D2bXu/L1Zm3lmxzP4OPnw22G/bdUclpjRawbHio+RVJTU4OsOkZGYcvMwFRUBdVs8S89p2yDFXETj5BOsEGfZlpLH1a9t4S9fHqZ/gBvf3Hsp/5g9FB9XB2uHJoQQQnRrO75cTm11FWObWOXLTC4iM7mYEdPCsbHgi9iCpe9jcHXFc951rY5r7fG1HMw7yP3R9+Nq335HO6ZFTMOgDKw93vAWT4d+pyt4nl7tc/NxoiSvCndp0C4sIEmfEEBaQQWL/ruHm97ZQWmVkTdvjubTu8cyKLj5imBCCCGEuDjFOdnsX/s1gydOxTcsotFxu749jqOrHYPGBzc7Z/Wx45SsWYvXjfOxcW1dslZlrOKVva8Q5R3FrL6zWjWHpXycfBgTOKbRLZ4Okecmfe6+jpQVVdf36pNzfaIpkvSJHq2ixsiL3x1h6kub2Hgkl4en9eeHhycxc2hQmx/SFkIIIUTDtn7+X5QycMn1Nzc65lRSIWkJhURPj8DOwabZOfPfeQfl4ID3ggWtjuu/Cf8lszyTP8T+AYNq/4/N03tNJ70snaOFRy94zTYoCIOzM9XJZ3r1OaHNGmONHQ7OLrLSJ5okLRtEj6S15qsDp3hmdSJZJVVcMyKYx64cSJCHk7VDE0IIIXqU7OMpJGzZyKhZc3H39WtwjNaaHV8dx9ndniGTGu7dd7aa9AyKv/oKr5tvwtbHp1Vx5Vfm8+6hd5kcNpnRQaNbNUdLjQ8ZD8DPp35mgPeAc15TSmF/VjGX+gqeeZW4+wdQLCt9ogmy0id6nEPpxVy3+Gfu/3Q/vm72rFg4jlfnj5SETwghhLCCzR8vw9HFldHXNH7uLj2xkFNJRcRc2Qs7ewtW+d59B2Uw4HPHHa2O6839b1JtrOahmIdaPUdLBboEEukZydZTWxt8/ewKnu4+dZ9bzlTwlKRPNEWSPtFj5JVV8+iKg8x6Ywup+eU8N3coX90zgdhe3tYOTQghhOiRju/fQ+rBfYydMx9Hl4bP3dWt8h3D1cuBwROaP8tXm51N8f9W4jFnDnaBga2KK7kwmRVJK7h+wPXt0qKhKeOCx7E3ey+VxsoLXnOIjMSUl4exsBBXbweUgtL8Kjz8AyjJzWm03YMQkvSJbq/WZObdzce47PmN/G9vOndN6M2Pf5jMDaPCMVjQ30cIIYQQbc9YW8uGZW/jFRTM8CtmNjou9XA+2cdLiJ3ZCxu75j+65r/7Htpsxuc3d7U6thf3vIiLrQsLhy9s9RytNT54PDXmGvZk77ngtTMVPGuSk7GxMeDi5UBJXiUe/gEYa6qpKC7q4GhFVyFJn+jWNh3NZcYrP/H3bxOI6eXFugcn8uerBuHu2LpePUIIIYRoG3u+/ZLCzFNMWfBbbBvpoWc2mdm2MgUPPycGXhLU7Jy1mZkUffopHtdcg31oaKvi2paxjS0ZW7h72N14OXq1ao6LER0Qjb3Bnm2ntl3w2gUVPM+0bfA73bZBtniKRkjSJ7qlE3nl3PWfXfz6/Z2YNby/IJZlt4+mr1/79dcRQgghhGVK8/PYvvJTIkeNpdeImEbHxW/NpDCznEvmRFrUly/vzbcA8Lvnd62Ky2Q28cKeFwhxDeGmqJtaNcfFcrJ1IiYghm0ZFyZ9toGBGFxdqU76pW3DmZU+kF59onGS9IlupazayLNrErni5Z/4OSWfP105kHUPTGTKwABrhyaEEEKI0zZ++B6YNZNv+02jY2oqjez8+hhBkR70HuHb7Jw1J05QtHIlnjfcgF1I8xU+G/Jl8pckFSbxYMyD2NvYt2qOtnBJ8CWkFKeQVZ51zvNKKRz69q1f6fMMcKaipAYnt7oKpdKrTzRGkj7RLWit+WJfOlNe2MjiTSnMGhHMhj9O5reT+mJvK/+ZCyGEEJ3FiYP7OPrzZkZfO69+haohe9amUllay4R5/SzqnZv72usoe3t8f3t3q+IqqynjtX2vMcJvBFdEXNGqOdrKJSGXAHWtG85n3++XCp5egS4AlBWZcfbwpDgn64LxQoAkfaIbOJhexNy3tvHgZwcI8nDki99dwgvzhuPv5mjt0IQQQghxlpqqSr5f8hpewaGMmjW30XEl+ZUc+CGN/mMC8I9wb3beqiNHKFm9Gu9bb8XWr+Fef81ZcnAJ+VX5PDr6UYuSzPbUz7Mffk5+jZ7rMxUUYCwowCvQGYDCrHJp2yCaJM3ZRZeVW1rN8+sSWb4nHR8XB56/bhhzo0OlIqcQQgjRSW355ANK8nKZ/9fnsLVvfPvk1uXJKAVjr+lr0bw5L72EwdUVnztb15fvRPEJPkz4kGsjr2WI75BWzdGWlFKMCx7HpvRNmMwmbAy/9CZ0iOwH1BVzcY+JxWBQFGZV4O4fQFbyEWuFLDo5WekTXU6tycx7W44z5YWNrNybwW8u7cOGP0xiXmyYJHxCCCFEJ5WRGM++dd8wcvrVhAwc1Oi4E4fyOLY/l9ireuHm3fyunbLNWyjf9BO+ixZh4+HRqthe2P0CDjYO3B99f6uubw9jg8ZSXF3M0cKj5zzvEFmXCFefbtvg4e9EUVYF3sGhFOfmUFtdZY1wRScnK32iS9mclMtTX8eTnFPGxP5+PHH1ICL9pSKnEEII0ZkZa2pY9/a/cff1Y8KNtzU6rrbGxObPjuIV6MyIy8ObnVcbjWQ/9yx24eF43XJzq2LbkrGFTembeCjmIXydmi8Y01FGBY4CYHf2bqJ8ouqftw0IwODqSs1ZxVwKs8rpOyIctKYgI52APpFWiVl0XrLSJ7qEk/kV3P3Bbm59byc1RjPv3hbLf24fJQmfEEII0QX89PFSCk+lM+03v8fe0anRcXvXplKSV8WkGwdgY0EhtqLly6lJTsH/j3/A0MR20cbUmmt5budzRLhHcEvULS2+vj0FugQS6hrK7qzd5zyvlMIhMrK+bYNXoAvFOZV4BtX1JczPSOvwWEXnJyt9olOrqDHy1sYU3v7pGDZK8cfpA7hzQm8c7Wyav1gIIYQQVnd83272rfmakVf+il7DoxsdV5hVzt51qfQfE0DIgOaboptKSsj992s4jxqF2+WXtyq2TxI+4UTJCV6f8jp2Ng03iLem2MBYNqRtwKzNGNQvSbBDv0hK1/8AgFeQM2azxsbWC4ONLflpqdYKV3RistInOiWtNd8cPMXlL27itR+TuXJIID/+YRL3XBYpCZ8QQgjRRVQUF7H2rVfwDYtg4k23NzrObNb8+EECdg42jJ/bz6K58954A1NREQF/eqxV1TbzK/N568BbjA8Zz8TQiS2+viPEBsRSXF1MUmHSOc87REZiKizEmJ+PV0Bd24bi3Bq8goJlpU80SFb6RKeTmFXCX7+KY/uxAgYFufPK/JGM7u1t7bCEEEII0QJaa9YtfpXqinKue/zvTVbrPLA+jaxjJUy7YxDO7s1v06xKSKDgw//iOW8ejoMaLwrTlNf2vUaVsYpHRj1i9RYNjYkNjAXqzvUN8B5Q/7x9ZN2ZveqkZDyHxwB1K6U+oeHkHE/p+EBFpycrfaLTKK6o5clVh5n56mYSs0r5+7VD+PreCZLwCSGEEF3Q7m++4NjeXUy8+Xb8wns1Oq4gs5wdXx2j93Bf+o1qvFn7GdpsJvOvf8XG0xP/hx9qVWzx+fGsTFrJjVE30sejT6vm6AghriEEuwSzJ3vPOc/Xt21IScbByRYXD3uKsirwCQ2nKCeL2ppqa4QrOjFZ6RNWZzJrPt+dxvPrjlBUUcPNYyJ4aFp/vFxafiBbCCG6G6VUGPABEAiYgSVa61etG5UQTUuLP8Tmj5fRb8wljJzxq0bHmU1mfvhP3bbOyTcPtGjFrejz5VQdOEjwv55rVYsGszbzzI5n8HL0YuHwhS2+vqPFBsayOX0zWuv6Px9bfz8M7u5Un6ngGehCYXYFoRPPquDZ27Ieh6JnkJU+YVV7Ugu59o2t/GnlISL9XPn63gk8fe0QSfiEEOIXRuBhrXUUMBa4RynVuv1sQnSAsoJ8vnnlOTwDg5m+8IEmE7k9a1PJOVHCxPn9LdrWaczLI+ell3AeMwb3XzWeTDZlZdJK9ufu58GYB3G3d2/VHB0pNiCWwupCUop+2bZ5poJnTX0FT2cKM8vxCQ0DID/9pFViFZ2XrPQJq8gpreLZNYms3JtBoLsjr84fwazhwZ12T70QQliL1joTyDz9+1KlVAIQAsRbNTAhGmAy1vL1K89RU1XJvL/8Awdn50bHnkoqZNc3x+k/JsCibZ0AWf/4B+bKSgKffKJVnxnyKvN4ac9LxAbEck3fa1p8vTXEBvxyri/S65f+ew59+1L6/fdorfEKdKGmyoS9sy8GGxtJ+sQFZKVPdKgao5l3fjrGlBc28c2BTBZN7ssPD0/imhEhkvAJIUQzlFK9gJHADiuHIsQFtNasf/dNTh2JZ/pv78M3LKLRsVVltXz/fjzuvk5MunFAo+POVrJmDaVr1uJ3zz049GndObzndz1PlbGKv4z7S5f53BHqFkqAcwC7s8/t1+fQLxJTURGm/Hy8AuuS65K8GryCQiTpExeQlT7RYTYn5fLXr+JIyS1nykB//nL1IHr7ulg7LCGE6BKUUq7A/4AHtNYlDbx+N3A3QHh4eAdHJwTs+up/HN7wPWPnzmfg+EmNjtNa8+OHCVSU1DD3kRjsHZv/OGrMyyPrqb/hOGQIPnfd2ar4tmVsY/Xx1SwavqhTF285n1KK2MBYfj718znn+hzOVPBMTsZrwAiA01s8w8k5IRU8xblkpU+0u7SCCn774W5ufW8nRrPmvV/H8v6CUZLwCSGEhZRSdtQlfB9prVc2NEZrvURrHau1jvXz8+vYAEWPl7RzG5s/+Q8Dxl3KJfNubnLs/vVpHD+Qx7jZffGPaP5MndaarKeewlxRQfCzz6BsW75mUVFbwd+2/41e7r24c2jrkkZrig2IpaCqgOMlx+ufO7ttg4unA3YONhRmV+ATGkZRtlTwFOeSlT7RbqpqTSzelMJbG1MwKMUfpw/gzgm9pbm6EEK0gKr7Wv89IEFr/ZK14xHifKeOJrD69RcJ6tuf6b9runBLWkIBP69Mpu9IP4ZPDbNo/pJvvqH0+/X4//EP9atbLfXK3lc4VXaKpTOW4mDj0Ko5rKn+XF/W7vpVSls/PwweHlQnJ6OUwivQmaKscgaMjgCtKTyVgX+vrrOiKdqXrPSJNqe1Zl1cFpe/tIlX1icxbVAAPzw8iXsui5SETwghWm48cCswRSm1//SvmdYOSgiAvJMn+OLZp3D18uaaPz6OnX3jCVVJXiXr3j2MV5ALU34dZdGZuprUVLL++hROMTF4L1jQqhh3Ze3ik8RPuDnqZmICYlo1h7VFuEfg5+R3zrm+MxU8z7Rt8Ap0If9UOb5hddu789NSrRKr6JxkpU+0qZTcMp76Op6fjuYyIMCNT34zlnF9fawdlhBCdFla6y1A16g4IXqU4pwsVvzzCWzt7bnuz3/HxdOr0bG11SZWLz4EGq5cONSic3y6poaMhx4GW1tCXngeZdPyL44rait4YusThLmFce/Ie1t8fWehlCI2IJbdWbvPPdfXL5KS1WvQWuMX7saRHVnYOXrXVfDMSLNy1KIzkaRPtImyaiOv/ZjE+1uO42hrw5O/GsStYyOwtZHFZCGEEKK7KS3IY8Xf/4KppoYb/vosHv6Nt1wwm8x89+5hCjLKuOqe4Xj6N97G4Ww5L71MVVwcoW+8jl1QUKvifGXvK6SXpbN0+lKc7Sy7b2cVGxjLmhNrOFl6kgj3usqojoMHU/TpZ9SePIl/L08A8jMq8QoKIS9NKniKX0jSJy6K1pqvDpzin6sTyC6pZl5MKI/MGIifW9fbLy+EEEKI5pXm5/H53/5ERUkRc//vaXzDezU6VmvN5s+TOHEon4nz+xMxxLLdP6U/bqBg2TK8br4Zt6lTWxXnT+k/8UniJ9wSdQuxgbGtmqMzOftc35mkz2noUAAqDx3G94oZKAU5qaX4hkWQmXzEarGKzkeWYUSrJWaVcMOS7dz/6X783RxZ+btLeH7ecEn4hBBCiG6qND+Pz5/6ExXFRcz9v78R3H9gk+P3f5/G4U0ZjJwWztDJoRbdo/rYcU498giOgwfj/8gfWxVnXmUef9n6F/p79eeBmAdaNUdn09ujN96O3uec63Po2xfl4EDVoUPY2dvgHexCTmoJwf0HUpKbQ2l+nhUjFp2JrPSJFiuurOXl74/y4fZU3B1t+efsodwwKgwbgxw5EUIIIbqrouwsVvzjcSpLSpj7f083m/Al/pzJtpXJRMb4M252X4vuYSotJf2ee1D29oS+/hoGh5Z/kay15i9b/0J5bTnvXfFel6zW2ZAz5/p2Ze2qP9en7OxwjIqi8vBhAPwi3DlxMI9RMwcBkHEknoGXTLRm2KKTkJU+YTGzWfP57jSmvLCRD34+wU2jw9nwh8ncNCZcEj4hhBCiG8tNPc6nT/yR6vJyrvtz8wlf8p4cfvwggdCBXkxdEIWy4HOCNps59cij1KSlEfLKy60+x/dx4sdsydjCw7EPE+nVuhYPnVVsYCzZFdmkl6XXP+c4dChV8fFoo5GACDeqympxdAvCzsGRU0cSrBit6Ew6JOlTSs1QSh1RSiUrpR5r4PWblVIHT//appQa3hFxCcsdSi9m7uJtPLLiIL18Xfjq9xN4+toheDrbWzs0IYQQQrSj9MQ4PvvrYygbG+Y/9RxB/QY0Of7EoTy+fy+OwD4ezFw0DFsL2zXlvvwyZRs2EPDoo7iMHt2qWA/kHuCF3S8wOXQy8wfMb9UcndnZ5/rOcBo6BF1ZSXXKMfxON7vPS6sgqF9/MhLjrRKn6HzaPelTStkAbwBXAoOAG5VSg84bdhyYpLUeBjwNLGnvuIRlCstr+L8vDjHrjS2kFVTy0vXDWbFwHENCPKwdmhBCCCHaWcLWTaz4++M4e3px49/+hU9oeJPjTxzKY+3bh/EJdeWq3w/HzsGyhK/wk0/If+ddPOffgNctN7cq1oKqAh7e+DABzgH8fcLfLeoD2O7K8yFjL1QUtMl0fT374ungec65PschdcVcqg4fwjfEFYONIvdkCcEDBpObepyayoo2ubfo2jriTN9oIFlrfQxAKfUpcA1Q/9WD1nrbWeO3A5ad9BXtxmTWfLrrJM+vO0JplZE7xvfm/sv74e5oZ+3QhBBCCNHOtNb8vOJjfl7xCaFRQ5j18P/h5Obe5DUp+3L47t04fEJcmXXfCBycLPuYWfrjBrKe/juukycT+PjjrUrWTGYTj/30GIVVhXw480M8HKz05bSpFuJXwb4PIeswVJxVSMU9BIJHwoQHIbR11UQNylDfr+8M+14RGFxdqTx0CM+5c/EJcSX7RCkjpkShtZlTSUfoNWzkxb4z0cV1RNIXApzdHTIdGNPE+DuBNQ29oJS6G7gbIDy86W+aROvtO1nIE6viOJRRzNg+3vztmiH0D3CzdlhCCCGE6AA1VZWsW/xvjv68mcGTL2fab+7BxrbpL32P7sxi/bIEAnq5cfXvh+PgbNmXxJX795Px8MM4DhpEyEsvomxb99H0jf1v8HPmzzw57kkG+Zy/oawDGKthx2LY8TaUZIB3HxhwJfgNBM9wKDxelwQe2wjvToWRt8LlfwUX3xbfKjYwlvUn13Oq7BTBrsEogwHHIUOoOlRXzMU/wo2k3TkE9o1BKQOnjsRL0ic6JOlr6Osa3eBApS6jLumb0NDrWuslnN76GRsb2+AcovXyy6p5bm0in+9OJ8DdgX/fOJJfDQvqHNsjhBBCCNHu8jPS+PqlZyjISOfSmxYwatbcZj8HHPgxjS3Lkwjp58nM3w3D3tGyj5eVcXGc/M3d2Pr5Ebb4LQzOrWue/u2xb3nn0DvM6TeHuf3mtmqOi3JyO3x1L+Qdhd6T4OqXIXIaGBo4RVVdCpv+BdvfhISv4eYVEDaqRberP9eXvZtZrrOAunN9+cv+g7mmBv8Id+I2n6KqXOEb0UvO9QmgY5K+dCDsrMehwKnzBymlhgHvAldqrfM7IC5xmtFk5qMdJ3nxuyNU1Jj47aQ+3DelHy4O0tFDCCGE6CmO/LyFdYtfxdbOjrl//hsRQ0c0OV6bNdtWJrN/fRp9Rvgx7Y5B2NpbdoavOimJtDvvwuDmSsSypdj6tnzFC+Bg7kGe2PoEMQExPD6mdVtDW622Cr7/C+x8BzzC4Jb/QeTlTV/j4AZXPA0jboZP5sN/58KCryHI8hqG/bz64W7vzu6s3czqW5f0OQ4eArW1VB85gn+vXgDkpJYQMmAQcRvXYzaZMNhY9u9GdE8d8al+F9BPKdUbyADmAzedPUApFQ6sBG7VWh/tgJjEabtPFPDEqjjiM0uYEOnLX2cNJtLf1dphCSGEEKKD1FZVseE/Szj043cERQ7g6gcfw93Xr8lrjDUmfvhPAsl7chg6OZQJ1/fDYGH7pupjx0m94w6UnR0RS5diFxzcqrizyrO478f78Hf25+XJL2Nn04F1BwpT4fNbIfMAjFkIU/4CDi34/OQ/EH79Fbx/JXw4GxasrnvOAgZlICYghl1Zu+qfcxo6BIDKQ4fwumEwNnYGclJLCRkQxf5135CbepyAPt2rfYVomXZP+rTWRqXU74F1gA3wvtY6Tim18PTri4EnAB/gzdPf0Bi11q074SoskltazbNrEvnf3nSCPRx58+ZorhwSKFs5hRBCiB4k+1gy3/77eQqzTjH6muu45Pqbmz2/V1pQxZrFh8hNK2XcnL6MnBZu8eeHqiNHOHn7HWAwEL7sfewjIloVd3F1MYvWL6LaVM1709/Dy9GrVfO0SvJ6+N9dYDbDjZ/Wnd1rDc/wusRv6ZXwwTXw203gFmjRpbEBsWxI20BWeRaBLoHYBgdj4+1N1aHDeN9kwDfUlZwTJQy77Jcm7ZL09Wwdsn9Pa70aWH3ec4vP+v1dwF0dEUtPZzSZ+XB7Ki99f5SqWhO/m9yX30+JxNletnIKIYQQPYWxtpYdKz9lx5fLcfH0Yt7j/yB8yLBmrzuVXMTatw9hrDVz1aJh9Bpm+bbMyoMHOfmbuzE4ORG+9H0cevduVeyVxkp+/8PvSS1JZfHli+nr2bdV87TKrvdg9R/AfxBc/wH4XOS9ffrCrV/AO1Pg6wfgxk/AggQ6NvCXc31X97kapRSOQ4dQefAgAMH9PDmwPg1HFy/cfP3ISIgj+spZFxer6NLkk34PsvN4AU+sOkxiVimX9qvbytnXT7ZyCiGEED1JVvJR1r71CvnpJxk0cQqTf/0bnFybrtKttWb/+jS2f5GCm68j1z40DO8gF4vvWb5tG+n33oeNtzfhS5diHxrSqthrzbX8cdMf65qwT3qB0UGta+LeYlrDj3+HzS9A/xlw3ftgb/n7b1LAYJj6BKz7PzjwCYy4qdlLBngNwM3Ojd1ZdUkfgHNsLLmbfqI2J4eIwT7s++4k6YmFhA8ZTvLOnzHW1mJrJ623eqp2b84urC+ntIqHPtvP9W//TGmVkcW3RPPBHaMl4RNCCCF6kKryMta/+yYfPf4w1RXlzH70Sa6856FmE76qslpWv3mQbf9LptdwX+Y9FtuihK9o5RecvPu32AUHE/HfD1ud8JnMJh7f8jib0jfx+NjHuaLXFa2ap+U3roVVv69L+KJvgxs+aruE74wxiyD8EljzGBRnNDvcxmBDdED0OU3aXS65BICK7dsJjPTAztGG1Lh8BoydQHVFOakH97ZtzKJLkZW+bsxoMvPBz6m8/P1Rqo1m7rmsL/dcJls5hRBCiJ7EbDYR/9MGNn+8jMqSEkbOuJrx19+Cg3PziUtafAE/fphARWkNl97Qn6GTQyw+v6e1Ju+NN8l7/XVcLhlHyKuvYuPWur6/JrOJx7c+zurjq7k/+n6uH3B9q+ZpseoyWL4Akr+HSY/B5Mcs2n7ZYgYDXPsGvDW+rv3DLf9r9j6xAbFsSt9EbkUufs5+OEZFYePpSfnWbXjMmkXYQG9OxuVz6Q2jcXRx5cjPW+gb01SrbNGdyaf/bmrXiQL+8uUvWzmfmjWYPrKyJ4QQQvQYWmtOHNjLTx8tJe/kCYIiBzDnsb9aVNCjttrEzyuTObQpA69AZ+YujME/wt3ie5vLyzn158cpXbsWj9mzCfrbU6hWbi00mU38Zetf+ObYN9w38j7uGtpBZSDKcuHjeXUVOn/1KsQsaN/7efepa9i+5hFI/Bairm5y+Nnn+q7sfSXKYMB53FjKt21Da03EEB+O7c+lOKeayNGXcHT7Zow1Ndja27fv+xCdkiR93UxuaTXPrElg5d4MQjydWHxLNNMHS1VOIYQQoifJPpbMTx8t5eThA3j4B3DV/Y8wYOwEVEMNw8+TnljAxo+OUJxbyfApYYy9to/F/fcAak6eJP2e31OdkoL/I4/gffuCVn8OqTXX8uctf2bN8TX8fsTv+c2w37RqnhbLT6nroVeaBfM/bn2FzpaKvRN2vQs/PFV3dtCm8Y/qA70H4mLnwu6suqQP6rZ4lq5ZS01KCuGDQwHqtnhecimHN3zH8QN76DdqXIe8FdG5SNLXTRhNZv67PZUXvztKldEkWzmFEEKIHqgoO4ttn/+XhC0bcXRz57Jf/4Zh02ZaVMCjqqyWrSuTSdyWibufE9c8OJLQAS1rhVC6fj2n/u/PKKUIe2cJruPHt/atUGms5KGND7ElYwsPRD/AnUPvbPVcLZKxFz6aB9oMv/4awkZ1zH2hLsmb+iR8djPs/2+Tq4u2BltG+o9kZ9bO+udcT5/rK9+2De/bbsMnxIWTcfmMmDocJzd3jmzbLElfDyUZQTew+0QBf1kVR0JmiWzlFEIIIXqg3NTj7Fy1giPbNmNja8voa+cx+prrLDq3ZzZr4recYseqY9RUGomeEcGomb1atLpnrq4m57l/UfjxxzgOHkzIKy9jHxbW6vdTXF3MvT/ey4HcAzw57kmu639dq+dqkaTv4fNfg4sP3LISfPt1zH3PNvAqCBsDG56BofOaLBozLmgcz2c8T1ppGmFuYdiFhGAXEU75tp/xvu02wgf7cOCHNIy1mn5jLiFh80Zqq6uwc3DsuPcjOgVJ+rqwvLK6Busr9qQTJA3WhRBCiB4nPTGOXatWcGzvLuwcnYi5+lpiZl6Dq7ePRddnHClk8/Ik8tPLCO7nycT5/fEJadkXx1VHjnDq0ceoTkzEe8EC/B96EHUR58ZOlpzknh/uIaMsg+cnPt9xVTr3fVRXRCVgENy8wuJG6W1OKZj2N3h/Omx/Cyb+odGhl4VdxvO7n2dT2iZuGXQLULfFs2TVV+ja2nNaNwwYdykH16/l+P499B/T+hVY0TVJ0tcFmcyaj3ak8sK6I1TWmlg0uS/3SoN1IYQQokcwGY0c27OTPau/JCMxHkc3dy65/mZGTL+62fYLZ2QfL2HHVymkJRTi6uXAFXcNJjLGv0VfHGujkfx33yX3jTexcXcndPFbuE2e3Mp3VWd31m4e2PgACsU7V7xDTEDMRc1nEa3r2jH8+HfoMxmu/xAcLS9a0y7Cx8KAq2DrqxB7Bzh7NzgszD2Mvh592Zi28Zykr+iTT6k8cIDAkdH1rRsm3TgEZw9PEjZvkKSvB5IsoYvZe7KQv3x5mLhTJYyP9OGpWUOI9JetnEIIIUR3V5qfx8Ef1nH4x3WUFRbg5uvHZQvuZuhlV2DnaNl2vbz0MnZ+fYzjB/JwdLVj/HWRDJkY0qKtnABVCQlkPvEkVYcO4T7zSgL+8hdsvVp2/u9sWmuWH13OMzufIcwtjDemvEGYe+u3h1rMVFtXLXP3+zDsBpj1Oth2kuqWU/4Mb10CO96Gy/7U6LDJYZNZFreM4upiPBw8cBkzBgwGyrf9jHNsLBGDfUjZm8Ol8/oxdMoV7PhyOfnpafiEdsCfr+g0JOnrIgrKa/jX2kQ+3ZVGgLsDr980kquGBslWTiGEEKIbM5tNnDy4n/3fr+HYnp1oNL2HRzP1rnvoMzIWg41lyVphVjm7vj1B0u5s7B1tGTOrN8OmhGHv2LKPgqaycvJee42CDz/ExsuLkJdfwv3Ki6tsWWms5O/b/85XKV8xPmQ8z136HB4OHhc1p2U3LqzrwXdsI4x/oK6AigXVTTtMwGAYeDXseAvG3dPo6uNl4Zfx3uH32JKxhav6XIWNuzuOQ4dQtmULfvfdy+BLg0nek0PynhyiZ17D3tVfseOLz5h5b+PbRkX3I0lfJ2c2az7ZdZJ/rT1CebWRuyf24b6p/XB1kH91QgjREyil3geuBnK01kOsHY9of1prco6nkLBlI0e2/URZYQFO7h6MmjWHoVNn4Blg2VkzrTUZRwo58EMaJw7lY2tvIHp6BCOnhePo0rKeedpspuSbb8h58SWM2dl43nAD/g89iI3HxSVnx4qP8cdNfySpMInfDf8dvx3+WwyqAxKv/BT4+HooTIVr3oCRt7T/PVvj0och8Zu6Ng6XPtTgkKG+Q/F29GZj2kau6nMVAO5XXEHO8y9Qffw4IQN64RXozKFNGQwcF8vwK2ay55svGXfdjXgFhXTgmxHWJJlDJ3YwvYi/fHmYA+nFjOntzdPXDqF/gGV79YUQQnQby4DXgQ+sHIdoZ4WZGSRu/YmErZsoPJWOwcaW3iNjmTxhMpGjxmBja1miZjKaSdqdzf71aeSnl+HkZseoq3szZGIIzu4t37pYsXcf2c8+S9XBg/WVOZ1HjmzxPGc7s53z+V3P42jryJuXv8mEkAkXNafFjm2Cz28Dgw38+iuIuKRj7tsaIdEQeTn8/DqM+W2DlTwNysDksMl8d+I7ak212NnY4f6rX5Hz4ksUf7kK/wcfYMikEDZ/lkROagmxV89m/9pv2PHFcmb87oGOf0/CKiTp64SKK2p5/rtEPtpxEh8XB165YQTXjAiWrZxCCNEDaa1/Ukr1snYcou1ps5msY0mk7N5Jyp4d5J08AUoRFjWE2Ktn03/MeBxdLT+3X5JfSeK2TOK2nKKiuAavIBcuu3Ug/UcHYGvXsjN7AFVHjpL76quU/fgjtn5+BD3zDB7XzLKowXtT8irzeGrbU2xM38glwZfw9Pin8Xf2v6g5Lbb7fVj9R/DpBzd9Cl69Oua+F2PiI/D+FbBnWd02zwZMDp3MyqSV7M7ezbjgcdj5++MyYTzFq1bhd9+9DBgbxM9fHuPQpgym3hbF0Muns3/dt4y7bj4e/laqUio6lCR9nYjZrFmxN51n1yRSVFHDr8f14qEr+uPu2LItGEIIIYTonGqqKkmLO8SxPTtJ2buT8sIClMFA6MDBTL7tLvqPnYCbj6/F8xlrTRzfn0f81lOkHykEICzKm6m3hRE2yLtVXxhXJSaSv+QdStasweDqit8D9+N9660YXJrv+dcUrTWrUlbx/K7nqTJW8eioR7kp6qaO2c5prIHv/gw7l0C/K2Due9av0Gmp8DHQ61LY+m+IvRPsLizaMzZ4LA42DmxM28i44Lrm656zZ5Px4EOUb9+O6/jxDBgdQOL2LMbPjWTUrLkc/H4NP6/4hBm/e7CD35CwBkn6Oon4UyU8seowu1MLiQ735Ok7RzM4uAMOMQshhOgWlFJ3A3cDhIeHWzkacYbJaCQrJYmTh/aTemg/mUlHMJuM2Ds50Wt4DH1jx9B7RAxObpYnIFprsk+UcHRHNkd3ZlFdYcTN25FRV/Vm4NhA3H2dWhyn1pqKnbvIf/ddyjdvxuDsjM+dd+Bz113YeHq2eL7zpZak8o/t/+DnzJ+J9o/mr5f8ld4evS96XosUpdUVbMnYDeN+X9cDz9DylU+rmvhH+GAW7P8vjLrrgpedbJ0YFzSO9SfX88ioR7Ax2OA6ZQoGd3eKv/gS1/HjGTIplLjNp0jYlsnIaeHEXHUtO1etYMC4S+k9MtYKb0p0JEn6rKy0qpaXv0/iPz+fwMPJjn/NHcZ1MaEYDLKVUwghhOW01kuAJQCxsbHayuH0WMaaGrKPJZNxJJ6MxDjSEw5TU1kJShHQuy8xV19LxJARhEQNxtbO8p08Wmuyj5eQvDeHlL05lBVUY2NroM8IX6LGBxM6wAvVis8O2mymdP168t99j6qDB7Hx9sbvgfvxuvHGiy7SAlBRW8GSg0v4IP4D7G3seXzM48wbMK9jVvcAkr6HlXfXtWa4/gMYdE3H3Let9Z4IoaNhyysw8rYG20pc3fdqNqZvZOuprUwMnYjBwQH3q2ZSvPILTKWl+Ia6EdzPk33fn2TQhGDGzbuZlD07Wff2v/n1C29Y3ONRdE2S9FmJ1pqvD2by92/iyS2r5sbR4TwyfQCezp2kN4wQQgghmqS1pqwwn+yUZE4dTSAjMZ7sY0mYjEYAvIJCGDh+EhFDRxA2eFiLVvMATCYzWcnFHD+YV5foFVZjsFGED/JmzKw+9B7mi4Nz646AGHNzKfriS4qWL6c2LQ27sDACn3wCj9mzMVjY86/J2M0mvj72Na/te42cihxm9Z3FgzEP4utk+dbVi1JbBeufhB2LIWBIXcLn07dj7t0elIJJj8BH18HBzyD61guGTAmbgrejN/87+j8mhk4E6rZ4Fn3yKSVr1+I1bx7jr4tk+bO72fn1MS69vj9X3vMQHz/+MBuWvi0tHLo5SfqsIDmnjCdWHWZbSj5DQtxZclssI8I8rR2WEEKITkgp9QkwGfBVSqUDT2qt37NuVD2P1pqygnxyTqSQlZJMzvFkslKSqCguAsBgY0tA30hGXjmL4AFRhPSPwtnDs8X3KS+qJjUun9TD+aQnFFBTZcJgqwiP8mbMNReX6GmzmfKt2yj6/HNKN2wAoxHnUaPwf/AB3K64AmV78R8LtdZsztjMK3tfIakwiSE+Q3hx0ouM8B9x0XNbLOswrPwN5MTDmIVw+VMNnoPrciIvh6DhsOUlGH4j2Jz778vOxo5r+l7DB/EfkFuRi5+zH45Dh2Lfty+Fn3yC59y5+Ee4M+TSEA5tSCfqkmAC+kQyZvb1/LziEyJHjaX/2A6qoCo6nCR9HaiyxsRrPybxzuZjONrZ8PQ1g7lpTAQ2spVTCCFEI7TWN1o7hp5Ea015USH5aSfJT08lL/3k6d+fpLqiHAClDHiHhNJreDQBffoR0CeSgN59sbVv+W6d6opaTiUXk3G0kPTEQvLTywBw8XQgMjaAiCE+hA70anET9XPukZxMyeo1FH/5JbWnTmHj5YX3bbfhed11OPRpm3N1Wmu2ntrKWwfe4mDuQcLdwnlh0gtcEXFFx1UfN9XC5pfgp+fByQtuXgH9pnXMvTuCUnVn+z67BeJWwrDrLxgyp98clsYtZVXKKu4aehdKKXx+cxeZj/2JkjVr8LjqKsZc04fkvTn89OkRZj8czZjZN3B8327WvPky7r7+BEb2t8KbE+1Nad01t/3Hxsbq3bt3WzsMi30fn81fv4ojo6iSOdEh/OnKKPzcHKwdlhBCdAlKqT1aa6k0YKGu9jPSGqrKyijKzqQ4J4uirLp/FpxKJz/tJFXlZfXjHN3c8Q0Nxyc0HJ+wcPwj+vD/7d13fFxXmfDx35leNOrNlmRLluUet7im98QJ6WQhtGzosCHAQpY3LBtCKAvL+wZYXiDkJWQpJoX0QBLSnAax48Rxl3uTbBVbdXq75/3jjpot2ZJtaSTN8+VzOfeee+bq3ImlM8/cU4orp2A/yS6QYX+Mxj0dHNzZzqEd7Ryu84MGq81C6ZRsJs0uYPKcAvInek8pWIrt30/n88/T+dfniO7cCUrhXb6M3JtuIuvii7GcRIDan4SR4OUDL/P7Lb9n05FNTPBO4NNnfJrrp16P3TqCs48fWg9P3wZNm2DOB2HFf4G3YOR+/kgxDLjvbDCS8MXV0M/yGbe+cCuNwUb+esNfsSgLOplk7w03YoRCVP/1LyiHg61/P8SqP2zj4ltmMmP5BILtbTz0H18nFolw83d/TF7pxDTcnDgZg20f5UnfMKtrDfGdZ7fwcm0z00qyeOSzy1g6ZRz+ERJCCCFGCcNIEmpvJ9DaQqCtlUBrC50th+loauwO9KLBYJ/XeHJyyS2dyLTl51BQPpnCCjPQ8+TknnTwlYwbHK7307S3M7V10HkkAvQEeYuvqqJsWi4lVdkntZZeF6010W3bCLz+Ov4XXyKydSsA7oULKfn3f8d3+WXYi0/fWnidsU6e2PEEf9r2JxqCDVT4Krhr+V1cV33dyAZ74XZY9X1Y+xvwFsGH/wQzrhq5nz/SLBY492vw+Kdg27P9Tkxz47QbufPNO1nTsIblE5ejrFaKv/516j7zGdoefoT8T3ycmcsnsO3tBl5/eAeFFVkUludxw5338PBdd/D4D+7i5nt+jDc3Lw03KIaLPOkbJtFEkt+8uZefv7oTi1J85ZIabj27Crt1hGarEkKIcUSe9A3NaG8jh0prTTwaIdzZSSTgJ+zvTG1+IoFOQh0d3cFdsK2FYHs7Wht9rmGx2sgpLiGnpJSc4lJyS0rJKSklt9hMHa6hL3PQWyyc4Eh9gMN1fo7UBzhS56f1UBAjaX7OyspzUlKZTUlVDiVV2RRX+k4pyAMwQiGCq1cTeO11Am+8QaKxEQDX3Llkr1hB9hWXY58w4ZR+xtH2d+5nZe1Kntr1FOFEmMWli/nYzI9xfvn5WEdyGQTDgA0Pwct3Q+gILP4MXPhNcOeOXB3SxUjCL5aA3Q2fe9Ps9tlLNBnlokcvYumEpdx7wb2A+Tt04JOfJFq7jeqXXsTq8xHsiPLn/3wXpeCmOxfjyXbQsHM7j373m/gKirjxzu+QU1ySjjsUQzDY9lGCvmHw1s4j3PX0ZvYcCXLlGaX8xwdmMSHn1BoTIYTIZBL0Dc2ZCxfqVX97ATA/7IFO7acKaI1O5aF1Tz6a7s8Fml6v090v1r0u1PszRO8yXcfJRAIjEScZj5NIJEjG4yRTx137sXCYWCRMPBLu3o+Fe++HiPg7u2fE7I8ry0dWXj5Z+QV48/Lx5Rek9gvwpfI8OTlYTkNQEgnGaWsM0dYYpK0hSFtTiLaGYPcTPAC3z05hhY+iiiyKK7MpqcwhK+/Uh3ToWIzw5s2E1qwhuOYdwuvWoWMxLF4v3rPPJuv888k671xsRUWn/LN688f8/G3f33hm9zO83/w+NouNK6uu5GMzP8bMgpmn9WcNyp7XzYXWGzdB2SK46v/AxPkjX490Wv8neOoLcPMjMP2KY07/5L2f8ODmB3n06keZkT8DgPCWLey78YPkf+qTlNxxBwCHD/h54sfvUViRxbVfXYDNbqW+djNP/fi72OwOrv9fd1NSNYZnPc0AEvSlQXNnhO/+tZZnNxxicoGH71wzmwumn76uFEIIkakk6BuaGVOq9GcWz053NQbNZndgd7txuN04XD2p3e3B4XLj9vlw+7Jx+Xy4s1KpL9vM82ZhsZ7eJ0yxSILOIxE6j4Txt0RobzKDvNbGEOHOWHc5q81CbombvFIvBWVZFFZkUVThw5PjOC2TlxiRCJEtWwitW0do9RpC69ahw2EAnDNm4F22jKzzz8Nz5pmo0zRGr0vSSLK6YTVP736aVw+8SjQZpSqnimuqr+G6qdeN3NILvR1cB6t+ALtegpxJcMm3YfYN/Y5rG/eScfj5QvAWw6dfPuZpX2eskyufuJKZ+TO5/9L7u/89HvrWt+h4/Akqfn0fWeeZyzrsXtfMC/dvpmpeIZd9ejY2u5WW+gM8/oNvEwkGuOr2r1N95tIRv0UxOBL0jaBE0uD3b+/n3pd2EEsafPGCaj5/fjWuU+y2IYQQwiRB39DMO2OO/uN/34si9UEw9YGvdyCilOrJ71UG1f0qc7/7Narnc6VSva7d6/r0uqYCq9WGxWbHardjs9ux2mxY7eaxNZVvd7qwnoalAgZLG5pwIE6wI0qwPUqgLUrnkTCdRyL4W8w0Eoz3eY3TYyOv1ENuqZe8Ug/5pV7yJnjwFbixnKYZuLXWxOvrCa/fQHj9esIbNhDZtg1STzgdU6vxLl2GZ+kSPIsXY8s7/eOt4sk4axrX8MqBV1h1YBUtkRayHdmsqFrBtdXXMqdwzsjNxNnbwXXwxv+G7X8Fdz6c8xVY8rnxsQzDqXj3t/CXr8LHn4LqC485/cetf+RHa3/Ery75FeeUmUsxGOEw+27+CPFDh6h6/DEcFRUAbFxVz5uP7GBiTS5XfnEuTreNQGsLT/7XPTTv3c2CFVdz3kduPakZasXwkqBvhLx/oI1vPbWZLYc6OW9aEfdcM5vKQm+6qyWEEOOKBH1DM1rayJGitSYWThAOxIkE4qk0RqgzRrAjRrDdDPCCHVFCHbHucXZdLDZFdoEbX4GL7AIX2YWp/UI32YUuXF77aQ12dCJBdPceottqiWytJVJbS2TbNozOTgCUx4N7zhzc8+fjnj8P97x52AqGZxK4YDzIWwff4pUDr/Bm/ZsE4gE8Ng/nlp/L5ZWXc375+Tisafigbxiw4wV4+//C/r+DMwfOus1cd881tEXux61EFH42H/Imw63PH/O0L56Mc+3T1+K0Onns6se6x1zG6urY+8GbsE+YQOVDf8LiNocg7XinkVf+p5a8iV6u/tI8vDlOEvE4b658kHXPP0PR5CpW3PY1iiZVjvCNiuORoG+YtYdi/OiF7Ty89gAlPhd3XT2LFXNK0/MNmBBCjHMS9A1NutvIoUgmDRIxg0QsmdoM4tEksXCCaDhBLJwgFjGPu7au/K4gLxKIYxj9f55xuG14cxx4c53mluPEm+vo3s/Kc+HNcaCGYc1cIxYjtncfsT27ie7aTXTPbmK79xDbtw8dM7uJKqcT5/TpuGbOxDVzJu7583BOnXpaFkrvT8JIsPnIZt5ueJvVh1az8fBGEjpBviufCyou4OJJF7N0wlKc1jQtKxUPmxO0vP0LaNkFORWw7Auw4OMS7PVn7QPw13+FDz4Ic2445vSL+17ka69/jbuX382N027szg+88QZ1n/s8WRdeSNlP7sXiNP97H9jSwvO/3oTdZeOSW2Yyabb5ZcOedWt54Vc/JRLws/DKaznrpo+c8uRH4vSQoG+YaK15fN1B/vO5WtrDcW49q5KvXDqNLKesfiGEEMNFgr6hmTt7vr7vnsf6TrrS9X8acxIX3Tuvb5mu153ofNe+kdSpzcAweu0nNMmu/V5leoI8Y8Bg7RgKHC4bDrcVh8uG02PD5bXj9jlwZdlxpzZXVq9jnwO7c3iHWuhYjHhDA7G6euL19cTqDhDbs5font3E6+rNJ1YASmEvL8c5ZQqO6mpcM8xAz1FVNWwBHphB3o62Hbzf/D7vNLzD2sa1+ON+FIpZBbNYPnE5Z088mwXFC0Z29s2jNW2F9SvNgC/UAhMXwPLbYNZ1YJXPWAMyknD/BRA8DLetBaevz2mtNbe8cAu72nbx0AceYnL25O5zrStX0vTd7+FZsoTyX/4Ca1YWAC0HA7z4wBZaDwWZd0kFy6+txmq3EPZ38uaf/odNr75IVkEh5918CzPOPh+ViWMqRxEJ+obBjiY/33pyM+/sa2XhpFy+f/0ZzJwg3zoJIcRwk6BvaGZOPUN/ecXPU2P1OGrsXs8Yvl7D9fqM0es7dq/XawZ4nbIqrFaFxaqwWC2p1NzvL9/msGJzWLA7rN37NocVeyq1Oa3dAZ7TbcPhsmF3WofladyJGKEQ8aYmEk3NJJqbiB882BPg1deRaGzqCewAZbfjqJyMY0o1zuopPWlVFZaTXNB9KPwxP5uObGJ983rWNa9j4+GNhBPm5C8TvRNZPnE5yyYuY1npMnJducNen+MKtcLmx81g79D7YLHBtCtg2Rdh8lnHdFcUA6h7Bx64FM66HS777jGnDwYO8uG/fJgCVwErr1qJ194zDKnj2Wc5dOc3cU6rYdL992MrNCfoScSS/OPxXWx6/SC5JR7O+acaJqee+h3cXsurv72P5n27KZpUyTkfuYWq+Yukt1uaSNB3GoViCX72yk4eeHMvWS4bd66YwU1nVpy2wdtCCCGOT4K+oRlL3TvTQRsGyY4Oki0tJFpbSba2kWxrJXH4sBngNTaZAV5Tc/c4u95sRUXYKyqwl5fhKK/AXlGBo7wMe0UFtuLiEXvy0RnrpLallq0tW7u3A/4DAFiUhWl505hfNJ8FxQtYULyACVmnd82+kxINmLNvbnkKtj8HyRiUnAELPgpn3ATeNMwKOh48/S+w4WH4wj+gaPoxp9c0rOFzL32O88vP5ycX/gSL6vk3GnjjDepv/zIWXxYTv//97lk9wezu+cYjO+hoDlM5t5Czbqgmr9SLNgy2vf0mf3/kD3Q0NVJaXcOS625i6qJl8uRvhEnQd5q8uKWR7zy7lYPtYT60qIJvrJhBvldmLhJCiJEkQd/QZErQpxMJjECApN+P4feT9Acw/J2p1G8Gdm2tJFrbzACvLRXgtbf3eTrXTSlshYXYSkqwlZZgLy7BVlKCvdRMbcUl2CdOGJEndr2F4iF2t+9mV/su9nTsYVf7Lna376Yh2NBdZqJ3IrMKZjGrYBazC2czt3AuWY6sEa3ngIJHYPvzsO0vsHsVJKPgKTCDvPkfhQlz013DsS94xFzCoWQO3PIs9NNVt2s2z0/N+RRfXvjlPk/mItu3c+iOfyO6Ywe5H/4QJXfcgcVrPhFMxg02vFrH2uf2kYwlqVlcwqIrK8kr9ZJMxNny2iusfeZx2psayC+rYMEVVzPr3AtwuD0jdvuZTIK+U1TXGuI7z27h5dpmZpT6+N51c1hUmT9sP08IIcTAJOgbmnQEfdow0IkExOPoeBydSPRN42ZKou95IxxBR8IY4QhGODzwfjiMEQqRDPgx/Gagp0OhE9bLmpODNT8fa34+tlRqzc/DlteVl2fm5Zn7ym4fgXfrWAkjQWOwkTp/HfWBeg50HjCDvPY9HAoe6i7nsDioyqmiOreaqblTuwO9PNfpX8LhpCUT0LAe9rwGu1+FA2+DNsy19WZ+AGZ8ACqWyli90+39lfD0F+Hcr8HFdx1zWmvNPavv4bEdj3Ft9bV8e/m3sVt7/r0b0SiHf/bftD74ILaiIopu/xI511+PSq2DGeqMsf6lA2x6vZ5k3KBqfhFzLyxnYk0u2jDYsfot1j77BM17d2N3uZl17gXMvWQFxZVTRuwtyEQS9J2kWMLg/725h5+/uhOLUnz1kmn889mV2K3yqFoIIdJFgr6hmV9VpR+fPRsMnZpwxdx0r30Mo3tCFwxjwDLdE7oMVEZrM9hLJk/fDdhsWNxuLC4Xyu3u2fe4sWb5sGT7zNTnw5rtw5LVk1p8WVizs7FkZWH1+YZ1kpShSBpJWiItNAWbaAw1Uu+vp95f3x3kNQQaSOhEd3m7xd4nuKvOqaY6t5pyXzk2y+i4p25aw5GdsPd1M9Db+yZEO8xzpWfAtBVmsFc6V8bpDbenb4P3/wAffghmXHnMaa019224j19u+CVLSpdw7wX3kuPM6VMm9P77NP/wR4Q3bMBZU0Phbbfhu+TiPsHfhlfq2PLWQaLBBIUVWcw6eyI1i0pwem007t7BhhefZ/s/3iARjzFh2gzmX3olNUvOwj7CT8kzgQR9J2H1nha+9dRmdjUHuGJ2KXddPYuJuTIdrRBCpJsEfUOzoKZGP3vlVeYHbNW1ELsldaxSE6KkZmyxWLone1G9ymDpWpi9q0z/10KBstlRdjvKZjNTuw269nufc/Qq05X2F9yl6WnbyTC0QUe0g5ZwC62RVloiLTSHmmkKNdEUbDLTUBOHQ4dJ6r6BcY4zh4qsCsp95VT4eqVZ5RR7itM7m+bxRP1w8D2oWwv170D9Wgi3medyJ8GUC2HKBVB1nozRG2nxCPz2MmjdB59dBQXV/RZ7dvez3PWPu8h35XPnkju5eNLFfbp7aq3x/+1Fmn9yL/H9B7BXVJD/iU+Qc/31WLPMbp/xWJIdaxrZ9NpBWg4GsFgVk+cUMH1ZKZVzColFg2x9/RU2vPQ8bQ0HsTmdTF20jJnnXsDkMxZgHSVfyIx1EvQNQUsgyg+e28bj6+opz3Nzz7WzuWhGyWm5thBCiFMnQd/QZMqYvtNNa004EaYz1mlu0U78MT9t0TYzoEsFdl3BXWu4lbZoG4Y+dnyg2+amxFNibt6SY/bLfGVkO8bADODhdmjaDA0boXGjmR6uNbtrAhROh4ol5lZ5LuRXpbW6AmjbB78+H7xF8ImnIaes32KbDm/i7rfvZkfbDs4vP59vLP4GFdkVfcroZBL/y6/Q+uCDhNevR7ndZF92GTk33IBn8aLuSVuO1PvZvrqRHe80EeqM4fTYmLqohOoFRUyoyaFxRy21b73GjtVvEQkGcPmymb7sHGacfR5l02fJ5C+nQIK+QTAMzSPv1vHD57cRiiX47HlTuO3CGtyOUfrNmhBCZCgJ+oYm04I+rTXRZJRQIkQwHiQUDxFKhAjFU8dH5ftj/u7Azh/z0xnt7D5OGIkBf47H5iHflU++O58CV4G578qnwN2zn+/Kp9hTTLYje2xNYR9qNRdDP7LD7Kp5ZKcZ7LXv7ymTVWJ215y40ByTV34muEfRWELRY/8/YOU/gTsXPv4kFNb0WyxhJFhZu5JfrP8F0WSUyysv51NzPsX0/GNnAA1v3Ej7Y4/T+dxzGIEA1qJCfBdehO+Si/EsW4bF4cBIGtRva2Pb6kb2rj9MIm7gcNuYPKeAqnmFlE3P5tD2jWx76zV2v/cOiVgUT04u1YuWUrN4ORVz5mEbQ0/6RwMJ+k6gtqGTf39yE+sOtLOkKp/vXzeHmhLfiV8ohBBixEnQNzSLFi3Sa9euxdDmuL2usXsGBlqncnqlXflAz2uOKtf1NKt3+a7rJnSChNFrO/o4lZc0ksSNeJ8ySSPZfRxJRIglY0SSEaLJqLklokSSvfIT0Z5zqS2SiBzTdXIgNmXD6/CS7cju2Zw9+z6Hr89xtiObPFceea483LYxPORDa4i0Q3sdtB+A1j3QsrMnwAsd6SlrsZvdAotnmuPwSueawZ5PekGNKQ0b4I83mk9lP/JnM0gfQHOomT9s/QOPbn+UUCLEktIlXFN9DZdOvhSPve8snEY4jP+VV/G//DLBN97ACIWweDx4zzsP38UX4T3rLGwFBcRjSeq2trJv4xH2bTpC2B/HYlWUTculYlYBpVVu2hu3sGvtava+/y7xaASH203V/EVMXbyMyfMW4s6Sz+YnIkHfAILRBD99eQe//fs+ctx2vnnlTG5cWDa2vo0TQogMI0Hf0FTMqtC538hNdzVOitPqxGl14rK6cFgduGyu7jyn1YnT5uy3jNfuxWPzmKndg9dmph67p885u8U+Ptv8RBQCTeBvgs6D0JEK7rqCvI46iB615qCn0HwCVFgDhdOgILWfO1lm1hwvWnbDH66DzkNw1pfg/G+AfeAvLzqiHfx5x595YucT1PnrcNvcnFd+HhdUXMC5ZeceM+mLEY0SWr0a/8uv4F+1iuQR88sDZ81UPEuW4lm6BM/ixVhycmna08HeDWYA2NZozrzrznZQMSOPiTU+0PXUb32X3e+tIdzZgVIWSmumUTXvTCrnL6R0So10A+2HBH1H0Vrz4tYm7n5mCw0dEW5eUsE3rphBrkfW3BNCiNFOgr6hmTlvpr79d7ebk6yk/mdRlu5gx6Is3flK9aRdCzZ3l+86r1S/5btSq7Jis9iwWWzYLXZsyobV0pNns9iwKdvxjy02HBbH+AzITobWEAtCuBVCLamtzXwi5280t0CjGeQFGnsmUunNmQ05FebkKrmpNKfC3M+rAo8sRZURQq3w4rdg/Urzv/vlP4BpV5gTNA1Aa836w+t5ZvczrDqwipZIC1ZlZXbBbM4sOZMzS85kfvH8PkGgNgwimzYRXPMOoTVrCK1bhw6HQSmc06fjXjAf97x5uOfNI5ZTSv22dupqW6nf1krYHwcgp8jNhKnZuLNaifp3c2jHRhr37AStcfuymXTGfCpmzaF85hnkl5XL3wsk6Oujvi3E3c/0rLn3/evncOZk+UMnhBBjhQR9Q5NpY/pGLSMJsQBEOs0ZL6O90t55oVRgF27t2Q+1mouY98diN8fX+UrNrWu/O2+CGeC5c0f0dsUot+d1+MtXoXU35FfDsi/A3A+B6/gTChnaYPORzbxW9xrvNb3HpiObiBtxFIqavBoWFi9kTuEcZuTPYErOlO61/3QsRnjzZjMAXLuW8IaNGMEgAJacHNxz5+KeOxfnrFmE8ippaLHSsKuDQ7vaiQbNsbVZeU6KJtmxWuoJtu2kee8Wgu2tAHhycimfMZuymXMonzmbgvJJGTkjqAR9QDxp8MBbe/nZyzsB+OqlNdx6dpWsuSeEEGOMBH1DI0HfEBhJiIdTW+iotL+8Ac7FAmYQ1zuYiwVO/POVxZwMxZ0PngLz6Zsnv9dxV15BT54777hPaYQYUDIOW5+G1b80l92wOszlNWZ8AKZeMuBMn71Fk1E2Hd7Ee03vsa55Heub1xNKmN017RY7U3OnMj1/OjPyZ1CVU0VVdhUl3hKUoYnt2UN4wwZzW7+B6K5d5lNtwJqTg3PWTJzTZxKdNItWexnNbTYa9nQSbDe/AFEWyC2O4XQ1E48coL1xF8E2s0upze6gqLKK0upplFbXUFJdQ/6EsnHfJTTjg761+1r51pOb2d7k57JZJXz7mtmUyZp7QggxJknQNzSLzjxTv/uP14HUAurHpAyQP9S0v+sYYCTMYMpIgk72Ok6kjpODOE4cda1UmowduyX6yUvGzXFuXfvJaCqNpfJT+0Z86G+wsoDdY46NsrvB5ganz9xc2an97NTWO88Hzpy+eXavBHAiPerfhS1PQu2zPbO05k6CSWfBxAVQMguKZ51wrcWkkeSA/wDbWrexrXUb21u3U9taS2uktbuM2+ZmcvZkKrMrqcyppDK7kqqcKibZirHsriOyrZZobS2R2m1Ed+xAx2IAKJcL59SpGFUz8RdOp91eQmvYzeGmOIlY11Ipnbi9rVitzcQjjQRa60jGzSDR5nRSUFZBQfmk7q2wYhLZhcXjJhjM2KCvLRjjh89v45F36yjLdXP3NbO5dJbMNiWEEGOZBH1Ds2jaRP3uR4LprsbwsTrNJxRWO9icZmp1pPJT+33yHcfP6w7gjk7d/Rx7zNfJWCIxXmgNTVtg35uw/++w/+2+s7l6i8yZXItnQf6UvuNEXTkDXFLTEmlhb8de9nbsZV/nPjPt2Meh4KE+a1vmu/Ipyyrr2dylVLRZKKrz49nbTHznTqJ79pJobOy5vs1OvGou4fI5BHMn02nJpz1sJxQCrQ200Qq6Gbu9DWglHm4mFu7ofr3Vbie7qITc4hJySkrJKS4lp7iE7MJivHn5eHJysFjGxhJuGRf0aa15Yt1Bvv9cLR3hOJ8+p4ovX1KDx5F5fXuFEGK8kaBvaBbNnaHf/dUXMGdyUYNIGWS5QVzHYgVlBYvN3D/m2JY67u9cf2UtfY8tNgm4hBhOWpszwTZvhebaXuk2iB/1ZZIzxwz+fKXgLYasrq3EDBazSsxjd1737200GaWus469nXvZ37mfen89BwMHORg4SEOggYTuWStToSjyFFHqKaXMkk9Vh4OyVihsiuJr6MBe14yuOwQJ8zUxu5eQdwKRomoiBZMJeUsIWnIJJJwYySjaaMVIHkEb7VgsfpTuJBFvw0hE+tyWslhwZ+fiyy/AV1CAJzsXl8+HO8uHy5eN25eN2+fDlZVKvVlpe3I4qoI+pdQVwM8AK/AbrfUPjzqvUuevBELAP2ut1x3vmr2Dvl3NAb711CZW72ll4aRcfnDDGcwoPf6gVCGEEGNHpgd9J2pHjyZj+oQQp53WEDySWgLkqOVAAk0QaDa3/rpMW2zmmFR3Xs/m6TrO7c5LunJpVpp6I8zBeCcHI60cCjXSHGqmOdRMU6iJ4FGBp8XQTA55qQx7KQs6Ke20UNhhkNMWw9sSwnnEj05qos48Iq58Is58M/UWE80qJubMIaIsxIww2giYmw6AEUQbAdABtA6jjTAwUNyksDmc2BwubA4Xdpcbu9OFw+3G4XbhdHtweNzYnQ7sTidWux2b3Y7VZsdqT22pfZvNhtXhwGqzoSwWLBZrKrWgLFasdju5JaU9P3mQ7eOwPwZTSlmBXwCXAvXAWqXUM1rrrb2KrQBqUttS4Fep9Lgi8SS/WLWL+17fjdtu5T9vOIMPLarAYpFvAIUQQowPg2xHhRBieCkFWUXmNtBC71qby4cED/cNBIPNZn7X1lEPjZvM/V5BnBWYkNoWd2Xa3H3GxQadhTTZnTTZHTRbLTRbNIdzk7TpODuNKGuMCG2JMO2JIAltoDTkBhS5gXZyg21kh3aTE4LckKIwbCOvzUJFUOMNgjXhxlBZxBw+YnafmTrKSdg8xG1O4jYHcasibrWQsGgSFo3WUTMo1HHisRixaBz8MbTuBN0COobWcdBxIAEYR79rQ2K1efnKykeG/LqR6Pu4BNiltd4DoJR6GLgW6N1YXQv8XpuPHVcrpXKVUhO01g0DXTQQSXD5T99gf0uI6xeU8e9XzaQwyzmc9yGEEEKkw2DaUSGESD+lemagLZo+uNckohBuN5cs6R0YRv39LnHijfqZEmxhSu+Zcvt5AqcBv0XRZrHSZrXQYbHgz7YQyLUQsFjwW+3UW23UWi0ELFb8Fgthi59E0g/RBiwxsEU1tqjGE6V7c8bBGdc4EuCIgjNux2F4sCXt2JN27IYDq2HHajiwagdK2QE3WtlBOdDKRsJqw1AWDKUwLApDgWFRaJSZpzSGAm2xoFFo1XOHKpk45l4HYySCvjKgrtdxPcc+xeuvTBkwYNC3tyXIBKVY+emlnD31+LMKCSGEEGPYYNrRPpo6I8z8jxeGtVJCjCYyzHM8sQKFqe3EFAZuYjiJ4SaKiyiurn2j1z4xXETxEMWh4mSRII8E9hNsVuLgTKBdCQyVIGHRxJW5JZWBoaIY+EkqMJQmqTRRdPdxQmmSgKEgrhQJIKmVOeGxYZ7QGpQB2lBgKFRqImQ0ZtDXa5Jkq7IB3x3yuzoSQV9/v4ZHh+ODKYNS6rPAZ1OH0dfuuHDzOXecYu0yRyFw5ISlRBd5v4ZG3q+hk/dsaAb5lfG4dFJtJN9bsXlYazV6ZfrvVibfv9x75sqo+7/je32ahUG1jyMR9NUDFb2Oy4FDJ1EGrfX9wP0ASql3M3lQ/1DJ+zU08n4NjbxfQyfv2dAopTJ5VhJpI4cgk+8dMvv+5d4z894hs+9/sO3jSMwtuhaoUUpVKaUcwIeBZ44q8wzwCWVaBnQcbzyfEEIIkUEG044KIYQQAxr2J31a64RS6jbgb5iddH+rtd6ilPp86vx9wHOYyzXswlyy4dbhrpcQQggxFgzUjqa5WkIIIcaQEVm5XGv9HGZg1zvvvl77GviXIV72/tNQtUwi79fQyPs1NPJ+DZ28Z0OT0e9Xf+3oCWTy+5XJ9w6Zff9y75krk+9/UPc+IouzCyGEEEIIIYRIj5EY0yeEEEIIIYQQIk3GZNCnlLpCKbVdKbVLKfW/0l2f0Uwp9VulVLNSKlOn7h4SpVSFUmqVUqpWKbVFKfXldNdpNFNKuZRS7yilNqTer++ku05jgVLKqpR6Xyn1l3TXZbRTSu1TSm1SSq3P8Bk8ByWT28dMbu8yve2Stihz25VMbyOUUrlKqceUUttSv//LByw71rp3KqWswA7gUsxprNcCN2utt6a1YqOUUuo8IAD8Xms9J931Ge2UUhOACVrrdUopH/AecJ38++qfUkoBXq11QCllB94Cvqy1Xp3mqo1qSql/BRYB2VrrD6S7PqOZUmofsEhrnTHrL52sTG8fM7m9y/S2S9qizG1XMr2NUEr9DnhTa/2b1OzOHq11e39lx+KTviXALq31Hq11DHgYuDbNdRq1tNZvAK3prsdYobVu0FqvS+37gVqgLL21Gr20KZA6tKe2sfVN0ghTSpUDVwG/SXddxLiT0e1jJrd3md52ZXpbJO1KZlJKZQPnAQ8AaK1jAwV8MDaDvjKgrtdxPRn0h02MHKVUJbAAWJPmqoxqqS4l64Fm4CWttbxfx/dT4N8AI831GCs08KJS6j2l1GfTXZlRTtpHkbFtV4a3RT8lc9uVTG4jpgCHgQdTXXt/o5TyDlR4LAZ9qp+8jPk2R4wMpVQW8DjwFa11Z7rrM5pprZNa6/lAObBEKZVR3aqGQin1AaBZa/1euusyhpyttV4IrAD+JdWFT/RP2scMl8ltV6a2RdKuZHQbYQMWAr/SWi8AgsCAY7nHYtBXD1T0Oi4HDqWpLmIcSo0HeBxYqbV+It31GStSXQpeA65Ib01GtbOBa1JjEB4GLlJK/TG9VRrdtNaHUmkz8CRmF0bRP2kfM5i0XaYMbIsyul3J8DaiHqjv9VT7McwgsF9jMehbC9QopapSAxY/DDyT5jqJcSI1GPwBoFZrfW+66zPaKaWKlFK5qX03cAmwLa2VGsW01ndqrcu11pWYf7te1Vp/LM3VGrWUUt7UpBSkuqxcBmTczIxDIO1jhsr0tiuT26JMblcyvY3QWjcCdUqp6amsi4EBJ2+yjUitTiOtdUIpdRvwN8AK/FZrvSXN1Rq1lFIPARcAhUqpeuDbWusH0lurUe1s4OPAptTYAIBvaq2fS1+VRrUJwO9SswZagEe11hk1XbQYViXAk+bnWWzAn7TWL6S3SqNXprePGd7eZXrbJW1RZpI2Ar4ErEx90bcHuHWggmNuyQYhhBBCCCGEEIM3Frt3CiGEEEIIIYQYJAn6hBBCCCGEEGIck6BPCCGEEEIIIcYxCfqEEEIIIYQQYhyToE8IIYQQQgghxjEJ+oQQQgghhBBiHJOgTwghhBBCCCHGMQn6hBillFK/VkqdfZzzlUqpzSNZJyGEEGKsUkp9Wim1SSk14ALWQoxXEvQJMXotBVYP18WVSf4GCCGEyBQ3AhcBN6W7IkKMNPnAJ8QIUEo9rJR6RCm1Rim1Xyl11QnKzwR2aK2TJ7i0TSn1O6XURqXUY0opT69rPKWUek8ptUUp9dlUXqVSqlYp9UtgHVBxqvcmhBBCpNOJesb0sgZoTqVCZBQJ+oQYGfOAPVrrpcBHgW+foPwK4IVBXHc6cL/Wei7QCXyx17lPaq3PBBYBtyulCnq95vda6wVa6/1DuQkhhBBipCmlrCcoMtieMVnAm0DOKVdKiDFGgj4hhplSyg0UAt9JZW0F8k7wsssZXNBXp7X+e2r/j8A5vc7drpTagNkQVgA1qfz9Wuth6zYqhBBCnCql1J+VUvcqpVYBdx6n3KB6xqSGM1wPfAK4fhCBpBDjii3dFRAiA8wBdmqtI6njhcCGgQqnumjmaq0PDeLaur9jpdQFwCXAcq11SCn1GuBKlQkOuuZCCCFEepwB1GqtLzxBucH2jLkI2Ki13pf6QvQi4KVTrKMQY4Y86RNi+M0DJimlXEopL+YTv58cp/yFwKpBXnuSUmp5av9m4K3Ufg7Qlgr4ZgDLTqLeQgghxIhTSrmAfOCeQRQfbM+YjwIPpfYfSh0LkTHkSZ8Qw28esBJ4DcgGftCrS2Z/VgCPDfLatcAtSqlfAzuBX6XyXwA+r5TaCGxnGGcBFUIIIU6z2cAarXXieIUG2zMmNcziWuBipdR/YT708Cml3Frr8OmqtBCjmQR9Qgy/ecBntNbfGGT5s4CvnqiQ1nofMGuAc1HM4LE/cwZZDyGEECIdzgA2DqLcYHvGXAM8r7XufrqnlPo9cDXw6EnVUIgxRrp3CjH8qjGfwg2K1nqh1jo+jPURQgghRrPBBn2DHc/3UeDJo/KeBD42xHoJMWYprY+eB0IIIYQQQojRTSm1DlgqX5QKcWIS9AkhhBBCCCHEOCbdO4UQQgghhBBiHJOgTwghhBBCCCHGMQn6hBBCCCGEEGIck6BPCCGEEEIIIcYxCfqEEEIIIYQQYhyToE8IIYQQQgghxjEJ+oQQQgghhBBiHJOgTwghhBBCCCHGsf8PxSJRRMgn9EIAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "import matplotlib.pyplot as plt\n", "\n", "f, ax = plt.subplots(1,2,figsize=(15,5))\n", "ax[0].plot(isotherm.pressure/BAR, isotherm.total_adsorption*NAV)\n", "ax[0].set_xlabel('$p~~/~~\\\\mathrm{bar}$')\n", "ax[0].set_ylabel('$N$')\n", "ax[0].set_xlim(0,5)\n", "ax[0].set_ylim(0,1)\n", "\n", "for profile in isotherm.profiles[::10]:\n", " ax[1].plot(profile.r/ANGSTROM, (profile.density/(KILO*MOL/METER**3)).T, label=f'{profile.bulk.pressure()/BAR:.2} bar')\n", "ax[1].set_xlabel('$r~~/~~\\\\mathrm{\\AA}$')\n", "ax[1].set_ylabel('$\\\\rho~~/~~\\\\mathrm{kmol/m³}$')\n", "ax[1].set_xlim(0,6)\n", "ax[1].set_ylim(0,8)\n", "ax[1].legend()" ] } ], "metadata": { "kernelspec": { "display_name": "feos", "language": "python", "name": "feos" }, "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.9.10" } }, "nbformat": 4, "nbformat_minor": 5 }