{ "cells": [ { "cell_type": "markdown", "id": "2dab4a3e", "metadata": {}, "source": [ "# Multi-fluid Parameter Fitting\n", "\n", "Here is an example of fitting the $\\beta_T$ and $\\gamma_T$ values for the binary pair of propane+$n$-dodecane with the multi-fluid model. It uses differential evolution to do the global optimization, which is probably overkill in this case as the problem is 2D and other algorithms like Nelder-Mead or even approximate Hessian methods would probably be fine.\n", "\n", "In any case, it takes a few seconds to run (when the actual optimization is uncommented), demonstrating how one can fit model parameters with existing tooling from the scientific python stack." ] }, { "cell_type": "code", "execution_count": 1, "id": "3f8d8bce", "metadata": { "execution": { "iopub.execute_input": "2025-10-15T23:08:32.326891Z", "iopub.status.busy": "2025-10-15T23:08:32.326781Z", "iopub.status.idle": "2025-10-15T23:08:32.943529Z", "shell.execute_reply": "2025-10-15T23:08:32.943025Z" } }, "outputs": [], "source": [ "import json \n", "import teqp, numpy as np, pandas, matplotlib.pyplot as plt\n", "import scipy.interpolate, scipy.optimize\n", "\n", "import pandas\n", "data = pandas.read_csv('VLE_data_propane_dodecane.csv')" ] }, { "cell_type": "code", "execution_count": 2, "id": "b97637e7", "metadata": { "execution": { "iopub.execute_input": "2025-10-15T23:08:32.945352Z", "iopub.status.busy": "2025-10-15T23:08:32.945199Z", "iopub.status.idle": "2025-10-15T23:08:32.952939Z", "shell.execute_reply": "2025-10-15T23:08:32.952248Z" } }, "outputs": [], "source": [ "def cost_function(parameters:np.ndarray, plot:bool=False):\n", "\n", " # Fitting some parameters and fixing the others\n", " betaV, gammaV = 1.0, 1.0\n", " betaT, gammaT = parameters\n", "\n", " # betaT, gammaT, betaV, gammaV = parameters\n", "\n", " BIP = [{\n", " 'function': '',\n", " 'BibTeX': 'thiswork',\n", " 'CAS1': '112-40-3',\n", " 'CAS2': '74-98-6',\n", " 'F': 0.0,\n", " 'Name1': 'n-Dodecane',\n", " 'Name2': 'n-Propane',\n", " 'betaT': betaT,\n", " 'betaV': betaV,\n", " 'gammaT': gammaT,\n", " 'gammaV': gammaV\n", " }]\n", " model = teqp.build_multifluid_model([\"n-Dodecane\", \"n-Propane\"], teqp.get_datapath(), \n", " BIPcollectionpath=json.dumps(BIP)\n", " )\n", " ancs = [model.build_ancillaries(ipure) for ipure in [0,1]]\n", "\n", " cost = 0.0\n", " \n", " # The 0-based index of the fluid to start from. At this temperature, only one fluid \n", " # is subcritical, so it has to be that one, but in general you could start \n", " # from either one.\n", " ipure = 0 \n", "\n", " for T in [419.15, 457.65]:\n", " # Subset the experimental data to match the isotherm \n", " # being fitted\n", " dfT = data[np.abs(data['T / K90'] - T) < 1e-3]\n", "\n", " if plot:\n", " plt.plot(1-dfT['x[0] / mole frac.'], dfT['p / Pa']/1e6, 'X')\n", " plt.plot(1-dfT['y[0] / mole frac.'], dfT['p / Pa']/1e6, 'X')\n", "\n", " try:\n", " # Get the molar concentrations of the pure fluid\n", " # at the starting point\n", " anc = ancs[ipure]\n", " rhoL0 = np.array([0, 0.0])\n", " rhoV0 = np.array([0, 0.0])\n", " rhoL0[ipure] = anc.rhoL(T)\n", " rhoV0[ipure] = anc.rhoV(T)\n", "\n", " # Now we do the trace and convert retuned JSON\n", " # data into a DataFrame\n", " df = pandas.DataFrame(model.trace_VLE_isotherm_binary(T, rhoL0, rhoV0))\n", " \n", " if plot:\n", " plt.plot(df['xL_0 / mole frac.'], df['pL / Pa']/1e6)\n", " plt.plot(df['xV_0 / mole frac.'], df['pL / Pa']/1e6)\n", "\n", " # Interpolate trace at experimental pressures along this \n", " # isotherm to get composition from the current model\n", " # The interpolators are set up to put in NaN for out\n", " # of range values\n", " x_interpolator = scipy.interpolate.interp1d(\n", " df['pL / Pa'], df['xL_0 / mole frac.'], \n", " fill_value=np.nan, bounds_error=False\n", " )\n", " y_interpolator = scipy.interpolate.interp1d(\n", " df['pL / Pa'], df['xV_0 / mole frac.'], \n", " fill_value=np.nan, bounds_error=False\n", " )\n", " # The interpolated values for the compositions \n", " # along the trace at experimental pressures\n", " x_model = x_interpolator(dfT['p / Pa'])\n", " y_model = y_interpolator(dfT['p / Pa'])\n", " if plot:\n", " plt.plot(x_model, dfT['p / Pa']/1e6, '.')\n", " \n", " # print(x_model, (1-dfT['x[0] (-)']))\n", " \n", " errTx = np.sum(np.abs(x_model-(1-dfT['x[0] / mole frac.'])))\n", " errTy = np.sum(np.abs(y_model-(1-dfT['y[0] / mole frac.'])))\n", " \n", " # If any point *cannot* be interpolated, throw out the model,\n", " # returning a large cost function value.\n", " #\n", " # Note: you might need to be more careful here,\n", " # if the points are close to the critical point, a good model might\n", " # (but not usually), undershoot the critical point of the \n", " # real mixture\n", " # \n", " # Also watch out for values of compositons in the data that are placeholders\n", " # with a value of nan, which will pollute the error calculation\n", " if not np.isfinite(errTx):\n", " return 1e6\n", " if not np.isfinite(errTy):\n", " return 1e6\n", " cost += errTx + errTy\n", "\n", " except BaseException as BE:\n", " print(BE)\n", " pass \n", " if plot:\n", " plt.title(f'dodecane(1) + propane(2)')\n", " plt.xlabel('$x_1$ / mole frac.'); plt.ylabel('$p$ / MPa')\n", " plt.savefig('n-Dodecane+propane.pdf')\n", " plt.show()\n", "\n", " return cost" ] }, { "cell_type": "code", "execution_count": 3, "id": "4ee30326", "metadata": { "execution": { "iopub.execute_input": "2025-10-15T23:08:32.953932Z", "iopub.status.busy": "2025-10-15T23:08:32.953822Z", "iopub.status.idle": "2025-10-15T23:08:33.512557Z", "shell.execute_reply": "2025-10-15T23:08:33.511850Z" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAi0AAAHKCAYAAADPQHvVAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjcsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvTLEjVAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAoWVJREFUeJzs3WV0VNfXgPHnzkzcSQhOCITgUkqLEyy4FggOQQptaan+W2jfutBSo1SRUtxdi6QQvLi7JlCcuM/MfT8kGZgkECHJRPZvrSy4vieBzJ5z9jlHUVVVRQghhBCigNNYOgAhhBBCiKyQpEUIIYQQhYIkLUIIIYQoFCRpEUIIIUShIEmLEEIIIQoFSVqEEEIIUShI0iKEEEKIQkGSFiGEEEIUCpK0CCGEEKJQkKRFiEd88sknKIqSq/esVKkSgYGBuXrPgiQ0NBRbW1t2796d7Wvv37+Pg4MDGzZsyIPIRE517tyZF198MUfX9u/fn4CAgFyOSIhkkrQIIZ7KZ599RqNGjWjWrJlp37lz53jzzTdp2rQptra2KIrC1atX013r7u7OqFGj+PDDD/MxYvEku3fvZvPmzbz33numfWfPnuXdd9+lfv36ODk5UaZMGbp06cLBgwfTXf/ee++xfPlyjh07lp9hi2JCkhYhRI7dvXuX2bNn89JLL5nt37t3L1OmTCEqKooaNWo88R4vvfQShw8f5p9//snLUEUWffvtt7Rt2xYfHx/TvhkzZjB9+nQaNmzI999/z1tvvcW5c+do3LgxW7duNbv+mWeeMZ0nRG6TpEUIkWPz5s1Dp9PRrVs3s/3du3cnPDycEydOMGjQoCfeo0aNGtSuXZtZs2blKIbAwEBatWqVo2ufRmxsbL4/M6/duXOH9evXp+veGTBgAKGhocyYMYPRo0fzv//9j3///ZcSJUrwySefpLtPQEAAK1asIDo6Op8iF8WFJC2i2Nq1axfPPfcctra2VKlShalTp2Z4nl6v5/PPP6dKlSrY2NhQqVIl3n//fRISEszOU1WVL774gvLly2Nvb0/r1q05depUhvcMDw/njTfeoEKFCtjY2ODj48M333yD0Wg0O89oNPLTTz9Rp04dbG1tKVmyJB07djRrlv/rr79o06YNnp6e2NjYULNmTX7//fd0z6xUqRJdu3Zl165dPP/889ja2lK5cmXmzJmT4/hWrVpFo0aNcHR0NNtfokQJnJycMnztGfH392ft2rXk56Lz27dvR1EUFi9ezPvvv0/p0qVxcHCge/fuhIaGmp3bqlUrateuzaFDh2jZsiX29va8//77QPIb/ciRIylVqhS2trbUq1eP2bNnm11/9epVFEXhu+++48cff8TLyws7Ozv8/Pw4efKk2bnHjx8nMDCQypUrY2trS+nSpRkxYgT37983Oy+1/urixYsEBgbi6uqKi4sLw4cPzzChmjdvHs8++yx2dnaUKFGC/v37p3ud69evR6/X065dO7P9zz77bLqfsbu7Oy1atODMmTPpnuXv709MTAxbtmzJ6FsvRI7pLB2AEJZw4sQJ2rdvT8mSJfnkk0/Q6/V8/PHHlCpVKt25o0aNYvbs2fTp04e3336bf//9l4kTJ3LmzBlWrlxpOu+jjz7iiy++oHPnznTu3JnDhw/Tvn17EhMTze4XGxuLn58fN27cYMyYMVSsWJE9e/YwYcIEbt68yeTJk03njhw5klmzZtGpUydGjRqFXq9n586d7Nu3j4YNGwLw+++/U6tWLbp3745Op2Pt2rW88sorGI1Gxo4da/bsixcv0qdPH0aOHMmwYcOYOXMmgYGBPPvss9SqVStb8SUlJXHgwAFefvnlp/55PPvss/z444+cOnWK2rVrP/X9suPLL79EURTee+897ty5w+TJk2nXrh1Hjx7Fzs7OdN79+/fp1KkT/fv3Z/DgwZQqVYq4uDhatWrFxYsXefXVV/H29mbp0qUEBgYSHh7O66+/bvasOXPmEBUVxdixY4mPj+enn36iTZs2nDhxwvRvb8uWLVy+fJnhw4dTunRpTp06xbRp0zh16hT79u1LVygeEBCAt7c3EydO5PDhw8yYMQNPT0+++eYbs9f44YcfEhAQwKhRo7h79y4///wzLVu25MiRI7i6ugKwZ88e3N3d8fLyytL37tatW3h4eKTbX7NmTezs7Ni9eze9evXK0r2EyBJViGKoZ8+eqq2trXrt2jXTvtOnT6tarVZ99L/F0aNHVUAdNWqU2fXvvPOOCqj//POPqqqqeufOHdXa2lrt0qWLajQaTee9//77KqAOGzbMtO/zzz9XHRwc1PPnz5vdc/z48apWq1VDQkJUVVXVf/75RwXUcePGpYv/0WfExsamO96hQwe1cuXKZvu8vLxUQN2xY4dp3507d1QbGxv17bffznZ8Fy9eVAH1559/Tvf8R3377bcqoF65cuWx5+zZs0cF1MWLFz/xXhkZNmyY6ufnl+3rtm3bpgJquXLl1MjISNP+JUuWqID6008/mfb5+fmpgPrHH3+Y3WPy5MkqoM6bN8+0LzExUW3SpInq6Ohouu+VK1dUQLWzs1OvX79uOvfff/9VAfXNN9807cvo57lw4cJ0P7uPP/5YBdQRI0aYndurVy/V3d3dtH316lVVq9WqX375pdl5J06cUHU6ndn+5s2bq88+++xjvmPmduzYoSqKon744YcZHvf19VU7deqUpXsJkVXSPSSKHYPBwKZNm+jZsycVK1Y07a9RowYdOnQwOzd1KO5bb71ltv/tt98GkpvTAbZu3UpiYiKvvfaa2SfhN954I93zly5dSosWLXBzc+PevXumr3bt2mEwGNixYwcAy5cvR1EUPv7443T3ePQZj7YGREREcO/ePfz8/Lh8+TIRERFm19WsWZMWLVqYtkuWLEm1atW4fPlytuNL7a5wc3NLF192pd7j3r17TzzPaDSaxXTv3j0SEhJISkpKtz8pKSlLzx46dKhZV1afPn0oU6ZMumHYNjY2DB8+3Gzfhg0bKF26NAMGDDDts7KyYty4cURHRxMcHGx2fs+ePSlXrpxp+/nnn6dRo0Zmz3r05xkfH8+9e/do3LgxAIcPH04Xf9oi6BYtWnD//n0iIyMBWLFiBUajkYCAALPvT+nSpalatSrbtm0zXXv//v0s/Tzv3LnDwIED8fb25t13383wnNR/P0LkJukeEsXO3bt3iYuLo2rVqumOVatWzewN5Nq1a2g0GrORFAClS5fG1dWVa9eumc4D0t2zZMmS6d4ELly4wPHjxylZsmSG8d25cweAS5cuUbZsWUqUKPHE17N7924+/vhj9u7dm66WISIiAhcXF9P2o0laKjc3N8LCwrIdXyo1F+pQUu+R2Rw5ISEheHt7Z3gsbbzbtm3LUoFu2p+Zoij4+PikG6Jdrlw5rK2tzfZdu3aNqlWrotGYf/5LHTGV+u/icc8C8PX1ZcmSJabtBw8e8Omnn7Jo0aJ03+u0SSik/5mm/nsLCwvD2dmZCxcuoKpqhs+G5CTrUZn9PGNiYujatStRUVHs2rUrXa3Lo/fJ7TmPhJCkRYgsyM1fvkajEX9//8d+QvX19c3yvS5dukTbtm2pXr06P/zwAxUqVMDa2poNGzbw448/piuc1Wq1Gd7n0TeqrMbn7u4OYJbw5FTqPTKqj3hU6dKl0xV3fvvtt9y6dSvdENt69eo9dVyPerQFJC8FBASwZ88e/ve//1G/fn0cHR0xGo107Ngx3c8TMv+ZGo1GFEVh48aNGZ77aNLh7u7+xJ9nYmIiL7zwAsePH2fTpk1PrD8KCwt7bKIkRE5J0iKKnZIlS2JnZ8eFCxfSHTt37pzZtpeXF0ajkQsXLpjNN3L79m3Cw8NNBYupf164cIHKlSubzrt79266N4EqVaoQHR2dboRGWlWqVGHTpk08ePDgsa0ta9euJSEhgTVr1ph94n60yT+7shpfxYoVsbOz48qVKzl+VqrUe2Q2p4utrW26uObNm0dCQkKm8T5O2n8Hqqpy8eJF6tatm+m1Xl5eHD9+HKPRaNbacvbsWdPxJz0L4Pz581SqVAlIfqMPCgri008/5aOPPnridVlVpUoVVFXF29s704S4evXqLF++PMNjRqORoUOHEhQUxJIlS/Dz83vsffR6PaGhoXTv3j3HcQuREalpEcWOVqulQ4cOrFq1ipCQENP+M2fOsGnTJrNzO3fuDGA2ogfghx9+AKBLly4AtGvXDisrK37++WezVou010HyJ+m9e/emexYkDzXW6/UA9O7dG1VV+fTTT9Odl/qM1E/Ojz4zIiKCv/76K+MXnwVZjc/KyoqGDRtmOCtqdh06dAgXFxfTCKb8lDqiJ9WyZcu4efMmnTp1yvTazp07c+vWLRYvXmzap9fr+fnnn3F0dEz3xr5q1Spu3Lhh2t6/fz///vuv6VkZ/Twh439HWfXCCy+g1Wr59NNP091XVVWzodRNmjQhLCzMrMYp1WuvvcbixYv57bffeOGFF574zNOnTxMfH0/Tpk1zHLcQGZGWFlEsffrpp/z999+0aNGCV155xfRGU6tWLY4fP246r169egwbNoxp06YRHh6On58f+/fvZ/bs2fTs2ZPWrVsDya0377zzDhMnTqRr16507tyZI0eOsHHjxnRdHv/73/9Ys2YNXbt2NQ03jomJ4cSJEyxbtoyrV6/i4eFB69atGTJkCFOmTOHChQum7oGdO3fSunVrXn31Vdq3b4+1tTXdunVjzJgxREdHM336dDw9Pbl582aOvjdZjQ+gR48efPDBB0RGRuLs7Gy6R0REBD///DOAaU2iX375BVdXV1xdXXn11VfNnrllyxa6detmkRqIEiVK0Lx5c4YPH87t27eZPHkyPj4+WVp7Z/To0UydOpXAwEAOHTpEpUqVWLZsGbt372by5Mnp5qrx8fGhefPmvPzyyyQkJDB58mTc3d1NXXHOzs60bNmSSZMmkZSURLly5di8efNTtWZVqVKFL774ggkTJnD16lV69uyJk5MTV65cYeXKlYwePZp33nkHSE7CdTodW7duZfTo0aZ7TJ48md9++40mTZpgb2/PvHnzzJ7Rq1cvHBwcTNtbtmzB3t4ef3//HMctRIbyf8CSEAVDcHCw+uyzz6rW1tZq5cqV1T/++MM0jPRRSUlJ6qeffqp6e3urVlZWaoUKFdQJEyao8fHxZucZDAb1008/VcuUKaPa2dmprVq1Uk+ePKl6eXmZDXlWVVWNiopSJ0yYoPr4+KjW1taqh4eH2rRpU/W7775TExMTTefp9Xr122+/VatXr65aW1urJUuWVDt16qQeOnTIdM6aNWvUunXrqra2tmqlSpXUb775Rp05c2a6YcZeXl5qly5d0n0f/Pz80g0Zzmp8t2/fVnU6nTp37lyz61OH+Gb05eXlZXbumTNnVEDdunVrutiy4mmHPC9cuFCdMGGC6unpqdrZ2aldunQxGwqvqsnfo1q1amV4n9u3b6vDhw9XPTw8VGtra7VOnTrqX3/9ZXZO6vfj22+/Vb///nu1QoUKqo2NjdqiRQv12LFjZudev35d7dWrl+rq6qq6uLioffv2Vf/77z8VUD/++GPTean/Vu/evWt2/V9//ZXhEPPly5erzZs3Vx0cHFQHBwe1evXq6tixY9Vz586Znde9e3e1bdu2ZvuGDRv22J9nRs9q1KiROnjw4Ay/X0I8DUVV83EKSiFEkTNy5EjOnz/Pzp07c3T9G2+8wY4dOzh06FC+trRs376d1q1bs3TpUvr06ZOnz7p69Sre3t58++23plaNgmrnzp20atWKs2fP5qiQ9ujRozRo0IDDhw9Tv3793A9QFGtS0yKEeCoff/wxBw4cMHUDZcf9+/eZMWMGX3zxhQyPLSBatGhB+/btmTRpUo6u//rrr+nTp48kLCJPSE2LEOKpVKxYkfj4+Bxd6+7uLovqFUAbN27M8bWLFi3KxUiEMCctLUIIIYQoFKSmRQghhBCFgrS0CCGEEKJQkKRFCCGEEIVCkSnENRqN/Pfffzg5OckoBCGEEKKQUFWVqKgoypYtm27x0bSKTNLy33//UaFCBUuHIYQQQogcCA0NpXz58k88p8gkLanTZYeGhppNJy6EEEKIgisyMpIKFSqkW/YiI0UmaUntEnJ2dpakRQghhChkslLaIYW4QgghhCgUJGkRQgghRKEgSYsQQgghCgVJWoQQQghRKEjSIoQQQohCQZIWIYQQQhQKkrQIIYQQolCQpEUIIYQQhYIkLUIIIYQoFCRpEUIIIUShIEmLEEIIIQoFSVqEEEIIUSgUmQUTiwNVNRAVdZqY2EtERV3jzOkjREa5UKf2cNzc3ChTpgzW1taWDlMIIYTIE5K0FAIxMZcJCZnO3XtbSEoKM+13dILYWG9WrVoFgEajoUyZMlSvXp3atWvj5uZmoYiFEEKI3KeoqqpaOojcEBkZiYuLCxERETg7O1s6nFyRkHCHCxe+5Pad9UDyj0mrdcTJqRY6XSkuXAjl7h0H7O3bcO/ePaKiosyur1ChAk2aNKF69epoNNITKIQQouDJzvu3tLQUUPfv7+DkqTfQ6yMA8PBoS4UKw3F1aYhGYwVAvbrm14SHh3Pp0iVOnjzJlStXCA0NJTQ0FA8PD5o3b06dOnXQarX5/VKEEEKIXCEtLQXQjRuLOHvuQ8CIk1MtalSfiJNTrWzdIzIykoMHD7J//37i4+MBKFGiBB07dsTX1zcPohZCCCGyLzvv35K0FDD//beEM2cnAFCmTF+qV/sUjcYmx/eLj4/n4MGD7N27l5iYGACqVatGx44dpeZFCCGExUnSUkiTlvsPdnHs2AhU1UDFCiPx8ZmAoii5cu+EhASCg4PZt28fRqMRnU5H8+bNad68OTqd9BIKIYSwDElaCmHSEh9/k3/3d0Gvj6B0qZ7UrPldriUsj7pz5w4bNmzg6tWrAJQuXZrevXtTsmTJXH+WEEIIkZnsvH/LkJICQFWNnD7zP/T6CJycalOjxld5krAAeHp6MmzYMHr37o2dnR23bt1i6tSpHDx4kCKSvwohhCiiJGkpAG7dWklY2F40Gjtq15r8VDUsWaEoCnXq1OHll1+mcuXK6PV61q1bx6JFi0x1L0IIIURBI0mLhSUlRXLh4tcAVPYeh729d74929nZmcGDB9O+fXs0Gg3nzp3j999/5/Lly/kWgxBCCJFVkrRY2LVrv5OU9AB7ex8qVAjM9+drNBqaNm3Kiy++iIeHB9HR0cydO5e9e/dKd5EQQogCRZIWC0pIvEfo9bkAVPV5D43GcusGlSlThtGjR1OvXj1UVWXTpk2sXr0avV5vsZiEEEKIR0nSYkEhITMwGuNwdq6Hu3trS4eDtbU1PXv2pEOHDiiKwtGjR5k1a1a65QGEEEIIS5CkxUL0+hj++28RAN6VXs2z0ULZpSgKTZo0YfDgwdja2nL9+nWmTZvGjRs3LB2aEEKIYk6SFgu5dWsVen0UdnZeuLu3snQ46VSpUsVU5xIVFcXMmTM5fvy4pcMSQghRjEnSYgGqqhJ6fTYAFcoPRVEK5o/B3d2dUaNG4evri8FgYMWKFQQFBUmBrhBCCIsomO+WRVxk5BFiYy+h0dhRpkxvS4fzRLa2tvTv358WLVoAsHPnTlavXo3BYLBwZEIIIYobSVos4OatVQB4enZAp3OybDBZoNFoaNu2Ld26dTMV6C5atIjExERLhyaEEKIYkaQlnxmNCdy+vQ6AMqVfsHA02fPss8/Sr18/dDodFy5cYM6cOcTGxlo6LCGEEMWEJC357MGDPej1EdhYl8LNrbGlw8m26tWrM3ToUNPIoj///JPw8HBLhyWEEKIYkKQln929twWAkiXboyhaC0eTMxUrVmTEiBE4Oztz//59/vzzT27fvm3psIQQQhRxkrTkI1U1cu9eEAAeJdtZOJqn4+npyciRIylZsqRpSPTVq1ctHZYQQogiTJKWfBQZeZTExHvodE64uT5v6XCemouLCyNGjKBixYokJCQwd+5cTp8+bemwhBBCFFGStOSju/f+AcDdvZVF1xnKTXZ2dgwZMoRq1aphMBhYunQpR44csXRYQgghiiBJWvLRgwe7AArkDLhPw8rKioCAABo0aICqqqxevZp///3X0mEJIYQoYgpM0mIwGPjwww/x9vbGzs6OKlWq8PnnnxeZ2VeTksKJijoJQAm3JhaOJvdptVq6detG48bJI6I2btzIzp07LRyVEEKIokRn6QBSffPNN/z+++/Mnj2bWrVqcfDgQYYPH46Liwvjxo2zdHhPLSzsX0DFwaEqNjalLB1OnlAUhQ4dOmBtbc2OHTsICgoiMTGRNm3aFJgFIYUQQhReBSZp2bNnDz169KBLly4AVKpUiYULF7J//34LR5Y7HoTtBsCtCLayPEpRFNq0aYO1tTVbt25l586dJCYm0rFjR0lchBBCPJUC0z3UtGlTgoKCOH/+PADHjh1j165ddOrUKcPzExISiIyMNPsqyMLDk5Ovotg1lJHmzZvTuXNnAP799182bNhQZLr6hBBCWEaBaWkZP348kZGRVK9eHa1Wi8Fg4Msvv2TQoEEZnj9x4kQ+/fTTfI4yZ5KSIomJuQCAi0sDC0eTf55//nl0Oh1r1qzhwIEDqKpK586d0WgKTK4shBCiECkw7x5Llixh/vz5LFiwgMOHDzN79my+++47Zs+eneH5EyZMICIiwvQVGhqazxFnXWTkUQDs7Cpibe1h2WDyWYMGDejZsycABw8eZN26dRiNRssGJYQQolAqMC0t//vf/xg/fjz9+/cHoE6dOly7do2JEycybNiwdOfb2NhgY2OT32HmSERE8rwlLs7Fp5XlUfXr10dRFFatWsXhw4dRVZVu3bpJi4sQQohsKTDvGrGxsenexLRabZH4VB4RmZK0uDxj4Ugsp169evTq1QtFUThy5Ahr1qwpEj9bIYQQ+afAtLR069aNL7/8kooVK1KrVi2OHDnCDz/8wIgRIywd2lNRVePDlpZinLQA1K1bF0VRWLFiBUePHkVVVXr06CEtLkIIIbKkwCQtP//8Mx9++CGvvPIKd+7coWzZsowZM4aPPvrI0qE9lbi4axgM0Wg0Njg4VLN0OBZXp04dFEVh+fLlHDt2DFVV6dmzpyQuQgghMlVgkhYnJycmT57M5MmTLR1KroqKSl5A0NGxOhpNgfl2W1Tt2rVNicvx48cxGo306tULrVZr6dCEEEIUYPLxNo9FRZ8BwNGxhoUjKVhq1apF37590Wg0nDx5khUrVmAwGCwdlhBCiAJMkpY8Fh11CgAnp1oWjqTgqVGjBgEBAWg0Gk6dOsXy5cslcRFCCPFYkrTksdSWFidpaclQ9erV6devH1qtltOnT7Ns2TJJXIQQQmRIkpY8lJBwl8TEu4CCo6MU4T5OtWrVTInLmTNnWLp0qSQuQggh0pGkJQ/FxCSvo2RvXwmt1t7C0RRsvr6+9O/fH61Wy9mzZ6WrSAghRDqStOShmNhLANjbV7FwJIVD1apVTYnL6dOnWblypUxAJ4QQwkSSljwUE5OctDg4+Fg4ksKjatWqpuLckydPsnr1aklchBBCAJK05KnYmIsAONhXtnAkhUu1atXo06cPiqJw7Ngx1q5dK4mLEEIISVrykql7SFpasq1mzZr07t3btFbRhg0bUFXV0mEJIYSwIEla8khSUmTKyCFpacmp2rVr07NnTwAOHjzI33//LYmLEEIUY5K05JHYlFYWG+tS6HROFo6m8KpXrx49evQA4N9//2Xz5s2SuAghRDElSUseedg1JK0sT+uZZ56ha9euAOzdu5egoCBJXIQQohiSpCWPxMWFAGBv723hSIqGhg0b0rlzZwB27drF9u3bLRuQEEKIfCdJSx5JTVrsbCtYOJKi4/nnn6dDhw4ABAcHs2PHDgtHJIQQIj9J0pJH4uJCAbCzq2jhSIqWJk2a4O/vD8A///zDrl27LByREEKI/CJJSx4xtbRI0pLrmjVrRps2bQDYunUre/futXBEQggh8oMkLXlAr48mKekBAHZ20j2UF1q2bImfnx8AmzZtYv/+/RaOSAghRF6TpCUPxMVfB8DKyk2GO+ehVq1a0bx5cwA2bNjAwYMHLRyREEKIvCRJSx6IlyLcfKEoCm3btqVJkyYArFu3jiNHjlg4KiGEEHlFkpY8kFqEaytdQ3lOURTat29Po0aNAFi9ejXHjh2zcFRCCCHygiQteSA+/j8AbG3LWjiS4kFRFDp27EjDhg0BWLVqFSdPnrRwVEIIIXKbJC15ID7hJgC2NmUsHEnxoSgKnTt3pkGDBqiqyvLlyzl9+rSlwxJCCJGLJGnJAwnxKUmLrSQt+Umj0dC1a1fq1auHqqosW7aMs2fPWjosIYQQuUSSljyQ2tJiIy0t+U6j0dCjRw/q1KmD0WhkyZIlnD9/3tJhCSGEyAWStOQyozGRxMR7gLS0WIpGo6Fnz57UrFkTo9HI4sWLuXTpkqXDEkII8ZQkacllCQm3ARWNxhorK3dLh1NsabVaevfuTfXq1TEYDCxcuJArV65YOiwhhBBPQZKWXJY6csjGpjSKolg4muJNq9XSp08fqlatil6vZ8GCBVy7ds3SYQkhhMghSVpyWULCLUDqWQoKnU5HQEAAVapUISkpifnz5xMaGmrpsIQQQuSAJC25LCHxDgA2NqUe7ouP4e+PRrDu/cEkxsdYKrRiy8rKiv79++Pt7U1iYiLz5s3j+vXrlg5LCCFENknSkssSE+4CYGNdkqRbt7j/55/cevMdvJbspcqKQ7y6cQw/HPyBbSHbiEmSBCa/WFlZMWDAALy8vEhISGDevHn8999/lg5LCCFENiiqqqqWDiI3REZG4uLiQkREBM7OzhaL4+SpN7l9ew2lrj6P9rvjYDSajkXbwovjtBi0ybUuNlob/Mr70ce3D43LNJYamHyQmrCEhoZiZ2fHsGHDKF26tKXDEkKIYis779+StOSyQ/v6Eh57GNdZWuz3a7Fv2BDH1q1ILFeSc6X03LHXc/r+afbf2k9o1MPaiuolqhNYK5COlTqi1WgtFn9xEB8fz9y5c7lx4wb29vYEBgbi6elp6bCEEKJYkqTFQkmLITyc3Zsbk+SRhOdcT7xHfodj82YZnquqKmcenGHVxVWsuriKOH0cAL5uvrzd8G2alm2an6EXO3FxccyZM4ebN2/i4OBAYGAgJUuWtHRYQghR7EjSYoGkRVVVQkeP4WL3IFQHeLbyXFwrZS3xCI8PZ/G5xcw+PZuoxCgAmpdrzjsN36GKa5W8DLtYi42NZc6cOdy6dQtHR0eGDx+Ou7vMrSOEEPkpO+/fUoibSyLXrSN67w5Uh+Rth3I1s3ytq60rY+qNYeMLGxlcYzA6jY5dN3bRZ00ffj7yM4mGxDyKunizt7dnyJAheHp6Eh0dzezZs3nw4IGlwxJCCPEYkrTkAmNCAncmfYshJUFUFGt0Opds38fFxoX3nn+P1T1W07pCa/SqnmnHp9FnbR+O3jmau0ELABwcHBg6dCgeHh5ERkYye/ZswsPDLR2WEEKIDEjSkgsiVqxAf/cuSpXkrgUba4+nGglU0bkiU9pM4YdWP+Bu686ViCsM3TiUif9OJDYpNrfCFikcHR0ZNmwY7u7uREREMHv2bCIiIiwdlhBCiDQkaXlKqtHI/Zl/AeDQvQ0A1tYeuXJvfy9/VvdcTU+fnqioLDi7gBfWvCCtLnnAycmJYcOG4ebmRlhYGLNnzyYyMtLSYQkhhHiEJC1PKfbgQZJCQ9E4OmL1XDUArKxL5Nr9XWxc+LzZ50z1n0pZh7LciL5B4N+B/H7sd/RGfa49R4CzszPDhg3D1dWVBw8eMHv2bKKioiwdlhBCiBSStDylyLVrAXDq2AG9kjzDrZWVW64/p2nZpizrvowulbtgUA38dvQ3RmwawY3oG7n+rOLM1dWVYcOG4eLiwv3795k9ezbR0dGWDksIIQSStDwVNTGRyL83AeDSvTtJifcBsLbKvZaWRzlZO/F1i6/5qvlXOFg5cOTOEfqs6cP6y+vz5HnFlZubG8OGDcPZ2Zl79+5J4iKEEAWEJC1PIfbIUYxRUWg9PLBv2JCkpDAArDJJWiKDQrg+fieRQSEZbmemW5VuLO22lHol6xGdFM34neOZsHMC0YnyxppbSpQowbBhw3BycuLu3bvMmTOH4MVz+b5fN/YtXwTA3uUL+b5fN/YuX2jhaIUQoniQpOUpxOzaBYBjs6YoGg2JSclzfFg/oaYlMiiEyC3Xkv++5Rp3Z5ww285q4lLBqQKzOs7i5Xovo1E0rLu8jgHrB3Ax7OLTvCTxCHd3dwIDA3F0dOTOnTsEHzlBpLMrC3ftZto3X7BnyXxAZc+S+ZK4CCFEPpCk5SnE7N4NgEOz5Kn6H7a0PL6mJTVBSZVwMfyJx59Ep9HxSv1XmNVxFp72nlyNvMrADQNZd3ldlu8hniw1cVH0SZyqVINpA95icfeRfNywM8erP2s6b8+SBRaMUgghigdJWnLIEB5O/OnTADg0TZ6uPykxJWl5QkuLs7/XE++b2fGMPOP5DEu7LaVxmcbE6eOYsHMCX+77UmbSzSUeHh5Ur1WTHb71UVPm31E1Gja37EGUQ/KMgs0CBlkyRCGEKBYkacmhuBMnAbD28kLnkTwvi6l76Ak1Lc5tK2Lj45rhMRsfV5zbVsxRPCVsS/BHuz8YXXc0AIvOLSLw70BuRt/M0f2EuXIdu5sSllSqRkOYizsV69Snce/+FopMCCGKD0lacijuxHEAbOvWBcBoTMBgSC6EfVIhbmRQSLouoVQJF8OzXNOSEa1Gy2vPvMavbX/F2dqZE/dOELAugD039uT4niJZ2La/UYxGs32KasQ1KoyQE0dNxblCCCHyjiQtORSf0tJiV6c2AElJqdO+a9DpnB57XWY1K9mpaXmcluVbsrjrYmqUqEF4QjgvbX2J34/9jlE1Zn6xyND5xXNov2O1KXFRVCMtzx9DV6IUqkbD7iXzLRyhEEIUfTpLB1AYqapK3IkTANjWrgNAkj45adHpnFGUx+eCzv5eZomJjY+rWctLTmpaMlLeqTxzO8/l6/1fs+z8Mn47+hvH7x5nYvOJuNq65sozipOmAQNhyXy8Qy8Q5uJOldKleRAdh8HeibgKVWn1XANLhyiEEEWetLTkgP7OXQz37oFGg23NGsn7UlparKycn3itc9uKpsTE2d+LkqPqmG3ntKYlIzZaGz5u8jGfN/scG60Nu27sImBdACfvncy1ZxQXTXoPoGnAIJxiohjQvBnj3vuAxr6VwaDHYO/EpZgkEhOl8FkIIfKSoqqqaukgckNkZCQuLi5ERETg7PzkxOFpxezdS8jwEVh7eVFl098A3L0XxPHjo3FyqsPzz63K0+fnxLkH53hz+5uERoVipbFi/PPj6evb96lWoxZw/fp15syZQ2JiIt7e3gwcOBArKytLhyWEEIVGdt6/paUlBxIuXwbAukoV0z5TS4vOxSIxZaZaiWos6rqI1hVak2RM4vN9n/PBrg+I08dZOrRCrXz58gwePBhra2uuXLnCokWLSEpKsnRYQghRJEnSkgOJl68AYFPZ27TPVNNiVTCTFgBna2d+av0Tbz77JhpFw9rLaxm4fiBXI65aOrRCrWLFigwaNAgrKysuXbrE4sWL0etlBW4hhMhtkrTkQMLlSwBYe1c27XtY01JwkxYARVEYUXsEM9rPwN3WnYvhF+m/vj9br221dGiFmpeXlylxuXjxoiQuQgiRByRpyYHUlhbrjFpaCmj3UFrPlX6Opd2W0sCzATFJMby5/U2+O/AdSUbp2sipSpUqMXDgQHQ6HRcuXGDp0qWSuAghRC6SpCWbjLGx6G/fBsCm8qMtLZEAWOnytgg4N5W0L8mMDjMIrBUIwOzTsxm1aRR3Y+9aNrBCzNvbmwEDBqDT6Th37hzLli3DYDBYOiwhhCgSJGnJpqSbydPia5yc0Lo8bFUpDDUtGbHSWPF2w7f5odUPOFg5cPjOYfqu7cuBWwcsHVqhVaVKFfr3749Wq+Xs2bOSuAghRC6RpCWbkv5LTlqsSpc226/XJ7e06ApRS8uj/L38WdRlET6uPtyPv8+Lm1/kr5N/UURGxOc7Hx8fU+Jy5swZSVyEECIXSNKSTUm3kpMWXdkyZvsNhpjk/VqHfI8pt1RyqcT8zvPpWrkrBtXAD4d+4I1tbxCVGGXp0AqlqlWr0q9fP1PiIjUuQgjxdCRpySZ9SveQVRnzpEWvT05atLrCm7QA2FvZ81Xzr/iw8YdYaaz4J/Qf+q/rz7kH5ywdWqHk6+tr1lUkiYsQQuScJC3ZlHTzFgBWpR/X0uKY7zHlNkVRCKgWwJxOcyjjUIaQqBAGbxjM6ourLR1aoVS1alUGDBiAVqvl3LlzLFmyRBIXIYTIAUlasim1ENcqXfdQNADaIpC0pKrtUZslXZfQrFwz4g3x/N/u/+PTvZ+SYEiwdGiFjo+Pj2lU0fnz52UeFyGEyAFJWrIpo+4hozEJozF5sTxdIe8eSsvV1pXf2v7GK/VfQUFh2fllDNkwhOtR1y0dWqHzaOJy4cIFFi9eLFP+CyFENkjSkg2qqpJ0K7l7SPdI0pLaNQSgLcSFuI+jUTS8XO9lfm/3O642rpx5cIZ+6/qx4/oOS4dW6FSpUsVsAjpJXIQQIuskackGY2QkakJy14jO09O0P7UIV6OxRqMpuiv8NivXjCVdl1DHow6RiZGMDRrLL0d+wWCUobzZUblyZbMp/2WRRSGEyBpJWrJB/+ABABpHRzTW1qb9RbGe5XHKOJZhVsdZ9KvWD4Cpx6fy8taXeRD/wMKRFS7e3t5miywuXLiQxMRES4clhBAFWoFKWm7cuMHgwYNxd3fHzs6OOnXqcPDgQUuHZWIICwNA6+Zmvj+le6godg1lxFprzf81/j++bvE1djo79t7cS8DaAI7dPWbp0AqVSpUqmRKXy5cvS+IihBCZKDBJS1hYGM2aNcPKyoqNGzdy+vRpvv/+e9zSJAiWZEhpadGWMI9Jr09uaSlqRbiZ6VK5Cws6L6CScyVux94m8O9AFpxZILPoZkOlSpUYPHgw1tbWXLlyhQULFkjiIoQQj1FgkpZvvvmGChUq8Ndff/H888/j7e1N+/btqVKliqVDM0ntHtK5lTDfX8xaWh7l4+bDwi4L8ffyR2/UM3H/RN7b8R6xSbGWDq3Q8PLyMiUuV69elcRFCCEeo8AkLWvWrKFhw4b07dsXT09PnnnmGaZPn/7Y8xMSEoiMjDT7ymuGByndQyXMkxaDvvBP4f80HK0d+d7ve9597l10io6NVzcyYP0ALoVfsnRohUbFihUZMmSIKXGZP38+CQkyH44QQjyqwCQtly9f5vfff6dq1aps2rSJl19+mXHjxjF79uwMz584cSIuLi6mrwoVKuR5jIawlJaWEmlrWlIKcXVFvxD3cRRFYUjNIczsOBNPO08uR1xmwPoBrL+83tKhFRoVKlRgyJAh2NjYcO3aNUlchBAijQKTtBiNRho0aMBXX33FM888w+jRo3nxxRf5448/Mjx/woQJREREmL5CQ0PzPEZ9akvLY7qHisIU/k/rGc9nWNxtMY1KNyJOH8f4neP5Yt8XJBqkuyMrHk1cQkJCmDdvniQuQgiRosAkLWXKlKFmzZpm+2rUqEFISEiG59vY2ODs7Gz2ldceFuJm3D1U2BdLzC0edh5M9Z/KmLpjUFBYfG4xQzbKLLpZVb58eYYOHYqNjQ2hoaHMmzeP+Ph4S4clhBAWV2CSlmbNmnHunPlKwufPn8fLy8tCEaWnf0z30JMKcSODQrg+fieRQSEZbhdVWo2WV595ld/a/YarjSun758mYF0A20O3Wzq0QqFcuXIMHToUW1tbs8Rl7/KFfN+vG/uWLwIwbe9dvtDCEQshRN4rMEnLm2++yb59+/jqq6+4ePEiCxYsYNq0aYwdO9bSoZk8thA3paYlbSFuZFAIkVuuJf99yzXuzjhhtl3UExeA5uWas7TbUuqWrEtUYhSv/fMaPx76Eb1RFgvMzKOJy/Xr1/n9p8n8vXEDIWUr8ff6NSz9/AP2LJkPqOxZMl8SFyFEkVdgkpbnnnuOlStXsnDhQmrXrs3nn3/O5MmTGTRokKVDM3k4udzjuofMa1pSE5RUCRfDn3i8qCrtUJpZHWYxuMZgAGaenMmozaO4G3vXwpEVfGXLlmXYsGHY2dmxz6UU0wa+xeLuI5k66B3W6XVm5+5ZssBCUQohRP4oMEkLQNeuXTlx4gTx8fGcOXOGF1980dIhmRgTE03rDmmdncyOPSzENW9pcfZ/ctdWZseLEiutFe89/x7f+X2Hg5UDh24fou/avuy/ud/SoRV4ZcqUoePAwezwrY+qKACoGg2bW/YgyuFhLVezgIKT4AshRF4oUElLQWaMijL9XeNgnpw8bhp/57YVsfFxzfB+Nj6uOLetmLtBFgIdKnVgUZdFVHWryv34+7y45UWmH5+OUTVaOrQCLcrR2ZSwpFI1GsJc3AGoWKc+jXv3t0RoQgiRbyRpySJjdHLdisbeHkWrNTtmMMQBoNXame2PDApJ1yWUKuFieLGoaclIJZdKzO88nx5VemBUjUw5MoVXg14lPD7c0qEVWJXtbFDSLI+gqEZcoiMACDlx1FScK4QQRZUkLVlkiEpJWpyc0h0zGpOHo2q0tmb7M6tZKS41LRmx09nxRfMv+KzpZ9hobdh5YycB6wI4cfeEpUMrkMraWtM+eBWKMblFSlGNtDx/DJ1HGYw6KwB2L5lvyRCFECLPSdKSRaaWFqf0E8gZDcm1LhqNjdn+tDUrabuKilNNy+P0qtqL+Z3nU9GpIjdjbjL076Gy6OJjvFS3OmPmf0e/NX/y4qIp1Ll1FaONHbFe1TDqrGkaMNDSIQohRJ6SpCWLjNHJNS1ah/RJiyGlpUWrMW9pcW5b0ZSYOPt7UXJUHbPt4ljTkpFqJaqxqOsis0UX393xLtGJ0ZYOrUBp0nsAHbp0o+J/V+ncoSNjX38DOysrVGtbjDWfpXqbjpYOUQgh8pSiFpGPtJGRkbi4uBAREZEns+OGr1zFzQkTcGjenIozzBdy3La9FkZjPE2bBGNnVz7Xn11cqKrKvDPz+OHgD+hVPV7OXnzn9x3VS1S3dGgFVkREBLNnz+bBgwc4OzszdOhQPDw8LB2WEEJkWXbev6WlJYtM3UOO5i0tqqo+tqZFZE/qoot/dfyL0g6luRZ5jUHrB7H47GLpLnoMFxcXhg8fjoeHB5GRkcyaNYs7d+5YOiwhhMgTkrRkkal7KE1Ni9H4cCFAbZqaFpEz9T3rs7TrUlqVb0WiMZEv/v2C/+34H1GJUZlfXAw5OTkRGBiIp6cn0dHRzJo1i1u3blk6LCGEyHWStGSRafSQo/noodRWFkhfiCtyztXWlSltpvBOw3fQKTo2Xd1Ev3X9OH3/tKVDK5AcHR0JDAykTJkyxMbGMmvWLG7cuGHpsIQQIldJ0pJFD7uHzCeQMxoTUv6mQVGs8jmqok1RFIbVGsbsTrMp61CW0KhQBm8YLKOLHsPe3p6hQ4dSrlw54uPjmTNnDqGhoZYOSwghco0kLVn0sHso45YWrdYWJc2MpSJ31C1ZlyXdltC6QmuSjElM3D+Rt4PfJjIx0tKhFTh2dnYMHTqUihUrkpCQwJw5c7h69aqlwxJCiFwhSUsWGVJbWtIMeTYYUopwpWsoT7nYuPBT659477n30Gl0bLm2hYC1AZy6d8rSoRU4NjY2DB48GG9vb5KSkpg3bx4XL160dFhCCPHUJGnJImNMLJB+3SHTyCFJWvKcoigMrjmYuZ3mUs6xHDeibzB442DmnZ4n3UVpWFtbM3DgQKpWrYper2fhwoWcOXPG0mEJIcRTkaQli9S45PWFNPbm6wuljh7SaGS4c36p7VGbJd2W0K5iO/RGPd8c+IY3tr1BREKEpUMrUKysrOjXrx81atTAYDCwZMkSjh8/bumwhBAixyRpySJjStKi2JonJ4ZHalpE/nG2duaHVj8w4fkJWGms+Cf0HwLWytpFael0Ovr06UO9evVQVZUVK1Zw6NAhS4clhBA5IklLFhnjU7qB7OzN90tNi8UoisLAGgOZ23kuFZwq8F/MfwzdOJTZp2ZLd9EjtFotPXr0oGHDhgCsXbuWvXv3WjgqIYTIPklassjUPWRn3qKSOuRZuocsp5Z7LRZ3XUx7r/boVT3fHfyOcf+Mk+6iR2g0Grp06ULTpk0B2LRpE8HBwZLcCSEKFUlasuhh91DampbUxRKlpcWSnKyd+M7vO/6v0f9hrbFm+/Xt9F3bl6N3jlo6tAJDURT8/f1p1aoVANu2bWPr1q2SuAghCg1JWrJANRpRE1JaVNIU4hpSW1qkpsXiFEWhX/V+zO8yHy9nL27G3CTw70BmnJiBUTVaOrwCQVEUWrVqRfv27QHYvXs3GzZswGiU748QouCTpCUL1PhHpupPU4grNS0FT/US1VncdTGdvDthUA38dPgnRm8ZzZ1YWUgwVdOmTenatSsABw4cYM2aNZK4CCEKPElasiC1awjSjx6SmpaCycHKgW9afMPnzT7HTmfHvzf/pc+aPuy4vsPSoRUYDRs2pFevXiiKwtGjR1m+fDl6vd7SYQkhxGNJ0pIFxrjk1hTF1hZFY/4tMw15lqSlwFEUhZ4+PVncdTHVS1QnLCGMsUFj+Wb/NyQaEjO/QTFQr149+vbti0aj4dSpUyxZsoSkpCRLhyWEEBmSpCUL1LiU2XBt0ycmD1tarPM1JpF13i7ezO88n8E1BgMw78w8Bm0YxJWIKxaOrGCoWbMmAwYMQKfTcf78eebPn09CQkLmFwohRD6TpCULUudoUdIU4QKoxuRPpYokLQWatdaa955/j1/b/oqbjRtnH5yl37p+rLywUkbPAFWrVmXw4MFYW1tz9epVZs+eTUxMjKXDEkIIM5K0ZEFqTYvGNn3SYlSTawA0ii5fYxI507J8S5Z1X0aj0o2I08fx0Z6PeG/ne0QlRlk6NIurVKkSw4YNw87Ojv/++49Zs2YRGSkraQshCg5JWrIgdfSQYpt+hJCqprS0SNJSaHjaezLVfyqvN3gdraJl45WN9F3bl+N3ZV2ecuXKMXz4cJycnLh79y4zZ87kwYMHlg5LCCEASVqyRE1MWRTROoOkxZjc0qJorPI1JvF0tBoto+qMYnan2aYVo4dtHCZzugCenp6MGDECNzc3wsPDmTlzJrdv37Z0WEIIIUlLVqQmLYp1+roVNaV7SFG0+RqTyB31StZjabeldKzUEb2qN83pcjf2rqVDsyg3NzdGjBiBp6cn0dHR/PXXX4SGhlo6LCFEMSdJSxYYn5C0GFO6h6SmpfBysnZiUstJfNb0M9OcLr3X9C72c7o4OTkxfPhwypcvT3x8PHPmzOHy5cuWDksIUYxJ0pIFasKTWloMyccU6R4qzBRFoVfVXizquohqbtVMc7pMOjCpWM/pYmdnx9ChQ6lcuTJJSUnMnz+fM2fOWDosIUQxJUlLFmSpe0gjLS1FQWWXyszvMp9BNQYBMPf0XAZvGMzViKuWDcyCrK2tGThwIDVq1MBgMLBkyRKOHj1q6bCEEMWQJC1Z8DBpSd+akjpPi3QPFR02WhvGPz+en9v8jKuNK2cenCFgXQCrLq4qtnO66HQ6+vTpQ/369VFVlVWrVrFv3z5LhyWEKGae6p329OnThISEkJho3nzevXv3pwqqoFGTnlTTklqIK91DRU2rCq1Y1m0Z7+96n/239vPh7g/Z898ePmr8EY7WjpYOL99ptVq6d++Ora0t+/bt4++//yY+Ph4/Pz8URbF0eEKIYiBHScvly5fp1asXJ06cQFEU06fP1F9cBoMh9yIsAB4OeZbuoeKmlEMppvlPY+bJmfx69Fc2XtnIibsn+KblN9QtWdfS4eU7jUZDhw4dsLOzY9u2bWzfvp3Y2Fg6duyIRiMNt0KIvJWj3zKvv/463t7e3LlzB3t7e06dOsWOHTto2LAh27dvz+UQLc/UPWSVJmkJnoR6bTeQ0j0UPAk+cU3+UxQZv267zFcLytLF43PKOpTlevR1Bq4fwvCVEzEYi1aCnhWKouDn50enTp0A2L9/PytWrJAVooUQeS5HScvevXv57LPP8PDwQKPRoNFoaN68ORMnTmTcuHG5HaPFZTjkOXgSbPsSY8p3UPlnImz7ElCJ3HKJ6+N3EBkUAkBkUAjXx+80bYvCY0rQBX7Ych4VmL9Dg1v4eyRF1kVRjByMXEDHxQO5dP0a18+FER0Wb+lw81WjRo3o3bs3Go2GkydPsnDhQhISEti7fCHf9+vGvuWLAEzbe5cvtHDEQojCLkd9GgaDAScnJwA8PDz477//qFatGl5eXpw7dy5XAywIMhw9tO2r5GMpXWLKrRMAROr7E6lPXk04css1Eq5EkHAx3LQN4Ny2Yn6ELXLBj1vOm23vuxgHDEAfXQ3bUqtxDXVm4xcXUNCgKNBqcHVqNitrmWAtoE6dOtjZ2bF48WIuXbrEb5N/JPziKcLLViJq/RpCT58g5OQxAPYsmQ9Ak94DLBmyEKIQy1FLS+3atTl2LPkXUaNGjZg0aRK7d+/ms88+o3LlyrkaYEGgJqasL/Ro0tL6/eRjKfWHSsrM75H6gWbXpiYsqVITF1E4vOnvm8FeBX3EsygX38Lvcj+UlP9Gqgrb550tdi0uPj4+poUW97mUYtqAt1jcfSRTB73DOr3556I9SxZYKEohRFGQo6Tl//7v/zAak9+lP/vsM65cuUKLFi3YsGEDU6ZMydUAC4IMW1r83gVvP4wpLS2alGJkZ92Tfyk7+3vlTZAiT4xrW5VmPu4ZHmvh4W1KWFKpKhw8X/wWXixfvjydhwxjh299U+ujqtGwuWUPohycTec1CxhkqRCFEEVAjrqHOnToYPq7j48PZ8+e5cGDB7i5uRXJoY8ZztMSPAmuBKOWdEs+ljJ9h7NuEQnG2iQY6wHm3wsbH1fpGipkpgRdYPfF+xke23ErjJrYmf2UjRh5//g7DLLpx5i6Y9AVo1FlEXYOpoQllarREObijlNMJBXr1Kdx7/4Wik4IURRkq6XFaDTyzTff0KxZM5577jnGjx9PXFwcACVKlCiSCQtkVtOSvKmktLRE6vuTYKxP2oQFkruKpBi3cElb0/KoaA1ssktESS3GViCi0VmirMP449gfDPt7GKFRxWeRwcp2Nqb/B6kU1YhTUgIAISeOmopzhRAiJ7KVtHz55Ze8//77ODo6Uq5cOX766SfGjh2bV7EVGBnO05JS02LqHir7HJC+piUtqWkpXNLWtDT38TDb9u9ahaFfNqXnm88w9Kum/N/wV/mmxTc4WTlx/O5x+qzpw5pLa4rFTLplba1pH7wKJaXrWDEaaXn+GDr30iQ5JbdI7k4pxhVCiJzIVtv1nDlz+O233xgzZgwAW7dupUuXLsyYMaNITyz12JoWQI3/PflYj1/h4AqctywwjR6C5C6hR4txpaalcBnXtiqQ3OLylr8vr7WtypSgC/y45Txv+vuajju62Zqu6Vy5M/U96zNh5wQO3znMB7s+YOf1nfxf4//DxcbFIq8jv7xUtzre878jzMUd18gHuFSsQqSiIb5cZdTbIbTyb2/pEIUQhZiiZuMjoI2NDRcvXqRChQqmfba2tly8eJHy5cvnSYBZFRkZiYuLCxERETg7O2d+QTZc6deP+GPHKf/rLzi1bWt2bHtwPQyGaJo0DsLevlJyLEEhRG65hrO/F85tK6bbFsWDwWhg5smZ/Hb0N/SqntIOpfmq+Vc8V/o5S4eWp/YuX8ieJQtoFjCI53sF8NcvUwh9EA6An58frVq1KrJdyUKI7MvO+3e2khatVsutW7coWbKkaZ+TkxPHjx/H29s75xHngrxMWi73eoGEM2eoMH0aji1amB3btr0mRmMCTZvswM6uXK4+VxQNJ+6eYPzO8YREhaCgEFg7kFfrv4q1Nv2yEEWRqqoEBwebZstu2LAhnTt3LtKts0KIrMvO+3e2uodUVSUwMBAbGxvTvvj4eF566SUcHBxM+1asWJHNkAs4fco8LVYZrPKspvTfK/ILWGSsTsk6LO22lG8OfMOKCyv46+Rf7Lmxh69bfI2Pm4+lw8tziqLQqlUr7O3t2bBhAwcPHiQ6OprevXtjlcH/KSGEeJxsvdMOGzYMT09PXFxcTF+DBw+mbNmyZvuKGlWfvL6MotVmcFSSFpE5eyt7Pm36KZNbT8bNxo1zYefot64fc0/PxZiS+BZ1zz//PH379kWr1XL27FnmzJlDbGyspcMSQhQi2Wpp+euvv/IqjgJNTV0ITpv+2/Wwd0366EXm2lZsS72S9fho90fsvLGTSQcmEXw9mC+afUFph9KWDi/P1apVCwcHBxYuXEhoaCgzZ85k8ODBuLq6Wjo0IUQhIM0DWaAakpMWxSqjHC8laZGWFpFFHnYe/Nr2Vz5s/CF2Ojv+vfkvL6x5gY1XNlo6tHxRqVIlRowYgZOTE/fu3ePPP//k9u3blg5LCFEIZKulZcSIEVk6b+bMmTkKpsB6YvdQctIi7SwiOxRFIaBaAM+Xfp4JOydw8v5J3t3xLttDt/NB4w9wts7dYvKCplSpUowaNYp58+Zx9+5dZs6cSf/+/S1e0C+EKNiy1Twwa9Ystm3bRnh4OGFhYY/9KmpUQ3LSkrZ7yHzglaQtIvsquVRiTuc5vFzvZbSKlg1XNvDC6hf49+a/lg4tz7m4uDBixAgqVqxIQkIC8+bN49SpU5YOSwhRgGWrpeXll19m4cKFXLlyheHDhzN48GBKlCiRV7EVGKk1Lem7hx4WUEohrsgpK40Vr9R/heblmjNh5wRCokIYtXkUQ2sOZVyDcdhobTK/SSFlZ2fHkCFDWLFiBWfOnGHp0qVERUXRuHFjS4cmhCiAsvVO++uvv3Lz5k3effdd1q5dS4UKFQgICGDTpk1Fe5ry1KQlTfeQtLSI3FS3ZF2WdltKH98+AMw5PYcB6wdw7sE5C0eWt6ysrOjbty/PPZc86d7ff//Nli1bTCvJCyFEqmw3D9jY2DBgwAC2bNnC6dOnqVWrFq+88gqVKlUiOjo6L2K0uMd1D5mKcAFJWkRusLey5+MmH/Nzm58pYVuCC2EXGLB+ALNOzsJgNFg6vDyj0Wjo3Lkzbdq0AWD37t2sWrUKferIPSGE4ClHD2k0GhRFQVVVDIai+ws1K91DkrSI3NSqQitWdF9Bq/KtSDIm8f2h7xmxaQTXo65bOrQ8oygKLVu2pEePHiiKwvHjx1m4cCEJCQmWDk0IUUBkO2lJSEhg4cKF+Pv74+vry4kTJ/jll18ICQnB0dExL2K0KFVVn9A99PDvspaKyG3udu5MaTOFT5p8gr3OnsN3DtN7TW+Wn19epLtjn3nmGQYOHIiVlRWXLl1i1qxZREVFWTosIUQBkK2k5ZVXXqFMmTJ8/fXXdO3aldDQUJYuXVq01xF5tF893ZBn6R4SeUtRFHr79mZ59+U08GxArD6WT/Z+wqv/vMrd2LuWDi/PVK1alcDAQOzt7bl58yZ//vknd+8W3dcrhMiabC2YqNFoqFixIs8888wTWxYssfZQXi2YaExM5FzdegD4HjyA9pHWJIMhlu3BdQDwa3kcnc4hw3sIkRsMRgPzzszjp8M/kWRMwsXGhQ8bf0iHSh0sHVqeuX//PvPmzSMsLAxbW1v69+9PpUqVLB2WECIXZef9O1vNI0OHDqV169a4urqarTWU9qtISUoy/fVJo4dkyLPIa1qNlmG1hrGk6xJqlKhBREIE7wS/w7s73iUiIcLS4eUJd3d3Ro0aRbly5YiPj2fu3LmcOHHC0mEJISwkWy0tBVletbQYIiM5/3wjAKofP4ZibW06ptdHE7wjuRWmld8ptFrbXHuuEE+SZEhi6vGpzDgxA4NqwNPOk8+afUazcs0sHVqeSExMZMWKFZw9exaAdu3a0axZM6klE6IIyLOWluJIfXRUlE6GPIuCwUprxavPvMrcTnOp5FyJO3F3eGnrS3y+93Nik4reysnW1tYEBATQqFHyB4itW7eyfv36Ij1qUQiRniQtmVBTu4c0GpQ0xcbm3UOStIj8V6dkHZZ0W8KgGoMAWHJ+CX3W9uHInSMWjiz3aTQaOnXqRIcOyTU8Bw8eZPHixSQmJlo4MiFEfpGkJTOGzBdLTCZJi7AMO50d458fz/T20yntUJrQqFAC/w7kx0M/kmgoem/oTZo0ISAgAJ1Ox/nz5/nrr79kSLQQxYQkLZkwdQ9ZWWVwVCaXEwVH4zKNWdF9Bd2rdMeoGpl5cib91/fn7IOzlg4t19WsWZNhw4aZhkTPmDFDhkQLUQxkK2n56KOPOHToUF7FUjClTiyXwTw00j0kChonaye+bP4lk1tPfrgMwLoB/H7sd5KMSZnfoBCpUKECI0eOpESJEkRERPDnn39y9epVS4clhMhD2Uparl+/TqdOnShfvjwvv/wyGzduzJP+5K+//hpFUXjjjTdy/d7ZpRpTEpMMJ897tHtIGq1EwdG2YltWdF9Bu4rt0Kt6fjv6G4M3DOZC2AVLh5ar3N3dGTlyJOXLl5ch0UIUA9l6p505cya3bt1i4cKFODk58cYbb+Dh4UHv3r2ZM2cODx48eOqADhw4wNSpU6lbt+5T3yt3JCcmGbWkqEhLiyi43O3c+aHVD3zT4hucrZ05ff80/db1Y8aJGeiNRWchQgcHB4YNG0aNGjUwGAwsX76cXbt2FemlDoQorrLdPKDRaGjRogWTJk3i3Llz/PvvvzRq1IipU6dStmxZWrZsyXfffceNGzeyHUx0dDSDBg1i+vTpuLm5Zfv6PJE6jX9GLS3yS1EUcIqi0LlyZ1b1WIVfeT+SjEn8dPgnhm0cxuWIy5YOL9dYWVnRt29fmjRpAiQPiV63bp0MiRaiiHnqPo0aNWrw7rvvsnv3bkJDQxk2bBg7d+5k4cKF2b7X2LFj6dKlC+3atcv03ISEBCIjI82+8oL6pKRFiEKipH1Jfm7zM180+wInKyeO3ztOwNoAZp+ajcFYNN7YNRoNHTp0oFOnTgAcOnSI+fPnExcXZ+HIhBC5JVffiUuWLMnIkSNZvXo177zzTrauXbRoEYcPH2bixIlZOn/ixIlmSwdUqFAhJyFnLrU1RZO++0fZ+9sjpxkheBJ84pr8pxAFzM//XOSNGTq6uX9Ps7LNSDAk8N3B7+iweAAhkSGWDi/XNGrUiAEDBmBlZcXly5eZOXMmYWFhlg5LCJELCkTzQWhoKK+//jrz58/H1jZrU+FPmDCBiIgI01doaGjeBJfS0qKkHdIcPAll14+mTXVOD9j2JaAm/ymJiyhApgRd4Ict51GBqdseEHltGPE3X0A1WHM78QzdV/Vi/pn5RD6I5fq5MKLD4i0d8lOpVq0aI0aMwMnJibt37zJ9+vS8+x0hhMg3BWLtoVWrVtGrVy+0j0zgZjAYUBQFjUZDQkKC2bGM5NXaQ/GnT3Plhd7oSpWiavD2hwc+cSVJCzuauQPQeuc9NGbfSQU+Cc+1OIR4Gt7j15PRf3RFF4Zt2WXoHC5R/XZj/C73Q0GDokCrwdWp2axsvseamyIjI1mwYAG3bt1Cq9XSq1cvateubemwhBCPKHRrD7Vt25YTJ05w9OhR01fDhg0ZNGgQR48ezTRhyUuPHfLc+n2zthc1be9R6w/yMiwhsuVNf98M96t6N+JCRuJn8zItUxIWSO4V3T7/bKFvcXF2dmb48OFUq1YNg8HAsmXLCA4OlpFFQhRST5W03LhxI0ejhNJycnKidu3aZl8ODg64u7tb/lORmtI9lHZIs9+7KF7NH56W8mekvj/X49cSqe+XvB0UwvXxO4kMKjo1A6LwGde2Ks183DM81tzHk/HPDECT5teBaoRL1wp/l4qNjQ39+vUzjSzatm0bq1atQq83H/a9d/lCvu/XjX3LF5lt712e/UEFQoi8kXbZ4izZvXs3gwcPJiQk+Y3Yw8ODwMBAPvjgg1ztmikQUkcPpU1agifBlZ1Q3sN0PFLfj0j9YAAit1wj4UoECRfDTdsAzm0r5kfUQpiZEnSB3RfvZ3hs18V7rCzljFYxH8VvxMgbh1/hJdtR9K3WF41SIBpmcyR1ZFGJEiXYsGEDx44dIzw8nJY9X+AmGsK2/c35JfMB2L1kHqGnTxBy8hgAe1L2N+k9wGLxCyGS5ei30JgxY6hRowYHDhzg3LlzfPvtt2zdupUGDRrkSssLwPbt25k8eXKu3OupqI/pHtr2VbruoUj9QLNTUhOWVKmJixD57cct5594/Ifdl2k1uDqpeYmiwNX6e7mnvcUX/37Bi5tfJDSq8Le6PPfccwwaNAgbGxv+ToDmRy7T5+glRrtU4Xj1Z03npSYsqfYsWZDfoQohMpCjpOXSpUtMnjyZBg0a4OPjw9ChQzl48CDPPPNMgZh6Pzc9rGlJ09LS+n2URz6VqpWa4ax78i82Z3+vXI5OiKxJW9PS3Mcj3fGazcoy9Mum9HzzGYZ+1ZRvx3zA+OfHY6ezY/+t/fRe05v5Z+ZjVI0UZj4+PnQZMowdvvVRU1pQVY2GzS17EOWQcUtxs4BB+RmiEOIxcpS01KhRgzt37pjtUxSFzz77jL///jtXAiswTDUtab5Vfu+iPFJsqwbMxtm/Cjaao5DBOA0bH1fpGhIWM65tVd7y90UB3vb3Zd6oRqbtt/x9Gde2KgCObraUq+aGo5stGkXDoBqDWN5tOc+Vfo44fRxf7/+aEZtGFPp5XSIdnEwJSypVoyHMJX3dT8U69Wncu39+hSaEeIIcJS2BgYG89tpr6eY9yO3hxgXCk2bE9XsX07dQNRCp70+CsT6kndOF5K4iKcYVljSubVWufN2F11ISlNTt1ITlcSo4V2BG+xl80OgD7HR2HLp9iN5rejPn1JxCO5tuZTubdL/8FNWIvc4q3UeOkBNHTcW5QgjLylHS8sYbb3Ds2DGqVq3KwIEDmTRpEhMnTmTkyJFMmlS0JlUzdQ89ZkFERUkejq2qhkxrVqSmRRRWGkVD/+r9WdljJY3KNCLeEM+3B78l8O9ArkRcsXR42VbW1prvqlUgdTIFxWik5flj2Dg4E1+uMmqaltXdKcW4QgjLytHooZs3b3L06FGOHTvG0aNHmTVrFhcuXEBRFCZNmsTGjRupW7cudevWpWPHjrkdc/5KKcRVMpjGH5K7jVQVVFXF2d/LLDGx8XE1K8aVmhZR2JVzLMd0/+ksu7CM7w9+z9G7R+m7ti9j649laM2haDWWm1MpuwaWdadVCSeuxCXwYNvfnDq8nfgyXuidS6B3dkN74TgafRIATQMGZnI3IUR+yLUZcePj400TxKUmMydPniQ8PDw3bp+pvJoRN2bvXkKGj8DG15fKa1anO75te22MxjiaNtmGnV1FIoNCiNxyDWd/L5zbpt8Woqi4GX2TT/d+yu7/dgNQ16MunzX7jCquVSwcWc7sXb6QHWtXo69SiySDAWudFt2Fk7To3lOGOwuRh7Lz/l0gpvHPDXmVtETv3k3oyFHYVK9O5VUr0x3fHlwPgyGaJo23Ym/vnWvPFaIwUFWVVRdX8e2Bb4lKisJKY8Ur9V8hsFYgOk2OGnItLiwsjIULF3Lnzh20Wi09evSgbt26lg5LiCKr0E3jX6BlWtOSOu154R4GKkROKIpCr6q9WNljJS3LtyTJmMRPh39i4PqBnHtwztLh5YibmxsjR440Tf2/YsUKtm7ditEo/8eFsDRJWjKVUtPymKQl9VuoIr/QRPFVyqEUv7T5hS+bf4mztTNnHpyh/7r+TDk8hQRDgqXDy7bUqf+bN09eqmPXrl0sWbKEhITC91qEKEokacnMk4Y888j8LdLSIoo5RVHoXqU7q3uuxt/LH72qZ/qJ6fRd25ejd45aOrxs02g0tGvXzrQC/dmzZ5k5c2a+1ekJIdKTpCUTaqZJy8Mhz0II8LDz4IdWP/Bjqx/xsPPgSsQVhm4cysR/JxKbFGvp8LKtXr16BAYG4uDgwO3bt5k+fbpp3TUhRP6SpCUz6pO7h0xJi3QPCWGmnVc7VvVYRU+fnqioLDi7gF6re7H7xm5Lh5ZtFSpU4MUXX6RUqVLExMQwe/Zsjh49aumwhCh2JGnJzONWeTZJ2S/dQ0Kk42LjwufNPmeq/1TKOZbjv5j/eGnrS3yw6wMiEiIsHV62uLq6MmLECKpXr47BYGDVqlVs3rxZCnSFyEeStGRCzWiV5+BJ8IkrBH/7sKVlepvk/UIIAKYEXcB7/Hp+DrpA07JN8XeeROL9ZoDCmktr6LGqB1uubbF0mNliY2NDQEAALVu2BGDPnj0sWrRICnSFyCeStGQm7SrPwZNg25eACtu+QIlKXjhSVdTk/ZK4CMGUoAv8sOU8KvD9lvMMmrGPn4NCSbjTjZirL+Gqq8D9+Pu8tf0t3tz2Jtf++4/r58KIDou3dOiZ0mg0tGnTht69e6PVajl//jwzZszgwYMHlg5NiCJPkpbMpF3ledtXZoeVlMJCNbX3KM1xIYqjH7ecN9veffG+6e/GOC+unxzDmLpj0Ck6rh+IZu1np1n94xHmvL+H07v/y+9wc6ROnToMHz4cR0dH7t69y7Rp07h06ZKlwxKiSJOkJTNpa1pav292WJNy2Gg6/kE+BSZEwfWmv+8Tj7/VriavPvMqs1rMx+9yP5TU+Y5U2D7vbKFocQEoX748o0ePply5csTHxzNv3jz27NlDEZloXIgCR5KWTJh++aQmJX7vgref6biSclzVAJVbgd//8jdAIQqgcW2r0szHPcNjzX08eK1tVQDcE8uYEpZUqgqL/12J3qjP8zhzg7OzM4GBgdSvXx9VVdm8eTMrV64kKSnJ0qEJUeRI0pJdwZPgSrBp06yl5fJ2CP7WMnEJUYBMCbpg1iX0qF0X7/Fz0AUAXD3t0g3MM2JkesivDFw/kNP3T+d1qLnCysqKHj160KlTJxRF4fjx48ycOZOIiMI1QkqIgk6SluxKU7OiSWlpMaZ+J7d9mc8BCVHwpK1pSeuHlOOObra0Glyd1JIxRQNu7ePROqmceXCGAesH8N2B7wrFpHSKotCoUSOGDh2KnZ0dN2/eZNq0aVy7ds3SoQlRZEjSkgnTpHKp3URpaloU2xLJhzUZ17wIURylrWlp7uPx2OM1m5Vl6JdN6fnmMwz9simDX+jK6p6r6eTdCaNqZPbp2byw5oVCMymdt7c3o0ePNpuI7sCBA5YOS4giQZKWzKRNWvzeTSm2VaD1/6Ep3xhI6R5q/UHycSGKuXFtq/KWvy8K8La/L/NGNTJtv+Xvy7iUmpZUjm62lKvmhqObLZC8FMCklpP4te2vlHEow43oG7y09SXe2/Ee9+My7nYqSFJXiq5VqxZGo5H169ezdu1a9PrCUacjREGlqEWkzD0yMhIXFxciIiJwdnbOvfv+/Tc33ngT+4YN8Zo3N93xEyde5c7djfj6fkKF8kNy7blCiGSxSbH8cvQX5p+Zj1E14mLjwv8a/o/uVbo/YfX1gkFVVXbt2kVQUBCQvBxAv379cHR0tHBkQhQc2Xn/lpaWp6TRWAOgGhMtHIkQRZO9lT3vPvcu8zvPp5pbNSISIvi/3f/Hi1teJCSyYC9cqCgKLVq0YODAgdjY2BAaGsrUqVO5ceOGpUMTolCSpCVTyZ/kVDJukFI0VgAYjTK8UYi8VNujNgu7LuTNZ9/ERmvDvzf/5YU1LzDjxAySCvj/P19fX1588UU8PDyIiopi5syZHDt2zNJhCVHoSNKSGVNNS8aHU1tajGrB/qUpRFFgpbFiRO0RrOy+ksZlGpNgSOCnwz/Rf11/Ttw9YenwnsjDw4NRo0bh6+uLwWBg5cqVbNq0CYPBYOnQhCg0JGnJqseU/ihKakuLLJgmRH6p4FyBaf7T+LL5l7jauHI+7DyDNgzim/3fEJMUY+nwHsvW1pb+/fubFlzcu3cvc+fOJTo62sKRCVE4SNKSmUzq/LSa5NEOkrQIkb8URaF7le6s7rmarpW7oqIy78w8eqzqQVBIkKXDe6zUBRcDAgKwtrbm6tWrTJs2jevXr1s6NCEKPElaMpN2yHMaGq0dAEZDXH5FJIR4RAnbEkxsMZE/2v1Becfy3I69zRvb3uC1f17jZvRNS4f3WDVr1jTVuURGRvLXX39x8OBBWbdIiCeQpCWrHvOLRKuxAcBgLBwLvAlRVDUr14yVPVbyYp0X0Wl0bA/dTo/VPZh9anaBXceoZMmSjBo1iho1amAwGFi3bh2rV6+WdYuEeAxJWjKRbkbcNB62tEjSIoSl2epsGddgHMu6LaOBZwPi9HF8d/A7+q/rz/G7xy0dXoZsbW0JCAigXbt2KIrC0aNHmTlzJuHh4ZYOTYgCR5KWzDxu8qrgSfCJK9oL2wEw3D0Fn7gm7xdCWNT6Q0Z2BPehpdtYXGxcOBd2joHrBxOw7B0iEyMtHV46iqLQvHlzBg8ebFq3aOrUqVy6dMnSoQlRoEjSkhklg3lagielLIyoojm9BgBj+EVATd4viYsQFjMl6AI/bDmPiob1eypQLvpjksIboCgqZ2I20W5xFzZe2YiqqkSHxXP9XBjRYQWjpbRKlSqMGTOGsmXLEhcXx7x589i5c6fUuQiRQmfpAAqNR39nPLLSs9aY/KdBo5gflzWIhLCItCtM/3sxEQggKeJZbEqvIs7mLu/ueJfgLceocOR5UJM/m7QaXJ2azcpaJuhHuLq6Mnz4cDZs2MCRI0cICgrixo0b9OzZE1tbW0uHJ4RFSUtLpjKoaXlkJWeNIXm/8dGkpfUH+RGYECIDaVeYTmWIrULslddp6Nwf16SSlDvc0PRhRFVh+/yzBabFxcrKih49etCtWze0Wi1nz55l+vTp3Llzx9KhCWFRkrRkJqOSFr93wdsPAK0x+beeqaWlcisi9f24Pn4nkUHJ66JEBoWYbQsh8s64tlVp5uOe4bHmVUrzV68PmNzgNzRpfv2pRoi4U7CmLnj22WcZPnw4zs7O3L9/n+nTp3Pq1KnHnr93+UK+79eNfcsXmW3vXb4wv0IWIk9J0pJVapqalivBAGhSkhajNvlQ5PnSRG65lvz3Lde4O+OE2bYkLkLkrSlBF9h98X6Gx3ZdvMfPQReo5l0pXY29ESPTrv1KWHxYPkSZdeXLl2f06NFUqlSJpKQkli5dyubNmwmNjWdXWBT/xScv1rp3+UL2LJkPqOxeMo+ln39g2t6zZL4kLqJIkKQlMxkNeX60piVl2ZDUlpZI/UCzyxMuhpttpyYwQoi8kbamJa0ftpzH0c2WVoOro6T8BlQVlR2VF7Pi5mK6rerG8vPLMarGfIg2axwdHRkyZAhNmzYF4M/LN3h+3xn6HL1Ew72nWfDfffYsWWB2TchJ8wUZ0x4XojCSpCUTija5CUU1PvIL7JGaFlP3kDY5aXHWPfkXg7O/Vy5HKIR4VNqaluY+Hhker9msLEO/bErPN58h8KtmfBA4Fl83XyISIvhk7ycM2TiEM/fP5FvcmdFqtbRv355WL/Rhh2991JQPVEbgf+dC8e039InXNwsYlA9RCpG3JGnJjCblW/ToSqx+76YU2ypom7wNgKpRMCoKzv5VsPFxzfBWNj6uOLetmLfxClHMjWtblbf8fVGAt/19mTeqkWn7LX9fxrWtajrX0c2WctXccHSzpb5nfRZ3Xcy7z72Lvc6e43eP0399fyb+O5GoxCiLvZ60dBW8TAlLKgPg1qoDFWvXy/CainXq07h3/3yIToi8pahFZAKAyMhIXFxciIiIwNnZOdfuG7NvHyGBw7H2qUKVdevSHVdVA/9sS/7k1qL5fuJ3xjyxC8jZ30sSFyEKuDuxd/juwHdsvLoRAHdbd9557h26eHd5OEu2hfwXn0jDvad5tPNKUY2MO7kLq71BKEZDhtc1CxgsiYsokLLz/i0tLZlI7R7CkHH/tqJo0WiSp/I3GJ6csIDUtAhRGHjaezLJbxLT20+nknMl7sffZ8LOCYzcPJJL4ZadpbasrTXfVatAym8mNKi0unichAcPiPGuicHWPsPrdi+Zn39BCpFHJGnJjKmmJeNPLwA6nQMAen1MupqVtF1FUtMiROHRuExjlndfzusNXsdWa8uBWwfos6YPPxz6gdikWIvFNbCsOwea1GR5/SocbFKLH7p1wM7aCtXahliv6iS6laRC7fpm1zQNGJjxzYQoRCRpyYSSWtOif3zSotWmJC2GaJzbVjQlJs7+XpQcVcdsW7qGhChcrLXWjKozilU9V9G6Qmv0qp6/Tv5F91Xd2Xptq8Wm2C9ra00zNyfK2lpTpkwZxr35Fp7OjqDRkFDaC02N+jzfewCg0DRgEE16D7BInELkJqlpyUTciRNc7RuArkwZqm77J8Nz9h/oTlTUKerVnYGHR+tce7YQouAJDg1m4v6J3Ii+AUCzcs14//n3qehs+Q8kqqqyd+9etmzZgqqqeHh4EBAQgKenp6VDE+KxpKYlFz2saXlSS4tjyikx+RGSEMKC/Cr4sarHKsbUHYOVxordN3bTa3Uvfjv6G/F6yy4DoCgKTZs2JTAwECcnJ+7du8f06dM5fvy4ReMSIrdI0pKZjOZpSUOX2j2kj86XkIQQlmWrs+XVZ15lZY+VNC3blERjIr8f+51eq3ux4/oOS4eHl5cXY8aMoXLlyiQlJbFixQrWrVtHUlKSpUMT4qlI0pKJhzUt+seeo9VJS4sQxZGXsxd/tPuD7/2+x9Pek+vR1xkbNJbX/nmN0KhQi8bm6OjI4MGD8fNLXift4MGDzJw5k7CwgrVMgRDZIUlLZqSlRQjxBIqi0L5Se9b0XENgrUB0io7todvpuaonvxz5hTi95RZh1Gg0tG7dmsGDB2NnZ8fNmzf5448/OH36tMViEuJpSNKSiazUtOh0yYVDen1kfoQkhCiAHKwceLvh2yzvvpzGZRqTaExk6vGp9FjVw6KjjAB8fHx46aWXKF++PAkJCSxZsoQNGzZId5EodCRpyUwWWlqsrFwBSNKH50NAQoiCrLJrZab5T+PHVj9SxqEMN2Nu8ub2Nxm9ZTSXwy9bLC4XFxeGDx9Os2bNANi/fz9//vkn9+9nvCK2EAWRJC2ZUDJaeygNU9KSFJ73AQkhCjxFUWjn1Y7VPVfzUr2XsNZYs+/mPnqv6c33B78nOtEyXclarRZ/f38GDRqEvb09t27dYurUqZw4ccIi8QiRXZK0ZEarA0B9UveQlQsASUkR+RKSEKJwsNPZMbb+WFb1XEWrCq3Qq3pmnZpFt1XdWHtprcW6jKpWrcpLL71ExYoVSUxMZPny5axZs4bExESLxCNEVknSkglFm4WWFp0bAHrpHhJCZKCCUwV+bvMzv7X9DS9nL+7F3eP9Xe8T+HcgZx+ctUhMzs7ODBs2jJYtWwJw+PBhZsyYwd27dy0SjxBZIUlLZlILcXl8XYvVY1paIoNCuD5+J5FBIRluCyGKlxblW7Ci+wpeb/A6djo7Dt85TL91/fhi3xdEJOR/S61Wq6VNmzYMGTIEBwcH7ty5w7Rp01g24w++79eNfcsXAbB3+UK+79eNvcsX5nuMQjxKpvHPhCEigvONGgNQ/cRxFCurdOckJNxm1+6mgIY2rc+hKBoig0LMVnS28XEl4WK4aVvWIRKieLsVc4sfDv7AxqsbAXC1ceX1Bq/Ty6cXWo02k6tzX1RUFCtWrODKlSsA6MLvYXsrBOdnnufUrdu4RdzHKSZS1jESuU6m8c9NKTUt8Pi6Fp3ONeVvRtNcLY8mLIBZwpLRcSFE8VLaoTST/CYxs8NMfFx9CE8I59O9nzJww0CO3T2W7/E4OTkxZMgQrO/cAFVF7+rBweZd+LhhZxZ3H8nUQe9wvPqz7FmyIN9jEyKVJC2ZMNW0AOpjVnrWam3QaOyAh3UtqSs7P05mx4UQxcNzpZ9jabeljH9+PI5Wjpy+f5rBGwYzYecEbsfcztdYNBoNrdu2we7aOWI0WoJrPo+aMoJS1WjY3LIH1foNy9eYhHiUJC2ZeLQ7SE16fGX9w7qWcACc21bExsc1w3NtfFyla0gIYaLT6BhUYxBre62ll08vFBTWXV5Ht1XdmH58OgmGhHyLpUnvAVSuUoX4iPuoimJ2TNVocGvdId9iESItSVoyoeh0DyeYS3xk9sjgSfCJKwR/C4BVQvKxxMPTgOSi27RdQqkSLoZLMa4QIh0POw8+a/YZC7sspH7J+sTp45hyZAo9VvVgy7UtZkOkpwRdwHv8en4OumC2PSVlO6f2Ll9IyMljuIfdQUkz+EBRjVzdsPqp7i/E05BC3Cw4+0wD1Lg4qmzdgnX58skJy7YvH57g7ccRhyM8KGFNjXNRlK35Dtc3Nsn0vuW/bpGrcQohig5VVdl4ZSM/HPqB27HJ3UTPlX6O9557j027jSzdfpUwrZFoDTTzcWf3xYcz277l78u4tlVz9Nzv+3UDkt8Wjld/ls0te6BqNCiqkZbnj1Hj5lXatG1L8+bN0Wjkc694etl5/5akJQvONWqMMSKCyhvWY1O5cnILC+bfttO+jtwsbUuVKzFUCo0nssVxGT0khHhqsUmx/HXqL/46+RcJhgRq3GlCi0v90KBgRGWzXRInbMzr7RTgytddcvS8vcsXsmfJfNO207NNOHXzFq5R4Vg5uaJ3cQegUqVK9OrVCxcXlxy/NiFARg/lOsU6ua5FTZ0tsvX76c6xTkxuRk2w1kDrD3BuW9FUbOvs70XJUXXMtiVhEUJkhb2VPWPrj2VNzzV09uxBi0sBaEiuNdGg0D7OCsc0U0i95e+b4+c16T2ApgGDAIVmAYMZ/e4H9G/eFOeocNo0a0LPnj2xsrLi6tWr/PHHH5w5c+YpXp0Q2SMtLVlwsU1bkv77j0pLFmNXt27yztnd4Uqw6ZzQsrac93HEM9aZOl2P5OrzhRAC4Pq5MFb/mP73yxKPi1zTlwOguY8H80Y1ytM47t+/z7Jly7h58yYADRs2pEOHDlhlMI+VEJmRlpZcplhbA4+0tARPMktY4JGWlqT7puJcIYTITa6edqQZ0IMRIwkV52BXfhaK1T12XbxnKs7NK+7u7owcOZKmTZsCcPDgQaZNm8bt2/k7RFsUP5K0ZEG6pGXbV+nOsUlJWhKtNeZFukIIkUsc3WxpNbg6xpSaOiMq/5Q6RLR1FDqnszhU+REbzw38EHQ8z2PR6XS0b9+ewYMH4+DgwN27d5k2bRr79++32EKQougrMEnLxIkTee6553BycsLT05OePXty7tw5S4cFPExajI+raancyqymRW01IT/DE0IUIzWblUXtUpZFDglMdY7nSEJtqus/RR/ti6IYsHbfQckaP7Ls/DIMxscv9JpbfHx8ePnll6latSoGg4ENGzawaNEiYmJi8vzZovgpMElLcHAwY8eOZd++fWzZsoWkpCTat29fIP7hp2tp8XsXWn8AKND6/2DoamwavwOAUatgaP6KhSIVQhQH47rVoG9nH2I08La/L8tG9eTl6l8TFxKIi64sccYIPt37KQHrAth3c1+ex+Po6MjAgQPp2LEjWq2Wc+fO8ccff5jWMRIitxTYQty7d+/i6elJcHCwaen0J8nLQtxrgcOJ3bePst99h0vXxw8j3B5cD4MhmsaNtuDgUDlXYxBCiKxIMiSx6Nwifj/2O1GJUQD4lffj7YZv4+3inefPv3nzJsuWLeP+/eR5Y5o3b06rVq3Q6XSZXCmKqyJRiBsRkbxMe4kSJTI8npCQQGRkpNlXXkk35PkxbGw8U2K7lWexCCHEk1hprRhScwgbem1gYPWBaBUtwdeDeWH1C0z8dyLh8eF5+vwyZcowZswYGjRoAMCuXbv4888/uXfvXp4+VxQPBTJpMRqNvPHGGzRr1ozatWtneM7EiRNxcXExfVWoUCHP4knXPfQYtrbJQw7j46/nWSxCCJEVrrauTGg0gRU9VtCqfCv0qp4FZxfQZWUX5p6eS5IhKfOb5JC1tTXdu3cnICAAW1tbbt68yR9//MGBAwekSFc8lQKZtIwdO5aTJ0+yaNGix54zYcIEIiIiTF+hoaF5Fo8mi0mLnV3yhHFxcXkXixBCZEdll8r83PZnpvlPw9fNl8jESCYdmESvNb34J+SfPE0iatasySuvvIK3tzd6vZ7169ezcOFCoqOj8+yZomgrcEnLq6++yrp169i2bRvly5d/7Hk2NjY4OzubfeUVxSolaXnCKs8AdrbJ8cZJS4sQooBpUrYJS7ou4ZMmn+Bu6861yGu8vu11Rm0exdkHZ/Psuc7OzgwZMoQOHTqg1Wo5f/48v//+O+fPn8+zZ4qiq8AkLaqq8uqrr7Jy5Ur++ecfvL3zvmAsq7LcPWSX3EX1aEtLZFAI18fvNK3qnHZbCCHyi1ajpbdvb9a/sJ5RdUZhrbFm/639BKwN4KPdH3E39m6ePFej0dCkSRNefPFFPD09iYmJYcGCBaxfv55dS+fzfb9u7Fue3LK+d/lCvu/Xjb3LF+ZJLKJwKzDl3GPHjmXBggWsXr0aJycnbt1KLmZ1cXHBzs7OorGlm6flMexSkpb4+OSkJTIoxLRoYuSWayRciTAtmpi6X9YgEkLkNwcrB15v8Dp9ffsy+dBkNl7dyMqLK/n76t+MrD2SYbWGYauzzfXnli5dmhdffJGgoCD27dvHgQMH0CTEYWtrx9/r1xB04xZJR/fjhGpatLFJ7wG5HocovApMS8vvv/9OREQErVq1okyZMqavxYsXWzq0LLe02NkmJy2JifcwGOLMVnkGzFZ5BtIdF0KI/FTWsSyT/CYxt9Nc6nrUJU4fxy9Hf6Hbqm6svbQWo2rM/CbZZGVlRceOHRkyZAhKUiJGGzsONe7A1MHv8HOdVkwd9A7Hqz8LwJ4lC3L9+aJwKzBJi6qqGX4FBgZaOjQUm5SkJeHJSYuVlQs6nROQ3EWUuqrz42R2XAgh8kN9z/rM6zyPb1p8QxmHMtyKucX7u96n/7r+7L+5P0+eWaVKFVrVrUl8YgI7qjVAVZLfjlSNhs0texDl4EyzgEF58mxReBWYpKUg09jZA2CMi830XDvb5O6e+PjrOLetiI2Pa4bn2fi4SteQEKLAUBSFzpU7s6bnGsY9Mw4HKwfOPDjDyM0jGRs0lkvhl3L9mX79BuNYqhRqmlUgVY0G6/qNaNy7f64/UxRukrRkgcY+JWmJfSRpCZ4En7g+XNE5Zds2OvmcuLgQIoNC0nUJpUq4GC7FuEKIAsdWZ8uLdV9kfa/19K/WH52iY8f1Hbyw5gX6LnuLyh8sMq0iPSXoAt7j1zMlh6tK712+EP2x/ShG824oRTUSefM/diye/9SvRxQtkrRkgcY+uRDYlLQET0pZyVmFbV/A7O6mbfuQUwDExFzMtGZFalqEEAWVu507HzT+gJU9VtK2YluMqpGzMVuw9/mWn4/8xoAZwUzbdJ7ySRqmbTqfo8Rlz5IFOMVE0n7HalPiohiNtDx3FFsbO7YdP8XFixdz+6WJQkySlixIbWlRY1KSlm1fmZ9wJdj0V8cYPQDR0WfT1ayk7SqSmhYhREFXyaUSk1tPJvbqSxjiKqBoErEpuZWku8sYE2lL/xgbxkTasmVd9ruPmgYMBKDu2UOMmf8d404EM2b+dzy7bxNKQjyqlTXz5s1j/fr1JGYyEEIUD5K0ZIGpeyguLnlH6/cfe65jTPJS8NExF3BqU8GUmDj7e1FyVB2zbalpEUIUFq83b0/s1VeIuz4Qu2gv/K71QENyLYoGhQ7x1kSHxWfrnk16D6BpwCBAoWOX7rw/7nU6dOmGNj6Wtg3q0KhRIwAOHDjAH3/8kaczn4vCocCu8pxdebnKc+yBA1wbMhRrb2+qbNyQvHN2d7MWllTGyn5sr3AeVU2iaZNg7OweP6uvEEIUJoNm7GP3xftU0Kv0j7ZPd7z+SDeaPfdMrj7z0qVLrF69msjISBRFoXnz5vj5+cmq0UVIkVjluSBR7NIU4gZPyjBhAdBcDsZBTf6mR8ecy5f4hBAir00JusDui/cBCNMoGDH/vGvEyNtHXmPCzgncjL6Za8+tUqUKL7/8MnXr1kVVVXbu3MmMGTO4c+dOrj1DFB6StGRButFDaWta0nC8/R+QXNcihBBFwY9bHq4VFK2BzXZJpsTFiMq2MnuIsYlg3eV1dF3ZlR8P/UhUYlSuPNvOzo4XXniBvn37Ymdnx61bt5g6dSp79uzBaMz9CfBEwSVJSxZoHB4mLaqqpq9pqdzKbNOxdEsAoqOlpUUIUTS86e9rtu1Sy42pzvEsckhgqnM8rZoPYFGXRTxX+jkSjYnMPDmTzis6M/f0XBINuVNEW6tWLV555RWqVq2KwWBg8+bNzJ49m7CwsFy5vyj4JGnJAo2DY/JfDAbUuDjwexdafwAo0Pr/YOjqR7Y/wKHuaABiYmQVUyFE0TCubVXe8vdFAd7292XeqEaM7uDLdSsjozv4Mq5tVWp51OLP9n/yS5tfqOxSmfCEcCYdmET3Vd1zbVkAJycnBg4cSLdu3bCysuLatWv8/vvvHD58mCJSoimeQApxs0BVVc7WrgMGAz7bt2FVuvQTz4+Pv8nuPc1RFC1+LY+j1eb+wmNCCFGQ6Y16Vl9czW9Hf+NOXHL9STW3arzx7Bs0K9sMJc0suDnx4MEDVq5caRpV5OvrS7du3XBycnrqe4v8I4W4uUxRFLQp30hDRGSm59vYlMbGuhSqaiAi8khehyeEEAWOTqOjt29v1r2wjtcbvI6TlRPnws7x8taXGbV5FCfvnXzqZ5QoUYLhw4fTrl07tFot58+f59dff+X48ePS6lJESdKSRalJizEyItNzFUXB1e15AMLD8maxMSGEKAzsdHaMqjOKDS9sYFjNYVhprNh/az8D1g/g7e1vcy3y6WYG12g0NG/enNGjR1OmTBni4+NZsWIFixcvJjo6OpdehSgoJGnJIo2LCwCGyMxbWgBcXZOTlrDwf/MsJiGEKCxcbV1557l3WNdrHd2rdEdBYfO1zfRY1YMv9n3Bvbh7T3X/UqVKMWrUKFq3bo1Go+Hs2bP89ttvnDp1KpdegSgIJGnJoux0DwG4uSbP5BgZeRSjMcG0PzIohOvjd5oWS0y7LYQQRVlZx7J82fxLlnVfRsvyLTGoBhafW0znFZ355cgvRCfmvHVEq9Xi5+fHiy++SKlSpYiNjWXp0qUsXbqU2NhY9i5fyPf9urFv+SIA0/be5Qtz6+WJPCaFuFl04+13iFy/Hs/x7+EeGJjp+aqqsnNXI5KS7vNsg8W4ujYkMijEbJFEGx9Xs1WgZWp/IURxc/DWQX48/CPH7x4HwM3GjTH1xtDXty/WWusc31ev17Njxw527tyJqqpY67TE3L5BtM4Kt4j71PL2JuTkMdP5TQMG0aT3gKd+PSL7pBA3D2hdUltaMq9pgeS6Frc0XURpV3V+NGHJ6LgQQhR1DUs3ZF6nefzY6kcqOVciLCGMr/d/TfdV3Vl3eV2Oh0nrdDratGnDqFGjKFmyJMc8yjOz1xgWdx/J1EHvsE5vvgzAniULcuPliDwmSUsWaVILcbPYPQSkK8bNbFVnWfVZCFEcKYpCO692rOyxko+afERJu5LciL7BhJ0T6LeuH7tv7M7xaKBy5crRbdhwdvjWR00ZZq1qNGxu2YMoh4ef6psFDMqV1yLyliQtWaR1zl4hLjysa4mIPIzRmIRz24rY+LhmeK6Nj6t0DQkhijWdRkdf376s67WOcc+Mw9HKkbMPzvLS1pcYtXkUx+4ey/wmGQhJMpgSllSqRsNt75qoGi0V69Snce/+ufESRB6TpCWLtCXcADA8uJ/laxwcqqLTuWIwxBIZeYzIoJB0XUKpEi6GSzGuEEIA9lb2vFj3RTa8sIEhNYeYhkkP3jCY1/55jfNh2ZttvLKdDUqalhpFNeKg1RFTuRaXL182FeeKgk2Sliyy8vQEICkbK4sqigYP91YA3Lq9NtOaFalpEUKIh9xs3Xj3uXdZ12sdvXx6oVE0bA/dTp81fXhvx3uERGbtg15ZW2vaB69CSVlcUTEaabfnb5yiw1GtrImr6MvWnbuIi4vLw1cjcoMkLVmkK1kSAP3d7M0lULp0TwBu316Ho39Zs2Npu4qkpkUIIdIr61iWz5p9xsoeK2nv1R4VlQ1XNtBjVQ8+2/sZt2NuZ3qPl+pWZ8z87+i35k/GzP+O7sTicPk0Vvdvg6qid/Xg119/5dw5Wei2IJOkJYt0KS0txogIjPHxWb6uRImmWFt7oteHk1T3sikxcfb3ouSoOmbbUtMihBCPV9mlMt+3+p7FXRfTvFxz9KqepeeX0mVlF74/+D1h8Y9f7blJ7wF06NKNiv9dpWOX7vT9vy9o1ncAtneu83yViri7uxMdHc3ChQtZvnw5sbGx+fjKRFbJPC1ZpKoq5+rVR01MpMrWLViXL5/lay9c+IqQ0D8pWbIjdev8muuxCSFEcXTo9iGmHJ7C4TuHAXCwcmBYzWEMrTUUByuHbN0rKSmJbdu2sXfvXlRVxcHBga5du1KjRo28CF08QuZpyQOKophaW/TZqGsBKF26FwD37v1DUlLW5nkRQgjxZM+WepZZHWfxW9vfqF6iOjFJMfx27Dc6Le/E7FOzSTAkZH6TFFZWVrRv356RI0fi4eFBTEwMixcvZunSpcTExOThqxDZIUlLNuQ0aXFyqoGjQzVUNZE7dzbmRWhCCFEsKYpCi/ItWNx1Md/6fWuaoO67g9/RZUUXlp1fRpIxKcv3K1++PGPGjKF58+YoisKpU6f49ddfOXnypKwcXQBI0pINpmLcO3ezfW1qQe6tW6tyMSIhhBAAGkVDx0odWdljJZ81/YzSDqW5HXubT/d+Sq/Vvdh4ZWOWZ9e1srKiXbt2vPjii3h6ehIbG8uyZctYsmSJrBxtYZK0ZENqS0vS7VvpDwZPgk9cIfjbNNuTAChVujugEB5xgLi46/kSrxBCFDc6jY5eVXuxrtc63nvuPUrYluBa5DXe3fEufdf2JTg0GFVVmRJ0Ae/x6/k56AKAaXtKyjZA2bJlGT16NH5+fmg0Gs6cOcOvv/7K8ePHpdXFQqQQNxsezJnD7a8m4uTvT/mfpzw8EDwJtn35cNvbD64EP9xu/QH4vcuRI0N5ELYb70rjqFz59TyJUQghxEMxSTHMOz2PWadmEZ2U3EpSyroaVy74YYitjKMRWpR2ZeetcKJTPsa/5e/LuLZVze5z8+ZNVq9eza1byR9afX196dKlCy4uLvn6eoqi7Lx/S9KSDdHBwYSOeQkbX18qr1n98MAnrsCTvo0KfBLOrdtrOXXqDXQ6J5o22YaVlVuexCmEEMJcREIEf578k4VnFhJvSJ62ompID1rfaI0GBSMqm+2SOGFjQAGufN0l3T0MBgO7du0iODgYo9GItbU1/v7+PPvss2g00nGRUzJ6KI9YeyXPqZIYEoJqfKRvtPX7T76w9QcAlPLsjKNjdfT6KK5e/T2vwhRCCJGGi40Lbz37FutfWE9Nh47Yx7vR+kYrNCSvSaRBoX2cFY7G5JaWjGi1Wvz8/HjppZcoX748iYmJrF+/ntmzZ3PvXvYmHhU5I0lLNliVKwdaLWp8vPkIIr93k7uEMlK5Ffj9DwBF0eJT5V0AQq/PldoWIYTIZ572nizu8y11Yt5Dk+YtUINCk9IaXkvTNZTuHp6ejBgxgo4dO2JlZcW1a9f4/fff2blzJwaDIS/DL/YkackGxcoKq/LlAEi8+sg6QcGTzGtYHnV5+8PiXKBEiZa4uTVBVRO5fPnHPIxWCCFERqYEXWDfLTCm6dY3YuS4w9f0X/4e9+Ke3HKi0Who3Lgxr7zyCpUrV8ZgMBAUFMT06dO5efNmXoZfrEnSkk2mLqKrVx/u3PbVky96pEhXURR8qrwHwK3bq4iKOpXbIQohhHiCH7ecJ1oDm+2STImLEZVt5bYRa/uAU9Eb6LyiMz8c+oHw+PAn3svNzY0hQ4bQs2dPbG1tuXXrFtOmTWPLli0kJWV9fhiRNZK0ZJONT3KzYfy5sw93pq1pqdzKfDvNcWfnOpQq1Q2AixcnPfZZkUEhXB+/k8igkAy3hRBCZN+bKTUrJ2wMTHWO55SvLVOd4zkc05TYa6PwtK5KnD6Ov07+RccVHfnlyC9EJkY+9n6KolC/fn1effVVatasiaqq7N69m99//52rj3zA3bt8Id/368a+5YvMtvcuX5inr7cokdFD2X3Ohg3ceOttbOvWxXvJ4ocHgiclt7i0/iC5hsW0/X5yzUsacXGh7N3nj6omUbfOH5Qs6W/+nKAQIrc87IKy8XEl4WK4aVsWWBRCiJybEnSBH7ec5y1/X15rW9W0/aa/L6+18WHH9R38evRXzjw4A4CTtROBtQIZVGNQpusanTlzhvXr15smomvYsCEuseEE/b2BMBd33CLuU8vbm5CTx0zXNA0YRJPeA/LuBRdgMuQ5D5OWxGvXuNShI4q1NdUOHUSxssrxvS5cnEhIyAx0Omeea7gKe3sv07Hr43dmen35r1vk+NlCCCGezKga+SfkH349+isXwy8C4GrjyojaI+hfvT92OrvHXhsXF8eWLVs4fDh5McezJcsRXP1ZVI0GxWik/Y7V1D176JErFN5evDYvX06BJUOe85BVxYponJxQExNJuHjxqe5VpfJbODs/g14fyYkTL2MwPFwK3dnf6wlXZn5cCCHE09EoGtp5tWNZt2V80+IbKjlXIjwhnB8O/UCn5Z2Yf2b+YxdltLOzo3v37gwbNgytZymCazRETZnLRdVo2NyyB1EOD9+gmwUMypfXVNhJ0pJNiqJgW6sWAHEnTz7VvTQaG+rU+QUrK3eiY85x5uwHpqmhndtWxMbHNcPrbHxcpWtICCHyiVajpXPlzqzssZLPm31OOcdy3I+/z9f7v6bLii4sObeEJEPGRbfe3t606NMfVVHM9qsaDQ9c3AGoWKc+jXv3z/PXURRI0pIDdvXqARB74MBT38vWpjR1av+Comi5fXsN16/PBpJrWh6tYXlUwsVwKcYVQoh8ptPo6OnTk7U91/JRk48oZV+K27G3+Xzf53Rb1Y2VF1aiN+rTXefr7ICSphJDUY3YOThj1FkTcuKoqThXPJkkLTng0KwpADG795jPjJtDbm7P4+MzAYALF78iLOxfsyLcjGR2XAghRN6w0lrR17cv619Yz/jnx+Nh58GN6Bt8tOcjeq7uyfrL6zEYH04yV9bWmvbBq1BS3i8Uo5GW545gp7MipkotEt082bVkvqVeTqEihbg5oCYmcq5xE9TYWLxXrsC2Ro2nv6eqcur0W9y+vQat1p7K+o/QbCltOp7R6CFITl5SRxKljjiSkUVCCJF/4vRxLD67mJknZxKWEAZAFZcqjH1mLG0rtkWjaNi7fCGb1q81jR7y8a7Mlag4DPZOAGjiomlStzb+gwLZu3whe5YsoGnAwGIxokhGD+Vx0gIQ+tLLRG/fjuc7b+M+atTjT8zGUGiDIY7jx1/iQdguFEVLJeM7WG+pkWFSAsiQaCGEKEBikmJYcGYBf536i6jEKACql6hORaUXK3Y587r7NQwH/0bTsCNT7lfitRJXSLh8iISS5UGrBdWIh0blbshlwp1L4BZxnw5duhX5xEWSlnxIWh7Mm8/tL77Atl5dvBcvzvik4Elms+Hi7Wc+3X/rD9IlLkZjImfOTODW7VUAVK78FpW8XkFJU8QlQ6KFEKJgikyMZO7pucw9PZeYpBgADHEVSLjrT6MyjTl+/gFuBg1hWiPDr/2OUWdFQqmK6J3dOFPaix2+9VEVxTQ0evann1r4FeUtGfKcD5w7dgCtlvhjx0m4fCXjk9JO7592faIMpv/XaKypWfM7vLxeAuDy5R84d/5jVNV8ES4ZEi2EEAWTs7UzY+uP5e8X/ibxnh+q0QqtXSj2FWcSf2s1YyJt6R9jw5hIW644NUKjT8LuxiWS7t8yJSyQPMJoi19P/otPtPArKjgkackhnYcHjs2bAxCxalXGJ6Wd3j/d8Q8y3J28PtH/8PX9GFC4cWM+R44OIzb2YXeQDIkWQoiCzdXWlbH1Xyfm4rsk3m+OfZw7ra+3Q0NyUqJBobpVM8pVrwNAtM4q3dBoo6Kw9eRpikinyFOTpOUpuPTqBUDE6tWoGS2M5fducpdQRiq3Sq5xeYIK5YdSp/YvaDQ2hIXt5d/9nbh6bSpGY5IMiRZCiEJgXNuqNPWuRMKdrthcfhVNmrddQ+w+bpw9AYBbxH3TCKNUimrk+Oa/mTt3Lg8ePMi3uAsqSVqegmOb1mhLlEB/+zYRa9akPyF4UvouoVSXt0Pwt5k+w9OzI42e34CbW1OMxgQuXZrEgYM9uf3vtideJ0OihRDC8qYEXWD3xfsAPMDOtKp0Kn38HtPfnWIiab9jtdnQ6DYHt+NiSOLy5cv89ttv7Nq1C4PBvFygOJFC3Kd0f+Zf3Jk0Caty5ajy90bztYg+cQWe9O1V4JPwLD1HVVVu3VrJhYtfkZQUBqqCa0g7PC71RKt3kNFDQghRAHmPX2/2LlAnQUv7OCs0KBhRuaBbitfd66bjCeXtSQzTmYZGd+jSDd9W7Vm3bh1XriTXT5YqVYru3btTrly5fH41eUMKcfOR24D+aD08SLpxg/DlK8wPZlbTgprcGpMFiqJQpswLNG60idKle4KiEu61hcst3ybcfx2Og0qYim8lYRFCiILhTX9fs22XWm5MdY5nkUMCU53j8en1NjW6d0YFDvuGs7DuGS6WC6XCf1c4bVuBg64NcXd3J7xcE3YleqPorLl9+zbTp09n48aNJCRkvPZRUSUtLbngwZw53P5qIhoXFyqvWY1VqVIPD6Yd9pyRDIY+Z+b+g11cuPAlMTHnAVAULZ4lO1Gx4kicnetmen3aiehkYjohhMgbU4Iu8OOW87zl78trbauatt/092Vc26qm8849OMevR39lW2hy97+qakiKaEBdhz6cOQ9uBg0x2gT6VAon8c5VAJydnenSpQvVqlUDME1M1yxgEI179y8UE9XJPC35nLSoSUlcHTCQ+JMncWjRggrTpprPq5KL3URmz1VVHjzYSUjIDB6E7Tbtd3VtRMWKI/Fwb4WiaNNdl5qgpJKuJSGEKDiqfDIV65Jb0DmeA6Da7ab4XQ4wdSlttkvix7E+rFu3jvDwcABq1qyJpzGB7Zs3mrqWanl7E3LymOm+TQMGFcjERZKWfE5aABIuXeJKrxdQExPxfO893IcHPjyYWWtL6//LdCRRZqKizhASOoPbt9ehqskLdllZuVOypD+eJTvg5tYYjcYakInphBCiIJsSdIEftpxHY3eNEi57GHY20GzUkapA4FdNsXbQEBwczJ49e1BVlTOe5dlRvQGqojFNTFf37KFH7qzw9uK1+f+CMiE1LRZgU6UKnu+8A8Cdb74h/NG5W55y6HNWODnVoFbN72naZDteFUej07mQlHSf//5bxNFjw9m5qxGnTr/N3bubcfAv9cR7ycR0QghhOePaVqWZjzvGOC/sQgelGyatqDB95xxijDH4+/szevRoHCp6saP6s6hK8rmqRsPmlj2IcniYBDQLGJSvryMvSNKSi9yGDMZt6BAAbn7wf0Ru2px8IBeGPmeVrW0ZfHzeo0Xzf6lffzblyg3E2toDvT6SW7dWcfzEyxzVvsCtptMILxdMgv1/qI90XcnEdEIIYVmPDpMO0xrTDZM2YmTxrTl0XN6RKYenYF/CnsY9eqebmE7VaHjg4g5AxTr1ady7f/68gDwk3UO5TDUaufnB/xGxciUoCh6vjsXjzgcoSu7XtGQ5JtVARMQR7tzdxN27m4iPv2F2XJPogF2ED3ZhVbELr4rnc364PlIcJoQQIv9kNkx6q1sINk03cubBGQCcrJzoUWMkk8MamCUuimpk8I61eFw7gy42mmYBgwtk4vL/7d17dBT12Qfw78zecyEJCQkJBDCEi9zkZvIGQsECogKip69QtDlIsVRFFHIqck8FBQ71tLSApVIVatFYrVhuRioW+wpYhRBFCMhV2pIEAiGXTfY287x/zO5kN9lcNpDdTPJ8zpkzM7/5zcyzP0L2m9nZXb6nJYShBQDI5ULxSy/hZq7yRYoRg7sjsfdR6M3uTzpMGadcYfFowbuHWlwbEc6v+TOq4o+hOuY72KIugHR1Ps1X1qFT9GBER41AVNQIRET0h8XS3e9NvYwxxm4vzz0tHpmpcSj4rlT9ksW5k/pi/g9T8enlT7H56804W3YWADD8wsPY/4OHQKJyT8vY08fQ/5ryR6rh5jWYSv6LX7zzYSgeUqM4tIQ4tHjc/GAnil98EWS3Q7QYEdevFJ1/Nh/C+MXum3PXKJ/lEqTA4uH97iESXJAH3EB5VT5qos+hJvosJFN5vX1E0YSwsN4ID0+tncL6wGLpAVHUB7V+xhhr75r7NmmZZOz/fj9eLXgVkUdvoM+VHuq7h/r26oWL12/CGRMPADDqdZjy4DQMHjzY9x2uIcahpY2EFgCwnTqFK8uXw35KuYxn6NkDXZ5+GpH33QfRZApZXQ19TkvkxB4wjhJw8+YxlFfko7y8ANXV5yDL/r9lVBAMCAu7wx1k+iDM0gtmSzeYTUkwmeL56gxjjAWBJEvYd3EfPnr790j5VsDxvuW4MtiAH9/4H9z84iLEfkNgtSu/x1NSUjBlyhR07tw5xFUrOLS0odACACRJKP/wQ1zdsAHStVIAgC46GtH/+yNEz5gBY3JyiCtsHJGEmpp/w2o9B2v1eVitZ5Vl6znIck2D+wmCHiZTAsymJJjNSTCZlbnZnKS26fURQXwkjDHWvjllJ3af340/fP0HXLFeAQB0De+KJwY+gfiSeHz+f8p3F+n1evzgBz/AqFGjoNeH9mo5h5Y2Flo8pCoryv78Fspy34WruFhpFARYhgxBxD3jEDFuHEz9+rWpy3aNIZJhsxXBWu0JMWdRU3MZNtsV2O3F6ufFNEav7wSzKRFGYxcYjbEwGDvDaKidG42dYTDEwmiMhU4XrpmxYYyxUHJKTuw8txN/+OYPuFp9FQDQLaIbZt8xGzUnanDp4iUAQJcuXTBlyhT07Bm6j7rg0NJGQ4sHuVyo+uwzlL39DqyHDvls03ftiogxmbAMHQrLkCEwpqRA0GnvJRYiCXbHNdhtV2DzTPYrsNmK1HWX62ZAxxRFozvAeIJNLIyGzjAYOkNv6AS9PhIGvTLX6zupkyiaOOwwxjoku2TH+9+9j63fbMV1m/I26p6RPfFIp0dQkl+C6upqAMDw4cMxceJEWCyWoNfIoaWNhxZvzuJiVH32T1QdPAjrkSMgm81nuxgWBvOgQbAMGQxTv/4w9uoF4x29oIvQ/ssqLpcVdrsSYhyO63A4r8PpuAGH8zocjutwOm+oc0mqbvF5BMHgFWQ8wcYTbiK9liOg04VDp7N4zcN8Jr5HhzGmRTWuGrx7+l288e0bKLOXAQD6hvfFePt4XD2rXIkJDw/Hfffdh0GDBgX1Dz0OLRoKLd5kmw3VX34J6xf/gu3ECdScPAmq9v9kresSB1PPXjDecQeMvXrBkJQIfUIC9PEJMMR3gWA0Brn61iVJ1XA4bqhBRgkznnkZXK5KOF0VcLkq3ZOyDMi3tQ5RNNULMj6TaIFOH67MdRaIohmiaIKoMyn7etZFM3Q6z7J7m9c6hyPGWGuwOq14u/BtbDu5DRWOCgDAMP0wDLw6EDXlyj2KvXv3xuTJk4N2oy6HFo2GlrpIkmA/fx62b75BzTcnYL9wHo6LlyBdv97kvrrYWOgT4mGIT1DCTGwsdNHR0MVEK3OvSYyIaJcvnxARJMmqBhhPmKkNNxXK5KyAS6qCy1UJSaqGJNW457XT7Q4/TREEg58wY4ZONEEQjRDdk9LPCNE9F8TadUE0+W1X1wX3Meq16yH4TDqIogGCoPNq4w/TZkzLKh2VeOvUW3jr1FuoclZBJBGZjkwklCRAlmTo9XqMHTsWGRkZrX6jLoeWdhJaGiJVVMDx/fdwXLwIx6VLcFy6BGdxCVzFxXBeuwY4nU0fxJteD11UlBJiIiMhhodDjIhQ5uHhECPCofMse7aFedbDIJpMEMxmCCYzRLMJQojvRL/diAiy7IAkWd2BxgpJrh9saieljyzZIMt2SLIyr123Q3a3Se42WbaDKMB/t5ASIQh6iKK+TpjRQRAMDQQddwAS9BDEum1efettqxOiRD/7eM4tGiB61SAIOkAQIUAHQRC91sXa7RDcy6L7cXn1hagu1+/r6ccBjmlXub0c205uw47CHahx1SDCGYExFWMQVhEGQLlRd+rUqQCu4L///T9ERaViwICHb2sNmg4tmzdvxq9+9SsUFxfjrrvuwsaNG5GWltbkfh0ptDSGZBlSWRlcJSVwlpTAVXIVrqslcN24AelmOaSbN32muvfQ3BYGgxpk1LnZXGfdBMFghGA0QjAYlMmz3Oy5AYJer9yorNND0ImATufbptcBoljb5pmLYpu7ukQk1QkyNiUsyTbIkmddaZNlJ2RygGSnsk7KnGQHZFK2k0+7092uTFSv3am2E7kgyy4QSe4g1aZ+RbRJSqDRQRAEr7DjJwhBBATvsFMbkJRtDe2nAwRBDU8+fdVlTzjz6gud1zZRaff0h6AGOKjLgte5Ba/2uvt61v1sV4/hntftD9F3u99jeerzV6/gU6/7X0A9prIs+Fn2+j/vfRy//UWvZfiOR73+dc+rTddrruPNb99E7plc2F12JFuTMaJsBHQu5aXqpKRC9E49Cod9KO6//6+39dyBPH+3qT+J3333XWRnZ2PLli1IT0/Hhg0bMGnSJJw5cwbx8fGhLk8TBFGEPjYW+thYmAcMaLK/bLMpAaa8HFLZTUiVFZCtVvdUDbmqymvdCtlaBcmzXGUF1dRAttlAdnvtQZ1OyE4nUFUFqRUf6y3T6XyCjLosioBeB0HnCTrKMnRivTYlHOnqt+l1EETvNuUXL0QRgij4XYbo/gXoWRZFQNT59NGJAnQ+ffSAYKztL4hNL4vuJwl12V2H3mvZq50EgAQCiTIgEEiQQYIMQJkTJBBkyJAA9zL5zN0TSSC46rTLXu2edZdXf6/9SHbPJa8+yrqstte2efqAZPe+EkBUe0xyPwaS1X0B8trmmTcd2pR9JbStPwFZaPgLS97ttSEHQJ3QI6rH8IQnT3hTWkX3cnPCk+C1ryf8+a9DaRYxThSQ2a8niqzFuFpzDpDPI9wZDrNkQUyM8jEdob7O0aautKSnp+Puu+/Gpk2bAACyLCM5ORnz58/H4sWLG92Xr7SEFhGB7HaQzaaEGJsNsme9xgay2yDb7MrcvU5OpzI5HD5z2eFQgk+ddnLU7w9JArknuFzKvE4bY/UIQr1JaKAdIEAnAKIykU7wXBxQlyECJLqXdQLI89ziXoZ7G7nzo9pX9GyH5/nNt81rm3c/EhrZJgIQqM463IHTs+zuA8/xyPt51N2vtk3djrrHUdoJ3vv72+59PPJd9/T3aiO1Nt/9yatmwKumusfwtNbt793mUwNrLltNBiZP/vNtPaYmr7Q4HA4cO3YMS5YsUdtEUcSECRNw5MiRev3tdjvsXn/dV1RUBKVO5p8gCBDMZsBsRlt63wsRAbLsE2Y84YZcEiC5QLJcv807/Hi1qe3ebS4JkJU5SS7A5ennAkmy2gaQci6ZlL/+Pcuy8te+uixLTfdpxrKyv3tZltWxaHy5meeQZIBInUgZ7HqT33Z3W0P7BOkHo965Aj1z3f5CM5ZZa7n1EfcNXvAKd57QhNuz3U9fn5DX0HZACaKBnNvn2E2d22s7AIhUux2ArBNRERYFV3gNMLnu6AVPmwktpaWlkCQJCQkJPu0JCQk4ffp0vf5r167Fiy++GKzymEYJglD78g3TDPIXZhoLR8oGP+21Ialeu/c+dfo0ejw/25psJ2rhPr7bGzpe7Xj560/qYbwGGD6N3gHOvUytsV/d/f3u572tgfMGup9PQPVto5DuV2es/G3Dre3nva3ev43fY9YfayJC2Z/eQgQqcP7ucoRSmwktgVqyZAmys7PV9YqKCiS38e/wYYw1j6C+NFOnPQS1MMaAopE98d/d7yFiyNCQ1tFmQktcXBx0Oh1KSkp82ktKStC1a9d6/U0mE0wh/JZkxhhjrKMYdu9jGHbvY6EuQ32FLOSMRiNGjBiBAwcOqG2yLOPAgQPIyMgIYWWMMcYYawvazJUWAMjOzsasWbMwcuRIpKWlYcOGDbBarZg9e3aoS2OMMcZYiLWp0DJjxgxcu3YNK1euRHFxMYYOHYq8vLx6N+cyxhhjrONpU5/Tciv4c1oYY4wx7Qnk+bvN3NPCGGOMMdYYDi2MMcYY0wQOLYwxxhjTBA4tjDHGGNMEDi2MMcYY0wQOLYwxxhjTBA4tjDHGGNMEDi2MMcYY0wQOLYwxxhjThDb1Mf63wvPBvhUVFSGuhDHGGGPN5Xnebs4H9Leb0FJZWQkASE5ODnEljDHGGAtUZWUloqKiGu3Tbr57SJZlXLlyBZGRkRAEocXHqaioQHJyMv7973/zdxi1Mh7r4OGxDi4e7+DhsQ6e1hprIkJlZSWSkpIgio3ftdJurrSIooju3bvftuN16tSJ/wMECY918PBYBxePd/DwWAdPa4x1U1dYPPhGXMYYY4xpAocWxhhjjGkCh5Y6TCYTcnJyYDKZQl1Ku8djHTw81sHF4x08PNbB0xbGut3ciMsYY4yx9o2vtDDGGGNMEzi0MMYYY0wTOLQwxhhjTBM4tDDGGGNMEzi0MMYYY0wTOmRo2bx5M3r16gWz2Yz09HR8+eWXjfZ/77330L9/f5jNZgwePBj79u0LUqXaF8hYb926FWPGjEFMTAxiYmIwYcKEJv9tWK1Af649cnNzIQgCHnroodYtsB0JdKxv3ryJefPmITExESaTCX379uXfI80U6Fhv2LAB/fr1g8ViQXJyMhYuXAibzRakarXrn//8J6ZOnYqkpCQIgoAPP/ywyX0OHjyI4cOHw2QyITU1Fdu2bWv1OkEdTG5uLhmNRnrjjTfo5MmT9LOf/Yyio6OppKTEb/9Dhw6RTqej9evX06lTp2j58uVkMBjoxIkTQa5cewId60cffZQ2b95Mx48fp8LCQnr88ccpKiqK/vOf/wS5cu0JdKw9Ll68SN26daMxY8bQtGnTglOsxgU61na7nUaOHEkPPPAAff7553Tx4kU6ePAgFRQUBLly7Ql0rHfs2EEmk4l27NhBFy9epI8//pgSExNp4cKFQa5ce/bt20fLli2jDz74gADQzp07G+1/4cIFCgsLo+zsbDp16hRt3LiRdDod5eXltWqdHS60pKWl0bx589R1SZIoKSmJ1q5d67f/9OnTafLkyT5t6enp9POf/7xV62wPAh3rulwuF0VGRtL27dtbq8R2oyVj7XK5aNSoUfTHP/6RZs2axaGlmQId69///veUkpJCDocjWCW2G4GO9bx58+iHP/yhT1t2djaNHj26Vetsb5oTWhYtWkQDBw70aZsxYwZNmjSpFSsj6lAvDzkcDhw7dgwTJkxQ20RRxIQJE3DkyBG/+xw5csSnPwBMmjSpwf5M0ZKxrqu6uhpOpxOdO3durTLbhZaO9apVqxAfH485c+YEo8x2oSVjvWvXLmRkZGDevHlISEjAoEGDsGbNGkiSFKyyNaklYz1q1CgcO3ZMfQnpwoUL2LdvHx544IGg1NyRhOq5sd18y3NzlJaWQpIkJCQk+LQnJCTg9OnTfvcpLi7227+4uLjV6mwPWjLWdb3wwgtISkqq9x+D+WrJWH/++ed4/fXXUVBQEIQK24+WjPWFCxfw6aef4rHHHsO+fftw7tw5PP3003A6ncjJyQlG2ZrUkrF+9NFHUVpaiszMTBARXC4XnnzySSxdujQYJXcoDT03VlRUoKamBhaLpVXO26GutDDtWLduHXJzc7Fz506YzeZQl9OuVFZWIisrC1u3bkVcXFyoy2n3ZFlGfHw8XnvtNYwYMQIzZszAsmXLsGXLllCX1u4cPHgQa9aswauvvor8/Hx88MEH2Lt3L1avXh3q0tht0qGutMTFxUGn06GkpMSnvaSkBF27dvW7T9euXQPqzxQtGWuPV155BevWrcMnn3yCIUOGtGaZ7UKgY33+/HlcunQJU6dOVdtkWQYA6PV6nDlzBr17927dojWqJT/XiYmJMBgM0Ol0atudd96J4uJiOBwOGI3GVq1Zq1oy1itWrEBWVhaeeOIJAMDgwYNhtVoxd+5cLFu2DKLIf6ffLg09N3bq1KnVrrIAHexKi9FoxIgRI3DgwAG1TZZlHDhwABkZGX73ycjI8OkPAH//+98b7M8ULRlrAFi/fj1Wr16NvLw8jBw5Mhilal6gY92/f3+cOHECBQUF6vTggw/innvuQUFBAZKTk4NZvqa05Od69OjROHfunBoMAeC7775DYmIiB5ZGtGSsq6ur6wUTT1gk/m7g2ypkz42teptvG5Sbm0smk4m2bdtGp06dorlz51J0dDQVFxcTEVFWVhYtXrxY7X/o0CHS6/X0yiuvUGFhIeXk5PBbnpsp0LFet24dGY1Gev/996moqEidKisrQ/UQNCPQsa6L3z3UfIGO9eXLlykyMpKeeeYZOnPmDO3Zs4fi4+PppZdeCtVD0IxAxzonJ4ciIyPpnXfeoQsXLtD+/fupd+/eNH369FA9BM2orKyk48eP0/HjxwkA/frXv6bjx4/T999/T0REixcvpqysLLW/5y3Pzz//PBUWFtLmzZv5Lc+tZePGjdSjRw8yGo2UlpZGX3zxhbpt7NixNGvWLJ/+f/nLX6hv375kNBpp4MCBtHfv3iBXrF2BjHXPnj0JQL0pJycn+IVrUKA/1944tAQm0LE+fPgwpaenk8lkopSUFHr55ZfJ5XIFuWptCmSsnU4n/fKXv6TevXuT2Wym5ORkevrpp6msrCz4hWvMP/7xD7+/fz3jO2vWLBo7dmy9fYYOHUpGo5FSUlLozTffbPU6BSK+ZsYYY4yxtq9D3dPCGGOMMe3i0MIYY4wxTeDQwhhjjDFN4NDCGGOMMU3g0MIYY4wxTeDQwhhjjDFN4NDCGGOMMU3g0MIYY4wxTeDQwhhjjDFN4NDCGNOccePGYcGCBbd0DCLC3Llz0blzZwiCgIKCgttSG2Os9XBoYYzdktmzZ2P58uWhLiNgeXl52LZtG/bs2YOioiIMGjQo1CUxxpqgD3UBjDHtkiQJe/bswd69e0NdSsDOnz+PxMREjBo1qsE+DocDRqMxiFUxxhrDV1oY60DeeecdWCwWFBUVqW2zZ8/GkCFDUF5eHvDxDh8+DIPBgLvvvtvv9nHjxmH+/PlYsGABYmJikJCQgK1bt8JqtWL27NmIjIxEamoqPvroI5/97HY7nn32WcTHx8NsNiMzMxNfffVVg3XIsoy1a9fijjvugMViwV133YX333+/wf6PP/445s+fj8uXL0MQBPTq1Uut95lnnsGCBQsQFxeHSZMmAVCuymRmZiI6OhqxsbGYMmUKzp8/73P+9evXIzU1FSaTCT169MDLL7/c3GFkjDUThxbGOpAf//jH6Nu3L9asWQMAyMnJwSeffIKPPvoIUVFRAR9v165dmDp1KgRBaLDP9u3bERcXhy+//BLz58/HU089hUceeQSjRo1Cfn4+7r33XmRlZaG6ulrdZ9GiRfjrX/+K7du3Iz8/H6mpqZg0aRJu3Ljh9xxr167Fn/70J2zZsgUnT57EwoUL8ZOf/ASfffaZ3/6//e1vsWrVKnTv3h1FRUU+gWj79u0wGo04dOgQtmzZAgCwWq3Izs7G0aNHceDAAYiiiIcffhiyLAMAlixZgnXr1mHFihU4deoU3n77bSQkJAQ8noyxJhBjrEPZvXs3mUwmeumllygmJoa+/fZbddtDDz1E0dHR9KMf/ahZx+rTpw/t2bOnwe1jx46lzMxMdd3lclF4eDhlZWWpbUVFRQSAjhw5QkREVVVVZDAYaMeOHWofh8NBSUlJtH79evW4zz33HBER2Ww2CgsLo8OHD/uce86cOTRz5swGa/vNb35DPXv2rFfvsGHDGn/QRHTt2jUCQCdOnKCKigoymUy0devWJvdjjN0avqeFsQ5mypQpGDBgAFatWoX9+/dj4MCB6rbnnnsOP/3pT7F9+/Ymj1NYWIgrV65g/PjxjfYbMmSIuqzT6RAbG4vBgwerbZ4rElevXgWg3GvidDoxevRotY/BYEBaWhoKCwvrHf/cuXOorq7GxIkTfdodDgeGDRvW5OOoa8SIEfXazp49i5UrV+Jf//oXSktL1Sssly9fRnV1Nex2e5PjwBi7dRxaGOtg8vLycPr0aUiSVO8ljHHjxuHgwYPNOs6uXbswceJEmM3mRvsZDAafdUEQfNo8Ly15gkCgqqqqAAB79+5Ft27dfLaZTKaAjxceHl6vberUqejZsye2bt2KpKQkyLKMQYMGweFwwGKxtKhuxljg+J4WxjqQ/Px8TJ8+Ha+//jrGjx+PFStWtPhYf/vb3zBt2rTbWJ2id+/e6j0lHk6nE1999RUGDBhQr/+AAQNgMplw+fJlpKam+kzJycm3XM/169dx5swZLF++HOPHj8edd96JsrIydXufPn1gsVhw4MCBWz4XY6xxfKWFsQ7i0qVLmDx5MpYuXYqZM2ciJSUFGRkZyM/Px/DhwwM61tWrV3H06FHs2rXrttcZHh6Op556Cs8//zw6d+6MHj16YP369aiursacOXPq9Y+MjMQvfvELLFy4ELIsIzMzE+Xl5Th06BA6deqEWbNm3VI9MTExiI2NxWuvvYbExERcvnwZixcvVrebzWa88MILWLRoEYxGI0aPHo1r167h5MmTar2bNm3Czp07Odgwdos4tDDWAdy4cQP33Xcfpk2bpj7hpqen4/7778fSpUuRl5cX0PF2796NtLQ0xMXFtUa5WLduHWRZRlZWFiorKzFy5Eh8/PHHiImJ8dt/9erV6NKlC9auXYsLFy4gOjoaw4cPx9KlS2+5FlEUkZubi2effRaDBg1Cv3798Lvf/Q7jxo1T+6xYsQJ6vR4rV67ElStXkJiYiCeffFLdXlpa6vMWacZYywhERKEugjHWdhw8eBCbNm1q9HNOHnzwQWRmZmLRokVBrIwx1tHxlRbGmGrChAn4+uuvYbVa0b17d7z33nvIyMio1y8zMxMzZ84MQYWMsY6Mr7QwxhhjTBP43UOMMcYY0wQOLYwxxhjTBA4tjDHGGNMEDi2MMcYY0wQOLYwxxhjTBA4tjDHGGNMEDi2MMcYY0wQOLYwxxhjTBA4tjDHGGNMEDi2MMcYY0wQOLYwxxhjThP8HDdBSoRww0bYAAAAASUVORK5CYII=", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "np.float64(0.47041664919729653)" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# The final parameter values, will be overwritten if \n", "# optimization call is uncommented\n", "x = [1.01778992, 1.17318854]\n", "\n", "# Here is the code used to do the optimization, uncomment to run it\n", "# Note: it is commented out because it takes too long to run on doc builder\n", "#\n", "# res = scipy.optimize.differential_evolution(\n", "# cost_function, \n", "# bounds=((0.9, 1.5), (0.75, 1.5)), \n", "# disp=True, \n", "# polish=False\n", "# )\n", "# print(res)\n", "# x = res.x\n", "\n", "cost_function(x, plot=True)" ] } ], "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 }