notebooks/02.04-Mixture-Design-Cold-Weather-Fuel.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "fIsduu71gNfY",
    "pycharm": {}
   },
   "source": [
    "# Design of a Cold Weather Fuel for a Camping Stove"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Imports"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "colab": {},
    "colab_type": "code",
    "id": "nm3gVnuzgUzZ",
    "pycharm": {}
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "\n",
    "import shutil\n",
    "import sys\n",
    "import os.path\n",
    "\n",
    "if not shutil.which(\"pyomo\"):\n",
    "    !pip install -q pyomo\n",
    "    assert(shutil.which(\"pyomo\"))\n",
    "\n",
    "if not (shutil.which(\"cbc\") or os.path.isfile(\"cbc\")):\n",
    "    if \"google.colab\" in sys.modules:\n",
    "        !apt-get install -y -qq coinor-cbc\n",
    "    else:\n",
    "        try:\n",
    "            !conda install -c conda-forge coincbc \n",
    "        except:\n",
    "            pass\n",
    "\n",
    "assert(shutil.which(\"cbc\") or os.path.isfile(\"cbc\"))\n",
    "\n",
    "import pyomo.environ as pyomo"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "fIsduu71gNfY",
    "pycharm": {}
   },
   "source": [
    "## Problem statement\n",
    "\n",
    "The venerable alcohol stove has been invaluable camping accessory for generations. They are simple, reliable, and in a pinch, can be made from aluminum soda cans. \n",
    "\n",
    "![alcohol-stove.jpeg](https://github.com/jckantor/ND-Pyomo-Cookbook/blob/master/notebooks/figures/alcohol-stove.jpeg?raw=1)\n",
    "\n",
    "Alcohol stoves are typically fueled with denatured alcohol. Denatured alcohol, sometimes called methylated spirits, is a generally a mixture of ethanol and other alcohols and compounds designed to make it unfit for human consumption. An MSDS description of one [manufacturer's product](https://www.korellis.com/wordpress/wp-content/uploads/2016/05/Alcohol-Denatured.pdf) describes a roughly fifity/fifty mixture of ethanol and methanol.\n",
    "\n",
    "The problem with alcohol stoves is they can be difficult to light in below freezing weather. The purpose of this notebook is to design of an alternative cold weather fuel that could be mixed from other materials commonly available from hardware or home improvement stores."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "UXtpYyp4gNfb",
    "pycharm": {}
   },
   "source": [
    "## Vapor pressure data\n",
    "\n",
    "The following data was collected for potential fuels commonly available at hardware and home improvement stores. The data consists of price (\\$/gal.) and parameters to predict vapor pressure using the Antoine equation,\n",
    "\n",
    "$$\n",
    "\\begin{align}\n",
    "\\log_{10}P^{vap}_{s}(T) & = A_s - \\frac{B_s}{T + C_s}\n",
    "\\end{align}\n",
    "$$\n",
    "\n",
    "where the subscript $s$ refers to species, temperature $T$ is in units of degrees Celcius, and pressure $P$ is in units of mmHg. The additional information for molecular weight and specific gravity will be needed to present the final results in volume fraction."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "colab": {},
    "colab_type": "code",
    "id": "Wy2WYE37gNfc",
    "pycharm": {}
   },
   "outputs": [],
   "source": [
    "data = {\n",
    "    'ethanol'          : {'MW':  46.07, 'SG': 0.791, 'A': 8.04494, 'B': 1554.3,  'C': 222.65},\n",
    "    'methanol'         : {'MW':  32.04, 'SG': 0.791, 'A': 7.89750, 'B': 1474.08, 'C': 229.13},\n",
    "    'isopropyl alcohol': {'MW':  60.10, 'SG': 0.785, 'A': 8.11778, 'B': 1580.92, 'C': 219.61},\n",
    "    'acetone'          : {'MW':  58.08, 'SG': 0.787, 'A': 7.02447, 'B': 1161.0,  'C': 224.0},\n",
    "    'xylene'           : {'MW': 106.16, 'SG': 0.870, 'A': 6.99052, 'B': 1453.43, 'C': 215.31},\n",
    "    'toluene'          : {'MW':  92.14, 'SG': 0.865, 'A': 6.95464, 'B': 1344.8,  'C': 219.48},\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "6Qy4mtDygNfg",
    "pycharm": {}
   },
   "source": [
    "The first step is to determine the vapor pressure of denatured alcohol over a typical range of operating temperatures. For this we assume denatured alcohol is a 40/60 (mole fraction) mixture of ethanol and methanol. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 312
    },
    "colab_type": "code",
    "executionInfo": {
     "elapsed": 5305,
     "status": "ok",
     "timestamp": 1555699125535,
     "user": {
      "displayName": "Jeffrey Kantor",
      "photoUrl": "https://lh5.googleusercontent.com/-8zK5aAW5RMQ/AAAAAAAAAAI/AAAAAAAAKB0/kssUQyz8DTQ/s64/photo.jpg",
      "userId": "09038942003589296665"
     },
     "user_tz": 240
    },
    "id": "CnvZREQYgNfi",
    "outputId": "6de811ab-a13c-4781-e583-94ba36b3da96",
    "pycharm": {}
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Vapor Pressure at 0C = 22.1 mmHg\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEWCAYAAAB8LwAVAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xl8FPX9x/HXJwkkQIDIFW7CjYAcBgFv4q1VUetRvI+CVlttq61H+/Nqba311nqheCEErGepolaDggKScCM3JBBuSAgJIeTYz++PmegaN8luks1sks/z8dhH5p73zsJ+duY7h6gqxhhjTEVRXgcwxhgTmaxAGGOMCcgKhDHGmICsQBhjjAnICoQxxpiArEAYY4wJyAqEMbUkIseLyHoRKRCRC4KY/jUR+Wt9ZKtvIjJHRH5Zw3lVRPrVcv2ZInKaV+tvbKxANAAi8omIPBhg+HgR2SkiMV7k8stxrYiUuV+QB0RkqYic62WmevYg8Kyqxqvq+16HKVebL2tjwApEQ/EacJWISIXhVwFvqWppfQWpohjNV9V4IAF4BZgpIu1CmL/ehCFDL2BVHS/Tc5HwWRlvWYFoGN4H2gEnlg8QkSOAc4E33P6ficgS9xf8VhG532/aJHf3eZKIbBeRHSJyu9/4WBF50h233e2OdceNE5FsEblTRHYCr1YVVFV9wBSgBdCnsvlF5Fx3T2O/iHwjIsP88twpIttEJF9E1orIqe7w0SKS7r7HXSLyuH9G/xz+hxpE5H4R+beITBWRA8C1IhIlIneJyEYR2SciAQua3/ImisgGEckRkQ9FpKs7fCPQB/iPuwcVG2DekSKy2H0/M4C4CuOr2haZInKHiCwXkTwRmSEice64I0RklojsEZFct7u7O+4hnH8vz7q5nvX7dxDjt/zv9zLcPcGvReQJEckB7neHXy8iq911fCIivfzmP11E1rjZngUq/ojxf5+jRWS++z53uJmaVzJtCxF5TESy3GXPE5EW7rjzRWSVu5w5InJkhdlHBNpeVX2OphKqaq8G8AImAy/79d8ILPXrHwcchVP0hwG7gAvccUmAAtOBVu50e4DT3PEPAguATkBH4BvgL37LLQX+AcQCLQJkuxaY53bHALcB+UDbQPMDRwO7gTFANHANkOmOHwhsBbr6Ze/rds8HrnK744GxfhmzK2TK9Ht/9wMlwAXu9mkB/NZ9z93d9b4ITK9k258C7HVzxwLPAF8FWleAeZsDWcDvgGbAxW6Wv7rjK90Wfsv+FuiK8yNhNXCTO6498HOgJdAaeBt432/dc4Bf+vWX/zuICTSN+zmWAr9xP8cW7jbbABzpDvsz8I07fQfggPuemrnvsdR/nRW2RTIw1l1Okvtefus3XoF+bve/3Gzd3O1ynLvtBwAHgdPddf7Rzdc8iO1V3ef4/frt5W4TrwPYK8gPCk4A8nC/oIGvgd9VMf2TwBNud/kXwyC/8Y8Ar7jdG4Fz/MadCWS63eOAYiCuinWVf7Hsd/8DLuCHL+efzA88j1uA/IatBU4G+uF8YZ4GNKswzVfAA0CHCsPHUX2B+KrC+NXAqX79XXC+uGMCvL9XgEf8+uPdaZMqrivAvCcB2wHxG/YNPxSISreF37KvrPC5vVDJukYAuX79cwi9QGypsMyPgRv8+qOAQpzDalcDC/zGCZBNJQUiQN7fAu/59av7+UcBh4DhAeb5P2BmhTzbgHHVba8gPkcrEBVedoipgVDVeTi/+seLSB/gGGBa+XgRGSMiae7hhjzgJpxfeP62+nVn4fzKwv2bVck4gD2qWlRNxAWqmqCqHVR1rKr+r4r5ewG3u4cI9ovIfqAHzl7DBpwvjvuB3SKS6ncY4AacX5BrRGSRhNYQvrVCfy/gPb/1rwbKgMQA8/5o+6hqAbAP59dtdboC29T9BnL5b+tKt4XfNDv9ugtxvtgQkZYi8qJ7GOYATgFNEJHoIHJVJtB2esovWw5OIejmZvx+evc9Vpz/eyIywD0MttPN+zd++m8Ud1gczg+Xiip+Fj53nf6fRcDtFWDeUD7HJskKRMPyBs6vtquAT1V1l9+4acCHQA9VbQu8wE+PB/fw6+6J88sW92+vSsaB88uqNirOvxV4yC0o5a+WqjodQFWnqeoJbibFOTyFqq5X1Qk4h8L+AfxbRFrhHHJoWb5w9wuyYxAZzq6QIU5VtwXI/6Pt466zPc4v1+rsALqJ/OgEg57Bbotq3I5zSG6MqrbB2VuBHz73iu/5oPu3pd+wzhWmCbSdbqyQr4WqfuO+t+//TbnvsQeVex5YA/R3895D4DaLvUAR0DfAuIqfRfk6g/ksavM5NklWIBqWN3AOvUwEXq8wrjWQo6pFIjIauDzA/P/n/uocAlwHzHCHTwf+LCIdRaQDcC8wNSzvwDEZuMnd6xERaSVOI3trERkoIqe4jb1FOIcaygBE5EoR6ej+atzvLqsMWAfEuctohnOc/CeNxRW8ADxU3uDqvvfxlUw7DbhOREa4uf4GLFTVzCDe63ycw2+3ikiMiFwEjA5mWwSx7NY422e/OA3s91UYvwunAR0AVd2D82V4pYhEi8j1BP4S9vcCcLf7bwYRaSsil7jj/gsMEZGL3IbvW/lpwamY9wBQICKDgF8Fmkh/ONHhcRHp6mY91t32M4Gficip7md9O3AY57BddWrzOTZJViAaEPcf8jc4Dc0fVhh9M/CgiOTjfMHPDLCIL3Ea9D4HHlXVT93hfwXSgeXACmCxOywsVDUdp8g9C+S6ma51R8cCD+P8ityJs7dwjzvuLGCViBQATwG/UNUiVc3Def8v43wBHsQ5Fl6Vp3C24afuNluA01AcKO/nOMe+38H51dwX+EWQ77UYuMh9f7nAZcC7fuOr2hbVeRKnIbm83Wd2hfFPAReLc/bR0+6wicAfcA6tDKGaL1ZVfQ9nby3VPSy0EjjbHbcXuATn89oH9MdpG6vMHTg/XPJxCuOMaqZdASzCOaz1DyBKVdcCV+I0MO8FzgPOc7dzlWrzOTZV8uNDo6YxEpEkYDNOo2+9XTNhjGnYbA/CGGNMQFYgjDHGBGSHmIwxxgRkexDGGGMCatA34+rQoYMmJSXVeP6DBw/SqlWrugtURyxXaCxXaCxXaBpjroyMjL2qWvFaoZ/y+lLu2rySk5O1NtLS0mo1f7hYrtBYrtBYrtA0xlxAutqtNowxxtSUFQhjjDEBWYEwxhgTkBUIY4wxAVmBMMYYE5AVCGOMMQFZgTDGGBOQFQhjjGlgXv8mk+/2lYV9PVYgjDGmAUnPzOHBWd8xN7sk7Otq0LfaMMaYpiSvsITbUpfSLaEFVw8J//psD8IYYxoAVeXOd5az60ARz0wYSYuYQI/zrltWIIwxpgF4a+EWZq/ayR/PGsjwHgn1sk4rEMYYE+HW7sznL7O+46QBHfnlCX3qbb1WIIwxJoIdKi7j19MW0zquGY9dMpyoqPAfWipnjdTGGBPBHpz1Het3F/DmDaPp2Dq2XtdtexDGGBOh/rt8B9O/3cJNJ/flxP7VP9+nrlmBMMaYCLQ1p5C73l3OiB4J3H7GAE8yhK1AiEgPEUkTkdUiskpEbnOHtxORz0Rkvfv3CHe4iMjTIrJBRJaLyNHhymaMMZGspMzHbalLQOGZCSNpFu3Nb/lwrrUUuF1VjwTGAreIyGDgLuBzVe0PfO72A5wN9Hdfk4Dnw5jNGGMi1uOfrWPxlv38/edH0aNdS89yhK1AqOoOVV3sducDq4FuwHjgdXey14EL3O7xwBvuI1MXAAki0iVc+YwxJhJ9vnoXz8/ZyOVjenLusK6eZhHn+dVhXolIEvAVMBTYoqoJfuNyVfUIEZkFPKyq89zhnwN3qmp6hWVNwtnDIDExMTk1NbXGuQoKCoiPj6/x/OFiuUJjuUJjuUJTn7n2FPq4f/4hOrSI4k9j4mgeXfkprbXJlZKSkqGqo6qdUFXD+gLigQzgIrd/f4Xxue7f/wIn+A3/HEiuatnJyclaG2lpabWaP1wsV2gsV2gsV2jqK1dRSame/8xcHXrfbM3ae7Da6WuTC0jXIL6/w9ryISLNgHeAt1T1XXfwrvJDR+7f3e7wbKCH3+zdge3hzGeMMZHib/9dzbLsPB69ZDg923vX7uAvnGcxCfAKsFpVH/cb9SFwjdt9DfCB3/Cr3bOZxgJ5qrojXPmMMSZS/GfZdl6fn8XEE3tz5pDOXsf5XjivpD4euApYISJL3WH3AA8DM0XkBmALcIk77iPgHGADUAhcF8ZsxhgTETbuKeCud5YzqtcR/PGsQV7H+ZGwFQh1Gpsra2E5NcD0CtwSrjzGGBNpDhWXcfPUxcQ2i+aZy7273qEydi8mY4zxgKry5/dXsm53Pm9cP5oubVt4HeknIqtcGWNMEzEzfSvvLM7m1lP6e3KfpWBYgTDGmHq2clse936wihP6deDWU/t7HadSViCMMaYe5Rws5sY3M2jXqjlP/mIE0fX4fIdQWRuEMcbUk9IyH7dOX8KegsO8feOxdIiv3+c7hMr2IIwxpp7889O1zNuwl79eMLTenitdG1YgjDGmHsxavp0Xv9zElWN7cumoHtXPEAGsQBhjTJit3ZnPH/+9nOReR3DvuUO8jhM0KxDGGBNGeYUlTHoznfjYGJ6/4miaxzScr11rpDbGmDDx+ZTfzljC9v2HSJ00lk5t4ryOFJKGU8qMMaaBefJ/60hbu4d7zxtCcq92XscJmRUIY4wJg09X7eTpLzZwSXJ3rhzT0+s4NWIFwhhj6tjanfn8bsZShnVvy18uGIrz9IOGxwqEMcbUoZyDxfzyjUW0io3hpatGEdcs2utINWaN1MYYU0eKS338amoGuw4cZuaNx9K5bcNqlK7I9iCMMaYOqCr3fbiShZtz+OfFwxjRAK6Uro4VCGOMqQOvf5PJ9G+3cktKX8aP6OZ1nDoRzmdSTxGR3SKy0m/YDBFZ6r4yyx9FKiJJInLIb9wL4cpljDF17at1e3hw1necPjiR208f6HWcOhPONojXgGeBN8oHqOpl5d0i8hiQ5zf9RlUdEcY8xhhT5zbuKeCWaYsZkNiaJy8bQVQE3747VOF8JvVXIpIUaJw453xdCpwSrvUbY0y45RWWMPH1dJpHR/HyNaNoFdu4zvsRVQ3fwp0CMUtVh1YYfhLwuKqO8ptuFbAOOAD8WVXnVrLMScAkgMTExOTU1NQa5ysoKCA+Pr7G84eL5QqN5QqN5QpNZbnKfMrjGUWsyfFx5+g4BhxRv6ez1mZ7paSkZJR//1ZJVcP2ApKAlQGGPw/c7tcfC7R3u5OBrUCb6pafnJystZGWllar+cPFcoXGcoXGcoUmUC6fz6d/fm+F9rpzls74dkv9h9LabS8gXYP4Dq/3s5hEJAa4CJhRPkxVD6vqPrc7A9gIDKjvbMYYE4xX5m3mzQVZ3HhSHy49pmE826EmvDjN9TRgjapmlw8QkY4iEu129wH6A5s8yGaMMVWavXIHD320mnOO6sydZw3yOk5YhfM01+nAfGCgiGSLyA3uqF8A0ytMfhKwXESWAf8GblLVnHBlM8aYmliyJZfbUpcyokcCj1/auM5YCiScZzFNqGT4tQGGvQO8E64sxhhTW1tzCvnl6+l0ahPL5Ksb9j2WgmVXUhtjTDXyCku49tVvKfUpr147mg7xsV5HqhdWIIwxpgqlPuXGqelsySnkxauS6dcp8k7FDZfGdVWHMcbUIVVlyspiFmwv5InLhjO2T3uvI9Ur24MwxphKPPX5er7ZXsrvThvAhSO7ex2n3lmBMMaYAGYu2sqT/1vP8V1juPXUfl7H8YQdYjLGmAo+X72Lu99bwYn9O3BN78IG+8jQ2rI9CGOM8ZORlcst0xYzuEsbnr8ymZhGfq1DVaxAGGOMa8PufG54fRGd28Tx6nXHEN/I7s4aKisQxhgD7Mwr4upXviUmKoo3rh/TZK51qIoVCGNMk5dXWMI1U77lQFEpr113DD3bt/Q6UkSwAmGMadKKSsqY+EY6m/YW8NJVyQzt1tbrSBGjaR9gM8Y0aWU+5bbUJSzKyuHpX4zkuH4dvI4UUWwPwhjTJKkq//fBSj5ZtYt7zx3MecO7eh0p4liBMMY0OarKw7PXMG3hFm4e15frju/tdaSIVO0hJhH5fYDBeUCGqi6t+0jGGBNez83ZyItfbuKqsb34w5kDvY4TsYLZgxgF3AR0c1+TgHHAZBH5Y/iiGWNM3Xv9m0z++claLhzZjQfOH9Jkr5IORjCN1O2Bo1W1AEBE7sN56ttJQAbwSPjiGWNM3XknI5v7PlzF6YMT+efFwxr9E+FqK5g9iJ5AsV9/CdBLVQ8BhyubSUSmiMhuEVnpN+x+EdkmIkvd1zl+4+4WkQ0islZEzqzBezHGmErNXrmTP/x7Gcf3a88zE0YSE21NsNUJZg9iGrBARD5w+88DpotIK+C7KuZ7DXgWeKPC8CdU9VH/ASIyGOdZ1UOArsD/RGSAqpYFkc8YY6r01bo93Dp9CSN6JPDSVU3jcaF1odoCoap/EZGPgBMAAW5S1XR39BVVzPeViCQFmWM8kKqqh4HNIrIBGA3MD3J+Y4wJKD0zh0lvptOnYytevXY0rZr4/ZVCIaoaeIRIu6pmVNWcahfuFIhZqjrU7b8fuBY4AKQDt6tqrog8CyxQ1anudK8AH6vqvwMscxJOQzmJiYnJqamp1cWoVEFBAfHxkff4QMsVGssVmqaUKzOvjH8sKqJNc+GeMS1oGxt6m0Nj3F4pKSkZqjqq2glVNeAL2Axscv8Wud3l/Zsqm6/CMpKAlX79iUA0TtvHQ8AUd/i/gCv9pnsF+Hl1y09OTtbaSEtLq9X84WK5QmO5QtNUcq3alqfDH/hEj/v755qdW1jj5TTG7QWkaxDf4ZXua6nq91eOiMgSVR0ZUokKvMxdfsucDMxye7OBHn6Tdge213Z9xpimac3OA1zx8gJaNItm+sSxdEto4XWkBinYZvzAx6FCJCJd/HovBMrPcPoQ+IWIxIpIb6A/8G1drNMY07Ss35XPFZMX0jwmiukTx9qdWWshbK01IjId54K6DiKSDdwHjBORETgFJxO4EUBVV4nITJyzokqBW9TOYDLGhGjD7gImTF5IVJQwbeJYkjq08jpSg1Zpgahwi41OFW+5oaqPV7VgVZ0QYPArVUz/EE67hDHGhGzz3oNcPnkBoEyfOJa+HSOvYbmhqWoPorVf9+QK/cYYEzGy9h1kwksLKPUpqZPG0q+TfV3VhaoaqR+ozyDGGFMTW3MKmfDSAopKy5g+cSwDEq041JVg7ubaG/gNzimr30+vqueHL5YxxlRva04hEyYv4GBxGW/9cgxHdmnjdaRGJZhG6vdx2g7+A/jCG8cYY4KT6bY5HCwuY+oNY+xRoWEQTIEoUtWnw57EGGOCtGF3AZdPdtocpk8cy+CutucQDsEUiKfcW3x/it/dW1V1cdhSGWNMJdbuzOeKlxcAQuoka3MIp2AKxFHAVcAp/HCISd1+Y4ypN6u253HlywtpFh3FtIlj6dfJTmUNp2AKxIVAH1UtrnZKY4wJk+XZ+7nqlW9p1TzaLoKrJ8HcamMZkBDuIMYYU5mMrFyumLyQ1nExzLjxWCsO9SSYPYhEYI2ILOLHbRB2mqsxJuwWbtrH9a8tomPrWKZNHEtXu/FevQmmQNwX9hTGGBNA2prd3DQ1g25HtGD6xLEktonzOlKTEswT5b4EEJE2wUxvjDF14cNl2/n9jKUM6tKa168bTfv4WK8jNTnBXEk9CfgLcAjnLCbBOYupT3ijGWOaqqkLsvi/D1ZyTFI7Xr5mFG3imnkdqUkKZo/gD8AQVd0b7jDGGPPcnA08MnstpwzqxHNXHE1cs2ivIzVZwRSIjUBhuIMYY5o2VeXvH6/mxS83MX5EVx69ZDjNooN9ppkJh2AKxN3ANyKykB+fxXRr2FIZY5qUMp/y2qpivszexFVje/HA+UOIihKvYzV5wRSIF4EvgBXYzfqMMXWsuNTH72Yu5cvsUn6d0o/bzxiAiBWHSBBMgShV1d9XP9mPicgU4Fxgt6oOdYf9EzgPKMY5dHWdqu4XkSRgNbDWnX2Bqt4U6jqNMQ1LflEJN03N4OsN+7hsYHPuOHOg15GMn2AO8KWJyCQR6SIi7cpfQcz3GnBWhWGfAUNVdRiwDufwVbmNqjrCfVlxMKaR23WgiEtfXMDCTTk8dslwzu5tZypFmmD2IC53//p/mVd7mquqfuXuGfgP+9SvdwFwcRDrN8Y0Mht253PNlEXsLyxmyrXHcNKAjsyZs8HrWKaCYC6U6x2mdV8PzPDr7y0iS4ADwJ9VdW6Y1muM8VB6Zg43vJ5Os+goZtx4rD3oJ4KJqoZv4c4exKzyNgi/4X8CRgEXqaqKSCwQr6r7RCQZ5yl2Q1T1QIBlTgImASQmJianpqbWOF9BQQHx8ZF3u2DLFRrLFRovc2XsKuWFZYdpHyfcPiqOji1/OMpt2ys0tcmVkpKSoaqjqp1QVcP2wnmO9coKw64B5gMtq5hvDjCquuUnJydrbaSlpdVq/nCxXKGxXKHxKtcb32zW3nfN0gv+NU/3FRz+yXjbXqGpTS4gXYP4Dq/0EJOIxKhqaY3KU+XLPAu4EzhZVQv9hncEclS1TET6AP2BTXW5bmOMN3w+5dFP1/LcnI2cdmQiz0wYSYvmdnV0Q1BVG8QCEckGZgOzVTUzlAWLyHRgHNDBXc59OA3dscBn7nnO5aezngQ8KCKlQBlwk6rmhPhejDERpqikjDveXsas5TuYMLonfxk/hBi7OrrBqLRAqOooEekFnA08KSLdgHnAx8CXqnq4snnd+ScEGPxKJdO+A7wTdGpjTMTbW3CYiW+ks3Trfu4+exCTTupjF8A1MFWexaSqWcALwAsi0gw4Eefahr+KyB5V/Vk9ZDTGNDDrduVz/WuL2FtwmOevOJqzhnbxOpKpgaCf76CqJTi33PgCwN2jMMaYH5m7fg83T11MXPNoZt54LMO62xOLG6oaPwBIVbfVZRBjTMM3beEW/u+DlfTvFM+Ua4+xx4M2cPaEOGNMrZX5lIc/Xs3kuZtJGdiRZy4/mvhY+3pp6IL+BEWklaoeDGcYY0zDk19Uwu9mLOV/q3dz7XFJ/PlnR9qZSo1EtZ+iiBwnIt/h3G0VERkuIs+FPZkxJuJt3nuQC5/7hrS1e3hw/BDuP99OY21MgtmDeAI4E/gQQFWXichJYU1ljIl4X63bw6+nLSY6Sph6wxiO7dve60imjgV1iElVt1Y4f7ksPHGMMZFOVXll3mb+9tFqBiS2ZvLVo+jRrqXXsUwYBFMgtorIcYCKSHPgVtzDTcaYpqWopIx73lvBu4u3cfbQzjx6yXBaWWN0oxXMJ3sT8BTQDcgGPgVuCWcoY0zk2ZlXxI1TM1i2dT+/P30Av07pZ8+NbuSqLBAiEg1cpapX1FMeY0wEysjK4VdTF3PwcCkvXpXMmUM6ex3J1IMqTzdQ1TJgfD1lMcZEGFXl1a83c9mLC2jRPJp3bz7eikMTEswhpq9F5Fmcp799fx2Eqi4OWypjjOcOHi7l7ndX8OGy7Zx2ZCKPXTqcti3sudFNSTAF4jj374N+wxQ4pe7jGGMiwaY9Bdw0NYMNuwv4w5kD+dXJfa29oQkK5pnUKfURxBgTGWav3Mkdby+jeUwUb1w/hhP6d/A6kvFItQVCRO4NNFxVHww03BjTMJWW+fjnp2t58ctNDO+RwPNXHG0322vigjnE5H//pTjgXOw6CGMald0Hirg1dQkLNuVw5die/N+5g4mNsceCNnXBHGJ6zL9fRB7Fve1GdURkCk5B2a2qQ91h7XAavJOATOBSVc0V51Ltp4BzgELgWmsINyb8vlq3h9/NWEphcRmPXTKcnyd39zqSiRA1uatWS6BPkNO+hvMEOn93AZ+ran/gc7cfnEeb9ndfk4Dna5DNGBOk0jIfj8xew9VTvqVDfCz/+c3xVhzMjwTTBrEC56wlgGigIz8+o6lSqvqViCRVGDweGOd2vw7MAe50h7+hqgosEJEEEemiqjuCWZcxJnjb9x/i1ulLSM/KZcLoHtx33hDimtkhJfNjwbRBnOvXXQrsUtXSWqwzsfxLX1V3iEgnd3g3YKvfdNnuMCsQxtShJbtL+e3Tcykp9fHUL0YwfoQ9PdgEJs4P9iomEOkLZKvqYREZBwzD+aW/P6gVOHsQs/zaIParaoLf+FxVPUJE/gv8XVXnucM/B/6oqhkVljcJ5xAUiYmJyampqUG90UAKCgqIj4+v8fzhYrlCY7mCU+pT3l5XzCeZpfRqE8WvhsfSuVXkPLsh0rZXucaYKyUlJUNVR1U7oapW+QKW4uxp9AM24jwf4qPq5vObPwlY6de/FujidncB1rrdLwITAk1X2Ss5OVlrIy0trVbzh4vlCo3lqt6G3fl67tNztdeds/T6f83WQ8WlXkf6iUjaXv4aYy4gXYP4/g7m54NPnUNKFwFPqurv3C/2mvoQuMbtvgb4wG/41eIYC+SptT8YUyuqyrSFWzj36XlszS3khSuTuWpwrLU3mKAE0wZRIiITgKuB89xhQd2QRUSm4zRIdxCRbOA+4GFgpojcAGwBLnEn/wjnFNcNOKe5XhfkezDGBLCv4DB3vbuCz77bxQn9OvDYpcNJbBPHnDlrvI5mGohgCsR1OM+EeEhVN4tIb2BqMAtX1QmVjDo1wLSKPWfCmDrx1bo93P72MvIKS/jzz47k+uN7272UTMiCuVDuO5ynyCEiRwCtVfXhcAczxoSuqKSMR2avZcrXmxmQGM/r141mcNc2XscyDVQw10HMAc53p10K7BGRL1X192HOZowJwXfbD/D7mUtZszOfa49L4q6zB1lbg6mVYA4xtVXVAyLyS+BVVb1PRJaHO5gxJjilZT6en7ORp79YT0LL5rx63TGkDOxU/YzGVCOYAhEjIl2AS4E/hTmPMSYE63flc/vby1iencf5w7vywPlDOKJVc69jmUYimALxIPAJ8LWqLhKRPsD68MYyxlSlzKe8Mm8Tj366jvjYGJ674mjOOao2Z58b81PBNFK/Dbzt178J+Hk4QxljKpe59yB3vL2M9KxczhicyEMXHkXH1rFexzKNUDBOGuJtAAAWIElEQVSN1ANw7qyaqKpDRWQYcL6q/jXs6Ywx3/P5lDfmZ/KP2WtpFi08cdlwLhjRDedO+cbUvWAOMU0G/oBzKwxUdbmITAOsQBhTTzbszufOd1aQkZXLyQM68o+fD6Nz2zivY5lGLpgC0VJVv63wK6U2d3M1xgSpuNTHi19u5JkvNtAyNprHLx3OhSNtr8HUj2AKxF73jq4KICIXY7fgNibslm3dz53vLGfNznzOHdaF+84bYm0Npl4FUyBuAV4CBonINmAzcEVYUxnThB0qLuPxz9byyrzNdGwdy+SrR3H64ESvY5kmqMoCISJRwChVPU1EWgFRqppfP9GMaXrmrd/LPe+tYEtOIZeP6cldZw+iTVxQ98Y0ps5VWSBU1ScivwZmqurBespkTJOzO7+Iv85azYfLtpPUviXTJ47l2L7tvY5lmrhgDjF9JiJ3ADOA74uEquaELZUxTUSZT3lrYRb/nL2Ww6U+bju1P78a19fuoWQiQjAF4nr3r/+tuBXoU/dxjGk6VmTn8af3V7A8O48T+nXgwfFD6NMx8h5taZquYK6k7l0fQYxpKg4UlfDYJ2t5c0EW7eNjeXrCSM4b1sVOXTURJ5grqeOAm4ETcPYc5gIvqGpRmLMZ06ioKu8v3cbfPlrD3oLDXD22F7efOdAaoU3ECuYQ0xtAPvCM2z8BeJMfHhVqjKnGym153PfhKjKychnevS2vXDOKYd0TvI5lTJWCKRADVXW4X3+aiCyr6QpFZCBOg3e5PsC9QAIwEdjjDr9HVT+q6XqMiQT7Cg7z6KfrSF20hfatmvPIxcO4+Oju9vhP0yAEUyCWiMhYVV0AICJjgK9rukJVXQuMcJcVDWwD3sN59vUTqvpoTZdtTKQoLfPxWVYJt86ZQ2FxGdcf35vbTutvh5NMgxJMgRgDXC0iW9z+nsBqEVkBqKoOq8X6TwU2qmqWNdCZxuKbDXt54D/fsXZXMSf068B95w2mf2Jrr2MZEzJR1aonEOlV1XhVzarxykWmAItV9VkRuR+4FjgApAO3q2pugHkmAZMAEhMTk1NTU2u6egoKCoiPj7zTCi1XaCIl1/YCHzPWFrNsTxkdWggX9vJxXK9WEXd2UqRsr4osV2hqkyslJSVDVUdVO6GqevICmgN7cZ4zAZAIRANRwEPAlOqWkZycrLWRlpZWq/nDxXKFxutce/KL9E/vLdc+d/9Xh947W59L26CHiks9z1UZyxWaxpgLSNcgvqeDOcQULmfj7D3sAij/CyAik4FZXgUzJhhFJWW8Mm8zz8/ZyKGSMq4Y05PbTu1P+3i746ppHLwsEBOA6eU9ItJFVctvI34hsNKTVMZUw+dzrmd49JO1bM8r4rQjE7n7nEH0taugTSPjSYEQkZbA6cCNfoMfEZEROBfjZVYYZ0xEmL9xH3/7aDUrtuUxtFsbHrt0hN1UzzRanhQIVS0E2lcYdpUXWYwJxorsPB75ZA1z1++la9s4nrhsOOOHd7PrGUyj5uUhJmMi3obd+Tz26To+XrmThJbNuOecQVx9bJLdbdU0CVYgjAlga04hT32+nncXZ9OiWTS3ndqfX57Ym9Z2oZtpQqxAGONnT/5h/pW2gbcWZiEiXH98b341rq+dmWSaJCsQxgC5B4t5ed4mpszLpLjMx6WjenDrqf3o0raF19GM8YwVCNOk5RwsZvLcTbzxTSYHi8s4f3hXfnf6AHp3aOV1NGM8ZwXCNEl7Cw4zee4m3pyfxaGSMn52VBduPbU/A+yeScZ8zwqEaVL25B/mpa82MnXBFg6XlnHe8K785pR+9OtkhcGYiqxAmCZh94EiXvhyE9O+zaK41McFI7pxyyn97OpnY6pgBcI0alv2FfLS3I28nZ5NqU+5cGQ3bknpZ20MxgTBCoRplFZuy+OFLzfy0YodxERFcdHR3fjVuL70am+FwZhgWYEwjYaqMn/jPp7/ciNz1+8lPjaGiSf14frje5PYJs7reMY0OFYgTIPnU+XjFTt4/suNLM/Oo0N8LH88ayBXjOlF2xZ25bMxNWUFwjRYhcWlvLt4G8/MPcSuwsUktW/J3y48iouO7mb3SjKmDliBMA3OjrxDvP5NFtO/3ULeoRKS2kTxr8tHctbQzkTb3VWNqTNWIEyDsWRLLlO+zuSjFTtQVc4c0pnrT+hNweZlpAzr4nU8YxodKxAmopWW+fh45U6mfL2ZJVv20zo2huuOS+Ka45Lo0a4lAHMyba/BmHCwAmEi0t6Cw8xM38rU+VlszyuiV/uW3H/eYC4e1YP4WPtna0x98Ox/mohkAvlAGVCqqqNEpB0wA0jCeezopaqa61VGU79UlfSsXKYuyOKjFTsoKVOO7dOeB8YP5ZRBnax9wZh65vVPsRRV3evXfxfwuao+LCJ3uf13ehPN1JeCw6W8v2QbUxdksWZnPq1jY7hiTC+uGNOT/nbzPGM843WBqGg8MM7tfh2YgxWIRmvtznymLsjivSXbKDhcyuAubfj7RUdx/vCutLLDSMZ4TlTVmxWLbAZyAQVeVNWXRGS/qib4TZOrqkdUmG8SMAkgMTExOTU1tcYZCgoKiI+PvJu1NeZch0uVb3eWMndbKetyfcREwejOMZzSM4a+baMQCf0wUmPeXuFguULTGHOlpKRkqOqoaidUVU9eQFf3bydgGXASsL/CNLlVLSM5OVlrIy0trVbzh0tjy+Xz+XTJlly9651lOuTe2drrzlk67p9p+sKcDbqv4LBnucLNcoXGcoWmNrmAdA3ie9qz/XhV3e7+3S0i7wGjgV0i0kVVd4hIF2C3V/lM7eUcLOa9JduYuWgra3flE9csinOO6sJlo3owune7Gu0tGGPqjycFQkRaAVGqmu92nwE8CHwIXAM87P79wIt8pubKfMq8DXuZmb6Vz1btorjMx/DubXnowqGcN7wrbeLs3kjGNBRe7UEkAu+5vyBjgGmqOltEFgEzReQGYAtwiUf5TIhW7zjAu4uz+WDpdnbnHyahZTMuH9OTy47pwZFd2ngdzxhTA54UCFXdBAwPMHwfcGr9JzI1setAER8s3ca7i7exZmc+MVHCuIEduXBkd049spPdMM+YBs7OJTQhKSwu5ZNVO3l38Ta+3rAXn8LwHgk8cP4Qzh3WhfbxsV5HNMbUESsQplpFJWV8uW4PU5YWcfMX/6OwuIxuCS24eVw/Ljy6mz3X2ZhGygqECai41MfXG/byn+Xb+WzVLvIPl9KqGYwf2YMLRnTjmKR2RNmtL4xp1KxAmO+VlvmYv2kfs5btYPaqneQdKqF1XAxnDu3MucO6ULptFaedMszrmMaYemIFookrLfPxbWYO/12+g9krd7LvYDGtmkdz+uBEzh3WlRMHdCA2xmlsnrPjO4/TGmPqkxWIJqiopIy56/fyyaqdfL56F7mFJbRoFs0pR3bivGFdGDfQzkAyxliBaDLyCkv4Yu0uPlm5iy/X7eFQSRmt42I4dVAnzhzSmZMHdqRlc/vnYIz5gX0jNGI784r47LudfLJqFws27aPUp3RqHcvPk7tx5pDOjOndnuYxUV7HNMZEKCsQjYjPp6zYlscXa3aTtnY3y7PzAOjToRW/PLEPZwxJZET3BDv7yBgTFCsQDdyBohLmrd/LF2t2M2ftHvYWHEYERvZI4I4zBnDmkM706xRvN8YzxoTMCkQDo6ps3FPAF2t288Wa3aRn5lLqU9rExXDywE6cMqgjJw/oRLtWzb2Oaoxp4KxANAB5h0qYv3Evc9c7ry05hQAM6tyaiSf14ZRBnRjZI4GYaGtPMMbUHSsQEajUpyzKzGHuuj3M3bCXZVv341No1TyaY/u2Z9JJfUgZ1IluCS28jmqMacSsQEQAVWXT3oPMW7+Xuev3MG9dIUWfzidKnBvh/TqlHyf078jIngk0s70EY0w9sQLhAVUlO/cQ8zftY8GmfSzYuI/teUUA9GzXkmO7xnDZycM4tk8H2ra0B+wYY7xhBaKebM0p/L4gLNyUw7b9hwBo16o5Y3q34+aUDpzYvwO92rdizpw5jBvaxePExpimzgpEGKgqW3MOsXDzPhZsymHBpn0/KQg3ntyHsX3a099OQTXGRKh6LxAi0gN4A+gM+ICXVPUpEbkfmAjscSe9R1U/qu98NVFS5mPV9gOkZ+aQkZVLelYue/IPA3BEy2aM7eM0LJcXBLtQzRjTEHixB1EK3K6qi0WkNZAhIp+5455Q1Uc9yBSSvMISFm/JJT0rh/TMXJZl76eoxAdAt4QWHNe3PaN6HcGopHYMTGxtBcEY0yDVe4FQ1R3ADrc7X0RWA93qO0ewfD5l094Clm7NIyMrl4ysHNbtKgAgOkoY0rUNE0b3ZFSvdiT3OoLObeM8TmyMMXVDVNW7lYskAV8BQ4HfA9cCB4B0nL2M3ADzTAImASQmJianpqbWeP0FBQXEx//4cZm5RT425fnYtN/HprwyMg/4OFTqjGsRA/0Soul/RBT9E6Lp0zaK2Ji63zsIlCsSWK7QWK7QWK7Q1CZXSkpKhqqOqm46zwqEiMQDXwIPqeq7IpII7AUU+AvQRVWvr2oZo0aN0vT09Bpn+O9nabRNOopl2ftZunU/y7P3s+uA03YQEyUc2aUNw3u0ZXj3BIb3SKBfx/ppP5gzZw7jxo0L+3pCZblCY7lCY7lCU5tcIhJUgfDkLCYRaQa8A7ylqu8CqOouv/GTgVnhWv+K7Dxum7GETXsKgYWAc8fT4/p2YHj3tgzrkcDgLm3soTnGmCbNi7OYBHgFWK2qj/sN7+K2TwBcCKwMV4bENrH06RDPiIRiLjxpJMO6JdgFacYYU4EXexDHA1cBK0RkqTvsHmCCiIzAOcSUCdwYrgCd2sTx8jWjmDNnDif27xiu1RhjTIPmxVlM84BAB/IbxDUPxhjTVNid34wxxgRkBcIYY0xAViCMMcYEZAXCGGNMQFYgjDHGBGQFwhhjTEBWIIwxxgTk6c36aktE9gBZtVhEB5z7P0UayxUayxUayxWaxpirl6pWe5Vwgy4QtSUi6cHcsKq+Wa7QWK7QWK7QNOVcdojJGGNMQFYgjDHGBNTUC8RLXgeohOUKjeUKjeUKTZPN1aTbIIwxxlSuqe9BGGOMqYQVCGOMMQE1yQIhImeJyFoR2SAid3mdp5yIZIrIChFZKiI1f9h27XNMEZHdIrLSb1g7EflMRNa7f4+IkFz3i8g2d5stFZFzPMjVQ0TSRGS1iKwSkdvc4Z5usypyebrNRCRORL4VkWVurgfc4b1FZKG7vWaISPMIyfWaiGz2214j6jOXX75oEVkiIrPc/vBvL1VtUi8gGtgI9AGaA8uAwV7ncrNlAh0iIMdJwNHASr9hjwB3ud13Af+IkFz3A3d4vL26AEe73a2BdcBgr7dZFbk83WY4DwyLd7ub4TwYfiwwE/iFO/wF4FcRkus14GIv/425mX4PTANmuf1h315NcQ9iNLBBVTepajGQCoz3OFNEUdWvgJwKg8cDr7vdrwMX1GsoKs3lOVXdoaqL3e58YDXQDY+3WRW5PKWOAre3mftS4BTg3+5wL7ZXZbk8JyLdgZ8BL7v9Qj1sr6ZYILoBW/36s4mA/zQuBT4VkQwRmeR1mAoSVXUHOF88QCeP8/j7tYgsdw9B1fuhL38ikgSMxPn1GTHbrEIu8HibuYdLlgK7gc9w9ur3q2qpO4kn/y8r5lLV8u31kLu9nhCR2PrOBTwJ/BHwuf3tqYft1RQLRKDnYUfErwTgeFU9GjgbuEVETvI6UAPwPNAXGAHsAB7zKoiIxAPvAL9V1QNe5agoQC7Pt5mqlqnqCKA7zl79kYEmq99UP80lIkOBu4FBwDFAO+DO+swkIucCu1U1w39wgEnrfHs1xQKRDfTw6+8ObPcoy4+o6nb3727gPZz/OJFil4h0AXD/7vY4DwCqusv9T+0DJuPRNhORZjhfwm+p6rvuYM+3WaBckbLN3Cz7gTk4x/oTRCTGHeXp/0u/XGe5h+pUVQ8Dr1L/2+t44HwRycQ5JH4Kzh5F2LdXUywQi4D+7hkAzYFfAB96nAkRaSUircu7gTOAlVXPVa8+BK5xu68BPvAwy/fKv4BdF+LBNnOPB78CrFbVx/1GebrNKsvl9TYTkY4ikuB2twBOw2kfSQMudifzYnsFyrXGr8gLznH+et1eqnq3qnZX1SSc76svVPUK6mN7ed0y78ULOAfnjI6NwJ+8zuNm6oNzRtUyYJWXuYDpOIceSnD2uG7AOeb5ObDe/dsuQnK9CawAluN8IXfxINcJOLv3y4Gl7uscr7dZFbk83WbAMGCJu/6VwL3u8D7At8AG4G0gNkJyfeFur5XAVNwznbx4AeP44SymsG8vu9WGMcaYgJriISZjjDFBsAJhjDEmICsQxhhjArICYYwxJiArEMYYYwKyAmEaHBFJEJGbvc4RDBH5rYi0DPM6ZotItwrDuorIFyLygXsldfnwq0VkpXu30u9E5I5wZjMNmxUI0xAlABFRIMRR1f+j3wIhFQi/q2ODmbYFzvUV2yqMuhX4Dc7N3a50pz3bzXOGqg7BuTNuXijZTNNiBcI0RA8Dfd178/8TQET+ICKL3Buqld/HP0lE1ojIy+6v5rdE5DQR+dq9h/5od7r7ReRN9xf3ehGZWL6iKpa7WkSeAxYDPUTkeRFJr/AcgVuBrkCaiKS5wwr8ln2xiLzmdr8mIo+70/3DvbJ+irvuJSJS2R2Hx+HcEqKiaJwbu/n44b49d+Pc5rv8li5Fqjo5pC1vmhavrgi0l71q+gKS+PEzIc7AeYC74PzomYXz7IgkoBQ4yh2eAUxxpxsPvO/Ofz/OFewtgA44d/vtWs1yfcBYvwzt3L/ROF/Yw9z+TPye8QEU+HVfDLzmdr/mLj/a7f8bcKXbnYBz5X+rANviaeCUAMN7AV8B/wFau8NygLZef372ajivoHdljYlgZ7ivJW5/PNAf2AJsVtUVACKyCvhcVVVEVuB80Zf7QFUPAYfcX/GjcW5VUdlys1R1gd/8l7q3aI/BeVDPYJxbNoTibVUt83tP5/u1EcQBPXHuWeTveOAn7QiqmoVTzIypMSsQpjEQ4O+q+uKPBjrPQDjsN8jn1+/jx//+K95zRqtZ7kG//t44X9LHqGque9gorpKs/uupOM1Bv24Bfq6qaytZDiLSB9iqzoOvgrEKSMa5t5Ax1bI2CNMQ5eM8QrPcJ8D15WfriEg3EQn14TzjxXkmcXuc4/qLQlhuG5wv9zwRScR5nkdlWXeJyJFuw/aFVeT5BPiNewdRRGRkgGnOBmYH9e4cfwceEZHO7jJj3XYSYwKyPQjT4KjqPreheSXwsar+QUSOBOa736cFOGfulFW1nAq+Bf6LcxjnL+o05G4PZrmqukxEluD8Qt8EfO03+iXgYxHZoaopOM+mnoXTzrES57BVIH/Buef/crdIZALnVpjmLJwzlYKiqh+5Bex/7jIVp03GmIDsbq6myROR+3Eajx/1OkuwxHns5deqOsrrLKbxsj0IYxogdZ5uZsXBhJXtQRhjjAnIGqmNMcYEZAXCGGNMQFYgjDHGBGQFwhhjTEBWIIwxxgT0/xW9XznE9VVTAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "def Pvap(T, s):\n",
    "    return 10**(data[s]['A'] - data[s]['B']/(T + data[s]['C']))\n",
    "\n",
    "def Pvap_denatured(T):\n",
    "    return 0.4*Pvap(T, 'ethanol') + 0.6*Pvap(T, 'methanol')\n",
    "\n",
    "T = np.linspace(0, 40, 200)\n",
    "\n",
    "plt.plot(T, Pvap_denatured(T))\n",
    "plt.title('Vapor Pressure of denatured alcohol')\n",
    "plt.xlabel('temperature / °C')\n",
    "plt.ylabel('pressure / mmHg')\n",
    "print(\"Vapor Pressure at 0C =\", round(Pvap_denatured(0),1), \"mmHg\")\n",
    "plt.grid(True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "_wMhAg6ugNfo",
    "pycharm": {}
   },
   "source": [
    "## Cold weather product requirements\n",
    "\n",
    "We seek a cold weather fuel with increased vapor pressure at 0°C and lower, and also provides safe and normal operation of the alcohol stove at higher operating temperatures. \n",
    "\n",
    "For this purpose, we seek a mixture of commonly available liquids with a vapor pressure of at least 22 mmHg at the lowest possible temperature, and no greater than the vapor pressure of denatured alcohol at temperatures 30°C and above."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 295
    },
    "colab_type": "code",
    "executionInfo": {
     "elapsed": 5673,
     "status": "ok",
     "timestamp": 1555699125915,
     "user": {
      "displayName": "Jeffrey Kantor",
      "photoUrl": "https://lh5.googleusercontent.com/-8zK5aAW5RMQ/AAAAAAAAAAI/AAAAAAAAKB0/kssUQyz8DTQ/s64/photo.jpg",
      "userId": "09038942003589296665"
     },
     "user_tz": 240
    },
    "id": "hfq_rUEagNfp",
    "outputId": "fa6b9c53-0416-4c31-cc7f-25dbefbf06ed",
    "pycharm": {}
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEWCAYAAAB8LwAVAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzsnXd8jdcfx98nESFBROwaMYpKJIigqBWKltpaVRWUVq3S2lKqVFVbqzVr71ltdfkZsWPWpkWFxExkSCLz5vz+uE9ubyLjBvfeSM779bqvPM8Zz/mc59483+es7xFSShQKhUKhSIuNtQUoFAqFImeiDIRCoVAo0kUZCIVCoVCkizIQCoVCoUgXZSAUCoVCkS7KQCgUCoUiXZSBUCiyQAjRWAhxRQgRLYTo9Iyu6SqEkEKIfM/iek+hY7IQYo01NeQUhBC+QoiD1taRk1AGIgcghPhTCDElnfCOQoi7OeAh4iuE0GkPyIdCiNNCiPbW1GRhpgDfSSkLSSm3W1uMMUKIQCFEK2vrUOROlIHIGawAegshRJrw3sBaKWWSpYRkYoyOSCkLAUWBpcAmIUSxbOS3GGbQUBG48IyvqVDkeJSByBlsB4oBr6QECCGcgfbAKu38dSHEX9obfJAQYrJR2pTuioFCiNtCiDtCiI+N4u2FELO1uNvasb0W11wIESyEGCOEuAssz0yolDIZWAYUBCpnlF8I0V5raUQIIQ4LITyM9IwRQtwSQkQJIf4WQvho4fWFECe0Ot4TQnxrrNFYh/Gbs9ZNskUIsUYI8RDwFULYCCHGCiGuCSEeCCHSNWhG1xsghLgqhAgTQvwshCirhV8DKgO/aC0o+3TyZlQfkzUIIZyEEEu17+6WEGKqEMI2jb5LWhkXhRB1hRCrgQpG2kZraRtq9zxCCHFGCNHc6DqVhBD7tOv8Dyie0T3R0nfUvseHWj3aauFltfsUpt23AUZ5JgshNmvfR5QQ4pwQopoQYpwQ4r72+33VKL2/EGK6EOKYECJSCPGT8X0SQrwhhLig1cdfCPGSUZwUQlQ1Ol8hhJiqHaf8Nj/Wyr0jhOhrlNZFq8NDIcQxoIpRnBBCzNLyRQohzgoh3DO7V7kSKaX65IAPsAT4wej8feC00XlzoBZ6o+4B3AM6aXGugATWA45auhCglRY/BQgASgIlgMPA50bXTQJmAPZAwXS0+QIHteN8wHAgCnBKLz9QF7gPNABsgT5AoBZfHQgCyhppr6IdHwF6a8eFgIZGGoPTaAo0qt9kIBHopN2fgsBHWp3LaeUuAtZncO9bAqGabntgHrA/vbLSyZtZfTLUYPSd5dPOt2vxjtr3dAx4X4vrDtwCvAEBVAUqpqcNeAF4ALym3YvW2nkJo3v8raanqfY9rsmgbvWBSO0aNtq1a2hx+4D5QAGgNvrfm4/R9xEHtEH/e1kFXAcmAHbAAOC6UTn+Wv3ctfpvTdEEVANiNA12wGjgKpBfi5dAVaNrrQCmpvltT9HyvgY8Apy1+A3AJq1Md01Dyu+8DXASfYtZAC8BZaz9nLD4c8naAtRH+yKgifbPWFA7PwSMyCT9bGCWdpzysKlhFP8VsFQ7vga8ZhTXBgjUjpsDCUCBTMry1f7RItA/SAP47+H8WH5gAZoBMgr7G2iG/uF2H2gF2KVJsx/4DCieJrw5WRuI/WniL6E9sLTzMuiNSL506rcU+MrovJCW1jVtWenkzaw+GWow+s7yAaWAeIyMM9AT2Ksd/wkMz6D8VNqAMcDqNGn+RG+kK2jfo6NR3DoyNhCLUn5jacLLAzqgsFHYdGCF0ffxP6O4DkA0YKudF9bqXlQ79we+NEpfU/tN2QJ+wCajOBv0D/Lm2nlWBiLW+DvXvquG2rUTSf0/8wX/GYiWwD9aWhtLPANy4kd1MeUQpJQH0b+FdRRCVEb/trguJV4I0UAIsVcIESKEiAQ+4PHugSCj4xtAWe24rHaeXhxAiJQyLguJAVLKolLK4lLKhlLKXZnkrwh8rHUJRAghItA/VMpKKa+if7OeDNwXQmxI6c4B+qN/Y7wshDgusjcQHpTmvCLwo1H5l9A/1EqlkzfV/ZFSRqN/634hq0KzqI+pGiqif8O9Y5R2EfqWBOjv3bWstBhdq3uae98EvXEqC4RLKWOM0t9I7yJZlFsWCJNSRqW5jvH9umd0HAuESil1RuegN8QppP3t2qH/faf9bpK1tFl+NxoPZOoxvEdauSXQG+e05aaUswf4DvgeuCeEWCyEKGJimbkGZSByFquAd9EPTu+UUhr/k60DfgbKSymdgIXom77GlDc6rgDc1o5vo39wpBcH+rewpyFt/iBgmmZQUj4OUsr1AFLKdVLKJpomib57CinlFSllT/QPxhnAFiGEI/ouBoeUi2t98yVM0NAujYYCUspb6ehPdX+0Ml3Qv6lmXfkM6pMNDUHoWxDFjdIVkVK6GcVXIX3Sq/fqNGU6Sim/BO4Azlr9UqiQSdUyKvc2UEwIUTjNdUy6XxmQ9rebiL61mva7EVralLIeYfTbAEqbWF4I+tZU2nINSCnnSim9ADf0Ly6jTLx2rkEZiJzFKvRdFQOAlWniCqN/a4sTQtQH3k4nv58QwkEI4Qb0BTZq4euBiUKIEkKI4sCngDnnvi8BPtBaPUII4Sj0g+yFhRDVhRAthX6wNw7926QOQAjxjhCihPaWGKFdS4e+qV9Au4YdMBF9H3pmLASmCSEqatcuIYTomEHadUBfIURtTdcXwFEpZWBWFc2sPqZqkFLeAXYC3wghigj94HYVIUQzLckPwCdCCC/tflZNuSb6N/XKRpdbA3QQQrQRQtgKIQpog7XlpJQ3gBPAZ0KI/EKIJui7fzJiqXZffDRNLwghakgpg9CPY03Xru+BvvW3Nqv7lQnvCCFqCiEc0I8ZbNFaHJuA1zUNdsDH6I3pYS3faeBtra5t0XdjZol27W3AZO1/pib6bjgAhBDe2u/XDv0LShz/fa95BmUgchDaA+kw+kGzn9NEfwhMEUJEoX/Ab0rnEvvQD+DtBr6WUu7UwqeifzCcBc4Bp7QwsyClPIHeyH0HhGuafLVoe+BL9G+Hd9G3FsZrcW2BC0KIaGAO8JaUMk5KGYm+/j+gf3OMAVLNakqHOejv4U7tngWgHzRPT+9u9H3dW9G/ZVcB3jKxupnVx2QN6FuO+YGL6O/ZFvTdQkgpNwPT0BuyKP6b9Qb6vv+JWnfSJ9rDu6OmIQR9K2AU//2vv61pCAMmoc2SSw8p5TH0Lxqz0I+P7eO/t/me6MdRbgM/ApOklP/L8C5lzWr04wd30Q98D9M0/A28g37iQCh6g9ZBSpmg5RuuhUUAvdDfG1MZgr676a5WtvEMviLoX3TC0Xc9PQC+znatnnOElGrDoOcdIYQr+lkidtKCayYUimeBEMIf/UD5D9bWokiNakEoFAqFIl2UgVAoFApFuqguJoVCoVCki2pBKBQKhSJdrO5Y7WkoXry4dHV1faK8MTExODo6Zp3Qwihd2UPpyj45VZvSlT2eRtfJkydDpZRp1xI9jrWXcj/Nx8vLSz4pe/fufeK85kTpyh5KV/bJqdqUruzxNLqAE1K52lAoFArFk6IMhEKhUCjSRRkIhUKhUKTLcz1InR6JiYkEBwcTF5e5c1InJycuXbpkIVWmo3RljxRdBQoUoFy5ctjZ2VlbkkKRa8h1BiI4OJjChQvj6uqKeGwHz/+IioqicOHCGcZbC6Ure0RFRVGoUCEePHhAcHAwlSpVsrYkhSLXkOu6mOLi4nBxccnUOChyF0IIXFxcsmw1KhSK7JHrDASgjEMeRH3nCsWzJ1caCIVCocityKQk7s34CpuwMLOXpQyEBfniiy8Mx4GBgbi7u5u1PEuUoVAoLIdMTOTWqFGELV+O/fnzZi9PGQgLYmwgFAqFIjvIxERuffwJUb//QclRo4ht2tTsZSoDYSbWrFlD/fr1qV27Nu+//z6jRo0iNjaW2rVr06tXLwB0Oh0DBgzAzc2NV199ldhY/V7uS5YswdvbG09PT7p27cqjR48A8PX1ZdiwYTRq1IjKlSuzZcsWQO8uZdSoUbi7u1OrVi02btyYviiFQvFcIhMSCB4xgqidOyk1biwu/ftZpNxcN83VmM9+ucDF2w/TjdPpdNja2mb7mjXLFmFSB7dM01y6dImNGzdy6NAh7Ozs+PDDD6lVqxYFCxbk9OnTgL7758qVK6xfv54lS5bQo0cPtm7dSseOHenSpQsDBgwAYOLEiSxdupShQ4cCcOfOHQ4ePMjly5d544036NatG9u2beP06dOcOXOG0NBQvL29aWqBtwuFQmF+khMSuDX8I6L37qXUxIkUe6eXxcrO1QbCWuzevZuTJ0/i7e0NQGxsLCVLlnwsXaVKlahduzYAXl5eBAYGAnD+/HkmTpxIREQE0dHRtGnTxpCnU6dO2NjYULNmTe7duwfAwYMH6dmzJ7a2tpQqVYpmzZpx/PhxPDw8zFxThUJhTpLj4wkeNoyYffspPelTnHv2tGj5udpAZPamb86FX1JK+vTpw/Tp01OFf/116j3P7e3tDce2traGLiZfX1+2b9+Op6cnK1aswN/fP908UtvsKeWvQqHIPSTHxRE8ZCgxBw9SespnOPfoYXENagzCDPj4+LBlyxbu378PQFhYGDdu3MDOzo7ExMQs80dFRVGmTBkSExNZu3ZtlumbNm3Kxo0b0el0hISEsH//furXr//U9VAoFNYhOTaW4A8/JObQIcpMm2oV4wDKQJiFmjVrMnXqVF599VU8PDxo3bo1d+7cYeDAgXh4eBgGqTPi888/p0GDBrRu3ZoaNWpkWV7nzp3x8PDA09OTli1b8tVXX1G6dOlnVR2FQmFBdNExBL3/ATFHAijzxRcU7drVemJM2TTiaT6ALfAXsEM7rwQcBa4AG4H8Wri9dn5Vi3fN6trpbRh08eJFkzbMePjwoUnpLI3SlT2MdZn63VuCnLrJjJQ5V5vSJWVSRIS83uNNebGmm4z4+edM0+aWDYOGA8ZuQGcAs6SULwLhQH8tvD8QLqWsCszS0ikUCkWeICksjBu+fYm7eJEXZs/CqUMHa0syr4EQQpQDXgd+0M4F0BLYoiVZCXTSjjtq52jxPkI52FEoFHmAxHv3uPFObxKuX6fc/PkUad3a2pIAENKMM2CEEFuA6UBh4BPAFwjQWgkIIcoDv0sp3YUQ54G2UspgLe4a0EBKGZrmmgOBgQClSpXy2rBhQ6oynZycqFq1apbannQdhLlRurKHsa6rV68SGRlpZUV6oqOjKVSokLVlpEtO1ZZXddmEhuI8ezY20TFEDP6QxBdfNLuuFi1anJRS1ssqndmmuQoh2gP3pZQnhRDNU4LTSSpNiPsvQMrFwGKAevXqyebNm6eKv3TpkknTV3Py/gZKl+kY6ypQoAB16tSxsiI9/v7+pP1t5hRyqra8qCv+33+5OWkyyYlJVFi9ioK1auUIXSmYcx1EY+ANIcRrQAGgCDAbKCqEyCelTALKAbe19MFAeSBYCJEPcALM765QoVAorEDc5cvc7NcfhKDiqpUUqF7d2pIew2xjEFLKcVLKclJKV+AtYI+UshewF+imJesD/KQd/6ydo8Xvkebs/1IoFAorEXvmDDfe7YPIn5+Ka1bnSOMA1lkHMQYYKYS4CrgAS7XwpYCLFj4SGGsFbVbn7Nmz/Pbbb4bzyZMnP7YC+1ljiTIUCoWemGPHuNm3H7ZFi1JxzRrsc/A2uRZxtSGl9Af8teN/gceW+Uop44DultCTkzl37hznz5/ntddes7YUhULxjInas4dbI0ZiV74cFZYuw67U4z7achJqJbUZCAwMpEaNGrz33nu4u7vTq1cvdu3aRePGjXnxxRc5duwYMTEx9OvXD29vb+rUqcNPP/1EQkIC06ZNY+PGjdSuXdvgtvvixYs0b96cypUrM3fuXEM5nTp1wsvLCzc3NxYvXmwIL1SoEBMmTMDT05OGDRsanPrduHEDHx8fPDw88PHx4ebNm5a9MQpFHiZi248EDx2GfbVqVFy9OscbB8jlzvr4fSzcPZduVEFdEtg+QfVL14J2X2aZ7OrVq2zevJnFixfj7e3NunXrOHjwID///DNffPEFNWvWpGXLlixbtoyIiAjq169Pq1atmDBhAufPn+e7774D9N0/ly9fZu/evURFRVG9enUGDRqEnZ0dy5Yto1ixYsTGxuLt7U3Xrl1xcXEhJiaGhg0bMm3aNEaPHs2SJUuYOHEiQ4YM4d1336VPnz4sW7aMYcOGsX379uzfA4VCkS0eLF3K/Zlf49ioEeXmzcXG0dHakkxCtSDMRKVKlahVqxY2Nja4ubnh4+ODEIJatWoRGBjIzp07+fLLL6lduzbNmzcnLi4uwzf6119/HXt7e4oXL07JkiUNLYK5c+caWglBQUFcuXIFgPz589O+fXsgtRvxI0eO8PbbbwPQu3dvDh48aOa7oFDkbaSU3Js5k/szv6Zwu7aUW7jguTEOkNtbEJm86ceaeV6/sVtuGxsbw7mNjQ1JSUnY2tqydetWqqeZvbBv375Mr2Vra0tSUhL+/v7s2rWLI0eO4ODgYDAyAHZ2dqQsQk9Jnx5qobpCYT5kUhJ3Pp1E5LZtFO35FqUnTkTkwMWmmaFaEFaiTZs2zJs3z7CXw19//QXoxw+ioqKyzB8ZGYmzszMODg5cvnyZgICALPM0atSIlJXna9eupUmTJk9RA4VCkRHJcXEEDxtO5LZtFB88mNKffvrcGQdQBsJq+Pn5kZiYiIeHB+7u7vj5+QHwyiuvcPHixVSD1OnRtm1bkpKS8PDwwM/Pj4YNG2ZZ5ty5c1m+fDkeHh6sXr2aOXPmPLP6KBQKPbqoKILeG2DYIrTE0CHPb2vdFJevOfWj3H1bjudBl3L3bRo5VVtu0JV4/7681rGTvOheS0bs2GE+UdIy7r5z9xiEQqFQWIiEGze4OWAgSSEhlF+wgEJNGltb0lOjDIRCoVA8JbFnzhD0wSCQkoorllPQ09Pakp4JagxCoVAonoKo3bu50ccXm8KFcd2wPtcYB1AGQqFQKJ6YsHXr9KujX3wR1/XryO/qam1JzxTVxaRQKBTZRCYnEzJrFg+W/EChFi144ZuvsXFwsLasZ44yEAqFQpENkhMSuDN+Ag937KDoW2/qF8Dly52PUtXFZAYaNWpkbQkmsWLFCoYMGWJyen9/f4MLj+yS3a0RlQtyRU5E9/AhQQMG8nDHDkqMHEnpSZNyrXEA1YIwC4cPHzbbtZOSksiXi3+QCkVOJfHOHYIGvk98YCBlv5qB0xtvWFuS2VEtCDOQ8rZ8584dmjZtSu3atXF3d+fAgQMArF+/nlq1auHu7s6YMWNS5Rs/fjx169bFx8eHkJAQAJo3b8748eNp1qwZc+bMydBtt6+vLx988AGvvPIK1apVY8eOHYB+dfbp06cN5TRu3JizZ89mqP/YsWM0atSIOnXq0KhRI/7+++/H0kRHR9O3b19q1aqFh4cHW7duzbRugHJBrnhuibt8mcC3epJ45w4VFi/KE8YBcnkLYsaxGVwOu5xunE6nw/YJfKPUKFaDMfXHZJ0QWLduHW3atGHChAnodDoePXrE7du3GTNmDCdPnsTZ2ZlXX32V7du306lTJ2JiYvD09GTevHlMmTKFzz77zOD2OyIiwuDIr0OHDhm67Q4MDGTfvn1cu3aNFi1acPXqVd577z1WrFjB7Nmz+eeff4iPj8fDw4NTp06lX8caNdi/fz/58uVj165djB8/nhUrVqRK8/nnn+Pk5MS5c3p36uHh4VnWTbkgVzyP5D97lsCRH2NbuDAV167JsduDmgPVgjAj3t7eLF++nMmTJ3Pu3DkKFy7M8ePHad68OSVKlCBfvnz06tWL/fv3A3pPr127dgXgnXfeSeWO+8033zQcZ+a2u0ePHtjY2PDiiy9SuXJlLl++TPfu3dmxYweJiYksW7YMX1/fTHVHRkbSvXt33N3dGTFiBBcuXHgsza5duxg8eLDh3NnZOdO6KRfkiucNKSVhK1dSdMFC7CtVwnXTpjxlHCCXtyAye9OPMrO7b4CmTZuyf/9+fv31V3r37s2oUaMoUqSIyfmNHXw5ZuJD3jhdWqdgQggcHBxo3bo1P/30E5s2beLEiROZluvn50eLFi348ccfCQwMpHnz5o+lkVI+VpbUPNOmh3JBrniekImJ3J02jYgNG4mvU5saS5fmymmsWaFaEGbkxo0blCxZkgEDBtC/f39OnTpFgwYN2LdvH6Ghoeh0OtavX0+zZs0ASE5ONnSvrFu3LkN33Jm57d68eTPJyclcu3aNf//917DfxHvvvcewYcPw9vamWLFimeqOjIzkhRdeAHisaymFV1991dD9BfoupszqlhHKBbkip6GLiiLog0FEbNiIy4D3iBwwIE8aB1AGwqz4+/tTu3Zt6tSpw9atWxk+fDhlypRh+vTptGjRAk9PT+rWrUvHjh0BfSvh0qVLeHl5sWfPHj799NN0r5uZ2+7q1avTrFkz2rVrx8KFCylQoACg79YpUqQIffv2zVL36NGjGTduHI0bN0an06WbZuLEiYSHh+Pu7o6npyd79+7NtG4ZoVyQK3ISCcHBBPbsSczRo5SZNpWSH38MNnn4MWmKy9ec+slt7r4dHR2fSlefPn3k5s2b0427deuWfPHFF6VOp3uia+fE+yWlcvf9JORUbdbWFXPylPz75Ubycv0GMjrgqCHc2roywhLuvvOwacw7rFq1igYNGjBt2jRs8vLbkEKRAZE7fuWmry82hQvhumE9jg3qW1tSjiBXD1I/b0RHR5u03WhGZDRe8O677/Luu+8+8XUVityKTE4m9LvvCZ0/H4d69Xhh3lzyOTtbW1aOQRkIhUKRJ9FFx3B77Biid+3GqUsXykyehMif39qychTKQCgUijxHQlAQwR8OJv7ffyk1fjzOvd9RU6zTQRkIhUKRp4gJCODW8I+QQIUli3F8TpxrWgM1YqlQKPIEUkrCVq/hZv/3yFeyBJU2b1LGIQuUgcihnD59mt9++83aMhSKXEFyQgJ3/Py4N20ahZo1o+L6DeSvUMHasnI8ykDkUJSBUCieDUmhodzs40vklq24DPqAct/Nw7ZQxq5rFP+hDISZ6NSpE15eXri5ubF48WIA/vjjD+rWrYunpyc+Pj4AxMTE0K9fP7y9valTpw6//vorCQkJfPrpp2zcuJHatWuzceNGwsLC6NSpEx4eHjRs2NDgrnvy5Mn069eP5s2bU7lyZebOnWvQsGbNGurXr0/t2rV5//33M1wVrVDkVmLPX+B6t+7EXbrEC7NnUXL4cIRaC2QyuXqQ+u4XXxB/KX1330k6HWFP4O7b/qUalB4/Pst0y5Yto1ixYsTGxuLt7U3Hjh0ZMGAA+/fvp1KlSoSFhQEwbdo0WrZsybJly4iIiKBevXp06NCBKVOmcOLECYO/o6FDh1KnTh22b9/Onj17ePfddw17PFy+fJm9e/cSFRVF9erVGTRoEFevXmXjxo0cOnQIOzs7PvzwQ9auXavWQyjyDBFbt3L3synYFnfBdd1aCtSsaW1Jzx252kBYk7lz5/Ljjz8CEBQUxOLFi2natCmVKlUCMDjM27lzJz///LNhe834+Ph0N805ePCgYVOeli1b8uDBAyIjIwF4/fXXsbe3x97enpIlS3Lv3j12797NyZMn8fb2BiA2NpaSJUuat9IKRQ4gOSGBe1OnEbFpEw4vN+SFb74hXxYOKhXpk6sNRGZv+uZ09+3v78+uXbs4cuQIDg4ONG/eHE9Pz3R3ZpNSsnXrVoPX1RRdR48efSxdWlLmbdvb2xvCUlxpSynp06cP06dPf5ZVUyhyNIl37hA8/CPizp7FZcAASnw0HPEEPQUKPaozzgxERkbi7OyMg4MDly9fJiAggPj4ePbt28f169cBDF1Mbdq0Yd68eQYDcObMGQAKFy6cyu1G06ZNWbt2LaA3QMWLF890bwkfHx+2bNnC/fv3DeXduHHj2VdWocghxAQEcL1rNxKuXeOFeXMp+fFIZRyeEmUgzEDbtm1JSkrCw8MDPz8/GjZsSIkSJVi8eDFdunTB09PTsEOcn58fiYmJeHh44O7uztSpUwFo0aIFFy9eNAxST548mRMnTuDh4cHYsWNZuXJlphpq1qzJ1KlTefXVV/Hw8KB169bcuXPH7HVXKCyNlJIHS5dys19/bJ2dcd28iSKtW1tbVq4gV3cxWQt7e3t+//33dOPatWuX6rxgwYIsWrTIcJ7SaihWrBjHjx9Plfann3567HqTJ09OdX7+/HnD8Ztvvplqq1KFIrehi47hzoQJRP35J4XbtKHMtGlqCuszRBkIhULxXBL/778EDx1GwvXrlBw1imL9+ip/Ss8Ys3UxCSEKCCGOCSHOCCEuCCE+08IrCSGOCiGuCCE2CiHya+H22vlVLd7VXNoUCsXzTeSOXwns1h1deDgVli3DpX8/ZRzMgDnHIOKBllJKT6A20FYI0RCYAcySUr4IhAP9tfT9gXApZVVglpZOoVAoDCTHx3Nn0mRuf/IJ9jVqUGnbVhwbNrC2rFyL2QyEtrNdtHZqp30k0BLYooWvBDppxx21c7R4H6FeCRQKhUbCjRsEvtWTiI0bcXmvPxVXrsCudGlry8rViPTm16dKIMTIdIIjgZNSytNZ5LUFTgJVge+BmUCA1kpACFEe+F1K6S6EOA+0lVIGa3HXgAZSytA01xwIDAQoVaqU14YNG1KV6eTkRNWqVTOtE4BOp8M2B06BU7qyh7Guq1evGhYPWpvo6GgKFSpkbRnpklO1ZabL/uRJiqxeA7a2RPbpQ4JHrRyhy5o8ja4WLVqclFLWyzJhVptWA+uAf4BvtM9lYDVwHBhtysbXQFFgL/AKcNUovDxwTju+AJQzirsGuGR2XS8vr8c24zZ143rjze5zEkpX9jDWZep3bwly6kb3UuZcbenp0sXHyzufTZEXq9eQ13u8KRNu3coRunICT6MLOCFNeHab0sXkAtSVUn4spfwYqAeUAJoCvqZYKyllBOAPNASKCiFSZk+VA25rx8GawUCLdwLCTLn+84q/vz/t27e3tgyFIkeSEBTEjZ4+Q2g2AAAgAElEQVRvE75uHcV8fam4ehV2ZctaW1aewhQDUQFIMDpPBCpKKWPRD0SnixCihBCiqHZcEGgFXELfkuimJesDpEzu/1k7R4vfo1k6hUKRx3i4cyfXu3QlISiIct9/R6mxY9R+0VbAFAOxDggQQkwSQkwCDgHrhRCOwMVM8pUB9gohzqLvjvqflHIHMAYYKYS4ir51slRLvxRw0cJHAmOfqEZW5vjx43h4eBAXF0dMTAxubm507tw51SK3Xr168fPPP6fKl+L2u1mzZtSpU8eQfsWKFXTp0oW2bdvy4osvMnr0aEOenTt38vLLL1O3bl26d+9OdHQ0CsXzTHJcHHcmT+bWsOHkd3Wl0ratFNZc4yssT5YL5aSUnwshfgOaAAL4QEp5QovulUm+s0CddML/BeqnEx4HdDdRt0kc2PQPoUHpPzSfdNC1ePlCvNKjWobx3t7evPHGG0ycOJHY2FjeeecdGjVqxKxZs+jYsSORkZEcPnyYlStXcvDgQUO+FLffc+bMQafTUb9+fVq1agXoNw/666+/sLe3p3r16gwdOpSCBQsydepUdu3ahaOjIzNmzODbb7/l008/zXadFIqcgO2tWwR27078lasU69uXEiM+wka1GqxKhgZCCGHsH/e69jHESSlz9fjA0/Dpp5/i7e1NgQIFmDt3Lra2tgwePJj79++zbds2unbtSr58qW99itvvr776ChsbG+Li4gxuv318fHBycgL0PpZu3LhBREQEFy9epHHjxgAkJCTw8ssvW7aiCsUzQEpJ+Pr1uHw5gyQnJ8ovWUKhV5pYW5aCzFsQJ9GvWxDou4tSBpOFFl7ZvNKensze9M3p7jssLIzo6GgSExOJi4vD0dGR3r17s3btWjZs2MCyZcseyyM1t99ly5ZNpevo0aMZuvNu3bo169evN0sdFApLkBQezp2JfkTv3k2CW02qLVpEvuLFrS1LoZHhGISUspKUsrKUshJwSTuunBJuQY3PHQMHDuTzzz+nV69ejBkzBgBfX19mz54NgJub22N50rr9/uuvvzIto2HDhhw6dIirV68C8OjRI/75559nWQ2FwqzEHD3G9U6did6/n5JjxxAxeLAyDjkMU1dSq9lEJrJq1Sry5cvH22+/zdixYzl+/Dh79uyhVKlSvPTSS/Tt2zfdfCluv19++WXc3d3x8/PLtJwSJUqwYsUKevbsadin+vLl9LdXVShyEjIxkfuzZnPT1xebggVx3bAeF19fUHtF5zxMWSwBnDIlnaU/z9NCuZiYGFm5cmUZERGRabrnYUFaTkItlMs+1tQWf+OGvN7jTXmxeg15a9x4qYuOzhG6MiM36sLEhXKZDVIbu9gomdblhpTyW/OYrNzHrl276NevHyNHjjQMNisUeQkpJRFbtnBv+pcIW1vKfvM1Tq+/bm1ZiizIbJDaeAR3SZpzRTZo1aqVYUaSQpHXSAoN5Y7fp0Tv3YtDw4aUnf4FdmXKWFuWwgQyNBBSys8sKUShUOQ+ovbs4c5EP5Kjoyk1bizOvXsj1FjDc0OWC+WEEJWAoYCrcXop5Rvmk6VQKJ5ndNEx3PtyOpFbtmL/0ku8sHIF9i++aG1Zimxiypaj29G7wfgFSDavHIVC8bzz6NQpbo8ZS+KtW7gMHEiJIYOVH6XnFFMMRJyUcq7ZlSgUiucamZBAyPfzebBkCXZly1Jx9SocvLysLSv3kZwMAfOxS6hg9qJM6Qycoznqe1kIUTflY3ZlzykRERHMnz8/0zSBgYG4u7tbSJFCYX5iL1zgerfuPFi0CKcunam0fbsyDuYgJhTWdoWdEyh1b6/ZizOlBVEL6I1+q9CULqaUrUMVaUgxEB9++KG1pSgUZkcmJBC6cBGhixeTz9mZcvPnU7hlC2vLyp3cOAxb+sGjMOgwh+CHFcl678ynw5QWRGegspSymZSyhfZRxiEDxo4dy7Vr16hduzajRo1i1KhRuLu7U6tWLTZu3PhY+hUrVjBkyBDDeffu3fH39wcyduft6urKpEmTqFu3LrVq1TKsoE5xGe7t7Z3KZbhCYQ7iLl7kevcehM6fj9Prr1N5xy/KOJiDZB3smwkr2oOdAwzYDV6+IITZizalBXEG/Zah982s5Zmzd8Vi7t/4N904XZIO23zZd/ddsmJlWvgOzDD+yy+/5Pz585w+fZqtW7eycOFCzpw5Q2hoKN7e3jRt2tSkckJDQzN15128eHFOnTrF/Pnz+frrr/nhhx8MLsOXLVtGRESEwWW4o6NjtuupUGSEajVYkIe3YdtACDwAtbrD699CgSIWK94UA1EKuCyEOI7RDnJqmmvWHDx4kJ49e2Jra0upUqVo1qyZYUOhrAgICMjUnXeXLl0A8PLyYtu2bcB/LsO//vprAIPL8JdeeulZV02RR4m7dInb48YTf/kyTh3foNS4cdgWLWptWbmTv/+A7YMgKR46LQDPnhZpNRhjioGYZHYVZiKzN31zuvtOQZqwY2q+fPlITv5v9nB8fLwhb2buvFNcgKe4/07Js3XrVqpXr/600hWKVMiEBEIXLSZ00SJsnYuqVoM5SYqH/30KRxdCaQ/otgyKW2cNSZZjEFLKfVLKfcBfwDmjjyIdChcuTFRUFABNmzZl48aN6HQ6QkJC2L9/P/Xrp95Mz9XVldOnT5OcnExQUBAnT54Ensydd3ZdhisUphB75gzXu3Yj9PvvKfJaO6r8osYazEboFfjBR28cGgyC93ZZzTiAaSupBwKfA7HoZzE9NxsGWQMXFxcaN26Mu7s77dq1w8PDA09PT4QQfPXVV5QuXZrAwEBD+saNG1OpUiVq1aqFu7s7np6eQGp33imtiqlTp1KtWsabIPn5+fHRRx/h4eGBlBJXV1d27Nhh1voqci/JMTHcnzOH8NVryFeqlGo1mBMp4fQ6+G0U5LOHnhuheltrqzKpi2kU4CalDDW3mNzCunXrUp3PnDkz1bmrqyvnz58HQAjB2rVrDXHGXV8tW7bk+PHjj13f2MDUq1fPMOupYMGCLFq06FlUQZHHiT5wgDuTJpF05y7OPXtSYuQIbAsVsras3EncQ/j1Yzi3CVxfgS6LoUhZa6sCTDMQ14BH5haiUCisT1JYGPemf8nDX34hf5UqVFy7Boe6al2s2Qg8BD9+AA+DocVEeGUk2GR/dqW5MMVAjAMOCyGOknoW0zCzqVIoFBZFSsnDX37h3hfT0cXEUHzwYFzeH4iN8qFkHpLiYc9UODwPilWCfjuhvLe1VT2GKQZiEbAH/cC0ctanUOQyEoJvcXfyZGIOHqSgpydlpn6uPK+ak3sX9Gsb7p0Hr77w6lSwz5ndd6YYiCQp5ciskykUiucJmZDAg+UrCF2wAGFjQ6mJE3Hu+RbCNud0ceQqkpMh4HvYPQUKFIW3N0G1NtZWlSmmGIi92kymX0jdxRRmNlUKhcKsxBw9xt0pU0i4do3CrVtTavw4tcubOYkI0i96CzwANdpDhzngWNzaqrLEFAPxtvZ3nFGYmuaqUDyHJD14QJHlK7h59Ch25cpRftFCCjVrZm1ZuRcp4ewm+O0TkMnQ8Xuo3cviK6KfFFMWylVK56OMg4lMnjzZ4PriWWGKS/GnpVA2pzS6uroSGpr9mdC+vr5s2bLF5PTKVfqTIZOTCd+wkWvtXqPAiRO4vP8+lX/5WRkHcxJ9Hzb1hh8HQik3GHQI6rzz3BgHMK0FochhPKlLcZ1Oh63qX85zxF28yJ3PPiPuzFkc6tfnZru21OzZ09qyci9Swvmt+lZDwiNo9Rk0Gpqjpq+aito93AxMmzaN6tWr06pVK/7++29D+LVr12jbti1eXl688sorBjfdvr6+DBs2jEaNGuHh4WF4o46OjsbHx8fg1jvFfXdal+L+/v60b9/eUM6QIUNYsWIFoH+znzJlCk2aNGHz5s0Zarh+/Tovv/wy3t7e+Pn5ZVi3Tp064eXlhZubG4sXL043zapVqwwryHv37g3AjRs38PHxwcPDAx8fH27evGlIv3//fho1akTlypUNdZdSZukqXZE5uocPufvFF1zv1p3E4FuU/WoGFVauQKfGGsxH1F3Y0Au29odiVeCDA9Dko+fSOAD6f8T0PkC+jOJyysfLy0um5eLFi4Zj9GMlZvlkxIkTJ6S7u7uMiYmRkZGRskqVKnLmzJlSSilbtmwp//nnHymllAEBAbJFixZSSin79Okju3XrJnU6nTx27JisUqWKlFLKxMREGRkZKaWUMiQkRFapUkUmJyfL69evSzc3N0OZe/fula+//rrhfPDgwXL58uVSSikrVqwoZ8yYYYjLSEOHDh3kypUrpZRSfvfdd9LR0TFVvR4+fCillPLBgwdSSikfPXok3dzcZGhoqKGckJAQef78eVmtWjUZEhKSKn379u3lihUrpJRSLl26VHbs2PGxul+4cMFQ9y1btshWrVrJpKQkeffuXVm+fHl5+/btx+qeoivtd29t9u7da7Wyk3U6Gb55s/z75UbyYo2X5O3Jk2VSRESO0JYZz7Wu5GQpT6+XcnoFKT8vKeXBOVLqkqyvKwOAE9KEZ2xmXUwBQohg4A/gDyll4LM0TLmVAwcO0LlzZxwcHAB44w29V/To6GgOHz5M9+7dDWlTfCyB/s3cxsaGGjVqcO/ePUBvvMePH8/+/fuxsbHh1q1bhrjs8Oabb2ap4dChQ2zduhWA3r17M2bMmHSvNXfuXH788UcAgoKCuHLlCi4uLob4PXv20K1bN4oX18/QKFasGABHjhwxuCXv3bs3o0ePfqzuNWvWNNTvaVyl52Viz5zh7tRpxJ07R8G6dSn9wxIK1KxpbVm5m4d3YMdH8M8fUL6BfiDaig72niUZGggpZT0hREWgHTBbCPECcBD4HdgnpYzPKG9eR6QzCJWcnEzRokU5ffp0unlS3HfDf27C165dS0hICCdPnsTOzg5XV1fi4uIey5vWZXjaNCkbBmWlIT3dxvj7+7Nr1y6OHDmCg4MDzZs3f6wsKWWW10lbVnp1T/mrMI2k0FDuf/MtkT/+SL4SJSg78yuKtG9v0neheEKkhDPr4Y+x+pXRbb6ABh88v91J6ZDpGISU8oaUcqGUshPQCP1aiFbAASHEr5YQ+DRk1nR6+PDhU3VvZUTTpk358ccfiY2NJSoqil9++QWAIkWKUKlSJTZv3mzQdubMmUz1R0ZGUrJkSezs7Ni7dy83btwAUrsUB6hYsSIXL14kPj6eyMhIdu/ene71MtPQuHFjNmzYAJDKeWBaPc7Ozjg4OHD58mUCAgIeS+Pj48OmTZt48OABAGFh+uUyjRo1SnX9Jk2aZFp3U1ylK0AmJvJgxQqutW1H5I4duAx4j8q//45Thw7KOJiT8EBY202/tqFkTRh0GF4enKuMA2RjkFpKmSil3COlHC2lrA9kvBtPHqZu3bq8+eab1K5dm65du/LKK68Y4tauXcvSpUvx9PTEzc0tyz2je/XqxYkTJ6hXrx5r166lRo0aQGqX4qNGjaJ8+fL06NEDDw8PevXqRZ06dTK8ZkYa5syZw/fff4+3tzeRkZHp5m3bti1JSUl4eHjg5+dHw4YNH0vj5ubGhAkTaNasGZ6enowcqV+EP3fuXJYvX46HhwerV69mzpw5mda9c+fOhoHuli1bGlylK/4j5vBh/u3UmftfzqBg3TpU/vknSn78MbaF1BazZkOXBIfmwvcN4cYRaDsDfH8DlyoWlXH16lWCg4PNX9DTvEVb+5PVIHVmGA9u5iSUruyRFwep4679K29+MEherF5DXmn9qny4Z49MTk7OEdqehhyvK/iklAsaSzmpiJRr35Qy/KbFtcTFxckpU6ZIe3t76eHhIXU63RNdh2cwSK1QKHIQSeHhhH4/n/ANG7Cxt6fkJx/j3Ls3NkZjOIpnj21SLPwxTr/Lm2NJ6LEKXnrD4gve9u7dy6BBgwxT58+ePcvKlSvp27ev2co02UAIIRyllDFmU6JQKNJFJiQQtm4dofMXkBwdTdEe3SkxdCj5jGaPKczE33/gfXwIxIdCvf7QahIUcLKohJCQED755BNWrVqVKrxatWpmn9VnypajjYAfgEJABSGEJ/C+lDJ7y3gtiDRxJo0i9yBz4awnKSVRu3Zx/+uvSbxxE8cmTSg1ZrRyxW0Jou7C72Pg4nZ0DuWh359Q4fExN3Oi0+lYsmQJ48ePJzw83BBeuHBhvvjiC1566SW8vLzMqsGUQepZQBvgAYCU8gzQNKtMQojyQoi9QohLQogLQojhWngxIcT/hBBXtL/OWrgQQswVQlwVQpwVQjzRNlYFChTgwYMHufKBoUgfKSUPHjygQIEC1pbyzIi9cIGb7/bh1tBhCDs7yi9ZTIUflijjYG50SRCwAL7zhr9/h5YTOVFvlsWNQ0BAAPXr12fQoEGpjEP37t25fPkyQ4YMsYjbHJO6mKSUQWneyHUmZEsCPpZSnhJCFAZOCiH+B/gCu6WUXwohxgJjgTHo11u8qH0aAAu0v9miXLlyBAcHExISkmm6uLi4HPlAUbqyR4quAgUKUK5cOWvLeWoSgoMJmTOXh7/8gm2xYpSePImi3boh8qnhQrNz4zD8+gncvwBVWkK7mVC8KlLb890S3L9/n7Fjx7J8+fJU4ZUrV+a7776jXbt2FtMCphmIIK2bSQoh8gPDgEtZZZJS3gHuaMdRQohLwAtAR6C5lmwl4I/eQHQEVmkj7AFCiKJCiDLadUzGzs6OSpUqZZnO398/0+mg1kLpyh45VVd2SQoPJ3TBAsLXb0DY2uIyYAAuAwdgW7iwtaXlfqLuwf8+hbMbwKk89FgNL3Ww6CB0UlISCxcuxM/Pj4iICEN4gQIFGDduHKNHj7bKC5rIqitGCFEcmIN+gZwAdgLDpZQPTC5ECFdgP+AO3JRSFjWKC5dSOgshdgBfSikPauG7gTFSyhNprjUQbQ1GqVKlvFIWX2WX6OjobLu0tgRKV/Z47nXFx+O4ew8OO3ci4uOJbdyYmNdfJ9m5aNZ5za3Nwlhal0jWUfb2b1S6vg6b5ASCynfiRsXuJNumfhCbW9e5c+eYM2cO165dSxXeuHFjBg8eTJkMnCs+ja4WLVqclFLWyzJhZnNgAVtghCnzZTO5RiHgJNBFO49IEx+u/f0VaGIUvhvwyuza6a2DMJUcP+c6h6F0ZY+sdCUnJMiw9Rvk302ayIvVa8ibgwfLuKtXc4Q2a2FRXYGHpZzfSL+mYVUnKUOuWFzXrVu3ZO/evR9zBFq1alX522+/ZZnf2s76kFLqhBAd0Q9UZxshhB2wFVgrpdymBd9L6ToSQpQB7mvhwUB5o+zlgNtPUq5CkVORUhK183+EzJpFQmAgBb28KDlnLg51n/9usueCh3dg12R9d1KRclZZ0xAbG8s333zDl19+SUzMfysHChYsyIQJE/jkk09S+SezJqaMQRwSQnwHbAQMtZFSnsosk9CPai8FLkkpvzWK+hnoA3yp/f3JKHyIEGID+sHpSJnN8QeFIqcipSTm4EFC5swl7vx58letQrn58ynUormakm0JEmPhyHdwYBYkJ0KTEdB0FOS3nFsSKSUbNmxgzJgxBAUFpYrr2rUr33zzDRUrVrSYHlMwxUA00v5OMQqTQMss8jUGegPnhBAp7kPHozcMm4QQ/YGbQIrv6d+A14CrwCPAfMsDFQoLEhNwlJA5c4j96y/sXniBMl98gVPHNxBqdz/zIyVc2Ab/mwSRQfrB59ZToJhld00OCAhgxIgRjzm4dHd359tvv6V169YW1WMqWRoIKWWLJ7mw1A82Z/Rq5JNOegkMfpKyFIqcyKNTfxEydy6PAgLIV6oUpSdPpmiXzoj8+a0tLW9w65TeRUZQAJSqBZ3mQ6Usl3A9U4KCghg7dizr1q1LFV6iRAk+//xz+vfvT74cPIXZlJXUn6YXLqWckl64QpHXiT1/gaLzvuPGhQvYurhQavw4ir75pvKZZCke3oHdU+DMOnAsAR3mQp13LOqKOzo6mpkzZzJz5kxiY2MN4fnz52f48OFMmDABJyfLuux4EkwxXcb+lwoA7TFhHYRCkdeI+/sfQr+bR9T/dmHn6EiJj0dSrFcvbLTdBRVmJuERHPkeDn4LyUnQ+CN45WMoUMRiEhITE/nhhx/47LPPHtv9sUuXLnz11VdUqWJZ1+BPgyldTN8YnwshvkY/oKxQKID4K1cIXbCQh7//jo2jI8WHDuFi5cq4WXjVa55FlwSn14L/dIi6Y5VxBiklW7duZfz48Vy5ciVVXJ06dZg1axbNmjWzmJ5nxZN0fjkAlh3hUShyIHGXLhG6YCFRO3ciHBz0q5/79cW2aFGLumfIs0ip95e0+zMIuQzlvKHbMqjYKOu8zxB/f3/GjBnDsWPHUoWXK1eOKVOm8O6771rEb5I5MGUM4hz6WUugXzhXgtQzmhSKPEXsuXOEzl9A9N692BQqhMsH71OsTx/yOTtbW1reIei43j3GzcPgUtUq7jHOnj3L2LFj+f3331OFFy1alHHjxjF06FAKFixoMT3mwJQWRHuj4yTgnpQyyUx6FIocy6NTpwidv4CYgwexcXKi+NAhFOvdG9siluvjzvOEXoXdk+HSL/rNe17/Fuq+C7Z2FpMQGBjIpEmTWL16dSqv0fb29gwbNoyxY8dSrFgxi+kxJ6YYiHxAsJQyXgjRHOgqhFglpYzIIp9C8dwjpeTRseOEzp/Po6NHsXV2psTIkTi/3RPbHOjPKNcSdRf2fQUnV4BdQWg+Hl4eDPaW+w5u3brFtGnT+OGHH0hMTDSECyHw9fXls88+o3z58plc4fnDFAOxFagnhKiKfmX0z8A69IvaFIpciX7l8yFCFy4k9uRJbEsUp+SYMTi/2UPNSrIkMQ/g0Cw49oN+BXS9ftBsNBQqaTEJYWFhjBgxggULFhAfH58qrn379kyfPh13d3eL6bEkphiIZCllkhCiCzBbSjlPCPGXuYUpFNZAJiXx8I8/ebB0KfGXLpGvdGlK+U2kaLduah2DJYmN0LvGCFgAiY/A4029YbDgzKQHDx4wc+ZM5syZQ1xcXKq4Jk2aMG3aNJo2tezCO0tjioFIFEL0BN4FOmhhluvwUygsQHJsLBFbtxG2fDmJt26Rv0oVykybhlOH9mrlsyWJj4ajC+HwXIiLBLfO0HwclKhuMQkRERF8++23zJ49m6ioqFRx9evX5/PPP6d169Z5woeWKQaiL/ABME1KeV0IUQlYY15ZCoVlSAoPJ3ztOsLXrEEXEUHBOnUoNWE8hZo3R9iYsiOv4pmQGEu5oJ/gWD949ACqtYMW46GMh8UkREZGMm/ePL755ptUm/YA1K5dmylTptC+ffs8YRhSMGWh3EX0u8ih7R9dWEr5pbmFKRTmJCH4FmErVhCxdSsyNpZCLVrgMuA9HOo+0VboiiclMRZOrYYD31A1+i5UbgEtJ0K5rPeyeVY8ePCA2bNnM2/ePCIjI1PF1axZkx49euDn54dNHnxhMGUdhD/whpb2NBAihNgnpRxpZm0KxTMn7uJFHixbzsPffwcbG5zat8elfz/sq1a1trS8RUIMnFiu70qKvgcVGvFX1aHU6TTEYhLu37/PN998w/z584mOjk4VV7VqVSZPnsxbb73FgQMH8qRxANO6mJyklA+FEO8By6WUk4QQZ80tTKF4Vkidjui9ewlbuYpHx49j4+BAsd69KebbB7vSpa0tL28RHwXHf4DD38GjUL131W7LwLUJkRZafX779m1mzpzJokWLUjnSA6hWrRoTJkzg7bffztFeVi2FSesgtJ3fegATzKxHoXhm6KJjiNy2jbA1a0i8eZN8ZctQcvRoinbrqha3WZq4SDi6GAK+h9hwqOKjn5VUoaHFJNy8eZMZM2awdOnSx6arurm5MXHiRLp37/7cusUwB6YYiCnAn8AhKeVxIURl4EoWeRQKq5EQfIvw1auJ2LqV5OhoCtapQ8mRIyncygeh3goty6Mw/aykgIUQHwnV2kLT0VDOy2ISzp8/z1dffcX69etJSkrtBKJOnTr4+fnRsWPHPNuNlBmmDFJvBjYbnf8LdDWnKIUiu0gpiT11irCVq4jatQtsbCjSpg3F+rxLQQ/LzYRRaEQGw5H5+pXPiTFQo71+i8+ytS1SvJSSAwcOMGPGDH777bfH4hs0aICfnx+vvfZanpqVlF1MGaSuBiwASkkp3YUQHsAbUsqpZlenUGRFQgIRP24nfO1a4s6fx8bJCZf+/XHu9bYaX7AG9y/DoTlwbpPe26p7V2jyEZRys0jxycnJ/PTTT8yYMYOjR48+Ft+sWTMmTJhAq1atlGEwAVPa20uAUcAiACnlWSHEOkAZCIXVSAgKInzDBkps2MidmBjyV65M6cmTcHrjDeUKwxrcDICDs+Gf38HOAbzf0/tKKlrBIsXHx8ezevVqZs6cyT///JMqTghBly5dGDVqFA0aNLCIntyCKQbCQUp5LI21Vd5cFRZH6nTEHDxI2Lp1xOw/ADY2JHh4UG34cBwa1FdvhJYmORmu/Kk3DEEBUNBZv+rZewA4ulhEQkhICAsXLmT+/PncvXs3VVz+/Pnx9fXl448/plq1ahbRk9swxUCECiGqoO0JIYToBtwxqyqFwoik8HAit20jfP0GEoODsS1RnOKDBlH0zR4cunQJx4bqrdCiJMbC2Y16P0khl8GpPLSdAXV7Q35Hi0g4f/48s2fPZs2aNY/NSHJycmLQoEEMHz6c0qqb8akwxUAMBhYDNYQQt4DrQC+zqlLkeaSUxJ07R/i69Tz87TdkQgIO9epR8uORFPbx+c8/0iW1PbrFiLqrX8NwYpneHUbpWtB5Mbh3sch+DMnJyfzxxx/MmjWLXbt2PRZftmxZRowYwcCBAymipjE/EzI1EEIIG6CelLKVEMIRsJFSRmWWR6F4GnQPHxL58y9EbN5M/N9/IxwccOraBee3elKguuomsAp3zkLAfH+MOgoAACAASURBVDi3BZKToHo7aPghuDaxyA5uMTExrFq1ijlz5vD3338/Fl+vXj1GjBhBt27dyK8cKz5TMjUQUspkIcQQYJOUMsZCmhR5DCklsSdPErF5Mw//+BMZH0+BmjUpPXkSRdq3VxvzWIPkZPjnD71hCDwAdo5Qry80+ABcqlhEwpUrV1i4cCHLly8nPDw8VZyNjQ2dO3dmxIgRNGrUKM+NPyUnS5KNdrMzF6Z0Mf1PCPEJsBEwGAkpZZjZVCnyBEnh4UT+uJ2ILVtI+PdfbBwdcerciaLdu1PQzTLTIhVpiA2H0+v0XUlh/0KRctB6in5bz4Lm33Nbp9Px66+/Mn/+fP7888/H4osUKcJ7773H0KFDcXV1NbuenMbDuES2ngxm9ZEbdKigo6WZyzPFQPTT/g42CpOA5XbuUOQaZHIyj44eJXzTJqJ27YbERArWrk2ZadMo0q6tmqJqLe6cgWNL9N1ISbFQvoHeq+pLHcHW/KvP7927x5o1a+jTpw83b958LL5y5coMHz6cvn37UrhwYbPryWn8cy+KVUcC2XbqFo8SdNSpUBR7W53ZyzVlJXUls6tQ5HoSbt4kcvtPRP70E4m3bmHj5IRzz7co2q0bBdQUROuQGAcXf4LjSyD4uH79gkcP/RoGC+zDIKXk8OHD/2/vzaPruO78zs99+4aHh30HARIgwX0DxUUkxUUUtZiStdgtybak2Iodu+2OMxOn3afnJE5nOnEyM0lnMtPOdGxHre52q1u2W4vVokhRAilSXMSdBDeAIEjsO/D2te78UfUeFoIkSOIRIFGfc+rUrVu36v1QeO/3rbv9Ln/+53/O22+/PWqdZ1DnLzz11FN873vfY/v27TMuFEYsofDxuS7+8mAzh5r6sZgMPL20mFfWzmJJqYe6exDccCIzqW3A94D1qDWHz4D/LqUM3/RCnRlPwu/Ht3Mng++8Q+joMRAC59q15P3wh2Q8tk1fwnOKsIa74eOfwPE31dFIOVXw+E9h6Utg96T9871eL7/61a/42c9+xunT1weGzs3N5Vvf+hbf+c53qKycee+nLf1B3vriGn9/tJUeX4QSj50fP1HDV2vLyHbe2074idQd3wR8wH/Tjl8C/gr4SrqM0rl/kYkEwcOHGfyHd/Dt3o0Mh7FUVpL3L/4Fmc88rYe/mCqUBFz+FL74OWsufQQCmPekWluYvSnto5GklBw4cICf//znvP322wSDwevKLFy4kB//+Me88MIL2Gy2tNoz3YglFPac7+ZXR67xWUMPAthSk89LD5WzaV4+RsPUdMJPRCDmSSmXjjj+VAhxKl0G6dyfRK5cSTUhxTs7MbjdZH75GTxf/jK2pUtn3CiTacPgNTjxN3Dir8HbCs48rpU/z6znfwKZpWn/+O7ubt58801+/vOfjztE1eFw8LWvfY3vfve7DA0NsWnTprTbNJ1o6Q/yd1+08PdHW+j2RSh02/iDLdX83qoyij32qTZvQgJxQgixRkp5CEAIsRo4kF6zdO4H4j09eD/8kKH3f0f4zBkwGHBuWE/BH/4rXFu26E1IU0U8qsZEOv4mNO5R8+Zsge1/CvOe5Mr+z5mVRnFIJBLs2rWLX/ziF7z77rvXhdgGWLx4Ma+//jqvvPIKHo/arHUv2tSnA/GEwp4L3fzq8DX2abWFTfPyefmhcjbNy8NknD59LRMRiNXAK0KI5NCCcuC8EOIMIKWUeizlGUTC58O3azfeD35H4NBhUBSsC+aT/6Mf4d7xJcz5+VNt4syl5xKceBNO/q26Wpu7RF2UZ/nX70nQvCtXrvDGG2/wy1/+ktbW1uvOu1wuXn75ZV5//XVqa2tnXK2yocvHr4+38g/H21K1hR9otYWSaVBbGI+JCMTjabdCZ1qjRCJYT5yg9bf/gL+uDhmNYi4rI/effQf3U09hnXNvJk7pjEM0oI5EOv4mXDsIBpM603nFq2qtwZDe1dGGhoZ4++23efPNN/nss8/GLfPwww/zrW99i6985Su4Ztikx6FgjPdOt/PrY62cahnEaBBsnpfP760qY/M0qy2Mx0SGuV69F4boTC9kPE7wyBGGPvgA367deHw+gjk5eL76VTJ3fAnbkiUz7g1w2qAo6uzmU2+p4hALqCORtv2JOhLJld5aXCwWY9euXbz55pu8++671wXLA3Uk0quvvsq3vvUt5s+fn1Z7phsJRfJZQw+/PtbKrnNdROMKNYUZ/G9PzeeZZSXkZdw/Ta/6+os6KZKi4N35Eb6PPybR34/B6SRj2zaay8tY++1v60t2TiU9F1VROP33aoez1Q2Ln1dFoXxtWkciSSk5ceIEb775Jr/61a/o6em5rozRaGT79u1885vfZMeOHTMuLlJjt59fH2vlH0600uWN4HGYefmhcl5YWcrCYvd9+UKl/9pnODIeJ3D4ML6kKAwMIBwOMjY9Qsb2x3E9shGDzcalujpdHKaCQC+c/Q2c+ltoPwHCCFVb4bE/UYepmtPbdn316lXeeust3nzzTc6dOzdumeXLl/PKK6/w0ksvUVBQkFZ7phvdvjC/O9XBu6faU01Im+bm8W+fLmVzTT5WU3qb+NKN/oufgchYjMDhI/g+2olv98ckBgc1UdhExuPbcW1URUFnioiF1YV4Tr0FDbvUCKqFS2D7v4dFL0BGep1wR0cHb7/9Nm+99RYHDx4ct0xxcTFf//rX+cY3vsGiRYvSas90wxuO8dHZTt471c6Bxl4UCQuL3fzxk/N5Znkx+Rnp/+14o16iSjTtn6MLxAxBRqMEDh/B+9FO/Ls/JjE0hMHhwLVlCxnbH8O1YYMuClNJIqZOZDv7G7jwAUR94CpUw2ovfTHtazr39vby29/+lrfeeou6ujrkOJFCHQ4Hzz//PK+88gqbN2/GaLy/345vh3AsQd3FHt492caeC91E4wrl2Q6+v7mKp5cVU5Wf/vhQwViQupY6Pmz+kANtB/hq1ld5jMfS+plpEwghxC+BLwHdUspFWl42alTYCqAZ+KqUckCojXP/FXgSCAKvSSmPp8u2mULC7yewbx++PZ/g37sXxe/H4HTi2rIF9+Pbca5fr89VmEqUBDTvV0Xh/HtqJFVbJix8BhY9D5WPpHUU0tDQEO+88w5/93d/x+7du8edr2A0Gtm2bRsvvfQSzz333IwahZRQJIea+vjFmQg/qPsYXzhOrsvCyw+V88yyYpaVedLerxBJRNjfup8Pmz9kb8tewokw+Y58Xqp5idKB9E90TGcN4g3g/0EN1ZHkx8AeKeVPhRA/1o7/EHgCqNa21cDPtL3ObRLr6sb/6Sf4Pt5D4PBhiMUwZmeT8fh2MrZsxfnwOl0UphJFUQPjnf0NnHsH/F3qWgs1T6qiMGcLmNL3/xkcHGTXrl382Z/9GR9++CHR6PXNFEIIHnnkEV588UWef/55cnNz02bPdCOeUDh8pZ8PznSwq76TXn8UmxGeWlrKM8uKWTcnJ+1DU2NKjEPth9jZvJM91/YQiAXItmXzTNUzPFH5BMvzl2MQhukRrO9OkVLuE0JUjMl+Btikpf8SqEMViGeAN6Varz0khPAIIYqklPra17dASkn08mV8H+/B98knhLXgZ+ZZ5WR/4xtkPLoV+9KliBnUHDDtUBRoOwbn34X6d2CoBYxWmPuYKgrV28GSvjDn3d3dvPPOO/z2t79lz54949YUANauXcuLL77ICy+8QHFxcdrsmW7EEgqHmvr4xzMdfFTfRX8gisNiZEtNPk8uLsLUfYHHti699Y3ugrgS51jXMXY27+Tjqx8zGBkkw5LBY7Me4/HKx3mo8CFMhnvfIyDGa2uctJurAvG7EU1Mg1JKz4jzA1LKLCHE74CfSin3a/l7gD+UUh4d557fBr4NUFBQsPKtt966I9v8fv+0rC5PyK5EAnNTE9bTp7GeOo2puxuAWEUF4WVLiSxZSqKocFKHPd7Xz2sKCPiGKIlfI7f3IHk9h7BG+1CEkYGsZXTnb6A3dzUJU3pFYd++fXz22WecPXsWRVHGLVddXc3mzZvZvHkzhVMcSPFe/i/jiuR8X4IvuhIc74rjj4HNCMvyjawqNLEo14jVKNJqV0zGuBi6yKngKc6EzhBQAliEhSX2JaxwrqDGXoNZ3Hit77uxa/PmzceklLW3KjddOqnH82TjKpeU8i+AvwCora2Vdxrcq66ubloGBruRXfHeXvyf7ce/dy+BAwdQfD4wm3GuWUPG976La/MWzAXpmyB1vz2vKSEehSt74dy7RM++iyXmBZMNqh6F+TswzN1Ojj2LnDR9/MWLF3nnnXf4zW9+wxdffHHDcvPmzeO1117jueeeY+40Wosj3f/LUDTBZw097DrXxe5zXQyFYrisJh5dWMyTi4vYODcPm/n6mvZk2hWMBdnftp+Pr33MvtZ9BGIBXGYXmyo2sa18G+tK1mE3TWzo8r347t9rgehKNh0JIYqAbi2/FSgbUa4UaL/Htk0rpKIQPnsW/959+PftU4PhAaa8PHXU0SOP4Fy7Vl+veaqJBuHyHjj3nrqGc8QLFhcDWcspeOR1VRys6fkfxeNxDhw4wHvvvcf7779PQ0PDuOWEEKxfv57nn3+eZ599lqampukjqmmmzx9hz4VudtV3sb+xh3BMIcNmYtv8Ap5cXMSGublpn6vgjXrZ27KXj69+zIH2A0QSEbKsWTxe8Thby7eypmgNZuONawpTyb0WiPeAV4Gfavt3R+R/XwjxFmrn9NBM7H9IeL1Yjx6j/cOd+D/7jER/PxgM2JcuJe+H/xzXxo1Y58+/L2dkPlB421UxuLhTrTHEw+p6zfOfhgVPQ+UjnD9wiIKFmyb9owcHB9m5cyfvv/8+H374IQMDA+OWM5lMbNmyheeff55nnnlm1AS2pqamSbdrOnGlN8Duc53sPtfFsasDKBJKPHZeXFXOtgUFPFSZjTnNHc2tvlb2tu6lrqWOo11HiStx8h35PFf9HNtmbWN5/vIp6VO4XdI5zPVvUTukc4UQrcC/QRWGvxdCfAu4xvCiQ/+IOsS1EXWY6z9Jl13TCRmPq7WEAwcIfH6Q0MmTeBIJ/JmZODduxLVxI871D2PKSv9i8To3QUp1zeZLO+Hih9BxUs33zIKVr6nB8WatT9vazZcvX+b999/n/fffZ9++fTfsZHY4HDz22GM8++yz7Nixg6wZ8r1JKJJTrYN8rDUdNXT7AVhQ5OYHW6p5bGEBC4rSG+oioSQ403uGupY69rbupXGwEYDKzEq+Mf8bbJ21lcW5izGI6R2cbyzpHMX00g1ObR2nrAR+P122TCeiLS0EDhwgcOBzAocOqX0JQmBbuJCcf/o6DW436159VR91NNXEwnBln7quwqWPwNsGCChdBVv/tRrmIq8mLfGPgsEgdXV1fPTRR+zcuZNLly7dsGxpaSk7duxgx44dbN68ecasxDYQiLKvoYdPL3Sz91IPA8EYRoNgdWU2X1tdzqMLCijNSt8gAIBALMDn7Z9T11LH/rb99If7MQojKwpW8KPaH7GpbBPl7vSHWU8n07+Oc5+T8HoJHDqkCsLnnxNraQHAVFykTlZbtw7HmjWpWsK5ujpdHKaK/iZ1gZ3Gj1VxiAXVOQpzNsPmP4bqx8CVN+kfK6Xk/Pnz7Ny5k507d7Jv375xI6Qmqa2tZceOHTz99NMsnSGr9SmK5FyHl08vdPPpxW5OtgyiSMhxWthck8/meflsrM4j05G+tnwpJc3eZg60HWB/236OdB4hpsTIsGSwvmQ9m0o38XDJw2RaM9Nmw71GF4hJRolECJ04SfDIYQIHPid05gwoCgaHA8eaNWS/+irOh9dhqaiYET/saU00oM5kbvxY3fq1tvmsClj2Msx9AirWg3ny38oHBwfZs2dPShTGW2Anid1u59FHH2XHjh089dRTM2aOwlAoxueNvXxyoZu6Sz30+FTRXFqayQ+2VLO5Jp8lJZkY0rhecyAW4HTwNJ8d/IwD7Qdo87cBUOGu4KWal9hUtoll+cswG6ZnJ/PdogvEXaJEo4RPnSJw+AjBw4cJnTqFjEbBYMC2eBG5/+w7ONetUyermR/ML9F9g5TQfV4dddT4MVz9HBJRMDugYgOs/q4aKTVn8hdAisViHDlyhD179rB7924OHjxIIpG4YfmFCxfy+OOP8/jjj7N+/foZ0XQUVySHm/rY39jLZw29nG5Vawlum4mNc/PUWsLcvLSupyCl5NLAJfa37edA+wFOdJ0gLuM4BhysLlrNNxd9k3XF6yjNSH+Yi+mALhC3iYzFCJ05q9YQDh8mdOIkMhwGIbDOryHr5ZdxrH4IR20txoz0B/DSuQX+HmjepwbCa9wDPm30dP4CWP0dmLNVXUthkmsJiqJw4sQJ9uzZwyeffMK+ffsIBAI3LJ+Zmcm2bdvYvn0727dvp6ys7IZlHxSklDR2+/msoZf9jb0caAgSSRzCIGBpmYfvb65ifXUeK8o9aQ1v0R3s5nDHYQ51HOJg+0F6QupaF/Oy5vHKwldw9bh4bdtr03YoajrRBeIWKNEo4bP1BI8eJXjkCMHjx5HBIADWefPwfPUrOFevxrFyJUaP5xZ300k7YS9cPaD2ITTthe56Nd+aCXM2qfMS5myFzJJJ/VgpJZcuXUoJwu7du/F6vTe9pra2NlVLWL16NaYZsN5Gjy/CAa2GcKCxl05vGICKHAcPF5v46iNLWTsnh0x7+pzxUGSIo51HOdRxiMOdh7kydAWATGsma4rWsL5kPeuK15HvUCee1tXVzUhxAF0griPh8xE6cYLgseOEjh0jdOYMUuswtFZX4fnyl3GsXo3joVX68NNpgCERgaa6YUFoPwEyoc5gLl8Di/+NGhW1aOmkDkOVUtLU1MS+ffv49NNP+eSTT2hra7vpNZWVlWzZsoUtW7bw6KOPkp+f3qVBpwN9/ghHrvRzqKmPw1f6udDpA8DjMPPwnFzWV+eyviqXsmyHOjN40eSH+wjFQ5zoOsGhzkMc6TjC+f7zKFLBbrKzomAFz1U9x+qi1czLnnffDUNNNzNeIGJdXYSOHSN47DjBY8eIXLyotlUbjdgWLiTrpZewr1yBY+VKTNnZU22uTjyqikDzPriyj/XNB+GzmLrSWmktbPhfoXIjlD00qVFRFUWhvr4+Fd9o3759dHTcfC5nQUEBW7ZsYevWrWzZsoXKyspJs2e60uePcPhKP4eb+jjU1M/FLlUQ7GYjtRVZ7FhazIbqXBYWZ2JMU+dyJBHhTM8ZjnYd5XDHYU71nCKmxDAJE0vylvCdJd9hddFqluQumbE1g4kyIwUicOgw7v/5Bo3/+58S00aPCIcDx7KlZPz+7+OoXYl9yRIMjvSOo9aZANEAtByBawfVTuXWL9SZywCFi2kreZKyjV+HWWvBOnl9PrFYjOPHj6fEYP/+/TectZwkMzOTTZs2sXXrVlwuF6+99toDP1Ktd0QN4VBTH5e61ElqDouRlbOyeHpZMWtm57C4JBOLKT1v54FYgJPdJznWdYxjXcc403uGmBIDoCa7hpdrXmZ10WpWFqzEYdZ/07fDjBSIaPMVLOfOYVuzhuxvfB37ipXYaubpo4ymA8F+uHZI7Ue4dhDaT6pNRsKgLrtZ+02YtQ7K14Ezh8t1dZTN3XTXHzs0NMSRI0c4ePAgn332GZ9//jlBra/pRrjdbtavX8/GjRvZsmULy5cvT/Uj1NXVPXDiIKWkuS/I0eZ+jl0d4OjVARq7hwWhtiKbLy8vSQlCusJZDIYHOd59PCUIF/ovkJAJjMLIgpwFfG3+11hZsJLl+csfqDkJU8GMFAjP889zqqCARZs3T7UpMxspYfCaWiu4+rkqCN3n1HNGC5TUwvofqoJQ+hDY3JPysYqicPHiRQ4ePJjazp07N+4ymyPJz89n48aNbNiwgY0bN7J48eIHetnNSDzB2bYhVQyaBzh+bYBev7rAkNtmYuWsLJ5bUcLa2TksSpMgSClp9bdyqudUqpaQDGNhMVhYkreE1xe/zoqCFSzLW6bXECaZGSkQwmxOS4gEnVsQC6n9By1HVFFo/UJdUQ3A4lL7DRY9p9YOSlZO2tBTr9fL4cOHU2Jw6NAhBgcHb3ldRUVFSgw2bNjA3LlzH7hawUj6A1GtZtDPseYBTrcNEY2r60hU5Dh4ZG4+tRVZrJyVRVWeKy0T1ELxEPW99ZzqOZXa+sP9ADhMDpbnL+fJyidZWbCSRbmLsBgtk26DzjAzUiB07gFSwuBVaNGEoPUIdJ4BRQs0l1Wpji4qe0jtXC5YPCmjjCKRCKdOneLo0aMcPXqUL774gvr6+lvWDoxGI0uWLGHt2rWsW7eODRs2UF5+f8fRuRnhWIL6di+nWgY52TLIqdZBrvapTWpmo2BRSSavrp3FylnZrJyVlZbJaVJK2gPtnOw+yameU+zv2E/HrzqIS/U7Mss9i/Ul61mat5SleUup8lRhTOMa3TrXowuEzuQQ7FdrB8mt5QgEtOU+zA61RrDuB2pTUemqSYlpFIvFaGhooKGhISUIZ86cIRaL3fLa3Nxc1q5dm9pWrVqF0+m8a5umI4oiaer1c7JliJMtA+w/F6J110fEFVU0izJtLC318OKqclbOymJJaea4C+fcLYPhQer76jnbe5b6vnrO9J6hN9QLgN1kp9RYymuLXmNZ3jKW5C0hy6YPI08S9A7R13qN/rYW+lpb6Gu9hqV8NsMrOKcHXSB0bp/wkBr+uv0EtB1X94NXh89nz4E5W9SaQdlDkL/wrmsHiUSC8+fPp4Tg6NGjnDx58qZB7ZKMrB2sXbuWNWvWMGfOnAeyuUhKSetAiPr2IU63DnGqdZDTLUP4IupbuctqoswJ3944m6VlHpaVeShwT34YD3/Uz/n+8ykxONt7NhXHCNRYRmuK1qRqB9VZ1ezft59NKzZNui33C1JK/P199Le10tfWQn/btZQYhHzDky7NNjs5pWWYb7CM7GSiC4TOzYn4ofM0pS3vwm/+WhWDvsbh855yKF4Otf8EileoE9Lsdzej3O/3c+bMGU6ePMnJkyc5deoUp0+fJhQKTej66upqamtrqa2tZeXKlaxcuXJarlt9tyiK5EpfgLNtQ9S3e1P7oZBagzIZBDVFGTyzvJilpaoYzMlzsW/fXjZtqpk0O4KxIJcGLlHfV099bz1n+87SPNSM1FYNLnYWszB3IV+Z+xUW5S5iQc4CMiwzNwxNLBymv6ONgfZW+tvbGOhoo7+tlYGONmKRcKqc1ekkp6ScqlVryCktJ6ekjOzScjJychFCUFdXl3ZbdYHQGcbfrfYTdJ6BrrPqvvcSSIUqAHeJKgZLX1T3RcvBeecrLEspaW9vHyUEJ0+epLGx8ZZ9BkkKCwtZv359ShBWrFjxQC6UE0soNHb7RwjBEOfavQSiasA/i8nA/MIMnlxcxKISN4uKM5lXmDHpTUU9wR4u9F/g4sBFLvRf4EL/Ba55r6XEINeey6KcRTxR+QSLchaxMHch2baZN8FUKgq+/l5VANpb6U+KQXsbvr6e4YJC4M7NJ7u4hJL5C8guKiWruISc0nKcnqwpr+XqAjETScSh//KwGCQFITmiCMBdCoWLYMGXoXg5nzeHWLf92Tv+yGAwyPnz56mvr+f06dMpUejr65vwPYqLi1m1atWo2kF9ff0Dt75yfyDKhQ4vFzp9XOz0cb5TTSdHFDksRhYUuflKbRkLi90sKsmkKt81qcNME0qCq76rXOxXhSC57wsP/79KXCXUZNfw1OynmJ89n5rsGgocBVPu1O4VUkqCQ4MMdLYz2Nmhbe1q7aCjjfiI5k+L3U5WUSml8xeSXVxKVnEp2cUleIqKMVvSF532btEF4kEn2A89F6CrHjpPQ+dZda5BcjaywQz5NWoAu8LFqigULALH6Le+aEfdhD4uEolw8eJFzp49S319fWrf1NQ04VqB0WikpqaGpUuXsmzZMpYtW8bSpUsfuNhF4ViCxm4/Fzt9XNBE4EKnL7XuAagL4swrzOC1dRUpMajIcU5qmIr+cD+NA400DDbQONjIpYFLNAw0EIqrTXomg4kqTxXrS9YzP2c+87LmMTd7Lm7L5MxLmc5IRSHq99FSf5qBzg4Gu1QRSArCyCYhYTDgzssnu6iEsgWLyS4uSYnBdKgN3Am6QDwohIeg+wL0nFf33edUYRhZK7BnqyKw6nVVBAoXQ+5cMN3+WPJYLEZjY+N1QtDQ0HDTdQ7GkpGRkRKC5H7hwoXY7fbbtmm6Ek8otAyENDEYFoIrvQES2kgii8nA3AIXG6vzmF+UwbzCDGoK3ZM6vNQX9XF58DIHfAc4ePgglwcv0zDYkJpnAGpE0ypPFc9XP09Ndg012TXMzpz9QMcsSsRjeHu6GerpZqirQxUCrTYw1NVJPBbljFbWYDSRWVBIVmERZQsW4ykswlNYjKewCHduPsYHLCLvg/XXzAQifui5OCwA3efVvXdEJFGzQ10vuepRdZ8/X13/wF182xMEe3t7uXDhAh988AEffPABFy9e5OLFi1y+fPm2hMBgMFBVVcXChQtZtGhRShAqKysxGB6MCJqReIIrvQEau/3sbojydvtxLnf7aeoJEE0Mjzgpy7Yzr8DNE4sKU0JQkeOYtDUP/FE/zd5mmoaaRtUMOgOdqTJ2r51qTzWbyjZR5alijmcO1Z5qcu259+Wb7s1QlAT+/j6GursY6u7C29OVSg/1dOHv71Pn7WiYzBYyCwrxFBZTsWwlPV4fD23YhKewiIzcXAwzaC6GLhDTEUUBbyv0NqgjhnoboK8BehvV/CQmm1oDqFivCcECtbkosxxuw+lGo1Gampq4cOFCSgCS6f7+/lvfYAyVlZUsWrQoJQYLFy6kpqbmgVkVzR+Jc7nbT2O3n8YePw1dfi73+LnaF0CrECCAsuwhqvJdPDI3jzn5LqryXVTnu8iw3f3buCIVugJdXBm6whXvFXWvbckFbwDMBjOzpIIG0gAAHKBJREFUM2ezsmAlVZ4qqj3V9F7s5dmtzz4woa2TfQFJh+/t7mKou1MTg268vT0oifjwBULgys4hM6+A8oVLcOcVkJmvbp6CIlxZ2YgRv5+6ujpmLVk2BX/Z1KMLxFQS9g47/r4G6G2g9upJ2N8F8RFDOq1uyKmCiochtxry5qu1gqwKmODbTDQa5cqVK1y+fDm1NTY2cunSJZqamm6rNpCkrKwsJQBJMZg/f/4DMeEsHEvQ3BeguTfAld4gV3r9NPcGaeoN0Osf7iMwGQSVuU5qCjPYsaQoJQRt54/z2Na7j/UVjAW55rtG81DzsAh4r3DVezXVRwCQYc6g0lPJuuJ1VGZWUpFZQWVmJeUZ5ZgMo3/mdZfr7itxiEXC+Pp68fX24u3rxtfbi6+vB29vj5bfQzw6ej6MI9NDZl4BBbOrmLvmYTLzC3FrIpCRk4dJD8w5IXSBSDcRH/RfgYErw/u+y2qtwD9c5UcYwDOLiDUX15IvqYKQWw051eDKn1DTkM/nGyUASRG4fPkyLS0tKHcwscbhcDBv3jyysrLYuHEj8+bNo6amhurq6vteCKJxhZaBoCYC6tbcF+BKT4D2ofCosrkuK5W5DrbU5FGR62R2rpOq/Axm5TjGHT3Uc2nizTTeqJcWbwvXfNe45r1Gi6+FFp96nJxpDCAQFLuKqcysZFXhKlUI3KoQ5Nhy7sumIUVJEBgYwNfXQ3/jBb7w9uHr61Wdf28Pvr6eUZPEABACpycLd04eeeUVzF5eS2Z+wXBNIK8A8wNSW51qdIG4W6RU5w+MFICR+2Dv6PKOHHWmcdXW0SKQXQkmK2fq6m44bDMSidDS0sLVq1dpbm5ObUkx6O7uvuM/o7y8nHnz5qUEILkvKSlJTcq534aTSikZDMa41h+kZSCo7vtDtGjHrQOhVCcxQKbdTGWuk9Wzc6jIcVKZ56Qyx0lFruOumoWklAxEBlSnrwnANd+1lCgMRkYHDsx35FOWUcaGkg2Uu8spzSil0l3JLPcsbKb7x/HFohH8/X3qNtA/nE4d9+Lv70MZUXu9AljsDty5eWTk5lFYVY07N5+M3DwycnJx5+bhys7BaNJrAPcCXSAmQtgLQy0w2KKGpx66polAs7qPjVyMXkBmqerwa55S91mVw/ubhKwOh8O0tLSwa9cumpubRwnB1atXaW9vn/BQ0bEIISgtLaWqqoo5c+YwZ84cqqqqqKqquq9rA9GEpLHbR0t/SBMATQgGVCHwR+Kjyuc4LZRmO1hcksnTS4tHCUGW884ig0op8Ua9tPnbaPe30+Zv44v+L/j1nl+n8oLx4bUlBIIiZxFl7jK2zdpGeUY5Ze4yyjNUMbCbpvcILqkoBL1D+Pv78PX3ERjoS6X9/X0ENDEIB/zXXWuyWsnIzsGVlUPJvAW48/LJyMklIzePi03NbHniSayO+/O7+CCiC4SUEBrQBODaCBFIHl+D8JjQ0Ear2v6fXQkVG0aLgKd83KUu4/E4HR0dtLbW09LSQmtra2praWmhubmZzs7O6667HcxmM5WVldeJwJw5c6ioqLjvOomllAwEY7QPhmgbDNGe3IbCat5AiG5fBHbvS11jMxsoy3JQnu1gdWU2ZdkOyrLslOc4KMty4LTe/lc+KQAdgY6Uw2/3t9Pqb02l/bHRztAmbJSbVIe/umg1Ja6SlBCUukqnZZjqeDRKcGiQwODAqC04pKUHBvAN9BEcHBj11g8ghAGnx4MrOwdPYREl8xepQqCJgSs7G1d2Dha744ZNYdeGAro4TDNmpkCc/x2LT/8XOPdHqgBEx7zpmJ2qo/eUqZFHk2nPLMgsA2feqFFC0WiUjo4OWi600Np67Drn39raSmdn5x31AYxECEFJSQkVFRXMmjUrtZ89ezZVVVWUlpbeVwvYhGMJOpLOfqQADIZpH1LT4djoZ2Y1GSjx2Cn22Hlkbh4JbxcbVyxUhSDbTp7Leltt8VJKfDEfXYEuOgOddAXV/ch0V7BrVIcwqGsTFLuKKXWVUltQS7GrmBJXSWp//PPjbJ4GC1IpSoKQ1zvK4XceP8anzZfU46EBAoODBIcGiAQC497DnuHG6cnC4cliVslSXNk5OLNUh5+RpYqAI9OD4T767ulMjJkpEOEhrJE+yFugLnDvKVcdf1IE7FkgBKFQiM7OTjo6Oui42EFn52E6Ot5Rjzs6Uue6urruuOlnJAaDgby8PGpqalICMFIMSktLsVim35vnWBRF0heI0uUNa1tkRDpMpzdCtzdMXyB63bX5GVaKPXbmF7rZMi+fYk0MVFGwke20jBKAuro6Ni0vGdcOKSVDkSG6Q930BHtuKAAjm38ADMJArj2XQkch1VnVbCjdQIGjgCJnESUZJZQ4S8i0Zt5UiNLVYSwVhZDfR8jrJegdJOQdIuj1ant1C2mbmvYi5fUvJj12u+r0M7PIK6/A6VmuiYAHpycLZ6aadrg9D9zkL52JM+P+81JKBiue5N28GOUF5aqzr++go2PvsBho29DQ0KR+dkFBAaWlpZSVlVFaWjpqq6iooKSkhP3790/bzuBkp2+vP0K3L0LnUJguX5iuIVUEOr1hur1hun2R1FoDSYSAHKeVAreV4kwby8o8FGfaRglAQaYVq+nWb6GKVOgP99Mb6qU+VE9/Qz89wR56Qj30hnpHpZOL16fsQKjO31nI7MzZrCteR6GzkAJHQWqf68jFbLg3naCJeJyw30fI5005+7GOf6TzD/t84zp8AJvThd2did2diaewiOK583F4slSHn3T8niyOnTrDlm3b7snfp3N/M+ME4oMPPmDHjh2Tek8hBIWFheM6/mRecXHxtHz7VxTJQDBKrz9Krz9Crz9Cjy+SOlbTEdr6gvh3fXid4wfIsJkodNsocNuYMyeXAreVAu04mc7LsN40mJwiFbwRL23+fvrCffSH+1NbX6gv5fR7Qj30h/pTq44BoA3eyrRmkmfPI9eeS627ljxHnnrsyCXfnk+hs5A8e15awkYoSoJIIKA5ex9DV5s4t09Rj7W8sH94Sx5HQ8Eb3tOW4caR4cbuziS7uJSSmgU4Mj3YMzJxuNV8hzsTR6YHmytjwm/6BvOFyfqzdR5wZpxAFBYWTrisyWSioKCAoqKiUVthYeF1afM0mXijKBJvOEZ/IMpAMEp/IMZAMMpAIEp/IEqPX3P+muPvC0RHDfVMYjYKcl1Wcl2qg88WQRZXz1LzMqzkZ1gpdNvId1txWMb/GoXiIc3Bt9LQrjn7Ec6/LzScHggPkJDXT9YTCLJsWeTac8lz5FGdVZ0SgXxHPi0XWnjs4cfItediNd5d3CIpJbFImEggQCTgJxwMqOmgduz3E/J7Cfv9hH3e0ccB/6hwDQCpVTOEwObKwO5yYXNl4Mj0kFNShi3Djc3lwu5yY8vIwOH2pBy/PcOtt+nrAKAkFKLhBJFgnGgoTiQUJxqME/XffbP2rZhxAlFUVITT6cTj8VBZWXlT55+TkzOlcYKklPgi8ZRzVx19THP80eF9IEa/JgIDwSjj+HtADQiX57KS67JQ7LGxpDRTEwELuRnWlCDkuay47aZUO7oiFXZ+spNFqxwMRgYZjLTQExmioXuQgWsDDEWGGIwMpvbJdDgRHtcOh8lBti2bbHs2xa5iFucuVo9t2eTYc1LpbFs2HqvnpusQ112po8RVknpe8Wgk5dTDgQCRoJ+I/3pnHwkERuT5U+fGjs657hnaHdgzMrC51C2zoBCby4XN5VYFIMON3ZXBuUsNrHvkEbWcwzkqdIPOzEB94UgQCydUpx6OEwuNSIcTRMOqw4+F4kS142gyHYoTCcaJRcb/ThatTP/EyBknEMXFxfj9/ns28UtKSSCaYCgUYygYU/ehGN7QcHrk1tIV4k+O1qWOx2vSATXEQ5bTQrbDQpbTzNwCF1kOC9lOy/B+xPkshxlhiOGP+fFFffiiPrzRwVS6I+rjQt8Qg+3XO/qh6BCKVKD1ejsMwoDb4sZj9ZBpzaTQWci87Hl4rB6ybFnk2HJGOf8sW9Z14/yloqhv7qEg0WCIqD9IpMdPf6iHjlCQaChINBgkGg4RDQWJBLW8cIi+7m4u//ZviISCRAKB0TF3xntuFitWpxOrw4nV6cThdpNVVIzV4cTmcqXyrQ4XVqcTW/LYqZ6baDPO1SE/WYXFEyqrM32QiiQWTaQceyySINAtaT7Tq+Zp+aMceTId0tJhNR0Lx8dWKsfFZDVisRmx2k2YbSYsNiOOTCsWuwmr3YTVYUqlLSOOT58/mvbnMeME4nZHl4RjCfyROP5wHF84ji8SwxdOHsfwR5L56t47jgDcyMkDGAS47WYytc1ugopiNx7tOOnwMx1GHNYEVksMsyWKQmiMs/fijXpTzt7X68PXPnzOF/Vd12E7FpvRRqY1E4/Vg8fqYW7W3JTj723tpXZhLR6rhwyTCxd2nNKKRTGRiESJhcNEIyFi4bCa9oWJRcLEwn6i4W7agiGuhIKqCISCREMhba+mJ4LRbMZid2C1OzDb7VjtDixOF4WlZVjs9pQTtzlHOnfXsCA4nJimYT+Qzp2RSCjEowrxpOOOJIhFVKcdG5kXHn1+ZF40Mrps/AZv682fnL4uz2QxYLGbsGhO3WI34ch0qGmbaficfeSxWs6cFASrEcMdRvE1XtZrEGnhVMsgHzRFORK+kHL+3nAcf9L5a87eH46PCtN8I6wmAxk2Ey6riUy7GbfdTGmWHbfNiMMWx26Nq47dHMVkiiGMEYQhgiLCJGSIYDyIP+YnGAvS3NFMJNPB5ZifQChA0Kudi9+4MzOJ2WDGbXGTYXaRaXLjES7KLQW4rA6cwoYDG3as2KUZi2LCkjBiSghMcYGIKSixGDG/6uBV5x4mGu4nFm4D7xAt71yiKRImEb/5W/pIhDBgcdhHOXab04U7Nx+L3YHFnjxnHz52ONS0zY41mbbbxw2vcD+GAHnQkVKiJCTxaIJ4VCEWTYxJK8S1t/RkOh5N0HlZYW/HxVTZeDRBbMT54bQqCspNXrzGYjAIzDajullVx2y2GslwmrHYjKljdTNp5dTt3IWz1K5egWVEvsV25479fmJGCsSRK/28fSmGsbEJl81Ahl3BYVVw2BK43QnyLXHMljhmUwKzKYbRGMNgjCIMMTBEkURQiBAnSlwJE0mECcVVRx+Kh7gWDxGKhwj5QuBTP1MoYFTE8JZQ96aEwCJNuFCduCUscVnj5GLBquRilSbMihGzYsCUMGBMgCEBhrhEJCTEEshYAiUWJxGNEotGiEejIL3AcJAzBfBr23gIgwGLzY7ZZsNstWG22bDY7Djcbsz5BYiBAcoqKrGkztsx26zq3mpT87VzFpsdk9WKxWbHaDbfl0HkHgSSjjoRU4jHFOKxBImYQmhA0tk0pObHFRJRhXhcPRePKiTiavnkdYmYQiKWGM6LK+OUHe3U5W047yTCCIGWbkwWA2arEZPFmHpLd7gtmCxGzBYDJu2c2WLQyoxw7mOcvcWmioHRdOfO/NqgoLAy846vv5+ZkQJhttcxt+wXoCgk4jHVWYcF8aAgqAgiisCYGHbohhGO3awYsUgTZmnELk2q45aGlNM3KE4MCSdCkYi4hISibrf1gwlqGyAEZosVk8WCyWrFZLEOHzvVvXps1c5bMFutmMzafmT+iPukxEBz8EaT6aaOXH9TvzkpZxxXUOKSWEji7Q2px1p+IqaQSJVRSMS1/PGONceeSAznK9q5RFyq6YRCIjZ8TTyaGO3AYwrc4GvX9NGxW/9RAkxmA0azAZPJgNFiVI9NBkwWdbM5TRjNRoxmoTlto3buxumR+5QQmAzs3beXTZs2TO4/RueumJEC4bwwwPZ9WRO/QAiMJhMmswWTxYLRbMZotmAymzGazZjMyTw1ncw3JsubzGpeMm0xp+5lsloxma2aM7dw7ORJHt6wAbPFlvqsB/0NXCoSRUpkQnWySkKSSChaWknljc5Xz/naJU0ne8Ypqzrj8e4xtpx63/Hzb3RupLNPisJYLr178M4figCjyYDRKDBqTtlgFOrepB4bTQKTxYDVYUo5cqPZgMls1PaG4b3m1I1mAxcunmfpsiVjrhm+Npk2GMUD/93TuTnTSiCEEI8D/xUwAj+XUv40HZ+zdevvsTeawdLly2/t2C1mDMbx366llEhFIiVj9hKpaG+VigSpzk8YmX+j65S4h8CgGanEkTJ+0/uNuu+oe2tOV1E/R0mMSI/cj81PSBQJMqFo++HzXV0KH5w7rd5/zHVKYuy9tc9PKEhFtTF5To4pn/w77oZr+87cupCGwSi0zTAirR4bxztnMmC2jsk3Cc1BGzAaDRjN2vVansEkuNzUyPwF80bkqQ5dLT/W4Y+4n1ZOGNLnnNsDF5i1KCct99Z5sJg2AiGEMAL/L7ANdUDlF0KI96SU5yb7szqvGOi9WMXnV5NOODgh5z32fLq4snsC1f9JQBgEBoNAGAUGgbo3iFR+8g3SYBREQ+AjfF2+wWjAZFb7MNR8zQmPuI8wjrjnmPxRxwZSTth4nRMfThtHpE+cOsGqVbW3LGcwptfpjmXAcJn56/Rhrjr3N9NGIICHgEYpZROAEOIt4Blg0gXC7jJjy4T8goyU0xAGtL3q5G6Yf8Oyw9cYDAKEOnLipvcaJ//s2TMsWbpkzP2A6+474rpR9xi+13iO2iA0hy1ub8iv2gfx0GT/K+6ai62CvPKMqTZDR+eBRExGFNLJQAjxAvC4lPJ17fgbwGop5ffHlPs28G2AgoKClW+99dYdfZ7f78flct2d0WlAt+v20O26faarbbpdt8fd2LV58+ZjUsraWxZMtV9P8QZ8BbXfIXn8DeC/3eyalStXyjvl008/veNr04lu1+2h23X7TFfbdLtuj7uxCzgqJ+CXp9NMj1agbMRxKdA+Rbbo6OjozHimk0B8AVQLISqFEBbgReC9KbZJR0dHZ8YybTqppZRxIcT3gY9Qh7n+UkpZP8Vm6ejo6MxYpo1AAEgp/xH4x6m2Q0dHR0dnejUx6ejo6OhMI3SB0NHR0dEZF10gdHR0dHTGZdpMlLsThBA9wNU7vDwX6J1EcyYL3a7bQ7fr9pmutul23R53Y9csKWXerQrd1wJxNwghjsqJzCS8x+h23R66XbfPdLVNt+v2uBd26U1MOjo6OjrjoguEjo6Ojs64zGSB+IupNuAG6HbdHrpdt890tU236/ZIu10ztg9CR0dHR+fmzOQahI6Ojo7OTdAFQkdHR0dnXGakQAghHhdCXBRCNAohfjzV9iQRQjQLIc4IIU4KIY5OoR2/FEJ0CyHOjsjLFkLsFkI0aPusaWLXT4QQbdozOymEeHIK7CoTQnwqhDgvhKgXQvxzLX9Kn9lN7JrSZyaEsAkhjgghTml2/Vstv1IIcVh7Xn+nRXWeDna9IYS4MuJ5LbuXdo2wzyiEOCGE+J12nP7nNZFFIx6kDTVS7GVgNmABTgELptouzbZmIHca2LERWAGcHZH3n4Afa+kfA/9xmtj1E+BfTvHzKgJWaOkM4BKwYKqf2U3smtJnBgjApaXNwGFgDfD3wIta/n8HvjtN7HoDeGEqv2OaTf8L8Cvgd9px2p/XTKxBpNa+llJGgeTa1zoaUsp9QP+Y7GeAv9TSfwl8+Z4axQ3tmnKklB1SyuNa2gecB0qY4md2E7umFKni1w7N2iaBLcCvtfypeF43smvKEUKUAk8BP9eOBffgec1EgSgBWkYctzINfjQaEtglhDimrb09nSiQUnaA6niA/Cm2ZyTfF0Kc1pqg7nnT10iEEBXActS3z2nzzMbYBVP8zLTmkpNAN7AbtVY/KKWMa0Wm5Hc51i4pZfJ5/an2vP6LEMJ6r+0C/gz4V4CiHedwD57XTBQIMU7etHhLAB6WUq4AngB+XwixcaoNug/4GTAHWAZ0AP/XVBkihHABvwF+KKX0TpUdYxnHril/ZlLKhJRyGerSwg8B88crdm+tut4uIcQi4I+AGmAVkA384b20SQjxJaBbSnlsZPY4RSf9ec1EgZi2a19LKdu1fTfwD6g/nOlClxCiCEDbd0+xPQBIKbu0H7UC/A+m6JkJIcyoTvhvpJS/1bKn/JmNZ9d0eWaaLYNAHWpbv0cIkVzEbEp/lyPselxrqpNSygjwP7n3z+th4GkhRDNqk/gW1BpF2p/XTBSIabn2tRDCKYTISKaBx4CzN7/qnvIe8KqWfhV4dwptSZF0wBrPMgXPTGsP/gVwXkr5n0ecmtJndiO7pvqZCSHyhBAeLW0HHkXtH/kUeEErNhXPazy7LowQeYHazn9Pn5eU8o+klKVSygpUf/WJlPJr3IvnNdU981OxAU+ijui4DPzxVNuj2TQbdUTVKaB+Ku0C/ha16SGGWuP6Fmqb5x6gQdtnTxO7/go4A5xGdchFU2DXetTq/WngpLY9OdXP7CZ2TekzA5YAJ7TPPwv8ay1/NnAEaATeBqzTxK5PtOd1FvhrtJFOU7EBmxgexZT256WH2tDR0dHRGZeZ2MSko6OjozMBdIHQ0dHR0RkXXSB0dHR0dMZFFwgdHR0dnXHRBUJHR0dHZ1x0gdC5rxBCeIQQ35tqOyaCEOKHQghHmj9jpxCiZExesRDiEyHEu9os6mT+K0KIs1qk0nNCiH+ZTtt07n90gdC53/AA00IghMrNfkM/BG5LIEbMjJ1IWTvq3Iq2Maf+APgBamC3r2tln9DseUxKuRA1Ku7Q7dimM/PQBULnfuOnwBwtLv//ASCE+JEQ4gstmFoyhn+FEOKCEOLn2lvz3wghHhVCHNDi5z+klfuJEOKvtDfuBiHEP01+0E3ue14I8efAcaBMCPEzIcTRMWsI/AFQDHwqhPhUy/OPuPcLQog3tPQbQoj/rJX7j9qs+l9qn31CCHGjaMObUMNBjMWIGtRNYThmzx+hhvhOhnMJSyn/x209eZ2Zx1TNCNQ3fbuTDahg9HoQj6Eu3i5QX3h+h7puRAUQBxZr+ceAX2rlngHe0a7/CersdTuQixrpt/gW91WANSNsyNb2RlSHvUQ7bmbE+h6Af0T6BeANLf2Gdn+jdvzvga9raQ/qrH/nOM/i/wa2jJM/C9gHvA9kaHn9QOZU///07f7aJlyd1dGZpjymbSe0YxdQDVwDrkgpzwAIIeqBPVJKKYQ4g+rok7wrpQwBIe0t/iHUMBU3uu9VKeWhEdd/VQvPbkJdpGcBariG2+FtKWVixN/09Ig+AhtQjhqvaCQPA9f1I0gpr6KKmY7OXaELhM79jgD+g5Ty/xuVqa5/EBmRpYw4Vhj93R8bb0be4r6BEceVqE56lZRyQGs2st3A1pGfM7ZMYERaAM9LKS/e4D4IIWYDLVJd9Goi1AMrUeMK6ehMCL0PQud+w4e6fGaSj4BvJkfrCCFKhBC3uzDPM0JdjzgHtV3/i9u4rxvVuQ8JIQpQ1/K4ka1dQoj5Wsf2szex5yPgB1r0UIQQy8cp8wSwc0J/ncp/AP6TEKJQu6dV6yfR0bkheg1C575CStmndTSfBT6UUv5ICDEfOKj5Uz/qyJ3Eze4zhiPAB6jNOP9Oqh257RO5r5TylBDiBOobehNwYMTpvwA+FEJ0SCk3o65L/TvUfo6zqM1W4/HvUOP9n9ZEohn40pgyj6OOVJoQUsp/1ATsY+2eErVPRkfnhujRXHVmNEKIn6B2Hv+fU23LRBHqkpcHpJS1U22LzoONXoPQ0bnPkOrKZro46KQdvQaho6OjozMueie1jo6Ojs646AKho6OjozMuukDo6Ojo6IyLLhA6Ojo6OuOiC4SOjo6Ozrj8/8IsOaLrKpXSAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "for s in data.keys():\n",
    "    plt.plot(T, Pvap(T,s))\n",
    "plt.plot(T, Pvap_denatured(T), 'k', lw=3)\n",
    "plt.legend(list(data.keys()) + ['denatured alcohol'])\n",
    "plt.title('Vapor Pressure of selected compounds')\n",
    "plt.xlabel('temperature / °C')\n",
    "plt.ylabel('pressure / mmHg')\n",
    "plt.grid(True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "Y_0VfdRGgNft",
    "pycharm": {}
   },
   "source": [
    "## Optimization model\n",
    "\n",
    "The first optimization model is to create a mixture that maximizes the vapor pressure at -10°C while having a vapor pressure less than or equal to denatured alcohol at 30°C and above. \n",
    "\n",
    "The decision variables in the optimization model correspond to $x_s$, the mole fraction of each species $s \\in S$ from the set of available species $S$. By definition, the mole fractions must satisfy\n",
    "\n",
    "$$\n",
    "\\begin{align}\n",
    "x_s & \\geq 0 & \\forall s\\in S \\\\\n",
    "\\sum_{s\\in S} x_s & = 1\n",
    "\\end{align}\n",
    "$$\n",
    "\n",
    "The objective is to maximize the vapor pressure at low temperatures, say -10°C, while maintaing a vapor pressure less than or equal to denatured alcohol at 30°C. Using Raoult's law for ideal mixtures,\n",
    "\n",
    "$$\n",
    "\\begin{align}\n",
    "\\max_{x_s} \\sum_{s\\in S} x_s P^{vap}_s(-10°C) \\\\\n",
    "\\end{align}\n",
    "$$\n",
    "subject to\n",
    "$$\n",
    "\\begin{align}\n",
    "\\sum_{s\\in S} x_s P^{vap}_s(30°C) & \\leq P^{vap}_{denatured\\ alcohol}(30°C) \\\\\n",
    "\\end{align}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "Y_0VfdRGgNft",
    "pycharm": {}
   },
   "source": [
    "## Pyomo implementation and solution\n",
    "\n",
    "This optimization model is implemented in Pyomo in the following cell."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 312
    },
    "colab_type": "code",
    "executionInfo": {
     "elapsed": 7130,
     "status": "ok",
     "timestamp": 1555699127382,
     "user": {
      "displayName": "Jeffrey Kantor",
      "photoUrl": "https://lh5.googleusercontent.com/-8zK5aAW5RMQ/AAAAAAAAAAI/AAAAAAAAKB0/kssUQyz8DTQ/s64/photo.jpg",
      "userId": "09038942003589296665"
     },
     "user_tz": 240
    },
    "id": "Dn7WgE9_gNfu",
    "outputId": "2ff5eb43-4f59-4e4d-bddf-98d0a0947de7",
    "pycharm": {}
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Vapor Pressure at -10°C = 17.48178543436185 mmHg\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEWCAYAAAB8LwAVAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzs3XmcTXUfwPHP14wha6Kkzdpmz6BsIUuUXdlJj5KiTVFPWhQtT3lahCckS8RYIm1SMokQYytE1pDsjLHPzPf545y5c2fcmbkzY+bO8n2/Xuflnt/Zvr97x/3ec37n/H6iqhhjjDGJ5Ql0AMYYY7ImSxDGGGN8sgRhjDHGJ0sQxhhjfLIEYYwxxidLEMYYY3yyBGFMCkSknoj8KSJRItLuEu2zjIioiARfiv2lI46hIjI1kDFkFSLSW0SWBjqOrMQSRBYgIt+JyGs+ytuKyD9Z4Eukt4jEuF+QkSKyTkRaBTKmTPYaMEpVC6nqvEAH401EdolI00DHYXImSxBZwySgp4hIovKewDRVjc6sQJJJRstVtRBwOTABmCkiV6Ri+0yTATGUBjZe4n0ak+VZgsga5gFXAA3iCkSkGNAKmOLO3ysia91f8HtEZKjXunGXK/qKyN8isl9EnvFank9E3neX/e2+zucuayQie0XkORH5B5iYXKCqGgt8AlwGlEtqexFp5Z5pHBeRX0Skqlc8z4nIPhE5KSJbRKSJW15bRFa7dTwgIu96x+gdh/cvZ/cyyWwRmSoikUBvEckjIs+LyHYROSIiPhOa1/4eFpFtInJUROaLyDVu+XagHPClewaVz8e2SdXH7xhEpKiITHA/u30iMlxEghLFt9k9xiYRqSEinwI3eMU22F33Dvc9Py4i60Wkkdd+yorIT+5+vgdKJPWeuOu3dT/HSLceLdzya9z36aj7vj3stc1QEZnlfh4nReQ3EblJRP4tIgfdv9/mXuuHi8ibIvKriJwQkS+83ycRaSMiG936hIvIrV7LVEQqeM1PEpHh7uu4v81n3OPuF5EHvdYt7tYhUkR+Bcp7LRMRec/d7oSIbBCRysm9VzmSqtqUBSZgPPCx1/wjwDqv+UZAFZykXhU4ALRzl5UBFJgOFHTXOwQ0dZe/BqwArgKuBH4BhnntNxr4D5APuMxHbL2Bpe7rYOBJ4CRQ1Nf2QA3gIHA7EAQ8AOxyl98M7AGu8Yq9vPt6OdDTfV0IuMMrxr2JYtrlVb+hwAWgnfv+XAY85db5Ove4Y4HpSbz3dwGH3bjzAR8CS3wdy8e2ydUnyRi8PrNgd36eu7yg+zn9CjziLrsf2AfUAgSoAJT2FRtwLXAEuMd9L5q581d6vcfvuvHc6X6OU5OoW23ghLuPPO6+b3GX/QSMAfID1XH+3pp4fR5ngbtx/l6mADuBIUBe4GFgp9dxwt36VXbrPycuJuAm4JQbQ15gMLANCHGXK1DBa1+TgOGJ/rZfc7e9BzgNFHOXzwBmuses7MYQ93d+NxCBc8YswK1AqUB/T2T691KgA7DJ/SCgvvuf8TJ3fhnwdDLrvw+8576O+7K5xWv528AE9/V24B6vZXcDu9zXjYDzQP5kjtXb/Y92HOeLdAXxX84XbQ/8DzcBeZVtARrifLkdBJoCeROtswR4FSiRqLwRKSeIJYmWb8b9wnLnS+EkkWAf9ZsAvO01X8hdt0ziY/nYNrn6JBmD12cWDJQEzuGVnIGuwGL39XfAk0kcP0FswHPAp4nW+Q4nSd/gfo4FvZZ9RtIJYmzc31ii8uuBGKCwV9mbwCSvz+N7r2WtgSggyJ0v7Nb9cnc+HHjLa/2K7t9UEPASMNNrWR6cL/JG7nxKCeKM92fuflZ3uPu+QML/M28QnyDuAra66+bJjO+ArDjZJaYsQlWX4vwKaysi5XB+LX4Wt1xEbheRxSJySEROAP24+PLAHq/Xu4Fr3NfXuPO+lgEcUtWzKYS4QlUvV9USqnqHqv6QzPalgWfcSwLHReQ4zpfKNaq6DeeX9VDgoIjMiLucA/TB+cX4h4isktQ1hO9JNF8amOt1/M04X2olfWyb4P1R1SicX93XpnTQFOrjbwylcX7h7vdadyzOmQQ47932lGLx2tf9id77+jjJ6RrgmKqe8lp/t6+dpHDca4Cjqnoy0X68368DXq/PAIdVNcZrHpxEHCfx325enL/vxJ9NrLtuip+N64gmbMM77R73SpzknPi4ccf5ERgFjAYOiMg4ESni5zFzDEsQWcsUoBdO4/RCVfX+T/YZMB+4XlWLAh/hnPp6u97r9Q3A3+7rv3G+OHwtA+dXWHok3n4P8LqbUOKmAqo6HUBVP1PV+m5MinN5ClX9U1W74nwx/geYLSIFcS4xFIjbuXtt/ko/YmiZKIb8qrrPR/wJ3h/3mMVxfqmmXPkk6pOKGPbgnEGU8FqviKpW8lpeHt981fvTRMcsqKpvAfuBYm794tyQTNWSOu7fwBUiUjjRfvx6v5KQ+G/3As7ZauLPRtx14451Gq+/DeBqP493COdsKvFxPVR1pKqGApVwfrgM8nPfOYYliKxlCs6lioeByYmWFcb51XZWRGoD3Xxs/5KIFBCRSsCDQJhbPh14UUSuFJESwMtARt77Ph7o5571iIgUFKeRvbCI3Cwid4nT2HsW59dkDICI9BCRK91ficfdfcXgnOrnd/eRF3gR5xp6cj4CXheR0u6+rxSRtkms+xnwoIhUd+N6A1ipqrtSqmhy9fE3BlXdDywE/isiRcRp3C4vIg3dVT4GnhWRUPf9rBC3T5xf6uW8djcVaC0id4tIkIjkdxtrr1PV3cBq4FURCRGR+jiXf5IywX1fmrgxXSsit6jqHpx2rDfd/VfFOfubltL7lYweIlJRRArgtBnMds84ZgL3ujHkBZ7BSaa/uNutA7q5dW2BcxkzRe6+PweGuv9nKuJchgNARGq5f795cX6gnCX+c801LEFkIe4X0i84jWbzEy1+DHhNRE7ifMHP9LGLn3Aa8BYBI1R1oVs+HOeLYQPwG7DGLcsQqroaJ8mNAo65MfV2F+cD3sL5dfgPztnCC+6yFsBGEYkCPgC6qOpZVT2BU/+PcX45ngIS3NXkwwc47+FC9z1bgdNo7iveRTjXuufg/MouD3Txs7rJ1cfvGHDOHEOATTjv2Wycy0Ko6izgdZxEdpL4u97Aufb/ons56Vn3y7utG8MhnLOAQcT/X+/mxnAUeAX3LjlfVPVXnB8a7+G0j/1E/K/5rjjtKH8Dc4FXVPX7JN+llH2K037wD07D9xNuDFuAHjg3DhzGSWitVfW8u92TbtlxoDvOe+OvATiXm/5xj+19B18RnB86x3AuPR0BRqS6VtmcqNqAQdmdiJTBuUskr2biMxPGXAoiEo7TUP5xoGMxCdkZhDHGGJ8sQRhjjPHJLjEZY4zxyc4gjDHG+BTwjtXSo0SJElqmTJk0b3/q1CkKFiyY8oo5RG6rL1idcwurc+pEREQcVtXEzxJdJFsniDJlyrB69eo0bx8eHk6jRo0uXUBZXG6rL1idcwurc+qISHJP0HvYJSZjjDE+WYIwxhjjkyUIY4wxPmXrNghfLly4wN69ezl7NqXOSaFo0aJs3rw5E6LKGnJbffPnz49cNEifMcZfOS5B7N27l8KFC1OmTJkUvxxOnjxJ4cKFk10nJ8lN9VVVjhw5kuvubDHmUspxl5jOnj1L8eLF7ZdjLiciFC9enKCgoJRXNsb4lOMSBGDJwQD2d2ByptjYWPr378+OHTsy/Fg57hKTMcbkZOPHj2fMmDHkyZOHHTt2MHx4hvXcnzPPILKSoUOHMmLEpe1G/vjx44wZM+aS7jOxQoUKpbySlzJlynD48OFUH6d3797Mnj3b7/V37dpF5cqVU30cY3KCgwcP8vzzzwPOmUSePBn7FW4JIhtKa4KIicl1A2IZk6MMHjyY48edwRavueYa/v3vf2fo8SxBZIDXX3+dm2++maZNm7JlyxZP+fbt22nRogWhoaE0aNCAP/74A3B+RT/xxBPUrVuXcuXKeX5RR0VF0aRJE2rUqEGVKlX44osvAHj++efZvn071atXZ9CgQYSHh9OqVSvPcQYMGMCkSZMA55f9a6+9Rv369Zk7d26SMezcuZM6depQq1YtXnrppSTr1q5dO0JDQ6lUqRLjxo3zuc6UKVOoWrUq1apVo2fPngDs3r2bJk2aULVqVZo0acJff/3lWX/JkiUX1V1VGTRoEJUrV6ZKlSqEhYX5PJYxucWSJUuYPDl+JOInnniCyy67LGMPqqrZdgoNDdXENm3a5HmNM6B7hkxJWb16tVauXFlPnTqlJ06c0PLly+s777yjqqp33XWXbt26VVVVV6xYoY0bN1ZV1QceeEDvu+8+jYmJ0Y0bN2r58uVVVfXChQt64sQJVVU9dOiQli9fXmNjY3Xnzp1aqVIlzzEXL16s9957r2e+f//+OnHiRFVVLV26tP7nP/9RVdXIyMgkY2jdurVOnjxZVVVHjRqlBQsW9Fm/I0eOqKrq6dOntVKlSnr48GHPcQ4dOqS///673nTTTXro0KEE67dq1UonTZqkqqoTJkzQtm3bJlv32bNna9OmTTU6Olr/+ecfvf766/Xvv/++qO4pWbNmjd/r5hSLFy8OdAiZLqfX+dy5c1qxYkXP90/Hjh3TVWdgtfrxHWuN1JfYzz//TPv27SlQoAAAbdq0AZyzgV9++YX777/fs+65c+c8r9u1a0eePHmoWLEiBw4cAJzk/cILL7BkyRLy5MnDvn37PMtSo3PnzinGsGzZMubMmQNAz549ee6553zua+TIkcydOxeAPXv28Oeff1K8eHHP8h9//JH77ruPEiVKAHDFFc7QycuXL+fzzz/37H/w4MHJ1n3p0qV07dqVoKAgSpYsScOGDVm1ahVVq1ZNdf2Nye7ee+89Nm3aBDjtg++//z7btm3L8ONagsgAvm6vjI2N5fLLL2fdunU+t8mXL5/ntbqDOE2bNo1Dhw4RERFB3rx5KVOmjM8nxIODg4mNjfXMJ14n7mGxlGJI6bbQ8PBwfvjhB5YvX06BAgVo1KjRRcdSVb9uL/Vex1fd4/41Jrfbvn07r776qmf+1Vdf5brrrsuUBJGj2yBSOn2KjIxM8+WtpNx5553MnTuXM2fOcPLkSb788ksAihQpQtmyZZk1a5YntvXr1ycb/4kTJ7jqqqvImzcvixcvZvdup4fewoULc/LkSc96pUuXZtOmTZw7d44TJ06waNEin/tLLoZ69eoxY8YMwElMScVTrFgxChQowB9//MGKFSsuWqdJkybMnDmTI0eOAHD06FEA6tatm2D/9evXT7bud955J2FhYcTExHDo0CGWLFlC7dq1k93GmJxGVXn00Uc5c+YMANWrV+eJJ57ItOPn6AQRCDVq1KBz585Ur16djh070qBBA8+yadOmMWHCBKpVq0alSpU8jc5J6d69O6tXr6ZmzZpMmzaNW265BYDixYtTr149KleuzKBBg7j++uvp1KkTVatWpXv37tx2221J7jOpGD744ANGjx5NrVq1OHHihM9tW7RoQXR0NFWrVuWll17ijjvuuGidSpUqMWTIEBo2bEi1atUYOHAg4FyamjhxIlWrVuXTTz/lgw8+SLbu7du39zR033XXXbz99ttcffXVyW5jTE4zbdo0vv/+ewDy5MnDuHHjCA7OxAs/af0F7ccv7OuBxcBmYCPwpFt+BfA98Kf7bzG3XICRwDZgA1AjpWOk1EidksjISL/XzQlyW31VrZE6t8iJdT58+LCWKFHC0zD91FNPJVieGY3UGXkGEQ08o6q3AncA/UWkIvA8sEhVbwQWufMALYEb3akv8L8MjM0YY7K0Z5991vPw6fXXX8+wYcMyPYYMSxCqul9V17ivT+KcSVwLtAXibuadDLRzX7cFprgJbgVwuYiUyqj4jDEmq/rxxx89zzIBjB49OtW9G1wKmdIGISJlgNuAlUBJVd0PThIBrnJXuxbY47XZXrfMGGNyjTNnzvDII4945u+//35at24dkFgyvLVDRAoBc4CnVDUymVsgfS246HYhEemLcwmKkiVLEh4enmB50aJFE9zhk5yYmBi/180Jclt9wWljS/w3ktNFRUVZnbOxCRMmeG5hLViwIJ07d/ZZt0ypsz8NFWmdgLzAd8BAr7ItQCn3dSlgi/t6LNDV13pJTdZInTq5rb6q1kidW+SUOv/2228aHBzsaZj+6KOPklw3WzdSi3OqMAHYrKrvei2aDzzgvn4A+MKrvJc47gBOqHspyhhjcrqYmBgeeughoqOjAefZpIcffjigMWVkG0Q9oCdwl4isc6d7gLeAZiLyJ9DMnQf4BtiBc5vreOCxDIwty0iqO/DM7Nb6/fff5/Tp0575S9kYlrgjQW9p7SLcl0A04BlzKb3//vusXLkSgLx58zJu3LgM7847JRnWBqGqS/HdrgDQxMf6CvTPqHhM0t5//3169Ojh6T8qPaKjozP3QR5jcoCtW7fy4osveuZfeuklKlasmPQGkZGZEJU9SZ0hUtvddZyIiAiqVatGnTp1GD16tM99P/bYY8yfPx9wnjb+17/+BTgNW3F/YFOnTqV27dpUr16dRx55xDMOxNNPP03NmjWpVKkSr7zyCuA84fz333/TuHFjGjdu7DnOkCFDqFatGnfccYenA71Dhw7RsWNHatWqRa1atVi2bBngnAX17duX5s2b06tXr4tijoyMpH379lSsWJF+/fol6DcqTlIxFypUyGcs/nZPbkxWFxsbS58+fTz9mlWvXt0zKJBPx49DrVpUGDkSLlzI2OD8aajIqlOKjdSQcVMSUtvd9SuvvOLpDrxKlSoaHh6uqqrPPvusz26tp0+frs8++6yqqtaqVUtvv/12VVXt3bu3LliwQDdt2qStWrXS8+fPq6rqo48+6unGe9euXaqqGh0drQ0bNtT169eranxX3fFvGzp//nxVVR00aJAOGzZMVVW7du2qP//8s6qq7t69W2+55RZPHWrUqKGnT5++KN7Fixdrvnz5dPv27RodHa1NmzbVWbNmJThucjEnFYu/3ZNbI3XukJ3rPHLkSE+jdHBwsK5duzbplaOjVe++O/57qFUr1djYVB+TQDdS51bJdXfdrVs3wOnueunSpQm2O3HiBMePH6dhw4aedXxp0KABP//8M5s2baJixYqULFmS/fv3s3z5curWrcuiRYuIiIigVq1aVK9enUWLFnkGN587dy41atTgtttuY+PGjZ7ugxMLCQnxtBuEhoaya9cuAH744QcGDBhA9erVadOmDZGRkZ7bZtu0aZPk4CW1a9emXLlyBAUF0bVr14vqnlzMScWybNkyunbtmux7ZUxWt2PHjgRnC88//zzVq1dPeoPBg+G77+Lnu3cHP3pPTiu7WHyJqaa+u+vUbHfttddy7NgxFixYwJ133snRo0eZOXMmhQoVonDhwqgqDzzwAG+++WaC7Xbu3MnIkSOJiIigWLFi9O7d22fX4eA0kMXFEhQU5LmrIjY2luXLl/tMBHFdivtTV1919xVzcrH42o8x2Ymq8tBDD3luEKlUqVKCdoiLTJgA78bfELq7e3dKd+mSoTHm7DOIFC4UnYyMTPtFpiSktbvryy+/nKJFi3p+XSfV5TZAnTp1eP/997nzzjtp0KABI0aM8PQa26RJE2bPns3Bgwc9x9+9ezeRkZEULFiQokWLcuDAAb799lvP/hJ3H56U5s2bM2rUKM98UuNKJPbrr7+yc+dOYmNjCQsLu6juScWcHH+6JzcmKxs3bhyLFy8GnJ5aJ06cmGBslATCw6Ffv/j5tm3Z6bY/ZqScnSACID3dXU+cOJH+/ftTp06dZMeabdCgAdHR0VSoUIEaNWpw9OhRT4KoWLEiw4cPp3nz5lStWpVmzZqxf/9+qlWrRtWqValUqRL/+te/qFevnmd/ffv2pWXLlgkaqX0ZOXIkq1evpmrVqlSsWJGPPvrIr/ekTp06PP/881SuXJmyZcvSvn37BMuTijk5/nRPbkxW9ddffzFo0CDP/LPPPkutWrV8r/znn9ChA8SdPVerBlOnQmbcAutPQ0VWnexJ6tTJbfVVtUbq3CI71Tk2NlbvvvtuT8P0zTff7PMGD1VVPXpU9aab4q9dXH216l9/qWo2f5LaGGPMxcaOHct3bkOziPDJJ5/4vmJw4QLcdx9s3erM588P8+fD9ddnWqyWIIwxJpNs27aNZ555xjM/cOBA6tate/GKqjBgAPz4Y3zZ5MmQ1GWoDJIjE4Qm04hscg/7OzBZSUxMDL169fLctRTX9ubTBx/AuHHx88OGQadOmRBlQjkuQeTPn58jR47Yl0Mup6ocOXLE80S2MYH2zjvvsHz5cgCCg4P59NNPyZ8//8Urfv01uDe3AM6zDkOGZFKUCeW45yCuu+469u7dy6FDh1Jc9+zZs74/oBwqt9U3f/78nDp1KtBhGMOGDRt4+eWXPfOvvPIKNWrUuHjF336DLl3ib6WvUwc+/jhDH4ZLTo5LEHnz5qVs2bJ+rRseHs5tt92WwRFlHbmtvkCKz1MYk9HOnTtHz549ueD2m1S7dm3ffS3t2wf33ANRUc58mTIwb57TOB0gOe4SkzHGZCVDhw5lw4YNgHNWO2XKlIt7PI6MdJLD3r3OfOHC8OWXcNVVBJIlCGOMySC//PILb7/9tmf+7bff5uabb064UtztrG4SITgY5syBTBoPJjmWIIwxJgOcPHmSXr16ebq3b9KkCf37JxryRhX69oXvv48vGz8emjXLxEiTlpFDjn4iIgdF5HevsjCv0eV2icg6t7yMiJzxWuZfHw7GGJNFPf7442zfvh2AIkWKMHHixItHiHvtNZg0KX5+6FDo3TuzQkxRRjZSTwJGAVPiClS1c9xrEfkv4N2JznZVTaafW2OMyR5mzJjB5MmTPfP/+9//uD7xE9ATJzoJIc6DD4LXnU5ZQUYOObpERMr4WiZOP82dgLsy6vjGGBMIu3fvpp9Xz6s9evTwjAXjsXChc2kpTvPmMHZswG5nTYpk5ANlboL4SlUrJyq/E3hXVWt6rbcR2ApEAi+q6s9J7LMv0BegZMmSoXFdPqdFVFRUrhrsPrfVF6zOuUVWqXNMTAxPP/00v/32GwClSpVi/PjxCcZLKbRtG9WfeILgM2cAiCpfnrUffEBMMmOq+JKeOjdu3Dgi7vs3Wf706JfWCSgD/O6j/H/AM17z+YDi7utQYA9QJKX9++rNNTWyUw+Ql0Juq6+q1Tm3yCp1fu211zy9tAYFBeny5csTrrB7t+o118T3znrddap796bpWDmyN1cRCQY6AGFxZap6TlWPuK8jgO3ATZkdmzHGpNXy5ct59dVXPfOvvPIKd9xxR/wKhw87l5L+/tuZL1IEvv0Wrr02kyP1XyBuc20K/KGqe+MKRORKEQlyX5cDbgR2BCA2Y4xJtcjISLp37+7p+6tBgwa88MIL8StERTkPwm3Z4syHhMDcuVniWYfkZORtrtOB5cDNIrJXRPq4i7oA0xOtfiewQUTWA7OBfqp6NKNiM8aYS6l///7s3LkTgKJFizJ16lSCgoKchefPOyPCrVrlzIs4I8LdlfXv0cnIu5i6JlHe20fZHGBORsVijDEZZdKkSUydOtUzP3bsWG644QZnJjYWevVK+CDc6NFw//2ZHGXa2JPUxhiTRhs3buSxxx7zzPfu3ZvOnd3HvVThySchLCx+g1dfhUcfzeQo084ShDHGpMGpU6fo1KkTZ9zbVW+99VZGjRoVv8Lrr4P3fP/+8NJLmRxl+liCMMaYNHj88cfZtGkTAJdddhmzZs2Kf95h7NiEyaBzZxg5Mss9CJcSSxDGGJNKU6ZMYeLEiZ75UaNGUalSJWdmzpyEl5GaNYMpUyBxP0zZQPaL2BhjAmjz5s086pUAevbsyYMPPujMLFgAXbvGjwhXq5aTMEJCAhBp+lmCMMYYP50+fZpOnTpx+vRpAG655RbGjBmDiMBPP0H79s74DgA33+yML124cAAjTh9LEMYY46cnn3yS3393RjDInz8/M2fOdPpD+vVXaNUKzp51Vixd2rm19corAxht+lmCMMYYP0yaNImPP/7YM//hhx9SpUoVWL8e7r47fizpUqXghx8gcffe2ZAlCGOMScG6desStDt069aNPn36wB9/OI3Qx487C4oXd5JDhQoBivTSsgRhjDHJOHbsGB06dOCse/moUqVKjB07Ftm1C5o2hUOHnBWLFHHGeahYMXDBXmKWIIwxJgmxsbH06NHD089S4cKF+fzzzyl04oSTHPbtc1YsUMDpmbVGjQBGe+ll5JCjxhiTrQ0fPpxvvvnGMz958mRuKlYM7rwTdrgdTufLB/PnQ926AYoy41iCMMYYHxYsWMBQrzGjBw8eTPsGDZxeWP/4wykMDobZs6FJk8AEmcHsEpMxxiSya9cuunXrFjfiJY0bN+b1gQOdROAOJ0qePDBtmnN7aw5lCcIYY7ycOXOGjh07cuzYMQCuvfZawsaMIbhFC9iwwVlJBCZPhk6dAhhpxrMEYYwxLlWlb9++rFmzBoC8efMyd8IEruzWDdatc1YSgUmToEePwAWaSVJsgxCRgT6KTwARqroume0+AVoBB1W1sls2FHgYcO8L4wVV/cZd9m+gDxADPKGq36WiHsYYk27vvvtugsF//vf669R64QVYu9YpEIGJE51BgHIBf84gagL9gGvdqS/QCBgvIoOT2W4S0MJH+XuqWt2d4pJDRZyhSCu524yJG6PaGGMyw3fffcfgwfFfaU/07Mm/wsLAPZtABCZMgAceCFCEmc+fBFEcqKGqz6jqMzgJ40qccaR7J7WRqi4B/B1Xui0wQ1XPqepOYBtQ289tjTEmXbZu3Urnzp2JjY0F4O7bb+e9TZuQiIj4lcaPh7heW3MJf25zvQE47zV/ASitqmdE5FwajjlARHoBq4FnVPUYzpnJCq919rplFxGRvjhnMZQsWZLw8PA0hOCIiopK1/bZTW6rL1idc4v01PnUqVM89thjnDhxAoDyxYvz2dGj5PnzT886W555hv3ly0MWel8z5XNW1WQn4CVgDfCKO60GXgYKAtNS2LYM8LvXfEkgCOfM5XWR1VifAAAgAElEQVTgE7d8NNDDa70JQMeUYgsNDdX0WLx4cbq2z25yW31Vrc65RVrrHB0dra1atVJAAb0mXz6NuuUWVWdEB2caO/bSBnuJpOdzBlZrCt+vqpryGYSqDhORb4D6gAD9VHW1u7h7KpPRgbjXIjIe+Mqd3Qt4d314HfB3avZtjDGp9fLLL/PVV87X0FXAxquuomDcQ3AAY8ZA376BCS4LSDJBiMgVXrM73cmzTFX9bV/w3mcpVd3vzrYHfndfzwc+E5F3gWuAG4FfU7t/Y4zx1/Tp03njjTcA50tn7RVXcPmePc5CEafNoU+fwAWYBSR3BhGBc9olQCnif9GLW14uuR2LyHScu51KiMhenMtTjUSkurv9LuARAFXdKCIzgU1ANNBfVWPSViVjjEnesmXL6N27N+A0sq647DKuOur+5g0Kch6C656qCyQ5UpIJQlXLxr0WkbWqeltqdqyqXX0UT0hm/ddx2iWMMSbD7Nixg3bt2nH+/HnKAz8HB1PqzBlnYXAwzJgBHTsGNMaswt/O+jRDozDGmExw/PhxWrVqxeHDh7kFWJwnD1dHRzsLQ0Kcjvdatw5ojFmJ9eZqjMkVLly4wP3338/mzZupAvwAXOU+98Bll8G8edC8eSBDzHKSa6T27mLjqsRdbqjquxkWlTHGXEKqSv/+/fnhhx+4A/gGKBa3sGBB+PpraNgwcAFmUcmdQRT2ej0+0bwxxmQb7777LuPHj6c58DnOQ1yAM0zot9/myMF+LoXkGqlfzcxAjDEmI8ydO5dBgwZxPzAVCIlbcOWVsGBBjhsm9FLypzfXssDjOE9Fe9ZX1TYZF5YxxqTfsmXL6NatGw+r8j+8Op+74Qb4/nu46aYARpf1+dNIPQ/n9tQvgdiMDccYYy6NP/74g9atWvH02bO84b3glltg4UK4/vqkNjUufxLEWVUdmeGRGGPMJbJ//35a3H03Q44f5xnvBTVrOm0OJUoEKrRsxZ8E8YGIvAIsBDy9t6rqmgyLyhhj0igyMpLWLVvyyl9/kaBz7rvucm5lLWz32/jLnwRRBegJ3EX8JSZ1540xJss4f/483du145X160nwuFv79vDZZ5A/f6BCy5b8SRDtgXKqej7FNY0xJkBUlad79ODFxYu53XvBv/4FY8c63WiYVPFnRLn1wOUZHYgxxqTH2488wlOzZiVMDs89Bx9/bMkhjfx510oCf4jIKhK2QdhtrsaYLOGXDz7goXnzuMqdjwXkww+RAQMCGVa250+CeCXDozDGmDRa8PjjPDlvnufp6HN58hAcFkae++4LaFw5gT8jyv0EICJF/FnfGGMyy6pHHqHpuHGeL6YTwcHkX7iQoMaNAxpXTpFiG4SI9BWRA8AGnPGoI9x/U9ruExE5KCK/e5W9IyJ/iMgGEZkrIpe75WVE5IyIrHOnj9JeJWNMjqfKth49qOWVHPaFhBC0YgX5LDlcMv40Ug8CKqlqGVUtp6plVTXZ0eRck4AWicq+ByqralVgK/Bvr2XbVbW6O/XzJ3hjTC509iyHWrSgwrRpnqLf8+YlX0QEhUJDAxhYzuNPgtgOnE7tjlV1CXA0UdlCVXVH52AFcF1q92uMycUOHeJUnTpcuXChpyg8f362jh9PicqVAxhYzuRPm8K/gV9EZCUJ72J6Ip3H/hcQ5jVfVkTWApHAi6r6czr3b4zJSTZt4vzdd1Nw715P0ZT8+blj9Wo4dCiAgeVcopr8aKIi8iuwFPgNr876VHVyijsXKQN8paqVE5UPAWoCHVRVRSQfUEhVj4hIKE4HgZVUNdLHPvsCfQFKliwZOmPGjJTCSFJUVBSFChVK8/bZTW6rL1idc4piERHc+vLLhJx2LmbEAi+EhFD+ww+58aabcmSdU5KeOjdu3DhCVWumuKKqJjsBv6S0TjLblgF+T1T2ALAcKJDMduFAzZT2HxoaqumxePHidG2f3eS2+qpanXOEjz7S2KAgVVAFjQK9P18+XbZsmWeVHFdnP6SnzsBq9eM73J82iMXunUylROSKuCkNSQsRaQE8B7RR1dNe5VeKSJD7uhxwI7AjLccwxuQQMTHwzDPQrx8SEwPAXqBpvnz0/+476toocBnOnzaIbu6/3nccKZDsnUwiMh1oBJQQkb04D9z9G8gHfC8iACvUuWPpTuA1EYkGYoB+qnrU546NMTlfZCR07w5ffeUpigDuDwlh3Fdf0dDGj84U/jwoVzYtO1bVrj6KJySx7hxgTlqOY4zJYbZuhXbtYPNmT9Fc4MHgYD6bO5emTZsGLrZcxp9LTMYYkzkWLIDatRMkh/8AnYOCmDRrFvfcc0/gYsuFLEEYYwJPFd55B+69F06cAOAM0B14IU8epn72Ge3atQtoiLlRkpeYRCRY4x9qM8aYjHHmDDz0kDOgj2sP0A5YlycPn376KZ06dQpYeLlZcm0QK9zG5QXAAlXdlTkhGWNyjT17nNHeIiI8RUuBjsCRoCA+mzaNzp07Byy83C7JS0zqPETxpDv7voisEpH3RKS5+2CbMcak3ZIlUKtWguQwDmcs46PBwYSFhVlyCLBk2yBUdbeqfqSq7YC6wJdAU+BnEfk6MwI0xuQwqvDf/8Jdd8GBAwBcAB4FHgHIm5dZs2bRsWPHAAZpIBXjO6jqBeBHd0JErs2ooIwxOdTJk84Y0bNne4qO5MlDh9hYlgAhISHMmTOHVq1aBS5G45HmAYBUdd+lDMQYk8Nt2gQdO8Iff3iKIoKDaRsdzT4gX758fP7553YraxZit7kaYzJeWJjzfINXchgfEkIdNzkUKFCAL7/80pJDFuN3ghCRgimvZYwxXi5cgKeegi5d4NQpAGLy5eNf+fLR9/x5LgBFixZl4cKFNGvWLLCxmov4M+RoXRHZBGx256uJyJgMj8wYk73t2QONG8MHH3iKokqVolZMDBPPOUPLXHnllSxevJh69eoFKkqTDH/OIN4D7gaOAKjqepzO9YwxxrevvoLq1WHZMk/RX6Gh3HDgAGujnedvr7vuOn7++Wduu+22QEVpUuDXJSZV3ZOoKCYDYjHGZHfnzztddLduDUfdDpmDglhyzz2UjojgWKwz5liFChVYunQpN998cwCDNSnx5y6mPSJSF1ARCQGewL3cZIwxHjt3Om0Nv/7qKdJrr2V0/fo8HhY/unCVKlVYuHAhV199dSCiNKngT4LoB3wAXIszXsdCoH9GBmWMyWbmzIE+fTwd7QHEtmzJo5ddxjiv5FC/fn3mz59PsWLFAhGlSaVkE4Q7yltPVe2eSfEYY7KTs2fh2Wdh9Oj4suBgzg4dSpvFi/n+2289xe3bt2fatGlcdtllAQjUpEVKXW3EAG3TunMR+UREDorI715lV4jI9yLyp/tvMbdcRGSkiGwTkQ0iUiOtxzXGZIKNG+H22xMmh9KlOfLFF9SdM4fvFy3yFD/66KPMmjXLkkM2408j9TIRGSUiDUSkRtzk5/4nAS0SlT0PLFLVG4FF7jxAS5yxqG8E+gL/8/MYxpjMpAqjRkHNmrBhQ3x5hw5sDQuj1oABrF271lM8fPhwRo8eTVBQUACCNenhTxtE3Mjgr3mVKU6ni8lS1SUiUiZRcVucsaoBJgPhwHNu+RRVVZyuxi8XkVKqut+PGI0xmeHAAacvpW++iS/Llw9GjGBxxYp0aNGC48ePAxAUFMTYsWPp06dPgII16SXO93EGHsBJEF+pamV3/riqXu61/JiqFhORr4C3VHWpW74IeE5VVyfaX1+cMwxKliwZOmPGjDTHFhUVRaFChdK8fXaT2+oLVudL6Yrly7nl7bcJcRMAQFS5cmx68UU+37KFESNGEBPj3AGfP39+Xn75ZerUqXPJ4/DFPufUady4cYQ7pEPyVDXZCXjZ15TSdl7blwF+95o/nmj5Mfffr4H6XuWLgNDk9h0aGqrpsXjx4nRtn93ktvqqWp0vidOnVfv3V3UuLsVPAwdqzOnTOmTIEMW5qqCAlipVSlevXn1pY0iBfc6pA6xWP76//bnEdMrrdX6gFel7DuJA3KUjESkFHHTL9wLXe613HfB3Oo5jjEmviAjo2RM2e/2XL1UKJk/mbIMGPPjgg3ifxVetWpWvvvqK66+/3sfOTHaTYoJQ1f96z4vICGB+Oo45H3gAeMv99wuv8gEiMgO4HTih1v5gTGBcuACvvw7Dh0OMV8cJ7dvDuHEcjI2lQ9OmLPPqSqNly5aEhYVRuHDhAARsMkJauvsuAJTzZ0URmQ4sB24Wkb0i0gcnMTQTkT+BZu48wDfADmAbMB54LA2xGWPS67ffnNtXX301PjkUKADjx8OcOazbu5datWolSA6PPvoo8+fPt+SQw6R4BiEiv+FcWwQIAq4k4R1NSVLVrkksauJjXcWe0DYmcGJiYMQIePllp0+lOPXrw8SJUKECs2bN4oEHHuDMmTMAiAgjRozg6aefRkQCFLjJKP60QXiP/RcNHFDV6AyKxxgTCFu3wgMPwIoV8WX58sEbb8CTTxIrwisvvcTw4cM9i4sUKcL06dNtkJ8czJ8EEQzsVdVzItII6CgiU1T1eArbGWOyupgY+PBDeOEFcM8KAKhVCyZPhltv5eTJk/Ts2ZMvvvjCs/jGG29k/vz53HLLLQEI2mQWf9og5gAxIlIBmACUBT7L0KiMMRlv0ybn8tHTT8cnh7x5nYbpX36BW29l27Zt1K1bN0FyaN68OStXrrTkkAv4kyBi3UtKHYD3VfVpoFTGhmWMyTDnz8OwYXDbbQkvKVWt6nTVPWQIBAfz5ZdfUrNmTX7/3dOVGgMHDuTrr7+23lhzCX8SxAUR6Qr0Ar5yy/JmXEjGmAyzapXTh5J3Q3TevDB0qLOsenViYmJ46aWXaNOmDSfc7rtDQkKYNGkS//3vfwkO9ufKtMkJ/PmkH8QZE+J1Vd0pImWBqRkbljHmkjp92kkK770H7qhugHM764QJUKkSAEeOHKFbt24sXLjQs8oNN9zA7NmzqVWrVmZHbQLMnwflNuGMIofbNXdhVX0r+a2MMVnGokXQty/s2BFfVqCA8yDc44+D28tqREQEHTt2ZPfu3Z7VmjVrxmeffUaJEiUyO2qTBaR4iUlEwkWkiIhcAawHJorIuxkfmjEmXf75B7p3h6ZNEyaHJk2ch+GeegqCglBVxo4dS7169RIkhyFDhvDtt99acsjF/LnEVFRVI0XkIWCiqr4iIhtS3MoYExixsfC//8G//51gCFAuvxzefRd69wb3obYTJ07Qt29fZs6c6VmtSJEifPrpp7Rp0yaTAzdZjV/PQbid6nUChmRwPMaY9Fi3jhoDBiTsXA+ga1cnOVx9tado1apVdOnShR1eZxdVqlRhzpw53HjjjZkVscnC/LmL6TXgO2C7qq4SkXLAnxkbljEmVU6ehIEDITSUIt7JoUIFWLgQPvvMkxxUlffee4969eolSA79+vVj5cqVlhyMhz+N1LOAWV7zO4COGRmUMcZPqjBzJjzzDOzbF18eEuJcYnr+ecif31N85MgRHnzwQb788ktPWZEiRfj444+5//77MzNykw3401nfTTjjQ5dU1coiUhVoo6rDU9jUGJORNmyAJ56An35KUHzsttsoNn063HxzgvKlS5fStWtX9u7d6ymrWbMmYWFhlCvnVwfNJpfx5xLTeODfwAUAVd0AdMnIoIwxyTh2zLk99bbbEiaHq66CqVNZ/9//JkgOFy5c4MUXX6Rhw4YJksPAgQNZtmyZJQeTJH8aqQuo6q+JuvK13lyNyWwxMfDJJ07HeocPx5cHBTlnEq+8AkWLQni4Z9HmzZvp2bMnERERnrIrrriCSZMm0bp160wM3mRH/iSIwyJSHndMCBG5D7CR3ozJTMuXO2cNXl/0gPNMw8iRULFigmJVZfTo0QwaNIizZ896yhs1asSUKVNsSFDjF38SRH9gHHCLiOwDdgLd03pAEbkZCPMqKge8DFwOPAwccstfUNVv0nocY3KEnTudhmav5xQAuOEG57bVDh08zzTEOXz4MC1atEjQXUZISAhvvvkmTz31FHnypGUgSZMbJZsgRCQPUFNVm4pIQSCPqp5MzwFVdQtQ3d1/ELAPmIvT59N7qjoiPfs3Jkc4ftzpCmPkyISju+XLB88950wFCly02ezZs+nTpw+RkZGesqpVqzJ16lSqVKmSGZGbHCTZnxKqGgsMcF+fSm9y8KEJzvMVu1Nc05jc4MIFZwCfChWc4T+9k0OnTs4DcK++elFyOHLkCN27d+f+++/3JAcRYdCgQfz666+WHEyaiDMUdDIriLwEnMG5LHQqrlxVj6b74CKfAGtUdZSIDAV6A5HAauAZVT3mY5u+QF+AkiVLhs6YMSPNx4+KiqJQoUJp3j67yW31hWxUZ1WKL1tG+XHjKLBnT4JFJypWZPujjxJZubLPTX/66Sc++OADjh2L/+9SsmRJnn/+eapXr56hYWcV2eZzvoTSU+fGjRtHqGrNFFdU1WQnnDaHxNOOlLbzY78hwGGc5ysASgJBOGc1rwOfpLSP0NBQTY/Fixena/vsJrfVVzWb1Hn1atWGDVWdx97ip7JlVcPCVGNjfW72zz//aMeOHRXnBhLP1KxZMz1+/Hjm1iHAssXnfImlp87AavXje9qfJ6nLpjI5+aslztnDAfc4B+IWiMh44gcnMiZn2rIFXnoJZs1KWF60KLz4onPXUr58F22mqkybNo0nn3ySo0fjT+SvueYaxo4dS6FChShatGhGR29yAX+6+84vIgNF5HMRmSMiT4lI/pS280NXYLrXcbyHMW0P/H7RFsbkBHv2wEMPOYP0eCeH4GAnKWzbBs8+6zM57N27l9atW9OzZ88EyeGhhx5i48aNtGrVKjNqYHIJf25znQKcBD5057sCnwJp7rhFRAoAzYBHvIrfFpHqOKfJuxItMyb7O3QI3nwTxoyBc+cSLuvQAd5446LuMeLExMQwbtw4nn/++QR3KJUuXZrx48fTrFmzjIzc5FL+JIibVbWa1/xiEVmfnoOq6mmgeKKynunZpzFZVmSk88zCf/8LUVEJlzVt6iSGZIbzXLt2Lf369ePXX39NUD5gwADefPPNXNc4azKPPwlirYjcoaorAETkdmBZxoZlTA5w5gx89JGTALy7xgCoXdspb9Ikyc1PnjzJyy+/zMiRI4n1Gkf6xhtv5OOPP+bOO+/MqMiNAfxLELcDvUTkL3f+BmCziPwGqKpWzbDojMmOTp+GsWPhP/+BAwcSLqtYEYYPh3btLnoCOo6q8vnnn/Pkk0+yz6sL75CQEF544QWee+458ue/FM2AxiTPnwTRIsOjMCYnOHXKOWN4552LE0Pp0s4Dbj16OJ3rJWHnzp0MGDCAb75J2MtMkyZNGDNmDDfddFNGRG6MT/7c5mpPORuTnFOnnDGg33kHDh5MuOzaa52+lB5+2OddSfG7OMV//vMf3nnnnQSd61111VW89957dO3aFUnijMOYjOLPGYQxxpeoKOeOpBEjnDuUvF13nTOiW58+ySYGVWXGjBkMHjw4wVgNIsIjjzzCG2+8QbFixTKqBsYkyxKEMal17BiMHg0ffHBx4/P11zvjNTz4YLKJAWDNmjU88cQTLFuW8J6PGjVqMHr0aO64445LHbkxqWIJwhh/7dsH773nNEAnvl31hhucxNC7d4qJ4eDBgwwZMoQJEybEdTsDOJeT3nzzTXr37m1dcpsswRKEMSnZssVpX5gyxelt1VuZMk5ieOABCAlJdjdnz57lww8/ZPjw4QkedsubNy9PPvkkL774onWRYbIUSxDGJGXVKnjrLZg71+k+z1vlys6YDJ07Q968ye4mJiaGadOm8eKLL7InUU+t9957L++++67dnWSyJEsQxniLjYXvvnMann/88eLl9es7dyXdc0+SzzF4W7hwIYMHD2b9+oSdD9x00028//77tGzZ8lJFbswlZwnCGHAebpsyxWl4/uOPi5e3auWcMdSv79fu1q5dy+DBg/nhhx8SlF955ZUMHTqUhx9+mLwpnHkYE2iWIEzutm8fjBoF48bB0URjYAUFQdeuMHgw+Dki27Zt2xg6dCjTpk1LUF6gQAGeeeYZBg0aROHChS9V9MZkKEsQJndavdq5I2nmTIiOTrisSBGnO+7HH3caof2we/duhg0bxqRJk4iJifGUBwUF0adPH4YOHUqpUqWS2YMxWY8lCJN7XLgA8+Y5l5GW+ehvsmxZePJJ5xmGIkX82uW+fft4/fXX+fjjj7mQ6A6ntm3b8uabb3LrrbdeiuiNyXSWIEzOt3cvjB/vXEb655+LlzdoAE8/DW3aJNtPkrd//vmHt956i48++ohzicZ2aNy4McOGDaNevXqXInpjAiZgCUJEduEMRBQDRKtqTRG5AggDyuAMGtRJVY8ltQ9jkqQKP/5IpVdfhV9+Aa/LPoAzeluXLvDUUxAa6vdu9+/fz7vvvsvo0aM5c+ZMgmX16tVj2LBhNG7c+FLUwJiAC/QZRGNV9e6r4Hlgkaq+JSLPu/PPBSY0ky0dP+7cjTRmDGzZwpWJl199NfTt60zXXuv3bnfv3s3bb7/NhAkTLjpjqFWrFsOGDaN58+bWoZ7JUQKdIBJrCzRyX08GwrEEYVKiCitXwscfw/Tpzi2riTVuDI89Bm3bpvhgm7ctW7bw1ltvMXXqVKITNWZXr16d1157jVatWlliMDlSIBOEAgtFRIGxqjoOKKmq+wFUdb+IXBXA+ExWd/gwfPopTJgAGzdevLxwYfY2bcp1w4c7A/Wkwrp163jjjTeYPXt2gv6SAGrWrMmQIUNo06aN9ZlkcjRJ/MefaQcWuUZV/3aTwPfA48B8Vb3ca51jqlos0XZ9gb4AJUuWDJ0xY0aaY4iKispV4/nmiPrGxlIsIoJSX39NiWXLyJP4FlUgqlw59rVty8FmzTgRE+N3nVWVNWvWEBYWxqpVqy5aXq1aNXr06EFoaGiWPmPIEZ9zKlmdU6dx48YRqlozxRVVNeATMBR4FtgClHLLSgFbktsuNDRU02Px4sXp2j67ydb13b1bdehQ1dKlVZ2LSgmnggVV+/RR/eUX1dhYz2b+1PncuXM6efJkrVq1quKc2SaYWrZsqUuXLs24ul1i2fpzTiOrc+oAq9WP7+aAXGISkYJAHlU96b5uDrwGzAceAN5y//0iEPGZLCIyEubMcS4jLV7se5077nAG5encGVL5hPKxY8cYN24cI0eO5O+//06wTETo0KEDL7zwAjVq1EhrDYzJ1gLVBlESmOuepgcDn6nqAhFZBcwUkT7AX8D9AYrPBEp0NPzwg3Mn0rx5kOhWUgCKF4eePZ3EULlyqg/x559/MmrUKCZMmMCpU6cSLCtQoAAPPvggTz31FBUqVEhrLYzJEQKSIFR1B1DNR/kRoEnmR2QCbv16JylMmwYHDly8PE8eaNrUSQpt26Y4KE9iMTExfPvtt4waNYrvvvvuouVXX301jz/+OP369eOKK65Iay2MyVGy2m2uJjfZvh3CwmDGDPjtN9/rVKkCvXpBt25wzTWpPkRkZCQjRoxgzJgx7Ny586LllStXZuDAgXTr1o18qUw6xuR0liBM5tqzx+kgb8YMp8M8X66+2kkIvXpBtYtONP2yZs0axowZw6effsr58+cTLBMRWrVqxYABA2jWrFmWviPJmECyBGEy3j//wOzZTlLw1UkewGWXQfv2TttC06ZOVxipFBkZyWeffcb48eNZs2bNRcuLFSvGQw89xKOPPkrZsmVTvX9jchtLECZj/P2308g8Zw6EhzsjtSWWNy/cfbfTJ1KbNqm+Cwmc27RXrFjB+PHjCQsL47SPp6irV6/O448/TpcuXShQoEAaKmNM7mQJwlw627c74zd//jksX+57naAguOsuJym0bw/FivleLwWHDh3ynC1s9PEUdf78+bnvvvu4/fbb6d+/v11GMiYNLEGYtFOF3393EsLnn8OGDb7XE3G61O7SBTp2hKvS1oPKuXPn+Oqrr5gyZQrffPPNRX0jAVSpUoWHH36YHj16UKxYMcLDwy05GJNGliBM6sTEwK+/OmcKc+fCtm2+1wsKgjvvhA4doF07uO66NB1OVVm+fDlTpkwhLCyM48ePX7ROwYIF6dKlCw8//DC1a9e2hGDMJWIJwqTsxAn47jv4+mv45hunkzxf8uWDZs2cpNC6NZQokeZD7tixg6lTpzJlyhS2b9/uc506derQu3dvunTpQhE/R4AzxvjPEoS5mCps3QpffeUkhZ9/vnjc5jiFCsG99zpJoWXLNDU0x9m3bx+zZs1i5syZLE+iDaNMmTL06tWLHj16cOONN6b5WMaYlFmCMI5z55xEEJcUkrp0BM5zCvfc4ySFJk0gf/40H3b//v3MmTOHsLAwli5d6nOdIkWK0KlTJ3r16kW9evWsi21jMokliNxKFbZsgYULnSk8HBL1S5RAaCi0auVMNWo4XV+k0cGDBz1JYcmSJReNtwAQFBTE3XffzQMPPEDr1q257LLL0nw8Y0zaWILIRYJPnHCeYo5LCnv2JL1ywYJOe0KrVs7ZQqlS6Tr2rl27+OKLL/jiiy/46aefiPXxXESePHm466676NSpEx06dKB48eLpOqYxJn0sQeRk5845Q3G6CaHe6tXOmUNSypVz2hNatYKGDVPdIZ43VWX9+vV88cUXzJs3j3Xr1vlcT0Ro2LAhnTt3pkOHDlyVxltgjTGXniWInOT8eVi1yrlctHgx/PJLgu6yL7r5s0gRpw2heXPnbKF8+XQdPjo6mqVLlzJv3jzmzZvH7t27k1y3fv36dO7cmY4dO1IqnWcnxpiMYQkiO7twASIinGQQHg5Ll4KPribiaJ48yO23OwmheXOoXTtNfR55++eff1iwYAHffvstCxcu9PmcAkBISAhNmjShbdu2tGnTxpKCMdmAJYjs5Nw5pwfUpUvjE0JUVPLbVKjgOUtYFhJC/Vat0nQuLB0AAA3hSURBVBVCdHQ0K1eu5Ntvv+Xbb7/12SlenKJFi3LvvffStm1bWrRoYc8qGJPNZHqCEJHrgSnA1UAsME5VPxCRocDDwCF31RdU9ZvMji9LOXLEuUy0bJmTDFatci4jJadcOWjUCBo3dv71eoI5Ojw8TWHs2bOHRYsWsWDBAhYuXMixY8eSXPe6666jTZs2tGvXjoYNGxISEpKmYxpjAi8QZxDRwDOqukZECgMRIvK9u+w9VR0RgJgCTxV27IhPBkuXwubNKW9XunR8MmjUyJlPp6NHjxIeHs4PP/zAokWL2Lp1a5LrBgcHU79+fVq2bEnLli2pXLmydXVhTA6R6QlCVfcD+93XJ0VkM3BtZscRcMeOOWcEK1c6fRv9+iscPJjydjfeCPXrO1PjxnAJxjU4ffo0y5YtY9GiRfzwww+sWbPG57MJca655hruueceWrZsSdOmTe3SkTE5lCT3RZDhBxcpAywBKgMDgd5AJLAa5yzjomsZItIX6AtQsmTJ0BkzZqT5+FFRURQqVCjN2/tLzp+n0I4dFNm8mcKbN1Pkjz8okNwzCK7YoCCibrqJE1WqcKJyZU5UrsyFNHaPDfH1jYqKYuPGjWzYsIENGzawZcsWLly4kOR2ISEhVKlShdDQUGrXrk25cuWyzVlCZn3GWYnVOXdIT50bN24coao1U1ovYAlCRAoBPwGvq+rnIlISOAwoMAwopar/Sm4fNWvW1NVJDVvph/DwcBo1apTm7X2KjnaeUF67Nv7MYO3alNsOAIoWhTp14s8QatWCSzDAzYEDB/j555+ZPn06O3bsYMOGDT4fVIuTJ08eatasSdOmTWnSpAl169Ylfzq60wikDPmMszirc+6QnjqLiF8JIiB3MYlIXmAOME1VPwdQ1QNey8cDXwUitlQ5d84ZD2HtWlizxvl3/foEzx4kKTjYGW/59tud201vvx1uuildXVgAxMbGsmXLFlauXMnSpUv5+eefk21DiHPrrbfSpEkTmjZtSsOGDbn88svTFYcxJvsLxF1MAkwANqvqu17lpdz2CYD2wO+ZHVuyoqKcL/+4RLBmDWzcmHQvp4mVK5cwGVSv7ozDnE6HDx9m5cqVrFixgpUr/9/e/cdWVd5xHH9/rLUIsjFAUH7jRemVX4UJLREIGuPAGRmbM5rJTGa2JVtUkgnRmSwuJnNui1uWzDnnDLqYubhfOjY0/sAYjQMULK3QlpZfCvhjQxgQKYN+98fz3HJbT0tLe3vpvd9XcnLPee655z7f3tv7Pec85zzPOtavX8/Bgwc7fc1ZZ51FRUUF8+fPZ8GCBcybN8/vYHbOfUo+jiAuB5YBNZIy/S98H7hJUgXhFNNO4Nt5qFv4wd+2DWpq2k7bt3d9G6NHw8yZoVO7TFLowdgIGUePHqW6upp169a1JoXtXahXWVkZc+bMYdy4cSxbtoy5c+d6w7Jz7pTycRXTayT0+gD07T0PZpR99BGsWXMyCdTWhktLm5u7vp1UKiSCTEKYOfO0h9TMdvjwYaqrq9m4cWPrtGXLlsRhNtsbMWIEVVVVVFZWMn/+fGbPns2AAQOK8jytc+70Fd+d1O+/DzfcADU1zO2gW4hEJSVQXh6SQCYRVFSEhuUeOnDgAJs2bWqTDOrr6zu91DSjrKyMWbNmUVlZSWVlJVVVVYwfP77fXGXknDtzFV+CGDoU3nij87aDMWNg2rQwTZ0aHtPpHvVuCtDc3ExdXR01NTXU1NRQW1tLTU0N73bhkteMSZMmtUkGM2bM8LuVnXM5UXwJ4pxzwpFAbS3HBw3i7IqKk8kgkxB6cK8BwIkTJ9i1a1drAshMDQ0NXTpFBKEhOZ1OM2vWrNapoqLC2w6cc32m+BIEwBNPwLBhvNbUxMIrrjjtzRw5coSGhgbq6upap61bt9LQ0EBzN9oxSktLmTp1aptkMH36dAb2wj0Qzjl3uoozQcycGR67cAVQS0sLe/bsobGxkfr6+tYkUFdXx+7du7v91hMmTGDatGltpksuuYTS0tJub8s553KpOBNEO8eOHWPnzp00NTXR1NREY2Nj6+OOHTu6dTSQMXLkSNLpdJtEMGXKFAYPHpyDCJxzrvcVXYI4dOgQDz/8ME1NTWzYsIH9+/eze/fuTrue6EhJSQmpVIry8nLS6TTl5eWUl5czefJkPtfDdgznnMu3oksQJSUlrFy5sluvGT58OKlUismTJ7cmgfLyclKplF9B5JwrWEWXIAYOHMioUaPYu3dvm/KxY8eSSqVIpVJMmjSpzaNfOeScK0ZFlyAAVqxYAYRxEJYuXcrEiRP7bW+lzjmXK0WZIJYvXw6E7nLT6XSea+Occ2emnvUt7ZxzrmB5gnDOOZfIE4RzzrlEniCcc84l8gThnHMukScI55xziTxBOOecS6SujFp2ppL0EbCrB5sYDvy7l6rTHxRbvOAxFwuPuXvGm9n5p1qpXyeInpL0ppldlu969JViixc85mLhMeeGn2JyzjmXyBOEc865RMWeIB7JdwX6WLHFCx5zsfCYc6Co2yCcc851rNiPIJxzznXAE4RzzrlERZcgJH1V0juSWiRd1u65uyU1SqqX9IV81TEXJC2KcTVKuivf9ckFSY9J+lBSbVbZUEkvSNoWHwtqsHBJYyWtlbQ1fq/viOUFG7ekAZLWS6qOMf8wlk+UtC7G/EdJBTUesKQSSZskrY7LOY+36BIEUAt8GXg1u1DSpcCNwBRgEfCQpJK+r17vi3H8ClgMXArcFOMtNKsIn122u4CXzOxi4KW4XEiOA98zszRQBXw3fraFHHczcKWZzQAqgEWSqoAHgJ/HmD8Gbs1jHXPhDmBr1nLO4y26BGFmW82sPuGpJcBTZtZsZjuARmBO39YuZ+YAjWa23cyOAU8R4i0oZvYqsL9d8RLg8Tj/OPClPq1UjpnZPjPbGOcPEX5ARlPAcVtwOC6WxsmAK4E/xfKCilnSGOCLwKNxWfRBvEWXIDoxGng3a/m9WFYICjm2UxlpZvsg/JgCI/Jcn5yRNAGYCayjwOOOp1veBj4EXgCagANmdjyuUmjf8V8AK4GWuDyMPoi3IMeklvQicEHCU/eY2TMdvSyhrFCuAS7k2Bwg6Tzgz8ByM/tv2MEsXGZ2AqiQNAT4K5A0uHxBfMclXQt8aGZvSVqYKU5YtdfjLcgEYWZXncbL3gPGZi2PAfb2To3yrpBjO5UPJF1oZvskXUjY4ywokkoJyeFJM/tLLC74uAHM7ICkVwjtL0MknR33qgvpO345cJ2ka4ABwGcIRxQ5j9dPMZ30LHCjpDJJE4GLgfV5rlNv2QBcHK96OIfQGP9snuvUV54FbonztwAdHUH2S/Fc9O+ArWb2YNZTBRu3pPPjkQOSzgWuIrS9rAWuj6sVTMxmdreZjTGzCYT/3ZfN7Gv0RbxmVlQTsJSwR90MfAA8n/XcPYRzmfXA4nzXtZfjvgZoiPHdk+/65CjGPwD7gP/Fz/hWwrnal4Bt8XFovuvZyzHPI5xa2Ay8HadrCjluYDqwKcZcC/wgll9E2KlrBJ4GyvJd1xzEvhBY3VfxelcbzjnnEvkpJuecc4k8QTjnnEvkCcI551wiTxDOOecSeYJwzjmXyBOE63ckDZH0nXzXoyskLZc0MMfv8Zyk0e3KRkl6WdIz8S7rTPnXJdXGXlC3SLozl3Vz/ZsnCNcfDQHOiAShoLP/o+VAtxKEpC73cBBvFBtqZnvaPXU7cBuhc7eb47qLY32uNrMpwCzgYHfq5oqLJwjXH/0YSEl6W9JPASStkLRB0uas8QEmSKqT9Gjca35S0lWSXo996M+J690r6fdxj3ubpG9m3qiT7W6V9BCwERgr6deS3mw3PsHtwChgraS1sexw1ravl7Qqzq+S9GBc7wFJgxTGt9gQxwDoqPfdhcArCeUlhI7dWjjZb8/dwJ1mthfAzI6a2W+79Zd3xSXfdwb65FN3J2ACUJu1fDVhAHcRdnpWAwvieseBabH8LeCxuN4S4G/x9fcC1cC5wHBCz7ejTrHdFqAqqw5D42MJ4Qd7elzeCQzPWu9w1vz1wKo4vypuvyQu/wi4Oc4PIdwFPyjhb/FLwtgI7cvHE8Y8+TswOJbtBz6b78/Pp/4zFWRnfa7oXB2nTXH5PEJfWruBHWZWAyDpHcIgOiaphvBDn/GMmX0CfBL34ucQurHoaLu7zOxfWa+/QdK3CB1gXkgYmGlzN+N42kIvpZmYrstqIxgAjKPtgDEQOnL7VDuCme0iJDPnTpsnCFcIBNxvZr9pUxjGR2jOKmrJWm6h7fe/fZ8zdortHslankj4kZ5tZh/H00YDOqhr9vu0X+dI1ryAr1jy4FaZ970IeNfCIFBd8Q7weeDlLq7vipy3Qbj+6BAwOGv5eeAbmat1JI2W1N0BcpYojHU8jHBef0M3tvsZwo/7QUkjCUO7dlTXDySlY8P20k7q8zxwW+ytFUkzE9ZZDDzXpeiC+4GfSLogbrMstpM4l8iPIFy/Y2b/iQ3NtcAaM1shKQ28EX9PDxOu3DnR2XbaWQ/8g3Aa5z4LDbl7u7JdM6uWtImwh74deD3r6UeANZL2mdkVhLGhVxPaOWoJp62S3Efo839zTBI7gWvbrbOIcKVSl5jZP2MCezFu0whtMs4l8t5cXdGTdC+h8fhn+a5LV0kqA143s8vyXRdXuPwIwrl+yMyaAU8OLqf8CMI551wib6R2zjmXyBOEc865RJ4gnHPOJfIE4ZxzLpEnCOecc4n+DxGApEUTfy69AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "m = pyomo.ConcreteModel()\n",
    "\n",
    "S = data.keys()\n",
    "m.x = pyomo.Var(S, domain=pyomo.NonNegativeReals)\n",
    "\n",
    "def Pmix(T):\n",
    "    return sum(m.x[s]*Pvap(T,s) for s in S)\n",
    "\n",
    "m.obj = pyomo.Objective(expr = Pmix(-10), sense=pyomo.maximize)\n",
    "\n",
    "m.cons = pyomo.ConstraintList()\n",
    "\n",
    "m.cons.add(sum(m.x[s] for s in S)==1)\n",
    "m.cons.add(Pmix(30) <= Pvap_denatured(30))\n",
    "m.cons.add(Pmix(40) <= Pvap_denatured(40))\n",
    "\n",
    "solver = pyomo.SolverFactory('cbc')\n",
    "solver.solve(m)\n",
    "\n",
    "print(\"Vapor Pressure at -10°C =\", m.obj(), \"mmHg\")\n",
    "\n",
    "T = np.linspace(-10,40,200)\n",
    "plt.plot(T, Pvap_denatured(T), 'k', lw=3)\n",
    "plt.plot(T, [Pmix(T)() for T in T], 'r', lw=3)\n",
    "plt.legend(['denatured alcohol'] + ['cold weather blend'])\n",
    "plt.title('Vapor Pressure of selected compounds')\n",
    "plt.xlabel('temperature / °C')\n",
    "plt.ylabel('pressure / mmHg')\n",
    "plt.grid(True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "5dX8AIH_gNfx",
    "pycharm": {}
   },
   "source": [
    "The Pandas library is useful for summarizing the solution in tabular form."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 235
    },
    "colab_type": "code",
    "executionInfo": {
     "elapsed": 7115,
     "status": "ok",
     "timestamp": 1555699127384,
     "user": {
      "displayName": "Jeffrey Kantor",
      "photoUrl": "https://lh5.googleusercontent.com/-8zK5aAW5RMQ/AAAAAAAAAAI/AAAAAAAAKB0/kssUQyz8DTQ/s64/photo.jpg",
      "userId": "09038942003589296665"
     },
     "user_tz": 240
    },
    "id": "X4ZEPlI5gNfy",
    "outputId": "087f4e2c-e502-4f80-8921-84448020b47e",
    "pycharm": {}
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>A</th>\n",
       "      <th>B</th>\n",
       "      <th>C</th>\n",
       "      <th>MW</th>\n",
       "      <th>SG</th>\n",
       "      <th>mole fraction</th>\n",
       "      <th>mass fraction</th>\n",
       "      <th>vol fraction</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>ethanol</th>\n",
       "      <td>8.04494</td>\n",
       "      <td>1554.30</td>\n",
       "      <td>222.65</td>\n",
       "      <td>46.07</td>\n",
       "      <td>0.791</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.0000</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>methanol</th>\n",
       "      <td>7.89750</td>\n",
       "      <td>1474.08</td>\n",
       "      <td>229.13</td>\n",
       "      <td>32.04</td>\n",
       "      <td>0.791</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.0000</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>isopropyl alcohol</th>\n",
       "      <td>8.11778</td>\n",
       "      <td>1580.92</td>\n",
       "      <td>219.61</td>\n",
       "      <td>60.10</td>\n",
       "      <td>0.785</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.0000</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>acetone</th>\n",
       "      <td>7.02447</td>\n",
       "      <td>1161.00</td>\n",
       "      <td>224.00</td>\n",
       "      <td>58.08</td>\n",
       "      <td>0.787</td>\n",
       "      <td>0.428164</td>\n",
       "      <td>0.2906</td>\n",
       "      <td>0.311695</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>xylene</th>\n",
       "      <td>6.99052</td>\n",
       "      <td>1453.43</td>\n",
       "      <td>215.31</td>\n",
       "      <td>106.16</td>\n",
       "      <td>0.870</td>\n",
       "      <td>0.571836</td>\n",
       "      <td>0.7094</td>\n",
       "      <td>0.688305</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>toluene</th>\n",
       "      <td>6.95464</td>\n",
       "      <td>1344.80</td>\n",
       "      <td>219.48</td>\n",
       "      <td>92.14</td>\n",
       "      <td>0.865</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.0000</td>\n",
       "      <td>0.000000</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                         A        B       C      MW     SG  mole fraction  \\\n",
       "ethanol            8.04494  1554.30  222.65   46.07  0.791       0.000000   \n",
       "methanol           7.89750  1474.08  229.13   32.04  0.791       0.000000   \n",
       "isopropyl alcohol  8.11778  1580.92  219.61   60.10  0.785       0.000000   \n",
       "acetone            7.02447  1161.00  224.00   58.08  0.787       0.428164   \n",
       "xylene             6.99052  1453.43  215.31  106.16  0.870       0.571836   \n",
       "toluene            6.95464  1344.80  219.48   92.14  0.865       0.000000   \n",
       "\n",
       "                   mass fraction  vol fraction  \n",
       "ethanol                   0.0000      0.000000  \n",
       "methanol                  0.0000      0.000000  \n",
       "isopropyl alcohol         0.0000      0.000000  \n",
       "acetone                   0.2906      0.311695  \n",
       "xylene                    0.7094      0.688305  \n",
       "toluene                   0.0000      0.000000  "
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "s = data.keys()\n",
    "results = pd.DataFrame.from_dict(data).T\n",
    "for s in S:\n",
    "    results.loc[s,'mole fraction'] = m.x[s]()\n",
    "    \n",
    "MW = sum(m.x[s]()*data[s]['MW'] for s in S)\n",
    "for s in S:\n",
    "    results.loc[s,'mass fraction'] = m.x[s]()*data[s]['MW']/MW\n",
    "    \n",
    "vol = sum(m.x[s]()*data[s]['MW']/data[s]['SG'] for s in S)\n",
    "for s in S:\n",
    "    results.loc[s,'vol fraction'] = m.x[s]()*data[s]['MW']/data[s]['SG']/vol\n",
    "\n",
    "results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {},
    "colab_type": "code",
    "id": "xU12DmVdgNf2",
    "pycharm": {}
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "colab": {
   "collapsed_sections": [],
   "name": "02.03-Mixture-Design-Cold-Weather-Fuel.ipynb",
   "provenance": [],
   "version": "0.3.2"
  },
  "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.9.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}