{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "06b263a1-b920-4b97-9bdf-238a87893e92",
   "metadata": {},
   "source": [
    "# Surface tension calculations using DFT for fluids experiencing quantum effects \n",
    "Classical density functional theory for interfacial properties of hydrogen, helium, deuterium, neon and their mixtures ([doi:10.1063/5.0137226](https://doi.org/10.1063/5.0137226))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "7f9082f3-e3e2-421b-b09a-a39160a83183",
   "metadata": {},
   "outputs": [],
   "source": [
    "from feos import *\n",
    "from feos.si import *\n",
    "from feos.dft import *\n",
    "from feos.saftvrqmie import *\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import seaborn as sns\n",
    "import json\n",
    "\n",
    "sns.set_context('talk')\n",
    "sns.set_palette('Dark2')\n",
    "sns.set_style('ticks')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0c514e92-9e5b-45a5-98b1-81f83eb4c104",
   "metadata": {},
   "source": [
    "# Pure fluid surface tension correlations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "b334936e-587b-44d0-a92e-6dacd5f08e73",
   "metadata": {},
   "outputs": [],
   "source": [
    "def surftens_mulero2012(fluid, tr):\n",
    "    \"\"\"\n",
    "    Calculate pure fluid surface tension using Mulero 2012 correlation (doi:10.1063/1.4768782)\n",
    "    Args:\n",
    "        fluid (str): Component name\n",
    "        tr (np.ndarray): Reduced temperature\n",
    "    Returns:\n",
    "        sigma (np.ndarray): Surface tension (mN/m)\n",
    "    \"\"\"\n",
    "    ff = open(\"mulero_2012_parameters.json\", \"r\")\n",
    "    complist = json.load(ff)\n",
    "    ff.close()\n",
    "    sigma = np.zeros_like(tr)\n",
    "    for i in range(len(complist[fluid][\"sigma\"])):\n",
    "        sigma[:] += complist[fluid][\"sigma\"][i] * \\\n",
    "            (1-tr[:])**complist[fluid][\"n\"][i]\n",
    "    return sigma * NEWTON / METER / (MILLI * NEWTON/  METER)    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "58b2d141-7ab7-4938-b888-bc950d942848",
   "metadata": {},
   "outputs": [],
   "source": [
    "def surftens(param, tr):\n",
    "    \"\"\"\n",
    "    Calculate pure fluid surface tension A(1-tr)**(B+C*tr+D*tr**2)\n",
    "    Args:\n",
    "        param (np.ndarray): \n",
    "        tr (np.ndarray): Reduced temperature\n",
    "    Returns:\n",
    "        sigma (np.ndarray): Surface tension (mN/m)\n",
    "    \"\"\"\n",
    "    A = param[0]\n",
    "    exponent = np.zeros_like(tr)\n",
    "    exponent[:] = param[1]\n",
    "    for i in range(2,len(param)):\n",
    "        exponent += param[i]*tr**(i-1)\n",
    "    sigma = A*(1-tr)**exponent\n",
    "    return sigma * NEWTON / METER / (MILLI * NEWTON/  METER)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f445dc40-3c4f-417d-a430-ed4b11d9216f",
   "metadata": {},
   "source": [
    "# Pure fluid parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "54513546-2b6a-4f44-9c2f-9f5f59fd9c1d",
   "metadata": {},
   "outputs": [],
   "source": [
    "fluid_param = {}\n",
    "fluid_param[\"hydrogen\"] = {\"Tc\": 33.145, \"Tnb\": 20.369, \"xlim\": [19.5, 34.0], \"ylim\": [-0.1, 2.2], \"param\": None} #Alternative correlation from Hammer et al. (2023): [0.00480413, 0.67152591, 0.47781899]\n",
    "fluid_param[\"para-hydrogen\"] = {\"Tc\": 32.938, \"Tnb\": 20.271, \"xlim\": [19.5, 34.0], \"ylim\": [-0.1, 2.2], \"param\": [0.00471037, 0.62668507, 0.51110266]} #Correlation parameters from Hammer et al. (2023)\n",
    "fluid_param[\"neon\"] = {\"Tc\": 44.49, \"Tnb\": 27.1, \"xlim\": [26.0, 46.0], \"ylim\": [-0.1, 5.5], \"param\": None}\n",
    "fluid_param[\"deuterium\"] = {\"Tc\": 38.34, \"Tnb\": 23.661, \"xlim\": [22.5, 39.5], \"ylim\": [-0.1, 3.1], \"param\": None}"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9ab9a113-6e49-48d4-acb5-cee8966087f2",
   "metadata": {},
   "source": [
    "# Select fluid and set up functional"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "aa462031-732f-479b-80c3-2fb66e48d6d2",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Calculating surface tension for \u001b[1mhydrogen\n"
     ]
    }
   ],
   "source": [
    "substances = [\"hydrogen\", \"para-hydrogen\", \"neon\", \"deuterium\"]\n",
    "substance = substances[0] # Setting fluid\n",
    "print(f\"Calculating surface tension for \\033[1m{substance}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "adfa3259-d508-482b-8314-290d3d1781f5",
   "metadata": {},
   "source": [
    "# Set up functional"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "4c8f318f-56e0-4a48-b21b-1ecab6bc10a3",
   "metadata": {},
   "outputs": [],
   "source": [
    "parameters = SaftVRQMieParameters.from_json([substance], \"../../parameters/saftvrqmie/hammer2023.json\")\n",
    "func = HelmholtzEnergyFunctional.saftvrqmie(parameters)\n",
    "Tc = fluid_param[substance][\"Tc\"]\n",
    "Tnb = fluid_param[substance][\"Tnb\"]\n",
    "state = State(func, temperature=Tnb * KELVIN, pressure=20.0 * BAR)\n",
    "model_tc = State.critical_point(func).temperature"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2ee15889-18da-4802-842f-5e042d673ecb",
   "metadata": {},
   "source": [
    "# Calculate and plot surface tension"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "c284f4b3-af36-4d50-8391-a7e44148b050",
   "metadata": {},
   "outputs": [],
   "source": [
    "dia = PhaseDiagram.pure(func, Tnb * KELVIN, npoints=50)\n",
    "sft_dia = SurfaceTensionDiagram(dia.states, n_grid=1024, l_grid=200 * ANGSTROM , critical_temperature=model_tc)\n",
    "surf_tens = sft_dia.surface_tension / (MILLI* NEWTON / METER)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "4e69ad65-f86a-432a-8c1b-b6a2378ff9ec",
   "metadata": {},
   "outputs": [],
   "source": [
    "tr = np.linspace(Tnb/Tc, 1.0, 100)\n",
    "if fluid_param[substance][\"param\"] is None:\n",
    "    s_corr = surftens_mulero2012(substance, tr)\n",
    "    corr_label = \"Mulero et al. (2012)\"\n",
    "else:\n",
    "    s_corr = surftens(fluid_param[substance][\"param\"], tr)\n",
    "    corr_label = \"Hammer et al. (2023)\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "0247cea8-f6cd-436f-9e98-7dc14130e5ec",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlsAAAHNCAYAAAApEr6yAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAACHKklEQVR4nOzdd1gU1xoG8HeBpUlXpDcLoqio2EWNLWqsQcVujLHHrknssWvUGBUjiXotsTeMLdHYY+8VsYCAivTel92d+weySAClLUt5f8+T517OnDnz7eDufsyc+Y5IEAQBRERERKQUaqoOgIiIiKg8Y7JFREREpERMtoiIiIiUiMkWERERkRIx2SIiIiJSIiZbRERERErEZIuIiIhIiTRUHQAB7dq1Q3R0NLS0tGBtba3qcIiIiCgf3r59i7S0NJiYmODChQt59mOyVQpER0cjNTUVqampiIuLU3U4REREVADR0dEf3c5kqxTQ0tJCamoqtLW1Ub16dVWHQ0RERPng7++P1NRUaGlpfbQfk61SwNraGnFxcahevTq8vb1VHQ4RERHlg7u7O3x8fD45BYgT5ImIiIiUiMkWERERkRIx2SIiIiJSIiZbRERERErEZIuIiIhIiZhsERERESkRky0iIiIiJWKyRURERKRETLaIiIiIlIjJFhEREZESMdkiIiIiUiImW0RERERKxGSLiIiISImYbBERUbnk6emJWrVqqTqMUmXmzJlo3759kcZYsGABvv7662KKqGR5eHhg5cqVJX5cJltERFTivL29UatWLdSqVQt37tzJsV0QBLRt2xa1atXCmDFjVBChal26dAmenp6qDiOHN2/e4NChQ9l+JyEhIdiwYQP69u2LJk2aoFmzZhg6dCiuXbuW6xjx8fGYN28emjdvjgYNGmDo0KHw8fHJ0e+vv/7CjBkz8Pnnn6NWrVoYOnRoruM9evQIixYtQrdu3dCgQQN89tlnmDx5MgICAnL0HTVqFPbs2YOIiIhCnoHCYbJFREQqo6WlhRMnTuRov3XrFkJDQ6GpqamCqFTv0qVL2LBhg6rDyOGPP/6AlZUVmjdvrmg7d+4cNm/eDDs7O0yZMgXjx49HUlISvv76axw+fDjb/nK5HKNHj8aJEycwZMgQfPfdd4iOjsbQoUMRGBiYre/evXtx7tw5mJubw9DQMM+YtmzZgn/++QctWrTAnDlz4OHhgTt37sDd3R0vXrzI1rdDhw7Q09PDnj17in4yCkCjRI9GRET0gbZt2+LUqVOYO3cuNDSyvpJOnDgBZ2dnxMbGqi64/0hJSYGOjo6qw1CZ9PR0HD9+HAMGDMjW3qxZM1y4cAEmJiaKtoEDB6JXr15Yv349+vTpo2g/deoU7t+/j3Xr1qFLly4AgK5du6Jz587w9PTEzz//rOi7cuVKmJmZQU1NDd27d88zruHDh2P16tXZEvMvvvgCPXr0wKZNm7B69WpFu5qaGjp37oyjR49i0qRJEIlEhT8hBcArW0REpDLdunVDbGwsrl69qmiTSCQ4ffo0evTokaP/zZs3UatWLdy8eTNb+9u3b1GrVi14e3t/8phHjx6Fu7s76tevj6ZNm2Lq1KkICQnJ1mfo0KHo3r07njx5gsGDB8PFxQVr1qwBAERFRWH27Nlo2bIl6tWrh549e+LIkSP5fs2XLl3CoEGD0KBBAzRs2BCjR4/Gy5cvFdtnzpyJ3bt3A4DiVuun5p6dPXsWo0ePhpubG+rWrYuOHTvi119/hUwmy3dcn3L37l3ExMSgZcuW2dpr1qyZLdECAE1NTbRt2xahoaFITExUtJ8+fRpVqlTB559/rmgzMTFB165dce7cOUgkEkW7hYUF1NQ+naY0atQoxxVQe3t71KxZE69evcrRv2XLlggODoavr+8nxy4uTLaIiEhlrKys0KBBA5w8eVLR9u+//yIhIQFffPFFsR/Py8sLP/zwA+zs7DBz5kwMGzYM169fx+DBgxEfH5+tb2xsLEaNGoXatWtj9uzZaNasGVJTUzF06FAcO3YMPXr0wPfffw99fX3MnDkTO3bs+OTx//zzT4wZMwa6urqYMWMGxo8fDz8/PwwaNAhv374FAPTv3x+tWrUCkHF1J/O/jzly5Ah0dXXx9ddfY86cOXB2dsb69euzXdUpqvv370MkEqFOnTr56h8REQEdHZ1sVwN9fX1Rp06dHElUvXr1kJKSkus8q8IQBAGRkZEwNjbOsa1u3boAgHv37hXLsfKDtxGJiMooeVoahLS0fPdX09fPcdtElpAACEK+9heJxVBTwm20Hj164Oeff0Zqaiq0tbVx/PhxNGnSBGZmZsV6nODgYHh6emLKlCkYO3asov3zzz/Hl19+iT179mRrj4iIwMKFC7PdNtuxYwf8/f2xatUq9OzZEwAwYMAADB06FGvXrkWfPn2gp6eX6/GTkpKwdOlS9OvXD4sXL1a0f/nll+jSpQt+//13LF68GA0bNoS9vT2uXr2KXr165eu1/fzzz9DW1lb8PHDgQMyfPx979+7F1KlTi2Xu26tXr2BoaJjn6/tQUFAQzpw5gy5dukBdXV3RHhERgcaNG+foX7VqVQBAeHh4sTxBeuzYMYSFhWHSpEk5tpmZmUEsFsPPz6/Ix8kvJltERGVU1KbNiPz113z3d7x1E+oGBtna/Np3gDwhIV/7G/buDcsVywsUY3507doVy5Ytw4ULF9C6dWtcvHgRc+fOLfbjnDlzBnK5HF27dkV0dLSivUqVKrCzs8PNmzezJVuamppwd3fPNsa///4LU1PTbHOIxGIxhg4dimnTpuH27dto165drse/du0a4uPj0a1bt2zHV1NTg4uLS45bowXxYaKVmJgIiUSCxo0bY//+/Xj16hWcnJwKPXam2NjYj05Uz5SSkoLJkydDW1sb06dPz7YtNTU118Qvsy2tAH885MXf3x+LFi1Cw4YN8eWXX+bax9DQEDExMUU+Vn4x2SIiIpUyMTFBixYtcOLECaSmpkImk6Fz587FfpzAwEAIgpBtvtCHPpygD2RcAflvYhAcHAw7O7sct8GqV68OAHj37t1Hjw8AX331Va7b83PFKC8vX77E2rVrcePGjWxzpAAgIZ/JdH4In7gKKpPJMHXqVPj5+WHz5s05rk5qa2tnm5eVKbNNS0urSPFFRERgzJgx0NfXx7p167JdVfuQIAglNjkeKKPJliAIuH//Ps6fP4+7d+/i1atXSExMhL6+PurUqYPevXujR48ehT6Rp0+fxq5du/Ds2TOkp6fDzs4OPXv2xLBhwyAWi4v51RARUffu3TFv3jxERkaiTZs2MPjPFbhMeX2uy+XyTx5DLpdDJBJh8+bNuX4J6+rqZvv5w6tFxSEzUVm5ciVMTU1zbM8rMfiU+Ph4DBkyBHp6epg0aRJsbW2hpaUFHx8frF69Ol/nJj+MjIxyzGv7r7lz5+LixYtYvXo1WrRokWO7qalprjWuwsPDAWTdTiyMhIQEjBo1CgkJCdi9e/dHb0PHx8fnOp9LWcpksnXjxg0MHz5c8bONjQ2srKwQHByMq1ev4urVqzh58iQ8PT0LfJ/6p59+wtatWwEAtra20NHRwcuXL7Fy5UpcuHABW7durbB1X4iodKk8ehRMvhqW7/5q+vo52mqcP1egOVvK0qlTJ/z444948OABfvnllzz7ZSZh/71aExwc/Mlj2NraQhAEWFtbw8HBoVBxWllZ4fnz55DL5dmubmU+9WZpaZnnvjY2NgCAypUr53ii778KcrHg1q1biI2NxYYNG9CkSRNFe+aE++JSrVo1HD9+HAkJCdDP5d/STz/9BG9vb8yePTvPUg1OTk64e/dujvP36NEj6OjoFPr3kpaWhrFjxyIwMBDbtm1DjRo18uwbFhaG9PR0xdXIklAmn0bMfLPMmTMH165dw9mzZ+Ht7Y2bN2/ip59+gqamJi5evIh169YVaNwzZ84okqmNGzfizJkzOHbsGI4fPw5ra2vcvn1b8egvEZGqqWlpQd3AIN//5fYFrq6vn+/9lTE5PlOlSpWwYMECTJw48aPLyVhZWUFdXR23b9/O1r53795PHuPzzz+Huro6NmzYkON2mCAI+ZrD06ZNG0REROCvv/5StEmlUuzcuRO6urrZkp3/at26NfT09PD7778jPT09x/YP53FlPsH3qStJABRJy4evSSKR5Ltw5+vXr/H69etP9mvQoAEEQcCTJ09ybNuyZQu2bt2KsWPH5nmbFAC6dOmCyMhI/PPPP4q26OhonDp1Cu3atSvUxQyZTIYpU6bgwYMHWLduHRo2bPjR/pnxf6pfcSqTV7bq16+PU6dO5XpLr3fv3ggNDcUvv/yCQ4cOYfr06fmq0wFAUa131KhR6NChg6K9evXqWLJkCYYPH47du3dj9OjROWqKEBFR0eQ1mflD+vr66NKlC3bt2gWRSAQbGxtcvHgRUVFRn9zX1tYWU6ZMwc8//4zg4GB07NgRlSpVwtu3b3H27Fl4eHjgm2+++egY/fv3x/79+zFz5kz4+PjAysoKp0+fxr179zB79uyPzrvS09PDggUL8P3338Pd3R1ffPEFTExM8O7dO1y6dAmNGjXC/PnzAQDOzs4AgCVLlsDNzQ3q6uro1q1bruM2bNgQhoaGmDlzJoYOHQqRSISjR49+cn5Vpsw7RefPn/9oP1dXVxgZGeH69evZbhGeOXMGq1atgr29PapVq4ajR49m269Vq1aoUqUKAKBz585o0KABZs2aBT8/PxgbG2Pv3r2QyWSYOHFitv1u376tSKqjo6ORnJyMjRs3AgCaNGmiSGxXrFiB8+fPo127doiNjc1x/P8+0Xnt2jVYWlrmu4RFcSiTydanJhG2adMGv/zyC2JjYxEdHa34JX9MYGAgnj17BiDjzfRfLVq0gJ2dHYKCgnDu3Dn069evcMETEVGRzJ07F1KpFPv27YOmpia6dOmC77///qNVxjONHj0a9vb22L59O359/ySnubk5WrVqla8FmrW1tbFz506sXr0aR44cQWJiIhwcHLB8+fIcTy7mpkePHqhatSo2bdqE//3vf5BIJDAzM0Pjxo2z7f/5559j6NChOHnyJI4dOwZBEPJMtoyNjfHbb7/hp59+wtq1a2FgYICePXuiRYsWn0weC0JTUxM9evTAqVOnMG3aNEV75ndnYGAgvv/++xz7/fHHH4rvYXV1dWzatAkrV67Ezp07kZaWhnr16mH58uWoVq1atv1u3LiRY8mizDtWEyZMUCRbmce/cOECLly4kOP4HyZbcrkcp0+fRt++fUt0grxIyG/qW4bcu3cPAwcOBJBR8TY/T3j8+eef+OGHH2BjY4OzZ8/m2mfWrFnw9vZGv379sGTJkmKL193dHT4+PnB2ds5X9WMiIiJVePPmDbp27YrNmzfnOgG+tDt79iymT5+OM2fOFGkyfqb8fn+XyTlbn5JZidjJySnfj9JmPpJra2ubZ5/MbcVV4ZaIiKgssbGxQZ8+fbBp0yZVh1IomzdvxuDBg4sl0SqIMnkb8WOePHmCffv2Aci4XJxfcXFxAPDRgm2Z2/IzYXHfvn04cOBAvo7t7++fr35ERESqtnDhQlWHUGj79+9XyXHLVbIVGRmJiRMnQiqVolOnTnne385NZtXaj9XRynxKIjU19ZPjRUREwMfHJ9/HJyIiovKp3CRbmcXM3r17B2dnZ6xYsaJA+2dWrc3tcdxMmRVu81PoztTUVPE0yaf4+/vnK4EjIiKisqdcJFtJSUkYOXIknj59ipo1a+J///tfgZc9yCyUl3k7MTeZ2/KqbPyhAQMGZFu89GMyJ9gRERFR+VPmJ8inpKRgzJgxePDgAezt7bFt27ZCleDPrFobFBSUZ5/Mom/29vaFipWIiIgqnjKdbKWlpWHcuHG4ffs2rKyssH379lzXm8oPFxcXABnLG4SFheXa5+7duwAyqugSERER5UeZTbbS09MxceJEXL9+HWZmZtixYwcsLCwKPZ6DgwMcHR0B5P60wvXr1xEUFASxWJytujwRERHRx5TJZEsmk2H69Om4dOkSTE1NsWPHDsUCn5/Svn17tG/fHqdOncqxbcKECQAy6nB8uGzBq1evMHfuXADAoEGDuFQPERER5VuZnCD/999/4/Tp0wAyyjHMnj07z77z5s3Ltv5R5srwycnJOfp27twZX331FXbs2IFx48bB1tYWurq6ePnyJWQyGVxdXTF9+vRifjVERERUnpXJZCuzBAOQkTxlJlC5SUhIKNDYs2fPRsOGDbFnzx74+voiPDwc1atXR8+ePTF8+PCP1uEiIiIi+q8ymWy5u7vna8HP3Dx//vyTfbp27YquXbsWanwiIiKiD5XJZIuIiMo+b29vzJo1S/GzpqYmDA0NUatWLbRt2xbu7u7ZaiZ6enpiw4YNuY61YMECbN68+aN3OjItX7680H+wExUGky0iIlKpSZMmwdraGlKpFJGRkbh16xaWLVuG7du3Y+PGjXBycsrWf8GCBdDV1c3W5uLiAlNTUyQlJSna/v33X5w4cQKzZs3KVn+xUaNGyn1BRP/BZKsCEwQBkoAAaFWrpupQiKgCa9OmDerVq6f4ecyYMbh+/TrGjh2L8ePH46+//sq2TFrnzp1zfSr8vwWnIyMjceLECXTs2BHW1tZKi5/oU8pk6QcqHlGbNiOgV2/EHvZWdShERNm0aNEC48ePR3BwMI4dO6bqcIiKhMlWBZVw8SIi1q6FkJ6OkDlzELp0GYSPLMJNRFTSevXqBQC4cuVKtva4uDhER0cr/vvYmrZEpQFvI1ZQ6vr6UDcxgSwqCgAQs3Mn0l68gNXaX6BRiLUliajkpUhkkMjkKo1BU10NOprqShnb3Nwc+vr6ePPmTbb2Ll26ZPvZysoqWyFqotKGyVYFpevqCodDB/H22wlIffoUAJB88yYC+/aD9a8boP2fCalEVLosPO6DHdcCIRdUG4eaCPiqpT1+7OGslPF1dXWzTXoHMp5K/PApRS0tLaUcm6i48DZiBSa2sIDdnt0w6NFD0ZYeHIzAgYMQn8tyRkRUevxxPUjliRYAyIWMWJQlOTkZlSpVytbWuHFjtGzZUvGfq6ur0o5PVByYbFVwatrasFz5E6p+9x2glvHPQUhJQfCUqQj/ZS0EuWpvURBR7oa1sIOaSNVRAOpqIgxrYaeUsUNDQ5GQkABbW1uljE9UUngbkSASiVD5mxHQcnRE8PTpkMfHAwCifv8dac+ewWrNz1D7z1+WRKRaP/Zwxvedncr1nK2jR48CANzc3JQyPlFJYbJFCnqt3eBwYD/efDsBEn9/AICQng4R50MQlUo6murQgXISHVW7fv06Nm7cCGtra/Ts2VPV4RAVCZMtykbT3h72+/fh3fc/IM3PD1ZrfoZIg/9MiEh5/v33X7x69QoymQyRkZG4efMmrl69CktLS3h5eXECPJV5/BalHNT19GC9wRPSyEioGxmpOhwiKufWr18PABCLxTAyMoKjoyNmz56dY21EorJKJAhCKXiepWJzd3eHj48PnJ2d4e1dMtXc9956jYdvYjGkuR3qWhnme7/EK1cR9+efsFi0EGr/WZuMiIioIsnv9zevbFVAz0MTMMv7MQDA+14w5nWvjSHN7SASffzRJklQEIKnTYM8Ph5pfn6w3rABmtZWJREyERFRmcXSDxWQgY4GdMQZk2olMjnmHfXBxL33kZgm/eh+Udu3K55UTHv2DIH9+iHp5i2lx0tERFSWMdmqgCwMdbB1eBNU0cuadHriUQh6el6Bb0h8nvuZz5kDk6++Uvwsi4nB6xEjEL1zF3g3moiIKHdMtiqoFtUr46/JbmhezUTR9ioyCb1/vYoDt9/kmjyJNDRgNmsmLFYsh0hTM6NRJkPY0qUImTsXcomkpMInIiIqM5hsVWBV9bWxe2RzTGxfA5nTtdKkcnx/+BFmHHyEZEnutxWNeveG3e5d0DAzU7TFHfZG0NChSA8LL4nQiYiIygwmWxWcupoI0z+vhe1fN4VJJU1F++F7b9Frw1X4hSfkup9OvXpwOHQQOo0aKdpSHz5CYN++SHnwQNlhExERlRlMtggA0NbRFCcnuaGxnbGi7WV4Inp4XsWf94Nz3UfD1BR227fByMND0SaNiMC7WbMhSD8+2Z6IiKiiYLJFChaGOtg7ujnGtKmmaEtJl2HK/geY5f0YqemyHPuINDVhsWghzBf8CGhoQKStzarzREREH+A3ImUjVlfDrC9qo7G9CaYfeID41IwrVHtvvcajt7HYOLgR7CrnXJTaeMAAaNWoAWlMDLRr1y7psImIiEotXtmiXHWqY4aTk1qjvnVWdXmfd/Ho7nkFp56E5rqPbuPGMOjUKUe75G0wUp8+VVqsREREpRmTLcqTjYkuDo5tga9a2CnaElKlGLvrLhafeAqJVP7JMeTJyXj77bcIHDQYcSdOKjNcIiKiUonJFn2UloY6FvaqC8+BDVFJU13R/r8rARiw6TrexaZ8dP+QefOR9vw5hNRUvJsxA+GrV0OQ5Zz7RUREVF4x2aJ86eFiiWMT3eBkrq9ou/c6Ft3WX8bF53nX1jJ0/xJqhlm3IqO2/A9vxoyFLC5OqfESERGVFky2KN+qm+rhyPhW6OdqrWiLSU7H19tv4+d/nkMmz1l1Xq9VKzgcOgitmjUVbUlXriCgnwfSXr4skbiJiIhUickWFYiOpjpW9XPByr71oS3O+OcjCIDneT8M/d9NRCSk5dhH08YG9vv2Qv/zzxVt6a9fI7D/AMSfOVNisRMREakCky0qFI/GNvjz21aoViWrDMQ1/yh0W38ZtwKic/RXq1QJVuvWwnTKZGSuDSRPTkbwxEmIWO8JQf7pyfZERERlEZMtKjQncwMcndAK3epbKNrCE9IwcPMN/HbJP8di1iKRCFXGjoX1xl+hpqenaI/cuBHvZswosbiJiIhKEpMtKhJ9bTE2DGyIBT3qQKyeccVKJhew4u9nGPXHXcQlp+fcp1072B/YD00HB0WbXrt2JRYzERFRSWKyRUUmEokwvJUDDoxpASsjHUX7Wd8wdN9wGY/f5nzyUKtaNdgf2A+9tm1h8vXXMOzRoyRDJiIiKjFMtqjYNLQ1xomJbvislqmi7U10Cvp4XcOuG0E5biuq6+vD2msjqk6flmMsQSLJ0Z+IiKgsYrJFxcq4kia2ftUE33WuBbWMu4qQyOSY++cTTN3/AMkSabb+IjW1HItWC4KAd3PmInjqNMiTkkoqdCIiIqVgskXFTk1NhG/b1cCukc1QRU9T0f7ng3foteEq/MITP7p/9LbtiD9+HAmnTiFw4CBI3rxRdshERERKw2SLlKZl9Sr4a1JrNHUwUbS9DE9Ezw1XcOzhu1z3EeRyJN24rvg57cULBPbth6Rr15QeLxERkTIw2SKlqmqgjT0jm2Fs2+qKtmSJDJP23sf8o0+QJs2+TqJITQ02GzfCZMQIRZssLg6vR45C1LbtnMdFRERlDpMtUjoNdTXM7OqEzcMaQ187a37WH9eD4PH7DbyNSc7WX6ShAbPvv4PlqlUQaWllNMrlCP/pJ7z7/gfIU1NLMnwiIqIiYbJFJaZTHTOcnNgada0MFG0P38Siu+cVXMhlMWvDHt1ht2c3NCyziqbGHz+OoEGDkf4u99uQREREpQ2TLSpRtpV1cWhsSwxqZqtoi01Ox9fbcl/MWsfZGQ6HDkG3aVNFW+rTpwjo2w9Jt26VWNxERESFxWSLSpy2WB3LvqyHNR4u0BGrK9ozF7OOTMy+mLWGiQls/7cFxsOGKtpk0dFIPHeuxGImIiIqLCZbpDLujawzFrM2zbmY9Z3A7ItZi8RimM+eDYtlyyDS1IRu06aoyvUUiYioDGCyRSpVy1wfxya4ofsHi1mHxadhwKYb2HL5VY6nD43cv4Td7l2wWvsLRGJxSYdLRERUYEy2SOX0tDTg+Z/FrKVyAUtO+mL87ntISM2+mLVOvXrQMDHJMU7MwYNIvnu3RGImIiLKLyZbVCpkLma9f0wLWBhqK9r/fhKKnhuuwjck/qP7J924idAFCxH01XDE7N3LelxERFRqMNmiUqXR+8WsW9esomgLiEzClxuv4tDdt7nuI09Lw7sffgBkMkAqRejCRQidPx9yiaSkwiYiIsoTky0qdSrraWH7100xuUNNiN4vZp2aLseMgw8x8/AjpKZnrzqvpqUFq7W/QMPUVNEWe/AQXg8dhvSwnPW7iIiIShKTLSqV1NVEmNrJEduGN4GxbtZE+H2336CP1zW8jspedV63YUPYHzoEHRcXRVvKw4cI7NsXyffvl1jcRERE/8Vki0q1z2pVxYlJrdHAxkjR5vMuHt09L+Ps07BsfcVmVWG78w8Y9euraJNGRCBo2FeIOXiwpEImIiLKhskWlXpWRjo4MKYFvmphp2iLT5Vi5B938NOpZ5DK5Ip2NU1NmC9aBPMFPwIa79dhTE9H6Lz5CFmwAALncRERUQljskVlgqaGGhb2qov1AxtCVzOr6rzXRX8M/d8tRCRkVZ0XiUQwHjAAdju2Q71K1kT72H37Ef7zzyUaNxEREZMtKlN6ulji6LetUKOqnqLt+quMqvO3/1N1XtfVFQ6HDkK7Xj0AgIa5OSqPGlWi8RIRETHZojKnppk+jn7bCj1cLBVt4Qm5V50Xm5vDbtdOGPXvD2tPT2h8cKWLiIioJDDZojKpkpYG1g9ogIU9nRVV52V5VJ1X09KCxcIF0KlXN8c4aa9ecR4XEREpFZMtKrNEIhG+ammP/WNawDKXqvPPQj9edT793TsEDRmKoK9HQBoRoexwiYiogmKyRWVeI1tjnJjUOkfV+d6/XoX3vTyqzqem4u2EiZBFRyPl7l0E9O2HlEePSipkIiKqQJhsUblgUkkz16rz0w48xJwjj5EmzV51HiIRtBwdFT9Kw8IQNHgIYg8fLsGoiYioImCyReXGh1XnjT6oOr/75mv0++063kRnVZ1X09KCxfJlMJszB1DPKCUhpKcjZM5chC5aDCE9Pcf4REREhcFki8qdz2pVxYmJbnCxNlS0PXobhx4bruDC86y1EkUiEUyGDoHttq1QNzFRtMfs2YOgr7+GNDKyROMmIqLyickWlUvWxro4MLYFhjS3VbTFJqdjxPbbWPPPc8jkWeUhKjVtCofDh6BdN+tpxZQ7dxHQpy/ncRERUZEx2aJyS0tDHUt618MaDxdoizP+qQsCsP68H4Zvu4WoxKyq82ILC9jt2gnD3r0VbdKwMAQNGYq4kydLOnQiIipHNFQdQGFFRETg6tWrePLkCR4/fgxfX1+kpaWhadOm2LlzZ6HG9PT0xIYNGz7aZ8GCBRg4cGChxifVcG9kjTqWBhi36x4CIpMAAJdfRqKH5xX8OrgRGtoaAwDUtLVhsXwZtJ2dEbZiBSCTAYIATRsbVYZPRERlXJlNtk6ePInly5crZezKlSvDzs4u122mpqZKOSYpl5O5AY5NaIXvDj7CKZ9QAMC7uFR4/H4d87vXwZDmdhCJRIp5XFq1HBE8ZSqqTp8Gnfr1VRw9ERGVZWU22dLT00PLli1Rr1491KtXD0+fPsXGjRuLZew2bdpgxYoVxTIWlR762mJ4DWmELZcDsOLUM8jkAtJlAuYd9cHdoBgsc68HXc2Mt0Slpk1R/e+/oG5omGMcQS6HSI134ImIKH/KbLLVt29f9O3bV/FzWFiYCqOhskIkEmFUm2pwsTHCt3vuISIhY97Wnw/e4WlIPLyGuKK6acYi17klWvK0NLwe8Q2M3L+EUZ8+JRo7ERGVTfzznCqkpg4mODnJDU0dsko+vAhLRE/PK/jrcUiu+wiCgNCFi5By9y5C5sxFyMKFXFeRiIg+iclWLp49e4bp06dj2LBhGDduHNauXYuXL1+qOiwqZlX1tbFnZDOMaVtN0ZYkkWH87ntYcuIp0mXybP2l4RFIvHRJ8XPs3n0IGv4111UkIqKPYrKVC19fX5w4cQI3b97E+fPn4eXlhR49emDZsmWQyWSfHoDKDA11NczqWhu/DXGFvlbWXfUtVwIwaPMNhMWnKtrEZlXhcOggtOvVU7Sl3LuXUY/rwYOSDJuIiMqQMjtnSxmqVq2KSZMmoXXr1rC2toaenh4CAgKwZ88e7Nu3Dzt27ICGhga+//77T461b98+HDhwIF/H9ff3L2roVERd6pqjlrk+xu26i2ehCQCA24Ex6Lb+CjYMaojm1SoDyKrHFbpwEeK8vQEA0vBwBA4dBvN5c2Hs4aGy10BERKUTk60P9O/fP0dbrVq1sHDhQlhbW2P16tXYsWMHBg0aBGtr64+OFRERAR8fH2WFSkrgUKUSjoxvhTl/Pob3vWAAQGRiGgZvuYnvO9fC6DbVIBKJMtZVXLoE2nWdEbZsOSCVAunpCJ3/I1Kf+MBs7hyoaWqq+NUQEVFpwWQrn0aMGIE//vgD4eHhOH/+PIYNG/bR/qampnB2ds7X2P7+/khNTf10R1I6HU11/NzPBa52xlh47CkkMjlkcgHL/36Gu0ExWO3hAgNtcUY9rkGDoF2rFt5OngLZ+3UUYw8cQNqLF7Batw5is6oqfjVERFQaMNnKJ3V1dbi4uODMmTMICgr6ZP8BAwZgwIAB+Rrb3d2dV8FKEZFIhMHN7FDX0hDjd99DcGwKAOCfp2Ho6XkFXkNcUdvCAACg6+oKh8OH8HbiJKS+X0cx5cEDJF27BqMve6vqJRARUSnCCfIFIBaLAQBSqVTFkVBJcLExwomJbmjjmLVqQGBUMr7ceBXe994q2sRmZhnrKvbNqLtl1L8/Ey0iIlJgslUAmeUfzM3NVRwJlRTjSprYNrwJpnSsCZEooy01XY5pBx5izpHHSJNmPJ2qpqkJi8WLYbV2LczmzFZhxEREVNow2cqnixcvKpKtVq1aqTgaKknqaiJM6eiIbcObwEhXrGjfffM1PH67rrjNKBKJYNClc47J8YIgIHLzZqSHhpZo3EREVDpUuGRr4MCBaN++PbZv356t/eXLl5g/fz6ePXuWrV0ul+PEiROYPn06AKBdu3aoz4WJK6TPalXFiYluqG+dtYzPw7dx6L7+Mi69yLuwaez+/Yj4eQ0C+vRF8u3bJREqERGVImV2gnxISAh69+6t+FnyftmUe/fuoVmzZor2kSNHYtSoUYqfw8LCEBwcjISEhGzjSaVS7N+/H/v374eRkREsLS2hrq6O169fIy4uDgDQuHFjrFy5Uomviko7a2NdHBzbAouOP8Xum68BADHJ6Ri+7RamdHDExPY1oKYmUvSXvHmD0KXLAACyqCgEfT0CZj/8AOMhgyESiXI9BhERlS9l9sqWTCZDbGys4r/k5GQAGUnTh+35LalgZWWFKVOm4LPPPoO+vj6CgoLw7NkziMVitGnTBitXrsQff/wBAwMDZb4sKgO0NNSx9Mt6+LmfC7Q0Mt5CggD8cvYFRuy4jdjkrPUSNW1sYPHjfIjeP1wBqRRhS5ciZOYsyFnug4ioQhAJgiCoOoiKLrP0g7OzM7zfVyWnsuHpu3iM230XQVHJijYrIx38NsQV9T643Zjy8CHeTpoMaViYok27Th1Ye66H2MqqRGMmIqLikd/v7zJ7ZYuoNKhjaYBjE9zQqY6Zoi04NgV9vK5h763XyPxbRsfFBQ6HD0GnsauiX+rTpwjo0xdJ16+XeNxERFRymGwRFZGhjhi/D3HFD12ckDldSyKTY5b3Y3x/6BFS0zPKQ2hUqQK7bdtgPHSoYl9ZbCxefzMSUf/bCl5kJiIqn5hsERUDNTURxn1WHbtGNkMVvazSDwfvvsWXG68hKCoJACASi2E+ZzYsf1oBkZZWRie5HOGrViH20CFVhE5ERErGZIuoGLWsXgUnJraGq52xos03JB7dPa/g7NOs+VqGvXrBfu8eiC0tAWTcZjTs1avE4yUiIuVjskVUzMwNtbFvdHN83cpe0ZaQKsXIP+5g1elnkMkzbhdq16kD+8OHYNCtG6zWr8tRDJWIiMoHJltESiBWV8OPPZzhObAhdDXVFe2/XvDHsK03EZmYBgDQMDaG1c+rITYzyzFGyhMfCHJ5icVMRETKUWxFTQMCAvDgwQOEh4cjOjoaaWlpMDIygomJCapXr45GjRpBR0enuA5HVCb0cLFEbQt9jNl5F/4RGfO2rvpFofv6K9g4pBEa2Rrnul/y/fsIGvYV9Nq2geWKFVDX0yvJsImIqBgVKdm6f/8+Dhw4gMuXLyMqKuqjfdXV1eHs7IwePXqgV69e0NfXL8qhicqMGlX1cXSCG344/AgnH4UAAELjU9H/9+uY260OhrWwy1ZNPj08HMGTJgPp6Ug8ew6B/Txg/esGaFWrpqqXQERERVCo24hHjx5Fjx49MGjQIBw5cgSRkZEQBAE6OjqwsLBA7dq10aBBAzg4OMDExAQikQhSqRQPHz7E0qVL0aZNG8ydOxchISHF/XqISiU9LQ1sGNgQ87vXgcb7+hDpMgE/HvPBlP0PkCyRKvqqaWlBq05txc+SgAAE9vNAwtmzJR43EREVXYGubN28eRM//fQTfH19IQgCDA0N0blzZzRu3BguLi6ws7PLdb+kpCQ8efIEDx8+xPnz5/HgwQMcOnQIx48fx7BhwzBmzBjo8TYJlXMikQgj3BxQ39oQ3+65h7D4jHlbRx+8g29IPLyGuKK6qR7UDQ1h4+WFyA0bELnRCwAgT0rC2wkTUXnMGJhOmgiRuvrHDkVERKVIgZbrcXJyAgC4ublhwIABaNu2LcSZa74VwJs3b3D06FHs2rULcXFxmDBhAr799tsCj1NecLmeiiciIQ0T997DjVfRijY9LQ2s6lsfXetZKNoSzp/Hu+++hzwpSdFWyc0NVqtXQd3IqCRDJiKi/1DKcj1ubm7Yv38/tmzZgo4dOxYq0QIAGxsbTJgwARcuXMC0adNgxC8NqmBM9bWw65tmGNu2uqItMU2KcbvvYenJp0iXZTyFqN++PewPHoRm9ax+SVeuIKBvP6Q+e1bicRMRUcEVKNnasmULXFxciu3gOjo6GDVqFAYPHlxsYxKVFRrqapjZ1Qm/D3WFvlbWHf3NlwMwePNNhMenAgC0qjnAfv9+6HfqpOiT/vYtAgcMROK//5Z43EREVDCss0WkYp2dzXFsohuczLOe0L0VGI1unldw81XGU77qepVgtX4dTKdOBd4/uaimpwetWrVUEjMREeUfky2iUsChSiUcGd8K7g2tFG0RCWkYtOUmNv/7CoIgQCQSocqY0bDZvBnqVarAeu0vuRZDJSKi0oXJFlEpoaOpjp89XLCkd11oqme8NWVyAUv/8sW4XfeQkJoOANBza4UaZ/6BbuPGOcaQSyQlGjMREX1asVWQBwCZTIbAwEDEx8dDKpV+tG+TJk2K89BE5YJIJMKQ5naoa2WIb3ffQ3BsCgDglE8onocl4Lchrqhlrg+1XFZjkEZEIHDAQFQePRpGHv2yFUolIiLVKZZkKzQ0FGvWrMGZM2eQmpr6yf4ikQhPnz4tjkMTlUsNbIxwfKIbJu+7j8svIwEAAZFJ6P3rVSx3r4feH9xuBABBIsHbKVORHhyM0B9/RMrjRzCfNw9qWlqqCJ+IiD5Q5NuIb968Qd++fXH8+HGkpKRAEIRP/ifn4rpEn2RSSRPbv26KSR1qKtpS0mWYsv8B5h99gjSpTNEuCQ6G5HWQ4ue4Q4cRNHgI0t+9K9GYiYgopyJf2VqzZg0iIyNhYmKC6dOnw83NDVWqVIE6K1wTFZm6mgjTOjmioY0Rpux/gLiUjHlbf1wPwqO3cdg4uBEsjXSg5eAAh0OHETxlClLu3wcApD55goA+fWH1yxpUat5clS+DiKhCK/KVrWvXrkEkEmHNmjXo06cPzMzMmGgRFbN2TlVxYqIb6lkZKtoevIlFt/WXcfllBABAbFYVdju2w2jgAEUfWUwMXo/4BlH/24oCLBZBRETFqMjJlkQigba2NprzL2cipbIx0cXBsS0wsKmtoi0mOR3Dtt7ChvMvIZcLEGlqwuLHH2GxbBlEmpoZneRyhK9aheBp07It+0NERCWjyMmWtbU1/2ImKiHaYnUsd6+HVX3rQ0sj4+0rCMDqf15g5B93EJeccZvRyP1L2O3ZAw3LD9ZZ/PsUAvr3R1pAgEpiJyKqqIqcbHXt2hVpaWm4fv16ccRDRPnQr7ENvMe3hK2JrqLt/LNwdPO8jCfBcQAAnbrOcDh8GLotsq46SwKDIIuKKvF4iYgqsiInWyNGjICTkxPmzZuHN2/eFEdMRJQPzpaGOD7RDR1rV1W0vY1JgbvXNey//RoAoGFsDNstW1B51CgAgNkPP+RaDJWIiJSnyE8jamtrY9u2bZg3bx569uyJzp07o169eqhUqdJH9+vdu3dRD01U4RnqiLFpaGP89q8/Vp9+DrkASKRy/HD4Me4GxWBRr7rQFquj6vRp0O/UEdr16qk6ZCKiCqdYipoGBwcjMjISKSkpOHr0KI4ePfrR/iKRiMkWUTFRUxNh/Gc10MDaCBP33kdUUsaSPQfuvMWT4Hj8NsQVtpV1oVO/fo59BYkE72bPQeVvRkC7du2SDp2IqEIocrL17NkzDB06FCkpGcuKiMViGBsbs/wDUQlrWaMKTk5qjfG77+Le61gAwNOQeHT3vIw1Hg3QsU7ORavDflqJ+BMnkHDmDCwWLYRhr14lHDURUflX5GTL09MTycnJsLGxweLFi9G0aVOoqXF9ayJVMDfUxr7RLbD8b19suxoIAIhPlWLkH3fwbbvqmNapFtTVMtZMTPPzQ8zevQAAIS0N736YiZSHj2A284esshFERFRkRc6K7t+/D5FIhF9++QXNmzdnokWkYpoaavixhzM8BzaErmbWFeZfL/hj2NabiExMAwBo1agB2y2boW5kpOgTs2cPgr4ajvSw8JIOm4io3CpyZpSSkgIdHR3UrVu3OOIhomLSw8USxya0QnXTrIdVrvpFofv6K7gbFA0AqNSyJRwOH4K2s7OiT8r9+wjo0wfJt2+XeMxEROVRkZMtOzs7SKVSyGSyT3cmohJVo6o+jk5wQ/f6WcVNQ+NT0f/3G9h2NQCCIEBsZQW7Pbth2Mdd0UcWGYmg4V8jescOFi0mIiqiIidbvXv3hkQiwfnz54sjHiIqZnpaGvAc2BA/9qgDjffztaRyAQuPP8XEvfeRlCaFmpYWLJYsgfnChRCJxRk7ymQIW74C76bPgDw5WYWvgIiobCtysjV06FA0b94c8+fPx/3794sjJiIqZiKRCF+3csD+Mc1hZqClaD/xKAS9fr0Kv/AEiEQiGPf3gN3uXdAwN1f0if/rLyScOaOKsImIyoUiP43o5eWFBg0a4OnTpxg0aBBcXV3zVdR0woQJRT00ERWQq50JTk5qjUl77+Oaf8ayPX7hiei54SpW9KmPni6W0KlfHw6HDyF42nQk37wJw149YdCzp4ojJyIqu4qcbG3YsAEiUcatCUEQcOfOHdy9e/eT+zHZIlKNKnpa2PlNM6w58xy/XvAHACRLZJi09z7uBcVg9he1oVm5Mmz/twXRu3bBuH9/xXuciIgKrsjJVpMmTYojDiIqQepqInzX2QkNbYwx9cADJKRKAQDbrwXi4dtYbBzcCBaGOqg8fHiu+8fs3Qv9Ll2gYWxcglETEZVNRU62du7cWRxxEJEKdKxjhpMTW2Psrrt4GhIPALj/Ohbd1l/B+gEN4VazSo59Yo/8idCFixC1eQus1q+HTl3nHH2IiChLgSfIX79+HXK5XBmxEJEK2FbWhff4lujf2EbRFp0kwdCtN7Hh/EvI5VmlHySvXyN0wQIAQPq7dwgaNAixhw+XdMhERGVKgZOtr7/+Gi1atMD333+PU6dOIZmPhBOVedpidfzUtz5+6lMPmhoZHwuCAKz+5wVG/nEHsckZi1uLbWxgOmki8H6lCEEiQcicuQiZNx/ytDSVxU9EVJoVONmqVasW4uLicOzYMUydOhXNmzfHmDFjsH//fkRERCgjRiIqIf2b2MJ7XEvYmugq2s4/C0d3zyt49DYWIpEIlb/5BrZb/wd1ExNFn9iDBxE0eAjS371TRdhERKVagZOto0eP4vz585gzZw6aNm0KuVyOS5cuYcGCBWjbti08PDywadMm+Pv7KyNeIlKyulaGOD7BDR1rV1W0vY1JQV+v69h9MwiCIKBS8+YZy/y41Ff0SX3yBAHufZB07ZoqwiYiKrVEQhHX4oiPj8fFixdx9uxZXLlyBcnJyYrHxG1tbdGhQwe0b98erq6ufHw8D+7u7vDx8YGzszO8vb1VHQ4RAEAuF/Dbv/5Yffo5Ppi2BfdGVljaux50NNUhl0gQtnw5Yvfuy+qgpgbTSZNQefQoiLgwPRGVY/n9/i5ysvUhiUSCGzdu4OzZs7hw4YLitqJIJIKxsTHatWuH9u3bw83NDVpaWp8YreJgskWl2TX/SEzaex+RiRJFm5O5PryGuMKhSkbx4tgjfyJ0wQIIH8zbMps7FyZDBpd4vEREJUUlydZ/PXr0CGfPnsW5c+cUtxVFIhG0tbXRsmVLdOjQAR06dIChoaGyQigTmGxRaRcWn4pvd9/DnaAYRZuelgZW96uPLnUzFrlO9fXF24mTkP72LbTq1Ib9nj1Q09ZWVchEREqX3+9vpV7jr1+/PqZNm4aTJ0/i9OnT+P7779GwYUOkpaXh3LlzmDNnDut0EZUBZgba2Du6OUa6OSjaEtOkGLvrHpaefIp0mRzatWvD4fAhGHTvDut165hoERG9V+SipvllZ2eHESNGYMSIEYiOjsb58+dx/vx56OjolFQIRFQEYnU1zO1eB652xvju0CMkpmVUnd98OQAP3sRiw6BGMDM0hNXqVbnun3zvPnTqOkOkqVmSYRMRqZxKZq+amJigb9++2LhxI7755htVhEBEhdS1ngWOTWgFJ3N9RdvtwBh0W38Z1/wjc90n5YkPXg8fjqBhXyE9LKykQiUiKhX4qBARFVg1Uz0cGd8K7o2sFG2RiRIM2XITv17wy1Z1XhoTg+BJkyBIJEh58CCjPMTNW6oIm4hIJQp8G/H27dvFcmAuYE1UtuloquPnfi5oYm+CH4/5QCKVQy4Aq04/x72gGKzxaABDXTHUdHWh27IF4g5lLOsji4rC6xEjUHXaVJiMGMGSMERU7hU42Ro6dGiRPxxFIhGePn1apDGISPVEIhEGNrVFXUtDjNt9F29jUgAA556Fo5vnZXgNdkU9a0NYLlkCHRcXhC1eAkEiAWQyhK9ajZQHD2GxfBnU9fRU/EqIiJSnULcRBUEo8n9EVH7UszbEyYmt0cEpe9X5Pl7XFFXnjfv1g93u3RBbWir6JJw5g8C+/ZD28qUqwiYiKhEFvrJ17ty5Ah8kOjoaGzduxMWLF5loEZVThrpibB7WGF6X/PHzPxlV5yUyOeYceYK7gTFY8mVd6NarC/vDh/Duu++RdOUKAEASGIgAj/6wWLwYht27qfhVEBEVvwInW1ZWVp/u9F5KSgq2bt2Kbdu2ISkpCYIgoFq1apg2bVpBD0tEZYCamgjftquBhrZG2arOe98PxpN3cfAa4orqpsaw+f03RP66EZEbNwIAhJQUvJsxA7KoSJh89ZUqXwIRUbFTytOIMpkMu3fvRqdOnbBhwwYkJibCzMwMS5YswYkTJ9CxY0dlHJaISomW1avg5KTWaGJvrGh7EZaInp5XcOLRO4jU1WE6aSJsfv8NagYGAAB1Q0Po87OBiMqhYk+2Tp48ia5du2LJkiWIjIyEgYEBZsyYgdOnT6Nv375Q48K0RBWCmYE29oxqjjFtqinakiQyTNhzHwveP72o17YtHLwPQ7tuXVj+/DPEBbhyTkRUVhRbBfkrV65gzZo18PX1hSAI0NbWxpAhQzB69GgYvP/LlYgqFrG6GmZ9URuN7Iwx4+BDJKRmVJ3ffi0QD97E4tfBjWBlbQ37A/shyuUPsfSwMGiYmua6jYiorCjyJ9jjx48xfPhwjBo1Ck+fPoWamhr69euH06dPY8aMGUy0iAidnc1xYqIb6lhkfR48eBOL7usv49KLiFyTKWlMDAIHDMTbCRMhi48vyXCJiIpVoZOtoKAgTJ48GR4eHrhx4wYEQUCnTp1w/PhxLF68GGZmZsUZJxGVcXaVK8F7fEsMbGqjaItJTsfwbbew5swLyD6oOi/IZHg3fQakISFIPH8eAX37IfXZM1WETURUZAVOtiIiIvDjjz+iW7duOH36NARBQJMmTbB//354enqiWrVqnx6EiCokbbE6lrvXx+p+LtAWZ3z8CAKw/txLfLX1FiIT0zLaUlMBcdYsh/TXrxHYfwBij/ypirCJiIqkwMlWp06dcODAAUilUjg5OWHTpk3YuXMnXFxclBEfEZVDfV2tcWR8KzhUqaRou+IXiW7rL+NOYDTUKlWCjZcXqkyaCLxfsUJIS0PIrFkImf8j5GlpqgqdiKjACjxBPjU1FSKRCCKRCPHx8Vi4cGGBDyoSiXD27NkC70dE5UdtCwMcm9AKPxx+hL8ehwIAwuLT0H/TDczs4oSRrR1gOn48dOq7ZNTgio0FAMQeOIBUHx9YrVsHTWs+vUhEpV+hnkbMrAL/7t27Qh2UC88SEQDoa4vx66BG2H4tEEtP+kIqFyCTC1j6ly/uBEVjZV8XGLq1goP3YbydMhWpjx4BAFJ9fBDQpw+sVv4EvbZtVfwqiIg+rsDJ1oQJE5QRR4FFRETg6tWrePLkCR4/fgxfX1+kpaWhadOm2LlzZ5HGvnHjBrZt24aHDx8iOTkZlpaW6NKlC0aPHg1dXd1iegVEBGT88fV1Kwe42Bhhwu57eBeXCgA47RMG35Ar2Di4EepaWcJu106Er1iBmD17AQDyuDi8GTMWVb//HpVHfK3Kl0BE9FFlNtk6efIkli9fXuzj7ty5E0uXLoUgCDA3N4eFhQX8/Pzg5eWFf/75B3v27IGRkVGxH5eoomtka4wTk1pjyv4H+PdFBADgdXQy3L2uYUEPZwxsagPz+fOh07AhQub/CCElBVBTg3adOiqOnIjo48pspUA9PT20bNkSY8aMwYYNGzB+/Pgij/nkyRMsW7YMALBo0SJcvHgRR44cwdmzZ+Hs7Ax/f3/MmzevyMchotyZVNLE9uFNML2TY+a8eEikcsw+8hjTDzxEskQKwx49YL9/HzTt7WE6eTIqNW+m2qCJiD6h2CrIl7S+ffuib9++ip/DwsKKPObGjRshl8vRu3dv9O/fX9FuZmaGNWvWoGvXrvjnn3/w7NkzODk5Ffl4RJSTmpoIEzvURENbY0zedx9RSdkXs944uBFqODrC4fAhiHR0cuwvT06GSEeHc0OJqNQo1mRLJpMhMDAQ8fHxkEqlH+3bpEmT4jx0kSUlJeHy5csAAA8Pjxzb7e3t0bx5c1y7dg2nTp1iskWkZG41Mxaznrj3Hm4HxgB4v5j1hqtY7l4PvRrkfBJRkErxZsxYaJiawmLxIqhVqpSjDxFRSSuWZCs0NBRr1qzBmTNnkJqa+sn+IpEIT58+LY5DFxtfX19IJBJoamqifv36ufZxdXXFtWvX8PDhwxKOjqhiMjfMWMx69enn+P3fVwCAZIkMk/c9wK2AaMzrXgfaYnVF/4h165B8+zYAIPX5c1ivXwet6tVVEjsRUaYiJ1tv3rzBwIEDERUVpSgJ8Sn57VeSAgICAACWlpYQi8W59rG1tc3W92P27duHAwcO5OvY/v7++YySqOLJXMy6sb0Jph94gPj3i1nvvvkaD9/GYuMgV9hW1oUglSLlyRPFfhJ/fwT084DFokUw7N5NVeETERU92VqzZg0iIyNhYmKC6dOnw83NDVWqVIG6uvqndy5F4uLiAACGhoZ59sncltn3YyIiIuDj41M8wREROtUxw8lJrfHtnnt49DbjPfgkOB7dPC9jdT8XdHY2h+3mzYhY74moTZsAAEJyMt7NmIGUe/dQdeYPUNPUVOVLIKIKqsjJ1rVr1yASibBmzRo0b968OGJSibT3y3/kdVULADTff1Cn5WOpEFNTUzg7O+fr2P7+/vm6/UpU0dmY6OLg2BZYcsIXO28EAQASUqUYs/MuRro54IeuTqg6bSp0GjTAux9+gDwhAQAQs2cPUp48gfXaXyC2tFTlSyCiCqjIyZZEIoG2tnaZTrQAQEtLCwCQnp6eZx+JRJKt78cMGDAAAwYMyNex3d3deRWMKJ+0NNSxuHddNHEwwazDj5AkkQEAtlwJwP03sfAc2BCW7dtlVJ2fPBlpT30BAKmPHiHAvQ8sV62CXms3Vb4EIqpgilxny9raulTOwSqo/NwizM+tRiIqGT1dLHFsohtqmekr2u4GxaDb+su4+DwcmjY2sN+7F0b9+im2y2Jj8Wb0aERs3KiKkImogipystW1a1ekpaXh+vXrxRGPytjb2wPIWO8xr6tbr1+/ztaXiFSruqke/vy2Ffq6WivaYpLTMXzbbaw+/RxyDTEsFi+CxfLlEGlrZ3QQBKjpcNktIio5RU62RowYAScnJ8ybNw9v3rwpjphUonbt2hCLxZBIJHj0frHb/7p79y4AoEGDBiUYGRF9jI6mOlb3c8HKvvWhLc76SNtwwQ9D/ncT4fGpMPqyd0bVeTs76HfqBJPhX6kwYiKqaIo8Z0tbWxvbtm3DvHnz0LNnT3Tu3Bn16tVDpU8UE+zdu3dRD12s9PT04ObmhgsXLuDAgQNwdXXNtj0wMBA3btwAAHTp0kUVIRLRR3g0tkF9a0OM330PryKSAAA3XkXji/VXsH5AA7SsVQv2hw8BQI7q8oJUCqirs+o8ESlFsRQ1DQ4ORmRkJFJSUnD06FEcPXr0o/1FIpHKkq2BAwciLCwMw4YNw/Dhw7NtGz9+PC5evIijR4+iUaNG8PDwgEgkQnh4OKZNmwa5XI6OHTuyejxRKeVkboBjE9ww2/sxjj18BwCITEzDkP/dxJSOjvi2XQ2oq+VMqMJXrUZ6SAgsli6Bur5+ju1EREVR5GTr2bNnGDp0KFJSUgBklE4wNjZWep2tkJCQbAlb5pOC9+7dQ7NmWQvTjhw5EqNGjVL8HBYWhuDgYCS8fyT8Q/Xr18fMmTOxYsUKzJ8/H15eXjA2Noafnx8kEgkcHBywePFi5b0oIioyPS0NrBvQAM2qmWDhsaeQyOSQC8CaMy9wOzAav/RvgCp6WU8Ux5/+B9E7dgAAUp8/g/W6ddDmH1REVIyKnGx5enoiOTkZNjY2WLx4MZo2bQo1tSJPBfskmUyG2NjYHO1SqTRbe0HrVw0fPhy1atXC1q1b8ejRI0RFRcHS0hJdunTB6NGjP3l7lIhUTyQSYXAzO7hYG2H87nt4HZ0MALj8MhLd1l+G58BGaOpgAgCI/+svxX7pQa8R2H8AzOfPg1GfPiqJnYjKH5FQxLoNLVu2RExMDA4ePIi6desWV1wVSmadLWdnZ3h7e6s6HKJyJT41Hd8ffIRTPqGKNnU1EaZ1csS4ttUhEuSI9PoNkb/+CnzwcWjo7g7zeXOhpqOjirCJqAzI7/d3kS9BpaSkQEdHh4kWEZVKBtpieA1phPnd60CsnjFfSyYXsOr0c4zYcRsxqTKYTvgWNls2Q93YWLFfnLc3AvsPQFo+1kIlIvqYIidbdnZ2kEqlkMlkxREPEVGxE4lEGOHmgINjW8LKKOtK1cXnEei2/jLuBEZDr1UrOBzxhk7DhortaS9eILBvP8SfOqWKsImonChystW7d29IJBKcP3++OOIhIlKaBjZGODnJDR1rmynaQuJS0X/TDfx2yR/qVc1g98cOmHzwpLI8KQnBU6YidOkyCPyjkogKocjJ1tChQ9G8eXPMnz8f9+/fL46YiIiUxkhXE5uHuWJut9rQUMu6rbji72cY+ccdxEoEmM38AVbr10FNT0+xnywmBiiBh3+IqPwp8tOIXl5eaNCgAZ4+fYpBgwbB1dU1X0VNJ0yYUNRDExEVikgkwsjW1dDQ1hgT99zDu7iMp5bPPwtHt/WXsWFwIzT6/HNo16qFt1OmQkhLg8XCBSx6SkSFUuRka8OGDYoPIEEQcOfOHcWyNh/DZIuIVM3Vzhh/TW6N6Qce4tyzcADAu7hUePx2HT90ccLI1g6w37sHsqgoqOXyB6Qgl0PEq11E9AlFTraaNGlSHHEQEalExm3Fxth8+RVWnn4OmVyAVC5g6V++uBkQhdX9XGBkZZVjv7hjxxB74CAsf/4ZYrOqKoiciMqKIidbO3fuLI44iIhURk1NhDFtq8PVzhgT995HyPvbimd9w/HFusvwHNQIrnZZZSHS/PwQ8uMCCCkpCHB3h9XqVajUooWqwieiUo7Xv4mI3mtsb4K/JrVGe6esK1Xv4lLR//fr+P2SP+TyjKKncSdOQHi/RJksKgqvR3yDiI0bIcjlKombiEo3JltERB8wrqSJLcMaY/YXTopFq6VyAcvfP60YkySB6eTJMF+wACJNzYydBAGR6z3xZtRoSKOjVRg9EZVGTLaIiP5DTU2E0W2q48CY5rA01Fa0n38Wji/WX8a91zEwHtAfdnv3QGxjo9iedPUqAr50R3I+HhIiooqjQMnWli1bCryw86c8fvwYly5dKtYxiYiKg6udCU5Oao0OH9xWDIlLhcfvGUVQtWrXgYP3Yeh36qTYLg0LQ9CwrxC1ZQtvKxIRgAImW6tXr0bHjh2xfft2xMfHF+nAd+7cwZgxY+Dh4YHHjx8XaSwiImUxrqSJLV81xpwvchZBHbHjNmJFmrBavw5ms2cBYnHGTjIZwlf/jLfjxkNWxM9KIir7CpRsjRkzBomJifjpp5/g5uaGSZMm4fTp04iKivrkvunp6Xj06BHWrl2Ljh07YujQobh06RLq1auHjh07FvoFEBEpm0gkwqg21XBgbIscayt+sf4ybgVEw2TYMNjv2gkNSwvF9vSwsKx5XURUYYkEQRAKskNYWBjWrFmDEydOQCaTKQqaWlhYoFatWjAxMYGhoSE0NTURHx+PuLg4vHnzBs+ePUN6ejqAjOKntra2mDx5Mrp161b8r6qMcXd3h4+PD5ydneHt7a3qcIjoI2KTJZhx8BHO+oYp2tREwNSOjhjfrgYQH4d3M2ch+fZtOBw+BE17e9UFS0RKld/v7wInW5nCwsJw4MABHD58GKGhoVkD5rKcReYhNDQ00LZtW/Tv3x+tW7fm0hfvMdkiKlsEQcC2q4FY/rcv0mVZH6GtalTGL/0bwLSSJiSBgdCqVi3nvjIZROrqJRkuESlJfr+/C13U1MzMDBMnTsTEiRPx4sUL3L59G48ePUJ4eDiio6ORlpYGIyMjmJiYoEaNGmjcuDEaNWoEvQ8WdiUiKotEIhFGuDmgsb0xJuy5j9fRyQCAq35R+GLdZazt3xBuNXMmWqlPnyJ4xnewXLEcOvXrl3TYRKQiRa4gDwCOjo5wdHTE4MGDi2M4IqIyob61EU5McsPMw4/w1+OMK/yRiRIM3XoT335WA1M61oSGesbUWFlCAt5OmYr0168ROHgIzL6bAeOhQ3mFn6gCYJ0tIqIiMNAW49dBjbCkd11oamR8pAoCsOGCHwZtvomQuIxK88m3biE9ODhjp/R0hC1bjuBJk/i0IlEFwGSLiKiIRCIRhjS3w5/jW6FalUqK9luB0fhi3WWc8w2DfocOsNuxHRpVs2p2JZw5iwD3Pkhh+Ruico3JFhFRMaljaYBjE93Qu4Gloi0mOR3f7LiDRcefQqNBIzgc8UalVq0U29PfvkXgoMGI/mMnCvm8EhGVcky2iIiKkZ6WBn7p3wAr+9aHjjjrqcOtVwPQx+sagqENm82bYDplMqD2/iM4PR1hy5bxtiJROcVki4iomIlEIng0tsHxia1Qy0xf0f44OA7d1l/B8cehqDJ2LGy3b4OGqalie+ZtxdTnz1URNhEpCZMtIiIlqVFVH0cntMKgZraKtsQ0KSbtvY9Z3o+g1sAVDn8eQaWWLRXb5QkJUNfXz204IiqjmGwRESmRtlgdy76sh18HNYK+Vla1nb233qDXr1fwSqoJmy2bYTp5EqCuDosVyyG2tPzIiERU1hQ52VqyZAkSEhKKIxYionKrW30LnJzUGi7Whoq2F2GJ6LnhCvbfeYvKY8ei+ulT0G/XLse+8rS0kgyViIpZkZOtXbt2oUuXLjh8+HBxxENEVG7ZVtbFwbEtMaq1g6ItNV2Omd6PMWHvfaRWMcuxjzQmBq+6dUfU1m18WpGojCpyslW3bl1ERUVh7ty56N+/P548eVIccRERlUuaGmqY060Otg5vDGNdsaL95KMQfLHuMu6/jlG0CXI5QmbOQvrbtwhfuRJvx42HNCYmt2GJqBQrcrJ18OBBLFq0CEZGRnj48CE8PDwwf/58xPADgYgoT+2dzPD35DZoXs1E0fY2JgX9frsOr4v+kMsFSMPCkPrsmWJ74sWLCHDvg+R791QRMhEVUpGTLZFIBA8PD5w+fRoDBw6ESCTCwYMH0aVLF+zZs4eXvYmI8mBuqI3dI5tjeidHqL1fIlEqF/DTqWf4atstxOqZZDyt2LaNYh9pSAiChg5D5KbNEORyFUVORAVRbE8jGhgY4Mcff8Thw4fRsGFDxMXFYfHixejTpw/u8a8wIqJcqauJMLFDTewf0wKWhtqK9ssvI9F13WVci5DCxssLVb/7DtB4/zSjTIaINWvwZtRoSKOiVBQ5EeVXsZd+cHJywp49e7BixQpUrlwZT58+xeDBg/HDDz8gMjKyuA9HRFQuNLE3wV+TW6Ozc9Yk+cjENAzbegsrTj+HwfDhsN+1M1tZiKSrV/Gqd28k3bihipCJKJ+UVmerd+/eOH36NIYPHw51dXUcO3YMXbp0wfbt2yGTyZR1WCKiMstIVxO/DXHF4t51oamR9fH8+6VX6PvbdYTb1ITDEW/od+qo2CaLiMTrr0cgavt2FURMRPmh1KKmlSpVwsyZM+Ht7Y1q1aohMTERP/30E7p3744VK1bg2LFj8PPz47wuIqL3RCIRhja3w9FvW6FGVT1F+8M3sei2/gpOBCTCav16mM2dC5E462lGbafaqgiXiPJB49Nd8i81NRUBAQHw9/fHq1evFP/7+vVrSCQSiEQiCIKAgIAABAYGKvbT1tZGrVq1sG/fvuIMh4iozKptYYBjE1ph0fGn2Hf7DYCMpX4m73uAf19EYmG//rBr2ADB06bBsHsPVGreTMURE1FeipxsxcfHY/Xq1bhy5QpCQ0OzXaX67xUrU1NTVKtWDXZ2doiMjISvry9CQkKQkpKChw8fFjUUIqJyRVdTAyv61IdbzSqY5f0YCalSAMDhe29xNygangMboc5hb6jpaOfYVxoTA/VKlSDS1CzpsInoP4qcbM2ePRvnzp1TJFaampqws7ODg4MDHBwcUK1aNVSrVg0ODg7Q09PLsX9MTAyePn2Kp0+fFjUUIqJyqXt9S7hYG2Hyvvu49zoWABAYlQx3r6v4vrMTvnFzgOiD/oJUiuCJkyCXSGD182po2tioJG4iylDkZOvG+6dgRo4cib59+8LOzg4ikegTe2UxNjZGq1at0KpVq6KGQkRUbtmY6OLAmBZYf+4lPC/4QRCAdJmApX/54t+XEfjZwwVV9TOucEV6/YbkO3cAAAFfusNi8SIYdO2qyvCJKrQiT5A3MDCAgYEBZsyYAXt7+wIlWkRElH8a6mqY9nkt7BnZHOYG2WtyfbHuMi4+D4cgCJBGZZXZkScmInjqNITM/xHylBRVhE1U4RU52Tp//jzOnz9fHLEQEVE+tKheGX9Pbo1OdT6sySXB8G23sfiEL0zmzoPVmp+hVqmSYnvsgQMI9PBA2suXqgiZqEIrltIPlT54QxMRkfIZV9LEpqGuWNzLOVtNrq1XA9D712sIb9waDke8oV23rmJb2ks/BPTzQMyBAyy5Q1SClFpni4iIlEckEmFoC3scn+CGWmb6inbfkHh097yCgyEC7HbvgsnXXyu2CampCJ3/I4KnTYMsIUEVYRNVOEy2iIjKuFrm+jg6oRWGt7RXtKWmyzHnyBOMO/AYmhOmwOb336BubKzYnvD3KQR69IcgkaggYqKKhckWEVE5oC1Wx4KezvjfV41hUimrttZpnzB0XXcZj6yd4fDnn9BtllX81Ki/B+twEZUAJltEROVIh9pmODW5NVrXrKJoC41PxeAtN/HL/WhYbN4M08mToNe+PUy++kqFkRJVHEy2iIjKmaoG2tjxdVPM7VYbYvWMcjyCAGy86I++m24iod8wWG/wzFGqRy6RIOnWLVWETFSuMdkiIiqH1NREGNm6Go6Mb4VqpllPjD96G4du6y/jwN23OZ5IjPj5Z7we9hXCVq7iXC6iYsRki4ioHKtrZYgTE90wsGnWkj3JEhl+OPwY43bdQ0xSRlKVcPEionf8AQCI3roVgQMHQRIYqIqQicodJltEROWcrqYGlrvXx29DXGGkK1a0n/IJRZd1/+KqXyTElpbQrFFdsS3VxwcB7n0Qd/SoKkImKleYbBERVRBd6prj9JQ2cKuRNXk+LD4NQ/53E7/4yWC5dz+M+vdXbJMnJ+PdDzMR/N33kCUmqiJkonKByRYRUQViZqCNP0ZkTJ7XVM/4ChAE4Pd/X6HvtntI+nY6rNatg5qBgWKf+OPHEfClO1IePVJV2ERlGpMtIqIKRjF5/tuWqFFVT9Hu8y6j8vxRA0c4HPGGTmNXxbb0N28QOGgwIjdthiCXqyJsojKLyRYRUQXlbGmI4xPcMLS5naItNV2OeUd9MOafYOh6/o4qEyYAau+/KqRSRKxZg6TLl1UUMVHZxGSLiKgC09FUx+LedfG/rxqj8geV5y88j0CXDddwv10f2O38AxqWFgAAg27dUKlNG1WFS1QmMdkiIqKMyvNT2qC9U1VFW3SSBKN33sXiIDHM9h+C8eDBMF/wY45iqET0cUy2iIgIAGCqr4X/fdUYS3rXhbY46+th76036Ln9IUK++hbq+vo59otY74nU589LMlSiMoXJFhERKYhEIgxpboeTk1qjnpWhoj0wKhl9f7uOdWdfQirLmiAfd/w4IjduRGA/D0T/sTNHVXoiYrJFRES5qG6qB+/xLTGhXQ2ovb9rKJML+OXsC/T7/ToCI5Mgi49H6KLFAABBIkHYsmV4M2YMpJGRKoycqPRhskVERLkSq6thRuda2D+mBayNdRTt91/H4ov1l7HfNxaWa9ZAvUpWkdSkfy/jVa/eSLx0SRUhE5VKTLaIiOijmtib4O/JrdGnkbWiLVkiw+wjjzHZXwv6uw9Ar21bxTZZVBTejBmL0CVLIU9LU0XIRKUKky0iIvokfW0xfvZwgdfgRjD+YH3F88/C8cWOx/CZ+CPM5s2FSEtLsS1m1y4E9u2H1OcvVBEyUanBZIuIiPKtaz0LnJ7SBp/VMlW0xSSnY+yue1imWReme/ZBy9FRsS3t5UsE9uuH6J27VBEuUamgoeoAiurGjRvYtm0bHj58iOTkZFhaWqJLly4YPXo0dHV1CzTWzJkzceTIkY/22bx5M9qwoB8RVWBVDbSxbXgT7L75GktP+iIlXQYAOHT3La776+DnVb/D/vBWxPyxE0DG5HlJUJAqQyZSqTKdbO3cuRNLly6FIAgwNzeHhYUF/Pz84OXlhX/++Qd79uyBkZFRgce1sLCAhYVFrtsMDQ1zbSciqkgyS0S0rF4ZUw88xMM3sQCA4NgUDNx+F6Nbf4kxLd0QOXcONIyNUHXGdNUGTKRCZTbZevLkCZYtWwYAWLRoETw8PCASiRAWFoZx48bBx8cH8+bNg6enZ4HH7tOnDyZOnFjcIRMRlTvVTPVweGwLbLzoj3XnXkImFyAIwO//vsIlc32s3LwHtbSlUNPWzrafIJdDSEuDmo5OHiMTlR9lds7Wxo0bIZfL0atXL/Tv31+xfISZmRnWrFkDNTU1/PPPP3j27JmKIyUiKt801NUwqUNNeI9riWqmlRTtz0IT0GfXY2wOkGYrhApkTJ4PcO+DlCc+JR0uUYkrk8lWUlISLr9fdd7DwyPHdnt7ezRv3hwAcOrUqRKNjYioonKxMcLJia3xVQs7RVu6TMDqf16g72/X8SoiEQCQ+vwFwlf/DElAAAIHDEDk75sgyGSqCptI6crkbURfX19IJBJoamqifv36ufZxdXXFtWvX8PDhwwKPf/PmTbx8+RKxsbEwMDCAs7MzevbsCSsrq6KGTkRUruloqmNhr7r43NkcMw4+REhcKgDgwZuMQqgzuzihe8BtCOnpGTtIpYj45RckXv4Xlit+gqY1P2ep/CmTV7YCAgIAAJaWlhCLxbn2sbW1zda3IG7fvo3Tp0/j5s2bOHPmDNauXYvOnTtj8+bNhQ+aiKgCaVWjCk5NaQP3RlnJU2q6HAuOP8Xk1OoQe/0PGh88iJRy5y4CevVC7J9/cn1FKnfK5JWtuLg4AB9/MjBzW2bf/LCzs8PMmTPRvHlzWFlZQVNTE8+fP8fWrVtx6tQprF69Grq6uhg8ePAnx9q3bx8OHDiQr+P6+/vnO0YiorLCUEeMNR4N0NnZHLO9HyMqSQIAuOoXhS/faGDe4t/R8sjvSDh5EgAgT0pCyMxZSLx4CRYLfoR6IZ4mJyqNymSylfZ++Ye8rmoBgKamZra++TFu3LgcbS4uLli3bh0WLlyIPXv2YO3atejduzcqVaqUywhZIiIi4OPDiZ9ERJ2dzeFqZ4w5Rx7jtE8YACAhTYrvT7xAx/oDMKvlZ5CuWAR5QkLGtlOnkHLvHiyWL4Neq1aqDJ2oWJTJ24ha75eDSM+8558LiUSSrW9RTZs2DWKxGPHx8bhx48Yn+5uamsLZ2Tlf/2n/55FoIqLypoqeFn4b4oqf+7lAXyvr7/yzvmHo46MJ31XboNukiaJdGh6ON9+MROKVq6oIl6hYlckrW/m5RZifW40Foa+vj5o1a+Lp06cIykcl5AEDBmDAgAH5Gtvd3Z1XwYio3BOJROjjao0W1Svj+0OPcMUvEgAQm5yOKacC8UXnyZjm9hDpG9YC6enQadgQlVo0V23QRMWgTF7Zsre3BwC8e/cuz6tbr1+/zta3OGTetpRKpcU2JhFRRWNppIOd3zTF4t51oauprmj/60koBkTYwG/1Vug0aADLlT9BpK7+kZGIyoYymWzVrl0bYrEYEokEjx49yrXP3bt3AQANGjQolmNKpVK8evUKAGBubl4sYxIRVVQikQhDm9vh1OQ2aOpgomiPTJRg4oUw/NL7eyRXNsuxX8yBA1xnkcqcMpls6enpwc3NDQByfeIvMDBQMa+qS5cuxXLM/fv3IyEhARoaGoqCqUREVDS2lXWxb1RzzO9eB9rirK+kow/e4fNf/sX5Z2GKtuS7dxG6YCFe9f4SMfv2s0QElRllMtkCgPHjx0MkEuHo0aPYvz/rTRceHo5p06ZBLpejY8eOcHJyyrZf+/bt0b59+xyV5a9evYpVq1YhMDAwW7tEIsHOnTuxfPlyABlzsapWraq8F0ZEVMGoqYkwws0Bf01qjUa2Ror28IQ0jNh+BzMOPkRMbCLe/TATkMshpKQgdMECvBk7Funh4aoLnCifyuQEeQCoX78+Zs6ciRUrVmD+/Pnw8vKCsbEx/Pz8IJFI4ODggMWLF+fYLzg4GACQnJycrT0lJQVbtmzBli1bUKVKFZiZZVy+DggIUPTt3LkzfvjhByW/MiKiiqmaqR4Ojm2JLZdf4eczLyCRZqyneOjuW1x+GYG5w6ehludCyN4/AJV06V8E9OgJ84ULYFBMdzGIlKHMXtkCgOHDh2Pbtm1o06YNUlJS4OfnB0tLS4wdOxaHDx+GiYnJpwd5z9nZGePHj0fLli2hpaWFgIAAvHjxAvr6+vj888/h5eWF9evXK+p3ERFR8VNXE2FM2+o4OdEN9a2zniYPi0/DxCdybBy/FrK27RXtsrg4BE+ZiuAZ3ymSMKLSRiTwprfKZZZ+cHZ2hre3t6rDISIqFaQyOTZdfoW1Z15CIpMr2k31tTDbNA61f1sG4YO7FBpVq8Ji6VLotXZTRbhUAeX3+7tMX9kiIqLyS0NdDeM/q4GTk9zgYmOkaI9ISMPUV9rY+O16SBpnPbAkDQ/Hm1GjELJgQdZC10SlAJMtIiIq1Wqa6ePw2BaY1dUJmhpZX1vHXsbiG6dBeDJ2NkQfLN8mjYwENMrslGQqh5hsERFRqaehroYxbavjr0mt0fCDJxYjEyX4LtQE6771RLJzA6gbG8Ni4UKIRCLVBUv0H0y2iIiozKhRVQ+HxrbE3G61ofXBVa6/AxMx0uVrPJy3Huq5PBwlefOmJMMkyobJFhERlSnqaiKMbF0Nf09ujcZ2xor2mOR0fPdvGL7ZcQchcSmK9oSzZ+HfpSvC166FIJGoImSq4JhsERFRmVTNVA8HxrTAgh51sq2xeP5ZOD5f8y/23HwNSVQUQn5cAMhkiPrtdwR49Efqs2eqC5oqJCZbRERUZqmpiTC8lQNOT2mD1jWrKNoT0qSYfeQxhux8iFAjC0V72rNnCOjbDxG//sonFqnEMNkiIqIyz8ZEF3+MaIqVfevDQDvrScRboSkY3XAE/u43GTK8nzQvlSLScwMC+vdH6vMXKoqYKhImW0REVC6IRCJ4NLbB2Wlt8XkdM0V7qlSO9ek2mDl8DV5Xr6doT3vqi4C+fRH522+8ykVKxWSLiIjKlaoG2vh9qCt+HdQIVfSyllh7EivDty7Dsdd9GtLU3l/9Sk9HxNp1769yPVdRxFTeMdkiIqJyRyQSoVt9C5yZ2hZfNrRStEvlAv6QW2LywJXwcWyiaE976ovYg4dUESpVAEy2iIio3DKupIlf+jfAjhFNYW2so2gPSpJjRp3+8PryOySIdSC2tkbVqVNUFyiVa0y2iIio3GvraIp/prbBSDcHqH1QXP6YYIZx7kvweMKPEOnqZttHkMtZl4uKBZMtIiKqEHQ1NTC3ex38+W0r1LYwULRHpQmYeiMOo/64m60YaszevQjo0xcpjx+rIlwqR5hsERFRhVLf2gjHJrTCD12csi35c9Y3DB1/voRtVwOQ8voNwn9eg7SXLxHYfwDCV6+GPDVVhVFTWcZki4iIKhyxuhrGfVYdp6e0QYtqlRXtSRIZFh5/Cvdt9/FC/H4pILkcUVv+h4BevZF8+7aKIqayjMkWERFVWPZVKmHPqGZY2bc+jHTFinbfFHVMaTcVm1sORrKGFgBAEhSEoKHDELpoEWSJiaoKmcogJltERFShZRZDPTetLdwbZZWJkAPwrtoQ43ouwDXLrGKoMXv24lWPnki4eLHkg6UyickWERERgMp6Wljj0QB7RjaDQ5VKivZwuRiLm36FJe2/RYSOIQBAGhKCt2PHIXjGd5BGR6sqZCojmGwRERF9oGWNKvh7cmtM6lATmupZX5NXDRwwpssceDu2g0yU0R5/4gTijh5TVahURjDZIiIi+g9tsTqmdXLEX5Nbo5mDiaI9RVDD5jrdMLnLTPiY2EOrTm2YDB2iwkipLGCyRURElIcaVfWwb3RzrOpbH8YfTKD31zLBjDYTsKHrZMSkybPtI5dIIEilJR0qlWJMtoiIiD5CJBKhX2MbnJv+Gfo3tsm2zds/Ee1/voi9t15DLhcAAJGengj06I8UHx9VhEulEJMtIiKifDCppImf+tbH4XEt4GSur2iPTU7HLO/HcPe6hntXHyJq6zakPn2KwH4eCFvxE+RJSSqMmkoDJltEREQF4GpnghMT3TCvex1U0lRXtD94E4u+x9/gtzrdkaShDcjliN6+Hf49eiDhwgUVRkyqxmSLiIiogDTU1fCNmwPOTf8M3etbKNrlEOFo9dYY1XkWztg0hhwiSN+F4O248Xg7aTLSw8JUGDWpCpMtIiKiQjI31MaGQY2w85um2WpzxYgrYY3rAMxoOwF+hhmFUhP++QevvuiG6D92cgJ9BcNki4iIqIha1zTFqSmtMb2TI7TFWV+tvsZ2mPTZFHi69EG8WBfypCSELVuWMYH+CSfQVxRMtoiIiIqBloY6JnaoibPT2qJrXXNFuyAS4S+HFhj5+SycsG8BGURIffoUkqBA1QVLJYrJFhERUTGyNtaF1xBX7PqmGaqbZt1aTBDr4NcGfTDpsyl41bY7DL74QoVRUkliskVERKQEbjWr4O/JbTDni9rQ09JQtL8yssK3xp9h6v4HCIlLUbRLIyOR+uKFCiIlZWOyRUREpCSaGmoY1aYazk9vC/eGVtm2/fngHdqvvoR1Z18iRSJD2MqVCHDvg7AVP0GWmKiiiEkZmGwREREpWVUDbazp3wCHxrZAHQsDRXtKugy/nH2Bdsv/wbG7byBIpRm1ubp2RdyxYxAEQYVRU3FhskVERFRCGtub4PhENyx3r4fKlTQV7aEpcvzUZAhmtP4WL4ysIYuIxLvvf0DQ0KFIffZMhRFTcWCyRUREVILU1UQY2NQWF777DKPbVINYXaTY9rSyAyZ/NgVrGvZHtJY+Uu7cRYB7H4QuWgRpTIwKo6aiYLJFRESkAgbaYsz+ojb+mdoWneqYZdt2xq4Jvuk0E3sdOyBVpI6YPXvxqktXxOzdC0EmU1HEVFhMtoiIiFTIoUolbB7WGLtHNkMts6wFrlM1tPBHna4Y2fEHnLFpDElcPEIXLsK7WbNUGC0VBpMtIiKiUqBVjSo4OckNi3vXhbGuWNEepWOENa4DMOmzKbhvWhPG/fqpMEoqDCZbREREpYSGuhqGNrfDxe/aYUybatBUz/qafmVkhdmtxmD8EwHPQxOy7SdIJCUdKhUAky0iIqJSxlBHjFlf1Ma56W3R08Uy27aLzyPQdd2/mHn4EcLjU5H64gX8OnTMmM/FBa5LJSZbREREpZSNiS7WD2yIo9+2QlMHE0W7XAD23X6Dz1ZfxIoNxxEXkzGf61Wv3ki8dIn1uUoZJltERESlnIuNEfaPbo5NQ11RrUrWeovJEhm26zphRKfZOFK9NRIDgvBmzFi8HjECKT4+KoyYPsRki4iIqAwQiUT43Nkcp6e2waJezjD5oChqvFYlbKrXCyM7/oB/bJsg4cYtBPbpi+DvvofkbbAKoyaAyRYREVGZIlZXw7AW9rj03WeY1KEmdDXVFdsidI3xS6P+GNduOq5a1EXc8eN41bUrwlb8xKKoKsRki4iIqAzS1xZjWidH/Pt9OwxvaZ+tEv0bAzMsaTYcU9tMwgND24z1Fjt2Qvq7dyqMuOJiskVERFSGVdHTwoKezjg//TO4N7KCKCvnwnMTW8x0G4fZLUfDz/UziC0t8x6IlIbJFhERUTlgY6KLNR4NcGpyG3SsnX35n/tVHfGtSTsM23oL919n3U6URkWxRlcJYLJFRERUjtQy18eWrxrj8LgW2cpFAMC/LyLw5cZrGL7tFh6+iUXIjz/Cv0tXxB46BCE9XUURl39MtoiIiMohVzsT7B/dHLu+aQZXO+Ns2y4+j0CvX69ieqIdfJNECJk7D/5fdMtIunilq9gx2SIiIiqnRCIR3GpWwaGxLbBjRFM0sDHKtv2WeR1MbDcVi5oNh2+8HCFz58GvSxfE7N0LOZOuYsNki4iIqJwTiURo62iKI+NbYtvwJqhvbZht+3WLupjYbirmtRiJB2naCF24CP6dPkfUtu2QJSapKOryg8kWERFRBSESidDOqSqOftsK//uqMZwtDbJtv2PmhO9af4sZbuNxXTBG2E8/wa99e4SvWwdBLldR1GWfhqoDICIiopIlEonQobYZ2jtVxZmnYfA874fHwXGK7T5VqmF+lWqoFhuM/i/O43MfH1RV4/WZwuKZIyIiqqAylwA6NqEVdn7TFM2rZX968ZWRFZY3HYrhtu7Yf/s1JNKsq1uy2FgueJ1PTLaIiIgqOJFIhNY1TbFvdAscHtcCHZyqZtselCjDD4cfw+2n89hw/iWikyQInj4Dr3r0QPTu3ZAlJqoo8rKByRYREREpuNqZ4H/Dm+Dvya3R08USah9UpA9PSMPqf16gxbKzWJ5oDr/QBIQtXoKXbdoiZP6PSPX1VV3gpRiTLSIiIsqhtoUB1g9siPPTP8PApjbQVM9KGdJkAv52aIExHb/HvBYjcaeSNWIOHEDAl+4IcO+D6J27uPD1BzhBnoiIiPJkX6USlrvXx7ROtbDrRhB23QhCVFJWDa47Zk64Y+YEu/hQ9Pb/F+2e3UPq0qUIW7kS+p+1hUHPntBv2xYiTU0VvgrV4pUtIiIi+iRTfS1M7eSIqzPbY2Wf+qhlpp9te5CBOdY19MCQLvPxe92eCNI2QcKZswieNh3y1FQVRV068MoWERER5Zu2WB0eTWzQr7E1rvpF4X9XXuHC8wjF9kRNXfxZow3+rNEGzpGv8KVOLBx0KkH9gzGkERGIP3MGeq1bQ9PGpuRfRAljskVEREQFlrkUkFvNKvALT8S2qwHwvheMlHSZoo9PlWrwAeC57BzcG1lhUFNb1DTTR+KlSwhbtBhhAMR2tqjUtCl0XFygXb8+tKpXh0hdPc/jlkVMtoiIiKhIalTVw9Iv62FmVyf8+eAd9tx8Dd+QeMX2uJR0bLsaiG1XA9HE3hgdXryCq1gbeumpSA96jdig14g9eAgAoFapErRqO0HT1g6atjbQtLWFfteuEImyHosUBCHbz6Udky0iIiIqFvraYgxtbochzWzx6G0c9t56jWMP3yFZknW163ZgDG5rOkPcdSEah/mi7dv7aBb6FNqydACAPCkJKXfuIuXOXQCAepUqMPjii2zHkfj5Ifj7H2DzmxfEZmYl9wILqcwnWzdu3MC2bdvw8OFDJCcnw9LSEl26dMHo0aOhq6tbqDFPnz6NXbt24dmzZ0hPT4ednR169uyJYcOGQSwWF/MrICIiKl9EIhFcbIzgYmOEOd1q4+j7q11PP7jala6mjusWdXHdoi50IEPLuAC0fnYZjUKfQSxkJWe5zemSp6QgzdcXIbNmwWbLFohK+VJCZTrZ2rlzJ5YuXQpBEGBubg4LCwv4+fnBy8sL//zzD/bs2QMjI6MCjfnTTz9h69atAABbW1vo6Ojg5cuXWLlyJS5cuICtW7dCswI/vkpERFQQ+tpiDGluh8Hvr3YdvPsGfz0ORfQH5SNSoI5zhjVwrlkNGGqK0EE/DW3SQlAv7Dl0a1TLNp40JgaSoNcAgKRr1xGzdy9MBg8u0ddUUKU7FfyIJ0+eYNmyZQCARYsW4eLFizhy5AjOnj0LZ2dn+Pv7Y968eQUa88yZM4pkauPGjThz5gyOHTuG48ePw9raGrdv38aaNWuU8XKIiIjKtcyrXUt618PN2R2w/esmcG9kBT2t7Nd94iQCvKM0MSXRDl0Nu2C6XnNs/vcVnocmQBAEqBsZIf7UKUX/8FWrIQkMLOFXUzBlNtnauHEj5HI5evXqhf79+ysmypmZmWHNmjVQU1PDP//8g2fPnuV7zA0bNgAARo0ahQ4dOijaq1evjiVLlgAAdu/ejejo6GJ8JURERBWLWF0Nn9WqijUeDXBnbkd4DW6ErnXNoamRPS2RSOW4/DISS//yRee1/6L58nP47tAj3B04EUmmlgAAITUV72bOgiCT5XaoUqFM3kZMSkrC5cuXAQAeHh45ttvb26N58+a4du0aTp06BScnp0+OGRgYqEjM+vfvn2N7ixYtYGdnh6CgIJw7dw79+vUr4qsgIiIibbE6utazQNd6FkhITcc/PmH4+0korvtHIkmSPYEKi0/DobtvcQiAqNU0WCaEwyYxHNYJEaj9y1406PcFqpvqwUi3dE33KZPJlq+vLyQSCTQ1NVG/fv1c+7i6uuLatWt4+PBhvsZ88OABAMDGxgZmeTzZ4OrqiqCgIDx8+JDJFhERUTHT1xajj6s1+rhaQyKV497rGFx+GYF/X0TicXBctr4CgGD9qgjWrwpYAIgE4HUdAFC5kiaqm+qhetVK6OlihRbVK5f4a/lQmUy2AgICAACWlpZ5Ph1oa2ubre+nBL6/35u5X1HH3LdvHw4cOJCvY/v7++erHxERUUWhqaGG5tUqo3m1yviuMxCVmIYrfpH490UkLr+MQHhCWp77RiVJEJUUjVuB0Thw5y0uffcZrI0LV6GgOJTJZCsuLiO7NTQ0zLNP5rbMvsU5Znx8fJ59MkVERMDHxydfxyYiIqKPq6ynhV4NrNCrgRUEQYB/RBL8whPhe8cHL/1D8K6qHfwjk5GQJs22nyAIKoo4S5lMttLSMrLZj9W8yizPkNm3OMdMzceCmqampnB2ds7Xsf39/fM1JhEREWU82Vijqh5qVNVDl7rminZBEBCRmAb/8CT4RyQiJC4FLapVUelVLaCMJltaWloAgPT09Dz7SCSSbH2Lc0xtbe1PjjdgwAAMGDAgX8d2d3fnVTAiIqIiEolEqKqvjar62iqfp/WhMln6IT+3CPNzW/BDBgYG+R4zsy8RERHRp5TJZMve3h4A8O7duzyvRL1+/Tpb309xcHAAAAQFBeXZp6BjEhEREZXJZKt27doQi8WQSCR49OhRrn3u3s1YwLJBgwb5GtPFxQUA8PbtW4SFhRXLmERERERlMtnS09ODm5sbAORaXiEwMBA3btwAAHTp0iVfYzo4OMDR0REAsH///hzbr1+/jqCgIIjF4mzV5YmIiIg+pkwmWwAwfvx4iEQiHD16FPv371c82hkeHo5p06ZBLpejY8eOOarHt2/fHu3bt8epD9ZVyjRhwgQAwObNm3H+/HlF+6tXrzB37lwAwKBBg2BiYqKsl0VERETlTJl8GhEA6tevj5kzZ2LFihWYP38+vLy8YGxsDD8/P0gkEjg4OGDx4sU59gsODgYAJCcn59jWuXNnfPXVV9ixYwfGjRsHW1tb6Orq4uXLl5DJZHB1dcX06dOV/tqIiIio/CizyRYADB8+HLVq1cLWrVvx6NEjREVFwdLSEl26dMHo0aNRqVKlAo85e/ZsNGzYEHv27IGvry/Cw8NRvXp19OzZE8OHD/9oHS4iIiKi/yrTyRaQsUB0ixYt8t3/+fPnn+zTtWtXdO3atShhEREREQEow3O2iIiIiMoCJltERERESsRki4iIiEiJmGwRERERKVGZnyBfHrx9+xYA4O/vD3d3dxVHQ0RERPnh7+8PIOt7PC9MtkqBtLQ0AEBqaip8fHxUHA0REREVROb3eF6YbJUCJiYmiI6OhpaWFqytrQs1hr+/P1JTU6GtrY3q1asXc4QVE8+pcvC8Fj+eU+XgeS1+5e2cvn37FmlpaZ9cWUYkZK5zQ2Wau7s7fHx84OzsDG9vb1WHUy7wnCoHz2vx4zlVDp7X4ldRzyknyBMREREpEZMtIiIiIiViskVERESkREy2iIiIiJSIyRYRERGREjHZIiIiIlIiJltERERESsRki4iIiEiJmGwRERERKRGTLSIiIiIl4tqI5YSHhwciIiJgamqq6lDKDZ5T5eB5LX48p8rB81r8Kuo55dqIRERERErE24hERERESsRki4iIiEiJmGwRERERKRGTLSIiIiIl4tOIJUgQBNy/fx/nz5/H3bt38erVKyQmJkJfXx916tRB79690aNHD4hEolz3T0pKwqZNm3D69Gm8e/cOurq6cHFxwYgRI9CsWbNCxxUVFQUvLy9cuHAB4eHhMDAwQJMmTTBmzBjUrl270OOWhMKe05iYGJw9exbXrl2Dj48PQkJCoKamBgsLC7i5uWH48OGwtrYucDxv375Fhw4dPtrHxcUFBw4cKPDYJako/1Zr1ar10bGrVKmCq1evFiouZb0HSkJhz+nMmTNx5MiRfB3j/PnzsLKyyndMyvxdlaS///5b8V4ODw9HbGwsxGIx7O3t0bZtW3z11VcwNjbOdV9+ruauMOeUn6t549OIJej69esYPny44mcbGxsYGBggODgYsbGxAIDPPvsMnp6e0NTUzLZvdHQ0Bg0ahICAAGhqaqJGjRqIjo5GaGgoRCIR5s2bh8GDBxc4pqCgIAwaNAiRkZHQ1dWFg4MDQkNDERUVBbFYjHXr1n3yH7kqFfacDho0CHfv3gUAVKpUCba2tkhNTcWbN28glUqhq6uLdevWoU2bNgWK58MPhUaNGuXap2bNmli0aFGBxi1pRfm3mvkFXrdu3RzbAMDIyAheXl4FjklZ74GSUthz+ttvv+HSpUt5jhsYGIjo6GhYWFjg/PnzUFPL/w0LZf2uSlqvXr3w7NkzaGpqwtTUFMbGxoiOjsa7d+8AAJUrV8bWrVvh5OSUbT9+ruatMOeUn6sfIVCJuXr1qtC+fXthx44dQmRkZLZtR44cEerWrSs4OjoKK1euzLHv2LFjBUdHR+HLL78UQkNDBUEQBLlcLuzbt09wdHQUateuLTx9+rRA8cjlcqFXr16Co6Oj8M033wjx8fGCIAhCenq6sHbtWsHR0VFo0KCBEBYWVshXrHyFPadDhgwRpk2bJty6dUuQSqWK9pCQEOHrr78WHB0dhUaNGgkREREFiufNmzeCo6Oj4OjoWPgXVQoU5d9q5ut/8+ZNscakjPdASSrKOc2LXC4XOnToIDg6Ogpr1qwpcEzK+l2VtP379wu3bt0SJBJJtvZnz54J3bt3FxwdHYUvvvgix378XM1bYc4pP1fzxmSrBCUkJOT4h/shLy8vwdHRUWjatKkgk8kU7T4+PoKjo6Pg5OQkBAYG5tjvu+++ExwdHYUJEyYUKJ4zZ84Ijo6OgqurqxAbG5tj++DBgwVHR0dh+fLlBRq3JBX2nEZHR+e5T3x8vNC8eXPB0dFR2LZtW4HiKQ8fCoJQ+PMqCMr5AlfWe6AkFeWc5uXmzZuK853befmU8pJsfczDhw8Vr9PPz0/Rzs/VwsvrnPJzNW+cIF+C9PT0IBaL89yeeWk1NjYW0dHRivbTp08DAJo3bw47O7sc+/Xv3x8AcOnSJSQnJ+c7nr///hsA0KVLFxgaGubY7uHhka1faVTYc5rX/A0A0NfXR4MGDQAAAQEBxRNoGVPY86osynoPlCRlnNPMuVyNGjXK9bwQUK1aNcX/T0lJUfx/fq4WXl7nlJ+reeME+VIkNTVV8f+1tbUV///BgwcAgMaNG+e6X/369aGpqYm0tDT4+vrC1dU1X8d7+PDhR8fNbA8NDUVYWBjMzMzyNW5pktc5/ZS0tDQAgI6OTqGPvWTJErx69QoikQhWVlZwc3NDx44dCzSnprTKz3nduHEjwsPDIZPJYGZmhubNm+OLL77IdW7QpyjrPVCaFPTfanJyMk6dOgUAcHd3L9Kxi/N3VdpkziHKnDuViZ+rhZfXOf2Uivy5ymSrFDl58iQAwMnJCXp6eor2wMBAAICtrW2u+4nFYlhYWCAoKAgBAQH5+lCQSCQIDg7+6LgWFhYQi8VIT0/Hq1evyuSHQl7n9GPCwsJw69YtAHl/YObHzp07s/28f/9+1K5dG56enrCxsSn0uKVBfs7r4cOHs/185MgRrF+/Hp6ennB2di7Q8ZTxHihtCvpv9fTp00hOToaOjg66du1apGMX5++qNJDL5YiIiMDVq1exevVqAMCMGTNQqVIlRR9+rhZMfs7px1T0z1UmW6XEkydPsG/fPgDA6NGjs22Li4sDgFwvSWfK3BYfH5+v4yUmJkIul390XJFIBAMDA0RFReV73NLkY+f0YxYvXoz09HTUqFED7dq1K9AxNTQ00LNnT3Tr1g01atRA1apVERMTg0uXLmHt2rXw9fXFN998A29v73wnf6XNp85rhw4d0KtXLzg5OcHc3BxJSUm4fv06fvnlF7x58wYjRozAn3/+CQsLi3wfUxnvgdKkMP9Wvb29AQCdOnUq9L8lZfyuVGn79u1Yvnx5trb69etjxYoVOZ6A4+dq/hTknH5MRf9cLd3X3SqIyMhITJw4EVKpFJ06dUK3bt2ybc+89Pqx+R6Zl/s/vBXxMZljfrhvcYxbWnzqnOZl06ZNOHPmDMRiMVasWAF1dfUCHdfc3ByrVq3CZ599Bmtra2hqasLMzAweHh7Yu3cv9PT0EBQUhD/++KMwL0vl8nNeN27ciM6dO8POzg5aWlowMTFBt27dcODAAVhaWiI2NhYbNmwo0HGV8R4oLQrzb/XNmze4ffs2gKLdQlTG70qVzMzM0KhRI7i4uMDU1BQikQi+vr44evRojsSGn6v5U5Bzmhd+rjLZUrmEhASMGjUK7969g7OzM1asWJGjj5aWFgAgPT09z3EkEgmA/M9Lyhzzw32LY9zSID/nNDdHjhzBmjVrIBKJsHTpUtSrV69Y47Kzs8PAgQMBAGfOnCnWsUtCYc9rJhMTE8VVm7Nnz0IoQIk/ZbwHSoPCntM///wTgiDAysoKzZs3L/a4ivK7UqWuXbti7969OHDgAK5cuYI///wTLi4uOHHiBIYNGwaZTKboy8/V/CnIOc0NP1czMNlSoaSkJIwcORJPnz5FzZo18b///S/XS6AGBgYAsi575yZzW2bfT9HT01NMKMxrXEEQFH+55HdcVcvvOf2vv/76C3PmzIEgCPjxxx/Rq1cvpcTXsGFDAFnzRcqKwp7X/8p8/bGxsYpCnvmhjPeAqhX2nAqCgD///BNARuHJvFacKKrC/q5KEycnJ/z+++8wNjaGr6+vYl4cwM/VwvrYOf0vfq5mYbKlIikpKRgzZgwePHgAe3t7bNu2Lc/HZu3t7QFkVCXOTXp6uqKqb2bfT9HU1ISlpSUA4PXr17n2CQkJUfzVV5AnTlSlIOf0Q2fOnMF3330HmUyGH374QfFXkjJk3rL41F+DpUlhz2tuPrxlU5BzoIz3gCoV5ZzeunULb9++hUgkKvJTiB9T2N9VaaOnp4emTZsCAHx8fBTt/FwtvLzO6Yf4uZodky0VSEtLw7hx43D79m1YWVlh+/btMDU1zbN/Zm2SzMdt/+vRo0dIT0+HlpZWgdbcyhz3zp07uW7PbDc3N4e5uXm+x1WFgp7TTJcuXcLUqVMhlUoxadIkjBgxQqlxvnz5EgBK/fnMVNjzmpfM16+lpQUjI6N876es94AqFPWcZtbWaty4sVKfvirs76o0kkqlALJ/GfNztWhyO6eZ+LmaE5OtEpaeno6JEyfi+vXrMDMzw44dOz75pE/nzp0BADdv3sz1r7D9+/cDyCiKmN/HcD8c99SpU7le8s5c1LNLly75HlMVCnNOgYy16iZOnIj09HSMGTMG3377rVLjTEpKwp49ewAArVq1UuqxikNhz2tepFIptm3bBiCjkKSGRv4fhlbWe6CkFfWcJiUlKYpxfvnll8oKs0i/q9ImNjZWUXLgw6SJn6uFl9c5Bfi5mifVFa+veKRSqTBx4kTB0dFRaNWqlfDq1at87ztmzBjFGl6Za2p9uIaXk5OT8OTJkxz7bdu2TWjXrp0wYMCAHNtkMpnQo0cPwdHRURg5cqRiDS+pVKpYw8vFxUWxZlhpVNhzeu/ePaFBgwaCo6OjsGzZsgId82PndO7cucLp06eFtLS0bO1+fn7CgAEDFOuilfblUQp7XletWiV4e3sLCQkJ2drfvXsnjBs3TnB0dBTq1KkjPHjwIMe+f//9t9CuXTuhXbt2uY5d2PdAaVGU93+mw4cPK/4NJSYmfrL/x85pUX5XpcnNmzeFX3/9Ndf31JMnTwR3d3fB0dFRaN26dY5zxs/V3BX2nPJzNW9l98+VMujvv/9W/FWqqamJ2bNn59l33rx5qFOnjuLnZcuWYeDAgfDx8UGHDh1Qo0YNxMTEICQkBCKRCLNnz861+GBCQoKiyN5/qampYd26dRg8eDD+/fdftGnTJsfq9KtWrSrVRfcKe05nzpyJ5ORkqKur49GjR3nOJ6hTpw7mzZuXre1j5/TRo0c4cOAAxGIxbG1toaenh5iYGMX8DUNDQ6xduxbW1tYFfq0lqbDn9dWrV9i8eTPmzJkDGxsbGBoaIiEhAQEBARAEAVpaWliyZAlcXFxyjJOcnJzneQUK/x4oLYry/s+UeQuxc+fO+bra8rFzWpTfVWkSHx+PdevWYd26dTA1NUXVqlWhrq6OkJAQREREAMgoX/D777/nOGf8XM1dYc8pP1fzxmSrBH34KHBwcPBHv1gSEhKy/WxiYoLDhw9j8+bNOHXqFPz8/KCrq4s2bdrgm2++KfTj3w4ODjh27Bi8vLxw4cIFvHjxAgYGBujcuTPGjh2b6wd+aVLYc5o5QVUmk+HevXt57lPQ2ydjxozB5cuX8eTJE0RGRiIoKAja2tpwdnZGmzZtMHjw4CLNeSophT2vAwcORJUqVfDkyROEh4cjODgYYrEYNWvWRIsWLTBkyJA8K2t/irLeAyWlKO9/IHttreK4hajM31VJatiwIWbNmoWbN2/Cz88PgYGBkEgkMDAwQLNmzdC+fXv07ds31yc9+bmau8KeU36u5k0kCGWkgAoRERFRGcQJ8kRERERKxGSLiIiISImYbBEREREpEZMtIiIiIiViskVERESkREy2iIiIiJSIyRYRERGREjHZIiIiIlIiJltERERESsRki4iIiEiJmGwRERERKRGTLSIiIiIlYrJFRKREv/zyC2rVqoVNmzYV+9jffPMNatWqhevXrxf72ERUfESCIAiqDoKIqCTNmjUL3t7ehdq3fv36OHjwYL76hoaGonPnztDV1cW5c+egq6ur2Obp6YkNGzYAAJ4/f57nGIcPH8a8efMgk8nQuHFj/Pbbb9DX1wcAPHjwAP3790edOnVw+PBhqKnx72ei0khD1QEQEZW0t2/fokqVKjnak5OTkZycDAC5bgeApk2b5vs4v/zyC1JTUzFp0qRsiVZ+bd++HStWrIAgCGjXrh3Wrl0LbW1txfYGDRrAzc0NV65cwbFjx9C7d+8CH4OIlI/JFhFVODt37sy1fd68eThw4AAsLCxw8eLFIh0jLCwMx48fh1gsRp8+fQq8/9q1a+Hl5QUA6NGjB1asWAENjZwf2QMHDsSVK1ewZcsWJltEpRSvORMRvefj4wMAqFOnTpHHOnDgAGQyGdq2bQsjI6N87ycIAhYtWqRItIYMGYJVq1blmmgBQJs2bWBkZISXL1/i7t27RY6biIofky0iIgBSqRQvX74EUPRkSxAEHDp0CADQvXv3AsXw3XffYffu3QCAb7/9FvPmzYNIJMpzH01NTXz++ecAMhI8Iip9mGwREQHw8/ODRCIBUPRk68WLFwgNDQUANG7cOF/7pKWlYcKECTh+/DhEIhHmzJmDSZMm5WvfzGNcuXKlcAETkVIx2SIiQtYtRABwdnYu0lh37twBAFhYWMDU1PST/RMTEzFy5EhcuHABGhoaWLFiBYYNG5bv47m4uAAAIiMj4e/vX7igiUhpmGwREQHw9fUFAFSuXBlmZmZFGuvhw4cAACcnp3z1HzZsGG7dugUtLS14enoWeKK7vb294mnHBw8eFGhfIlI+JltERMi6slW7du0ijxUeHg4AMDY2LtCx+/bti/bt2xfqmJnHyjw2EZUeTLaIqMKTy+V49uwZgKLfQgSA6OhoAMj3U4gNGzYEAOzevRt//PFHoY5paGiY7dhEVHow2SKiCi8gIEBRzLQ4yj6kpaUByHhSMD+2bNmCRo0a/b+9u1dprQnDMHwrQVd+PAKt1NomiIiCpYWdjZV4AornIEhEFBGDhmVEi3UGtmLhCQgiYm3soyA2QtzFJtkI+yeJa/zky33VM/NOFR5m3swCYHNzk7Ozs45rDg4Ofqgt6fswbEnqec1+LUjnZKt5pff8/NzW+EKhQLVapVgsAlAqlTg9Pe2oZrNWu1eXkr6OYUtSz2v2TA0NDTEyMvLp9ToNWwD5fJ5qtcrk5CQAW1tbnJyctD3fsCV9X4YtST3v7u4O+Nkc/7cHRNs1Pj4OQK1W62heLpcjjuPW9xe3t7c5Pj7+57yXlxfq9ToAY2NjHe5WUmiGLUk9L83mePj1yOj9/X3rodR2NQPX1NQUADs7O8Rx/Nc5t7e3NBoNMplMq/dL0vdh2JLU0x4fH3l6egLSefYBoFgskslkeHt7+9AP1q5sNkscx0xPTwOwu7tLpVL54/ibmxvgZ3N/Pp/vbtOSgjFsSeppzStESO9kq1AoMDc3B8Dl5WVXa0RRRKVSYWZmBoC9vT0ODw9/O7ZZo5PvMEr6OoYtST2tGbay2Syjo6Oprbu0tATA+fk57+/vXa0RRRFHR0fMzs4CsL+/T7lc/jCmVqtxfX1NFEUdvzwv6Wv0vXf7KyBJ+qNGo8H8/DwPDw8kSdL6l2HayuUyBwcHLC4uUiqVgtSQ9DmebElSAP39/ayvrwP8s8G9W6+vryRJwsDAAKurq0FqSPo8w5YkBbKwsMDExARXV1etJvY0JUlCvV5neXmZ4eHh1NeXlI7Mf70BSfq/6uvrY2Njg4uLiyDfLMzlcqytrbGyspL62pLSY8+WJElSQF4jSpIkBWTYkiRJCsiwJUmSFJBhS5IkKSDDliRJUkCGLUmSpIAMW5IkSQEZtiRJkgIybEmSJAX0A4C39Q7+RCZZAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, ax = plt.subplots(1, 1)\n",
    "ax.plot(tr * Tc, s_corr, color='tab:red', label=corr_label, ls=\"--\")\n",
    "ax.plot(dia.liquid.temperature/KELVIN, surf_tens, color='tab:blue', label='DFT')\n",
    "xlim = fluid_param[substance][\"xlim\"]\n",
    "ylim = fluid_param[substance][\"ylim\"]\n",
    "ax.set_xlim(xlim)\n",
    "ax.set_ylim(ylim)\n",
    "ax.set_xlabel(r'$T$ (K)')\n",
    "ax.set_ylabel(r'$\\gamma$ (mN/m)')\n",
    "legend = ax.legend(loc='best', frameon=False, fontsize=12, ncol=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4a3385ce-d12a-4402-b139-10d852df2e2f",
   "metadata": {},
   "source": [
    "# Calculate mean absolute deviation (MAD)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "dc17b6c6-47e5-414e-a432-cec8d01b493a",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "MAD=4.4 %\n"
     ]
    }
   ],
   "source": [
    "tr_mad = dia.liquid.temperature / Tc / KELVIN\n",
    "tr_filer_mad = np.less_equal(tr_mad, 0.95)\n",
    "trx_mad = tr_mad[tr_filer_mad]\n",
    "if fluid_param[substance][\"param\"] is None:\n",
    "    s_corr_mad = surftens_mulero2012(substance, trx_mad)\n",
    "else:\n",
    "    print(\"Correlation\")\n",
    "    s_corr_mad = surftens(fluid_param[substance][\"param\"], trx_mad)\n",
    "mad = 100*np.sum(np.abs((surf_tens[tr_filer_mad]-s_corr_mad)/s_corr_mad))/s_corr_mad.shape[0]\n",
    "print(f\"MAD={mad:.1f} %\")"
   ]
  }
 ],
 "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.10.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}