Trabalho Final da Disciplina COQ712 - Termodinâmica de Soluções - Parte 2 (linha de inversão térmica).ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Cálculo da linha de inversão térmica do SH2 com CPA EOS"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from scipy import optimize as opt\n",
    "from sympy import *"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "R = 8.31447"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def load_param():\n",
    "    a = 3.47972 \n",
    "    b = 0.0285 \n",
    "    c = 0.41107\n",
    "    Tc = 373.3 #http://webbook.nist.gov/cgi/cbook.cgi?ID=C7783064&Mask=4\n",
    "    Pc = 89.70e5 \n",
    "    #conversão de unidades\n",
    "    a = 0.1*a\n",
    "    b= b/1000\n",
    "    \n",
    "    association_type = '2B'\n",
    "    eps = 80.8884\n",
    "    beta = 0.08581\n",
    "    eps = 100*eps\n",
    "    \n",
    "    return a,b,c,Tc,Pc,association_type, eps, beta"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "def T_calc(V):\n",
    "    a,b,c,Tc,Pc,association_type, eps, beta = load_param()\n",
    "    def eta_zero(T):\n",
    "        eta_zero=-T*(-4.157235*T*(1.0 + 0.475*b/(V*(1.0 - 0.475*b/V)))*(-0.060136124130582*V*eps*(1.0 - 0.475*b/V)*(sqrt(1.0 + 4.0*b*beta*(exp(0.120272248261164*eps/T) - 1.0)/(V*(1.0 - 0.475*b/V))) - 1.0)*exp(0.120272248261164*eps/T)/(T**2*b*beta*(exp(0.120272248261164*eps/T) - 1.0)**2) + 0.120272248261164*eps*exp(0.120272248261164*eps/T)/(T**2*sqrt(1.0 + 4.0*b*beta*(exp(0.120272248261164*eps/T) - 1.0)/(V*(1.0 - 0.475*b/V)))*(exp(0.120272248261164*eps/T) - 1.0)))/V + 8.31447/(V - b) - 4.157235*(1.0 + 0.475*b/(V*(1.0 - 0.475*b/V)))*(-0.5*V*(1.0 - 0.475*b/V)*(sqrt(1.0 + 4.0*b*beta*(exp(0.120272248261164*eps/T) - 1.0)/(V*(1.0 - 0.475*b/V))) - 1.0)/(b*beta*(exp(0.120272248261164*eps/T) - 1.0)) + 1)/V + a*c*sqrt(T/Tc)*(c*(-sqrt(T/Tc) + 1) + 1)/(T*V*(V + b)))/(-8.31447*T/(V - b)**2 - 4.157235*T*(1.0 + 0.475*b/(V*(1.0 - 0.475*b/V)))*(-0.5*V*(1.0 - 0.475*b/V)*(-2.0*b*beta*(exp(0.120272248261164*eps/T) - 1.0)/(V**2*(1.0 - 0.475*b/V)) - 0.95*b**2*beta*(exp(0.120272248261164*eps/T) - 1.0)/(V**3*(1.0 - 0.475*b/V)**2))/(b*beta*sqrt(1.0 + 4.0*b*beta*(exp(0.120272248261164*eps/T) - 1.0)/(V*(1.0 - 0.475*b/V)))*(exp(0.120272248261164*eps/T) - 1.0)) - 0.5*(1.0 - 0.475*b/V)*(sqrt(1.0 + 4.0*b*beta*(exp(0.120272248261164*eps/T) - 1.0)/(V*(1.0 - 0.475*b/V))) - 1.0)/(b*beta*(exp(0.120272248261164*eps/T) - 1.0)) - 0.2375*(sqrt(1.0 + 4.0*b*beta*(exp(0.120272248261164*eps/T) - 1.0)/(V*(1.0 - 0.475*b/V))) - 1.0)/(V*beta*(exp(0.120272248261164*eps/T) - 1.0)))/V - 4.157235*T*(-0.475*b/(V**2*(1.0 - 0.475*b/V)) - 0.225625*b**2/(V**3*(1.0 - 0.475*b/V)**2))*(-0.5*V*(1.0 - 0.475*b/V)*(sqrt(1.0 + 4.0*b*beta*(exp(0.120272248261164*eps/T) - 1.0)/(V*(1.0 - 0.475*b/V))) - 1.0)/(b*beta*(exp(0.120272248261164*eps/T) - 1.0)) + 1)/V + 4.157235*T*(1.0 + 0.475*b/(V*(1.0 - 0.475*b/V)))*(-0.5*V*(1.0 - 0.475*b/V)*(sqrt(1.0 + 4.0*b*beta*(exp(0.120272248261164*eps/T) - 1.0)/(V*(1.0 - 0.475*b/V))) - 1.0)/(b*beta*(exp(0.120272248261164*eps/T) - 1.0)) + 1)/V**2 + a*(c*(-sqrt(T/Tc) + 1) + 1)**2/(V*(V + b)**2) + a*(c*(-sqrt(T/Tc) + 1) + 1)**2/(V**2*(V + b)))-V\n",
    "        return eta_zero\n",
    "    \n",
    "    T0 = Tc\n",
    "    try:\n",
    "        Temp = opt.newton(eta_zero, T0, fprime=None, args=(), tol=1.48e-10, maxiter=10000, fprime2=None)\n",
    "    except: ValueError\n",
    " \n",
    "\n",
    "    return Temp"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.02461723168023686"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "V0 = R*300/101325 #300K y 1 atm=101325 Pa\n",
    "V0 #m3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY4AAAETCAYAAAAh/OHhAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzt3Xl8XHW9//HXJ3vTfQ3dW6AU2soaW0TgFwShqJd6vaLUDQGtesHld69X4Xp/4vJA/elPvSIut0pZrkqtC1K1gMgltChLWyhtWigNpUu60CVp2jRplsnn98eclGk6k8xJMjPJzPv5YB6Z+ZzvOefTL0k+Oed7zveYuyMiIpKsvEwnICIiA4sKh4iIhKLCISIioahwiIhIKCocIiISigqHiIiEosIhIiKhqHCIiEgoKhwiIhJKQao2bGZLgHcB+9x9ThD7NTAzaDICOOTu5wbLbgNuAiLAZ9z90SA+H/gBkA/83N2/1d2+x4wZ49OmTetx7kePHmXw4ME9Xj/XqL/CUX+Fo/4Kpzf9tXbt2gPuPrbbhu6ekhdwKXA+UJVg+XeBLwfvZwEvAsXAdOBVooUiP3h/KlAUtJnV3b4vuOAC740nnniiV+vnGvVXOOqvcNRf4fSmv4A1nsTv95Qdcbj7SjObFm+ZmRnwPuBtQWgBsNTdm4HXzKwamBssq3b3rcF6S4O2m1KVt4iIdC1TYxyXAK+7+5bg80RgZ8zymiCWKC4iIhmSsiOObiwEHoj5bHHaOPELW9zpfM1sEbAIoKysjMrKyh4n19DQ0Kv1c436Kxz1Vzjqr3DS0V9pLxxmVgC8B7ggJlwDTI75PAnYHbxPFD+Buy8GFgOUl5d7RUVFj3OsrKykN+vnGvVXOOqvcNRf4aSjvzJxquoK4GV3r4mJLQeuM7NiM5sOzACeA1YDM8xsupkVAdcFbUVEJENSVjjM7AHgaWCmmdWY2U3Bous48TQV7r4RWEZ00PsR4GZ3j7h7G3AL8CjwErAsaCsiIhmSyquqFiaIfzRB/A7gjjjxFcCKPk1ORER6LFOD4yKShdyddodIu9Pu0Vf0PbQfj0XbRfyNuDvH23csd6KxHYcjVO2qh+BzR1sHPKZtovd0fMaDrzHLggbHlzkx67yxrY5WxHyO3XbHNjih/YntOrYZv99i23vceKL2nfdRs6+Nivir9RkVDpEBwN1pbmunqSVCU2uEY60RjrW209Qaobk1wrG2CM2t7TS3tdPS1k5zW4Tmtjc+t0Taae34GonGWyNOWyT42h6Nd8Qi7U5rJPpLv629nbZ2py34HAmKQVuknXaHtvZo+44CkRJ/fypFG84+pw7P43Mp3ocKh0gKuTtHWyLUN7VyuKmV+qZWnn+9jf1rdnL4WBtHjrVytLmNhuY2jhxrO/6+oTlCU0sbR1siNLVEaGxp69Uv5cJ8oyg/j8KCPArz86Lv843C/DwK8vMoyjcK8vMoyDNKiwooyDcK8oyCvDzyg/f5eR1f845/PuFlRl7wNT+P4+/zgniecfy9wfF17HgcjDeW51l0nU2bNjJnzhzyLIgH7YL/onGLxiyIEWzTgOj9xhxfZkGDNz4Hy4+3saANMeu9sY0TvnaOc+LymK2csL2To2/kcXI8/v9T63QXQ0e7Nc89E3+FPqTCIRJSc1uEAw0t7D/SfMLrQEMztY0t1B1tofZoC3WNLdQdbaUl0n7yRl5Yf/xtaVE+g4sLGFpcwODiAoYUFzBxRCGlRQWUFuUzqCifwUUFDCrKZ1Bh/vGvJYV5FBd2vM+nuCAv+irMpyg/j+LCaIEoLsg74ZfSQDPo4GYqZp+S6TQGjFdLUn+xrAqHSAx3Z39DMzsONlJT18Tu+ib2HDrGnvom9tQfY0/9MWqPtsRdd0RpIaMGFzGqtIjJo0o5Z9IIRg4uYmRpISNKCxlWUsjwQYW8svFFLr/kLQwrKWRISQH5eQP3l7rkJhUOyUm1R1vY8voRtuxrYNuBo2yvbWTHwUZ21DbS1Bo5oe2wkgImjBjE+OElnDN5BOOHlTBuWDFjhxYzZkj06+jBxRQVJPeXXktNPpNHlabinyWSFiocktWaWiJs2nOYTbvreeX1BrbsO8KW1xs4GHPUUFyQx5RRpUwdXcpbTx/D1NGlTBldyuSRgxg/fBCDi/VjIhJLPxGSNY61Ri/bXF9TT9Xueqp21VO9r+H4oPLQkgLOKBvK22eVcfq4IZw+bggzyoYyflgJeTpdJJI0FQ4ZsOqbWnl+ex3Pbatl9Wu1rK+pPz4QPXZoMW+aOJz5s09hzsThzJk4nPHDSwb0ILFIf6HCIQNGS1s7a7bVUvnKflZtOcDLew/jDgV5xpyJw/noW6dRPnUk50weQdmwkkynK5K1VDikX9tT30Tl5v088fI+/lZ9gKMtEQrzjfKpo/js5TOYO20U504ZQWmRvpVF0kU/bdLv7D7UxJ/W72b5i7up2nUYgIkjBrHgvIlcNnMcF502WgPWIhmknz7pFw42NLNiwx6Wv7ib1dvqADhn0nBuu/pMLjtzHDPGDdH4hEg/ocIhGRNpdyo37+MXz2xn5ZYDRNqdM8qG8Pkrz+BdZ09g2pjBmU5RROJQ4ZC0qz3awrI1O/nFM9upqWti3NBiPnHpqVxz7gTOPGVYptMTkW6ocEjavLjzEPc/vZ0/rt9NS1s786aP4rarz+LK2WUU5mfiYZQi0hMqHJJyz++o4/uPvcKqLQcYXJTP+8sn8+G3TOWMsqGZTk1EekCFQ1Jm3c5DfP+xV3jylf2MGlzEbVefyQfmTWFoSWGmUxORXlDhkD734s5D/OdfX+GJzfsZWVrIrVefyYcvnKpLaEWyhH6Spc/sO3yMO1a8xEPrdjOitJAvzJ/J9W+ZpoIhkmVSNiJpZkvMbJ+ZVXWKf9rMNpvZRjP7dkz8NjOrDpZdFROfH8SqzezWVOUrPRdpd5Y89RqXf/dJHt6wl8+87XSe+uLb+OeK01U0RLJQKn+q7wXuAu7vCJjZZcAC4Gx3bzazcUF8FnAdMBuYAPzVzM4IVvsR8HagBlhtZsvdfVMK85YQ1m6v5StPH2PnkU1cesZYvnrNbKbr/guRrJaywuHuK81sWqfwp4BvuXtz0GZfEF8ALA3ir5lZNTA3WFbt7lsBzGxp0FaFI8PqjrbwzYdfYtmaGkaVGD/54PnMn3OK7u4WyQHpPo9wBnCJmd0BHAM+7+6rgYlA7BPWa4IYwM5O8XnxNmxmi4BFAGVlZVRWVvY4yYaGhl6tn+1eOhjhv9Y3c6TFecf0Qi4/pYVBBzfz5JObM53agKDvr3DUX+Gko7/SXTgKgJHAhcCbgWVmdioQ789UJ/4YjMfbsLsvBhYDlJeXe0VFRY+TrKyspDfrZ6u2SDs/eHwLd62pZvqYwfxy4XnMnjBc/RWS+isc9Vc46eivdBeOGuD37u7Ac2bWDowJ4pNj2k0CdgfvE8UljWrqGvns0nWs3V7HtRdM4qsLZmsqc5Ecle6f/D8AbwMqg8HvIuAAsBz4lZl9j+jg+AzgOaJHIjPMbDqwi+gA+gfSnHPOe3jDHr74u/W0O/zgunNZcO7E7lcSkayVssJhZg8AFcAYM6sBbgeWAEuCS3RbgOuDo4+NZraM6KB3G3Czu0eC7dwCPArkA0vcfWOqcpYTtUXa+fqfNnHf09s5Z9Jw7lx4HlNH64opkVyXyquqFiZY9KEE7e8A7ogTXwGs6MPUJAmNLW18+lcv8PjL+7jp4ul8cf6ZFBVoIkIR0Z3jEsf+I83cdN9qqnbV8/V3z+HDF07NdEoi0o+ocMgJXt3fwEfveY4DR1pY/OFyrphVlumURKSfUeGQ49Zsq+Vj968h34yliy7knMkjMp2SiPRDKhwCRK+c+uyv1zFpxCDuvWEuU0aXZjolEemnVDiEP6/fw6cfeJ7zpozk5x8pZ+TgokynJCL9mApHjnvi5X187tcvcP6Ukdx/01zd1Cci3dL1lTns6VcP8slfrGXmKUNZcsObVTREJCkqHDlq3c5DfOy+1UwZVcr9N85jmB7nKiJJUuHIQS/tOcz1S55j9JBifvGxeYzSmIaIhKDCkWO27m/gw3c/x6DCfH75sXmUDSvJdEoiMsCocOSQfYeP8aGfP0u7O7/42Dwmj9IltyISngpHjmhpa+eff/k8dY2t3H/jXE4fNyTTKYnIAKXLaHLEHX/exJrtdfxw4XnMmTg80+mIyACmI44c8Pvna7jv6e187OLp/MM5EzKdjogMcCocWa5qVz23/X4DF546iluvPjPT6YhIFlDhyGKHGlv45C/WMrK0iLs+cD4F+frfLSK9pzGOLBVpdz6zdB37Djfz609cyJghxZlOSUSyhApHlvrPv77Cylf28833vInzpozMdDoikkV07iILPbF5Hz/8n2reXz6ZhXOnZDodEckyKSscZrbEzPaZWVVM7CtmtsvM1gWvd8Qsu83Mqs1ss5ldFROfH8SqzezWVOWbLeqbWrn1d+uZWTaUry6Ynel0RCQLpfKI415gfpz499393OC1AsDMZgHXAbODdX5sZvlmlg/8CLgamAUsDNpKAt9c8RL7jzTznWvPpqQwP9PpiEgWSlnhcPeVQG2SzRcAS9292d1fA6qBucGr2t23unsLsDRoK3Gs2rKfpat3sujS0zh7kh77KiKpkYkxjlvMbH1wKqtj1HYisDOmTU0QSxSXThqa27j1dxs4dexgPnfFjEynIyJZLN1XVf0E+DrgwdfvAjcCFqetE7+webwNm9kiYBFAWVkZlZWVPU6yoaGhV+tnwn9vamb3oTb+fV4Jz/xtVVr3PRD7K5PUX+Gov8JJR3+ltXC4++sd783sZ8Cfgo81wOSYppOA3cH7RPHO214MLAYoLy/3ioqKHudZWVlJb9ZPt2e3HuTxR57hhrdO4+P/kP4B8YHWX5mm/gpH/RVOOvorraeqzGx8zMd/BDquuFoOXGdmxWY2HZgBPAesBmaY2XQzKyI6gL48nTn3d00tEb74u/VMGVXKv101M9PpiEgOSNkRh5k9AFQAY8ysBrgdqDCzc4mebtoGfALA3Tea2TJgE9AG3OzukWA7twCPAvnAEnffmKqcB6Lv/mUz2w428quPz9Mzw0UkLVL2m8bdF8YJ391F+zuAO+LEVwAr+jC1rPH8jjru/ttrfHDeFC46bUym0xGRHKE7xweo9nbn9oc2csqwEs16KyJppcIxQP1x/W427Krn366aydCSwkynIyI5RIVjAGpui/DtRzZz1vhhvPtc3dYiIumlwjEA/ffT29l1qIl/f8eZ5OXFuwVGRCR1uh0cN7PRwEXABKCJ6CW0L7h73BvxJLXqG1v54f9Uc+kZY7lkxthMpyMiOShh4TCzS4DbgFOAdcA+oITovRRTzWwp0QkLG9KRqET9uLKaw8dauXW+BsRFJDO6OuJ4D3CLu2/tvCC4Ge8aojPZ/jZFuUknNXWN3PP3bbznvEnMmjAs0+mISI7qqnB82d2PJFg2x91VMNLsu395BQP+9cozMp2KiOSwrgbH/2pmwzsHzexy4KHUpSTxVO2q58EXdnHjxdOZMGJQptMRkRzWVeG4B3giGBwHwMzeR/Tu72tSnZi8wd355sMvMbK0kE9VnJbpdEQkxyU8VeXuPzWzZqLF40rgvcBngLfFG/eQ1Fm55QB/qz7Il981i2G62U9EMqzLy3Hd/R4zOwY8T3Q687e6+/60ZCZA9Gjj24+8zJRRpXzowqmZTkdEpMvLcV8gOoutAYOBYcCjZmaAu/v56Ukxtz1VfYCNuw/z7X86m6IC3a8pIpnX1RHHe9OWhST0X09uZdzQYhacNyHTqYiIAF2PcbyazkTkZFW76nmq+gBfmD+T4oL8TKcjIgJ0cVWVmT1hZp8yswmd4gVmdqmZ3W1mN6Q+xdz1s1VbGVyUzwfnaWxDRPqPrk5VvRP4GPCgmU0EaolOOVICPA78yN3XpD7F3FRT18if1u/hhoumMXyQrqQSkf6jq1NVjcCdwJ1mVgyMA5rc/UC6kstlS57ahgE3Xjw906mIiJwgqUfHunszsDPFuUigvrGVpat38A/nTNBd4iLS7+j6zn7oF89up7ElwscvOTXTqYiInCRlhcPMlpjZPjOrirPs82bmZjYm+GxmdqeZVZvZejM7P6bt9Wa2JXhdn6p8+4tjrRHu+ds2Lj1jrGbAFZF+KZVHHPcSnXb9BGY2GXg7sCMmfDUwI3gtAn4StB0F3A7MA+YCt5vZyBTmnHF/eGEXBxqa+cSlOtoQkf6p28JhZm82s2fMrN7MjplZs5kd7m49d19J9Eqszr4PfIHoXekdFgD3e9QzwAgzGw9cBTzm7rXuXgc8RpxilC3a253Fq7Yye8IwLjptdPcriIhkQDKD4z8GPgQsJfpX/0eByT3ZmZldA+xy9xejM5ccN5ETB99rgliieLxtLyJ6tEJZWRmVlZU9SRGAhoaGXq3fU8+/3sbW/c188uxinnzyybTvv6cy1V8DlforHPVXOOnor2QKR567bzazAndvBX5mZn8HvhxmR2ZWCnwJuDLe4jgx7yJ+ctB9MbAYoLy83CsqKsKkd4LKykp6s35P3fWTvzNxRB6ff38FBfkD57qFTPXXQKX+Ckf9FU46+iuZ305Hg0fFvmhm3zCzTwNDerCv04DpwXa2AZOA583sFKJHErFHMZOIzsabKJ51Nu0+zJrtddx48fQBVTREJPck8xvqo0G7W4AI0QHs0BMguvsGdx/n7tPcfRrRonC+u+8FlgMfCa6uuhCod/c9wKPAlWY2MhgUvzKIZZ3frN1JUX4e7zkv7pk4EZF+o8tTVWaWD9zu7tcDx4D/k+yGzewBoAIYY2Y1wXbuTtB8BfAOoBpoBG4AcPdaM/s6sDpo9zV3jzfgPqA1t0V48IVdvH12GSMHF2U6HRGRLnX3IKeImY03s8JgfCNp7r6wm+XTYt47cHOCdkuAJWH2PdD8ddM+DjW28r7yHl1zICKSVskMjm8FVpnZQ8DRjqC735myrHLMsjU7mTC8hItPH5PpVEREupVM4dhP9P6J0uAlfWj3oSZWbtnPpy87nfy8eBeRiYj0L90WDndPelxDwvvd2hrc4b0X6DSViAwM3RYOM3uMOPdOuHu8+zEkhPZ25zdra3jLqaOZMloHcyIyMCRzquo/Yt6XAP8ENKcmndzy7Gu17Kht5H+/fUamUxERSVoyp6qe7RR60swGznwY/dhv1uxkaHEB82ePz3QqIiJJS+ZUVezc3nnABYB+0/XS4WOtrKjawz+dP4lBRfmZTkdEJGnJnKrayBvzRrUBrwEfT2VSueCPL+7mWGu77t0QkQEnmcJxaueb/8wsqUfOSmLL1tQws2woZ08anulURERCSWauqs5jHADP9XUiuWTz3iO8uPMQ15ZPotP08iIi/V7CIwczG0d0LGOQmb2JN6Y4H4ZuBOyVZWt2Uphv/KMmNBSRAairU07vBG4kOpX5j2PiRwgx2aGcqKWtnQdf2MUVZ5UxekhxptMREQktYeFw93uAe8zsfe6+LI05ZbWVr+yn9mgL771gUqZTERHpkWTu41hmZlcBs4neANgR/0YqE8tWD1ftZWhJAZfMGJvpVEREeiSZ+zh+DIwALgXuIXrn+DMpzisrtUba+etLr3PFWWUUFegpfyIyMCXz2+tid/8AcDCY8HAe0XEPCemZrQepb2pl/pxTMp2KiEiPJVM4jnV8DZ4PfgyYlrKMstgjVXsZVJjPpTpNJSIDWDI38q0wsxHA/wPWEX3u+H0pzSoLRdqdRze+zmVnjtUUIyIyoHX3zPE84GF3PwT8xsz+BAzKxud+p9rzO+o40NDMVbN1mkpEBrYuT1W5ezvwg5jPTckWDTNbYmb7zKwqJvZ1M1tvZuvM7C9mNiGIm5ndaWbVwfLzY9a53sy2BK/rQ/8L+4mHN+ylKD+Pt505LtOpiIj0SjJjHI+Z2YIebPteYH6n2Hfc/Wx3Pxf4E/DlIH41MCN4LQJ+AmBmo4DbiQ7IzwVuN7ORPcglo9ydRzfu5eIZYxhaUpjpdEREeiWZwnEL8KCZNZlZrZnVmVm3Rx3uvhKo7RQ7HPNxMG88WXABcL9HPQOMMLPxwFXAY+5e6+51RJ993rkY9XsbdtWz61CTrqYSkayQzOD4mL7coZndAXwEqAcuC8ITgZ0xzWqCWKL4gPJI1V7y84y3n1WW6VRERHotmTvHI2Z2HdHp1b9hZpOAMmBtT3bo7l8CvmRmtxE9mrmdNyZQPKFpF/GTmNkioqe5KCsro7KysifpAdDQ0NCr9WO5O79/romZI40XV/+9T7bZ3/Rlf+UC9Vc46q9w0tFfydw5fhdQSPTO8W8AjcBPgTf3ct+/Av5MtHDUALFPNJoE7A7iFZ3ilfE25u6LgcUA5eXlXlFREa9ZUiorK+nN+rFeef0Iex9dyc1Xzqbiwql9ss3+pi/7Kxeov8JRf4WTjv5KZozjInf/BMGNgMFVVUU92ZmZzYj5eA3wcvB+OfCR4OqqC4F6d98DPApcaWYjg0HxK4PYgPHwhr2YwVWzdJpKRLJDMmMcrcH9HA5gZqOB9u5WMrMHiB4tjDGzGqJHFu8ws5nB+tuBTwbNVwDvAKqJHtHcANEiZWZfB1YH7b420O4heWTjXi6YMpJxw0q6bywiMgAkUzh+BPwOGGtmXwXeB3y1u5XcfWGc8N0J2jpwc4JlS4AlSeTZ72w/eJSX9hzmP955VqZTERHpM8kMjt9vZmuBK4LQte5e1dU6EvVI1V4A3S0uIlklmSMOgHyglejpKs0HnqSHq/YyZ+IwJo/Sk3ZFJHt0WwTM7EvAA8AEolc1/Sq4lFa6sKe+iXU7D3H1nPGZTkVEpE8lc8TxIeACd2+E4zfwrQW+mcrEBrrHNr0OwFWzdTWViGSXZE47befEAlMAbE1NOtlj5Sv7mTKqlNPHDc10KiIifSqZI45GYKOZPUp0jONK4Ckz+x6Au/9LCvMbkFoj7Tz96kHefd6Amx1FRKRbyRSOPwevDnreeDde2HGIoy0RLpnRp9N8iYj0C8lcjhv33gtJbNWW/eQZvOU0FQ4RyT7JXFU138xWBw9lSnpa9Vy2assBzpk8guGD9OwNEck+yQyO3wV8guh05mOJTrM+NpVJDWSHGltYX3OIS2aoi0QkOyUzxlEDrAseIyvd+PurB2l3uFTjGyKSpZIpHF8A/mhmlUBzR9Dd70xVUgPZqi0HGFJcwDmTR2Q6FRGRlEimcHyV6HQjI0hiVtxc5u6s2rKft5w2msJ8zcwiItkpmcIxzt0vSHkmWWD7wUZq6ppYdOmpmU5FRCRlkvmz+HEze1vKM8kCq7bsB9DAuIhktWQKx8eBv5pZgy7H7drKLQeYNHIQ00ZrNlwRyV7JnKrS5UFJaI2088yrB3nXOeMxs0ynIyKSMt0ecbh7BLgW+GLwfjxwbqoTG2he3HmII81tOk0lIlkvmTvH7wIuAz4chBqBn6YyqYFo1ZYDmMFFp43OdCoiIimVzKmqi9z9fDN7AcDda82sKMV5DTirtuzn7EkjGFGqrhGR7JbM4HirmeURnVIdMxtNEvdzmNmSYH6rqpjYd8zsZTNbb2YPmtmImGW3mVm1mW02s6ti4vODWLWZ3RrqX5cm9U2tvFhTzyWnazhIRLJfwsJhZh1HIz8CfgeMNbOvAk8B/zeJbd8LzO8UewyY4+5nA68AtwX7mgVcB8wO1vmxmeWbWX6w/6uBWcDCoG2/8vSrB4m0u6ZRF5Gc0NWpqueA8939fjNbC1wBGHCtu1d1sR4A7r7SzKZ1iv0l5uMzwHuD9wuApe7eDLxmZtXA3GBZtbtvBTCzpUHbTd3tP51WbdnP4KJ8zpsyMtOpiIikXFeF4/g1pe6+EdjYx/u+Efh18H4iJz4gqiaIAezsFJ8Xb2NmtghYBFBWVkZlZWWPE2toaAi1/mMbGjl9eB5/f2plj/c5kIXtr1yn/gpH/RVOOvqrq8Ix1swSPhbW3b/X052a2ZeANuCXHaF4uyD+qTRPkM9iYDFAeXm5V1RU9DQ9KisrSXb9HQcb2ffIE/zz5TOpeOv0Hu9zIAvTX6L+Ckv9FU46+qurwpEPDCH+L/UeM7PrgXcBl7t7RxGoASbHNJsE7A7eJ4r3C6uqo9OMXKz7N0QkR3RVOPa4+9f6cmdmNh/4IvC/3L0xZtFy4Fdm9j1gAjCD6BiLATPMbDqwi+gA+gf6Mqfe+lv1AcYPL+G0sYMznYqISFokNcbRE2b2AFABjDGzGuB2oldRFQOPBdNyPOPun3T3jWa2jOigdxtwc3CXOmZ2C/Ao0SOgJcF4S7/g7qzeVsdFp43WNCMikjO6KhyX92bD7r4wTvjuLtrfAdwRJ74CWNGbXFJlZ20T+480Uz5VV1OJSO5IeB+Hu2sG3G6s2R7togumjspwJiIi6aPH1PXCmu11DC0uYOYpQzOdiohI2qhw9MLabXWcO2UE+Xka3xCR3KHC0UP1Ta28su8I5TpNJSI5RoWjh57fUYc7lE/TwLiI5BYVjh5au62O/Dzj3Mkjum8sIpJFVDh6aM32Ws4aP5TBxck80kREJHuocPRAa6SddTsPaXxDRHKSCkcPbNp9mGOt7RrfEJGcpMLRA2u21wHoiENEcpIKRw+s3V7LxBGDOGV4SaZTERFJOxWOkNydNdvqdJpKRHKWCkdINXVN7NPEhiKSw1Q4QtLEhiKS61Q4QlqzTRMbikhuU+EIae12TWwoIrlNhSOE+qZWNr+uiQ1FJLepcITwgiY2FBFR4Qhj7XZNbCgikrLCYWZLzGyfmVXFxK41s41m1m5m5Z3a32Zm1Wa22cyuionPD2LVZnZrqvJNxpptdZrYUERyXiqPOO4F5neKVQHvAVbGBs1sFnAdMDtY58dmlm9m+cCPgKuBWcDCoG3aaWJDEZGolP3p7O4rzWxap9hLAGYnXZG0AFjq7s3Aa2ZWDcwNllW7+9ZgvaVB202pyjuRl/Ycpqk1wgW68U9Eclx/GeOYCOyM+VwTxBLF027NtmBiQw2Mi0iO6y8n6+PdFOHEL2wedwNmi4BFAGVlZVRWVvY4mYaGhpPWf3jdMUaXGJtfeJbNPd5ydorXX5KY+isc9Vc46eiv/lLoKrafAAAKRklEQVQ4aoDJMZ8nAbuD94niJ3D3xcBigPLycq+oqOhxMpWVlXRe/7anH+eimaOoqDivx9vNVvH6SxJTf4Wj/gonHf3VX05VLQeuM7NiM5sOzACeA1YDM8xsupkVER1AX57u5PYfaWZP/THOmTQ83bsWEel3UnbEYWYPABXAGDOrAW4HaoEfAmOBP5vZOne/yt03mtkyooPebcDN7h4JtnML8CiQDyxx942pyjmRjbvrAZg9QYVDRCSVV1UtTLDowQTt7wDuiBNfAazow9RCq9oVFI6JwzKZhohIv9BfTlX1a1W7DjNtdCnDSgoznYqISMapcCRhw6565kzUaSoREVDh6Fbd0RZ2HWpS4RARCahwdGPj7sMAzNHAuIgIoMLRrQ3BwPgcDYyLiAAqHN2q2l3PpJGDGFFalOlURET6BRWObmzcVc+bNL4hInKcCkcXDh9rZdvBRg2Mi4jEUOHowsZd0YHx2RM0viEi0kGFowsdU43oiENE5A0qHF3YsKue8cNLGDOkONOpiIj0GyocXajSHeMiIidR4UigobmNrQeO6sY/EZFOVDgSeGnPYdx145+ISGcqHAl0TKWuezhERE6kwpHAhl31jB1azLhhJZlORUSkX1HhSGDjrsM62hARiUOFI47miLNl3xHm6MY/EZGTqHDEsfNIO+0Os3XEISJyEhWOOLYfbgc0MC4iEk/KCoeZLTGzfWZWFRMbZWaPmdmW4OvIIG5mdqeZVZvZejM7P2ad64P2W8zs+lTlG2v74XZGDS5i/HANjIuIdJbKI457gfmdYrcCj7v7DODx4DPA1cCM4LUI+AlECw1wOzAPmAvc3lFsUmlbfTtzJg7HzFK9KxGRASdlhcPdVwK1ncILgPuC9/cB746J3+9RzwAjzGw8cBXwmLvXunsd8BgnF6M+1dwWYVdDuwbGRUQSKEjz/srcfQ+Au+8xs3FBfCKwM6ZdTRBLFD+JmS0ierRCWVkZlZWVPUrwtfoIEQc7VENl5d4ebSPXNDQ09Li/c5H6Kxz1Vzjp6K90F45E4p0T8i7iJwfdFwOLAcrLy72ioqJHiex+dgewgeuuvIjJo0p7tI1cU1lZSU/7Oxepv8JRf4WTjv5K91VVrwenoAi+7gviNcDkmHaTgN1dxFOmanc9gwth0shBqdyNiMiAle7CsRzouDLqeuChmPhHgqurLgTqg1NajwJXmtnIYFD8yiCWMlW76pk6LE8D4yIiCaTsVJWZPQBUAGPMrIbo1VHfApaZ2U3ADuDaoPkK4B1ANdAI3ADg7rVm9nVgddDua+7eecC9z7RG2nl5zxEun5Kfql2IiAx4KSsc7r4wwaLL47R14OYE21kCLOnD1BKqPdrCWROGcdrwpnTsTkRkQNKd4zHKhpXw0M1vpfyU/nLNgIhI/6PCISIioahwiIhIKCocIiISigqHiIiEosIhIiKhqHCIiEgoKhwiIhKKCoeIiIRi0Zu2s4uZ7Qe2x4SGA/UhPo8BDqQovc776qt1umuTaHm8eC70V3ft1F/h2vWmvzrH1F/hY7Gfe9NfU919bLet3D3rX8DikJ/XpCuXvlqnuzaJlseL50J/dddO/ZW+/uocU3/17nsulf3V8cqVU1V/DPk5lXqyr2TW6a5NouXx4rnQX921U3+Fa9eb/uocU3+Fj6Wzz7LzVFVvmdkady/PdB4DhforHPVXOOqvcNLRX7lyxBHW4kwnMMCov8JRf4Wj/gon5f2lIw4REQlFRxwiIhKKCoeIiISiwiEiIqGocIRkZoPNbK2ZvSvTufR3ZnaWmf3UzH5rZp/KdD4DgZm928x+ZmYPmdmVmc6nvzOzU83sbjP7baZz6a+C31n3Bd9XH+yLbeZM4TCzJWa2z8yqOsXnm9lmM6s2s1uT2NQXgWWpybL/6Iv+cveX3P2TwPuArL+cso/67A/u/nHgo8D7U5huxvVRf21195tSm2n/E7Lv3gP8Nvi+uqYv9p8zhQO4F5gfGzCzfOBHwNXALGChmc0yszeZ2Z86vcaZ2RXAJuD1dCefAffSy/4K1rkGeAp4PL3pZ8S99EGfBf4jWC+b3Uvf9VeuuZck+w6YBOwMmkX6YucFfbGRgcDdV5rZtE7huUC1u28FMLOlwAJ3/yZw0qkoM7sMGEz0f0qTma1w9/aUJp4hfdFfwXaWA8vN7M/Ar1KXceb10feYAd8CHnb351ObcWb11fdYLgrTd0AN0eKxjj46WMiZwpHARN6oxBDt4HmJGrv7lwDM7KPAgWwtGl0I1V9mVkH0MLkYWJHSzPqvUH0GfBq4AhhuZqe7+09TmVw/FPZ7bDRwB3Cemd0WFJhclajv7gTuMrN30kdTk+R64bA4sW7viHT3e/s+lQEhVH+5eyVQmapkBoiwfXYn0R/0XBW2vw4Cn0xdOgNK3L5z96PADX25o1wa44inBpgc83kSsDtDuQwE6q/w1GfhqL96Lm19l+uFYzUww8ymm1kRcB2wPMM59Wfqr/DUZ+Gov3oubX2XM4XDzB4AngZmmlmNmd3k7m3ALcCjwEvAMnffmMk8+wv1V3jqs3DUXz2X6b7TJIciIhJKzhxxiIhI31DhEBGRUFQ4REQkFBUOEREJRYVDRERCUeEQEZFQcn3KEREAzGwVMBRoc/esnwJepDd0H4dIHzKzs4DPAmOAx939J2nY56nAl4Dh7v7eVO9PRIVDJAXMLA/4WTofMmRmv1XhkHTQGIfkNDM73cw2dIoVm9lrZjbbzJ4MHpATZpt9/vAqM2sIvib1QCMzKzKzlWam09HS53TEITktODKoBUZ1PF/FzD4NXEB00rgCd/9BD7f9Z3d/Zx/l2eDuQ7ppc8IRh5ndTvTBPr/sixxEOuiIQ3JaUCx2ANMAzGwQ8K/AV4APAg8F8Wlm9rKZ/dzMqszsl2Z2hZn9zcy2mNncoF2Fmd1pZv9F8PCqZNcN2v5L0KbKzD6XzL/BzEab2U8JHmYUs+gPwb9BpE/piENynpn9GrjP3VeY2eeBKcDngR3ufkrQZhpQDZwHbCR6NPIicBNwDXCDu787wfaTWtfMLiD6LOkLiT6U51ngQ+7+QjJHHHH2mw/sdfexYdYT6Y6OOESiU1DPNLMhwM1EH0U6BjjUqd1r7r4hOErZSPSqKQc2EByxdCGZdS8GHnT3o+7eAPweuKSn/yh3jwAtZja0p9sQiUeFQyRaOM4kehntL939daAJKOnUrjnmfXvM53a6vycqmXXjPfqzt4qBYynYruQwFQ6RaOGYC9wIfAfA3euAfDPrXDxSaSXwbjMrNbPBwD8Cq3q6MTMbDex399a+SlAEVDhEADYDbwIWu3t9TPwvRE8fpYW7P090jOM5ouMbP3f3F3qxycsIBuhF+pIGx0USMLPzgH9x9w9nOpeeMLPfA7e5++ZM5yLZRUccIgkEf+0/EfYGwP7AzIqAP6hoSCroiENERELREYeIiISiwiEiIqGocIiISCgqHCIiEooKh4iIhKLCISIioahwiIhIKCocIiISyv8H0ZxDvswXSy8AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "V_grid = np.logspace(-4,0, num=50) #-4 e zero sao os exponentes do volume\n",
    "T_grid = np.zeros(len(V_grid))\n",
    "\n",
    "for j in range (len(V_grid)):\n",
    "    T_grid[j]=T_calc(V_grid[j])\n",
    "\n",
    "plt.semilogx(V_grid,T_grid)\n",
    "\n",
    "plt.ylabel('Temperatura (K)') \n",
    "plt.xlabel(r'$\\bar{V} (\\mathrm{m^{3}mol^{-1}})$')\n",
    "plt.grid(True)\n",
    "               \n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 983.36670922, 1094.13355026, 1188.71501482, 1268.91369617,\n",
       "       1336.56157268, 1393.39882219, 1441.01117914, 1480.80537253,\n",
       "       1514.00699987, 1541.67053658, 1564.69514603, 1583.8425696 ,\n",
       "       1599.75501061, 1612.97191798, 1623.94516259, 1633.05243759,\n",
       "       1640.60890196, 1646.87718382, 1652.07590456, 1656.38689864,\n",
       "       1659.96130056, 1662.92465892, 1665.38122138, 1667.41751699,\n",
       "       1669.10534661, 1670.50427533, 1671.66370746, 1672.62461323,\n",
       "       1673.42096506, 1674.08092669, 1674.6278485 , 1675.08108431,\n",
       "       1675.45667737, 1675.7679262 , 1676.02585044, 1676.23958441,\n",
       "       1676.41669852, 1676.56346285, 1676.685084  , 1676.78587278,\n",
       "       1676.86937214, 1676.93857745, 1676.9959371 , 1677.04344317,\n",
       "       1677.08279453, 1677.11537839, 1677.14234966, 1677.16517715,\n",
       "       1677.18349508, 1677.19866584])"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "T_grid"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "def pressure_total(T,V):\n",
    "    a,b,c,Tc,Pc,association_type, eps, beta = load_param()\n",
    "    \n",
    "    P_cubic=8.31447*T/(V - b) - a*(c*(-sqrt(T/Tc) + 1) + 1)**2/(V*(V + b))\n",
    "    \n",
    "    P_association=-4.157235*T*(1.0 + 0.475*b/(V*(1.0 - 0.475*b/V)))*(-0.5*V*(1.0 - 0.475*b/V)*(sqrt(1.0 + 4.0*b*beta*(exp(0.120272248261164*eps/T) - 1.0)/(V*(1.0 - 0.475*b/V))) - 1.0)/(b*beta*(exp(0.120272248261164*eps/T) - 1.0)) + 1)/V\n",
    "    \n",
    "    P = P_cubic-P_association\n",
    "    return P"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY4AAAEaCAYAAAAG87ApAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzt3Xl4VeW5///3JyOEMAQCgQARkOEIDkgpWOuQKlVsHTrp0WpFa4tttf46nVY7aLXV46+nw1Gr9VBLra2KqG3FltY6RZxwFgsqEpEhgAwyhiGQ5P7+sVZgsdlJdjbZe2cn9+u61pW9nudZa9159nCvecnMcM455xKVk+kAnHPOZRdPHM4559rEE4dzzrk28cThnHOuTTxxOOecaxNPHM4559rEE4dLG0nLJE1JYrrjJS1ORUwHQ9LVklZIOkLSk5mOx7l08cTRDEnDJJmkPEn/kFQbDnsk7Y6M397G+Z4s6W1JOyQ9KemQg4xzkqS5kjZL2ijpRUkXH8w8Oxoze9rMxmQ6jjiOBE4Cfgk8neFY0kLSJeHnd5uktZL+LqlnWHenpJ/GtN/7PQrH/0vSwnD69yT9Vyb+j9ZIOlXSvDDO9ZKeknRmWHeRpIbw+79V0uuSTo+ZfrikRkm3tbKcpv6pjRn+M9LmWElPhLFskfSwpLEx8/l+2J+1kmok3dee/RHLE0cCzOw0Mys2s2LgbuBnTeNm9pVE5yOpFPgz8COgL/AykPQbLOkjwBPAU8BIoB/wVeC0JOeXm2wsnY0CLX4/zOxzZlZtZh83s6vTFVtLUvkeSjoRuAE4z8x6AocBs9s6G+BCoASYClwu6dx2DfQgSfoccD9wFzAEKAOuBs6INHs+/D3oA/wOmC2pb6T+QmATcK6kwgQW2yfym1JsZveFsXwE+BfwEFAODAcWAM9KGhG2mQZ8AZgSxjQReDy5/z5BZtZlB8CAkZHxO4Gfhq+HhfV5MdPsbZPE8qYDz0XGewA7gf9Icn7PALe2UH8R8Exz/3P4v/wGmAtsB34IvA/kRtp/GngjfD0JeB7YDKwBfg0UtLD8LwDLgQ+AHwDLCD7cEKy0XAm8G9bPBvo2M59KoCYyvgz4DvAGsIUg+XYL694CTo+0zQM2ABPC8WOA58L/YQFQGWlbBVwPPBu+LyPDPlwKbAPeA84P2x5KkLQ/COd/N8GXv2leh4Xz2wwsAs5soZ+qgP8GXgz/n4eifUHwI/Z+WDcPGBfzeYy+h1OATwKvAVuBlcCPW/kcfRmoBjYCc4DyZtp9B/hrC/O5k5jvBs18jyL1NwO3tDDP4yLv10rgorC8N8EP+/rwM/ZDICfyuX8W+FU43VLg2LB8JbAOmNbM8gSsAP4r0e8VwffYgImRsncJVuLWAp9rYV6t9c/TwG1xyv8B3BW+/jXwv8n8hiQ7+BZHO5BUEe4qam74fNh0HMGPFQBmtp3gAzYuiWUWAR8BHjjI8D9P8GPZE/g5wY/PSTH194SvG4BvAqXhsk8GvtZMfGMJftC+QLCm1I9g7a3JFcCngBPD+k3ArW2I+xyCNdbhBLuMLgrL7wXOi7Q7FdhgZq9KGgz8HfgpwRbfd4AHJfWPtP8CQYLvSfCjdDNwmgVr2McCrzf9iwQ/9uUESWIo8OPwf88HHiZYUxwAfB24W1JLu9suBL4Yzq8+XG6TfwCjwnm9SpCkoqLv4TME7+GFBGvDnwS+KulT8RYq6aTw/zgHGETwIzyrmRhfAE6VdK2kjya4Jt0sSQKOJ0is8eorCP73W4D+wHj29f8tBMljBMFn6EIguot2MsGKRT+Cz+8s4MMEKwMXAL+WVBxnsWMI3suEvlfhLrgvAbXAkrDseILP+iyCFaILE5lXnHkXEXzm7o9TPRv4ePh6PnBhuBtwYlr2HKQzS3W0gfRvcfwOuDGm7FnCtag2zmtwGF+zWysktsVxV0z9T4GZ4eueBD9ChzQz/28Af2mm7mpgVmS8B7CbfVscbwEnR+oHAXti+zusq+TALY4LIuM/A24PX48k2DooCsfvBq4OX38P+GPMvB8hXPskWPO/LibmzcBnge6tvB+fAl4LXx9PsIWQE6m/l2bW/MPl3hgZHxv2VW6ctn3C97B3c+9hnGn+F/hVC5/Jn0XGi8P3YVgz7U8jSIqbCX4sf9kUZxjLrrCuadhKM2vUwLUEK1KFzSzrqnifLyAXqAPGRsouBaoin/slkbojwhjKImUfAOPjzPujYdturXyv6sP/bwPBD/eUSP0dhFtmBCtYe4ABzcxrWLi8zTHDYQTJJ+53nGClaU9k/HzgMYLv6wfAlS19Jg528C2O9KoFesWU9SL4oduPpPMjB8r+EWdem4BGgh/cg7EyZvwe4DPh2uRngFfNbHkY02hJf5P0vqStBPu7S5uZb3l03hZsXX0QqT8E+EvTVhlBImkg2J+ciPcjr3cQ/OBhZtXhvM4I19jOZN8W0yHA2dGtQYJdIdE+jI35P4GvAGvCA8H/EfbFAEmzJK0K++JPkb4oB1aaWWNkvssJkn1zou/DciAfKJWUK+lGSe+Gy1kWtiltZlokTVZw4sV6SVvC+Ft6n5ZH/udagvcpbqxm9g8zO4Ngi+0sgh/RL0Wa/NzM+jQNBFuDB5B0OcGa+CfNrK6Z2IYSbJHHKgUKonFzYP+ujbzeGcYeWxZvi6PpM9ra92p++D+WmtkxZvYYgKTuwNmEW4Vm9jzBrq/PNz8rAEqj/WZmb9Hyd3wQQdIiXM7dZjaFYMXiK8B1kk5tZZlJ88QB3SKv+yQzg3BXVexZEdHh/LDpIuCoyHQ9CPaVH7CpHn4Qmg6UHXCw28x2EBxv+GwLoW0HiiLLGxinzX63RzazNwm+hKex/24qCHY9vQ2MMrNewPcJdtnEs4bgi9+07CKC3QZNVhLsAop+WbqZ2aoW/p9ENe2uOgt4M0wmTcv8Y8wye5jZjZFpY/vjETP7OMEX9W3gt2HVf4dtjwz74gL29cVqYGjMwfUKoKX/bWjkdQXBWuoGgvfgLIJjF70J1lBh/36PvcX1PQTHKoaaWW/gdpp/n1YTJNRgpsFnsl8rsWJmjWb2OMFxnsNbahtL0hcJjm+dbGY1LTRdSfD9iLWBoH+iZyS21r+JWhwut6XvVUs+TbAyeFu4gvU+QUJr8+6qcMXleYJEFOsc4hwAN7M9ZnY/wW66Nr0vbeGJAy4O1+rGE+yz7xnuo06Yma2w/c+IiB2a9kn/BThc0mcldSPYnfOGmb2dZOzfBS4K9232A5B0lKSmfdQLgHGSxofL+3GC872H4BjECey/f7Unwa6H2nDN+6stzOMB4HRJx0kqAK5j/8/b7cD1Ck9HltRf0lkJxteaWcApYXzRxPcngi2RU8P3vJukSklD4s1EUpmkM8Mf0zqCLcaGsLpnOL45PHYSPa30BYKk/V1J+ZIqCc7Iae7YAcAFksaGCfY64AEzawiXU0ewJlxEsJXXmp7ARjPbJWkSLa/t3kPwHRgfbmXeALxgZstiG0o6S9K5kkoUmERwfGF+AjE1zeP8cBkfN7OlrTS/G5gi6RwFp8X3kzQ+7JfZBJ+fnuFn6FsE7+9BsWC/z7eAH0m6WFIvSTnh53hGArOYBswk2D02Phw+CoyXdEQSIV0JTJN0Rfi/lig45fkjBLv6mk4P/mRYnyPpNILjpi8ksbyEeOIIvoxrCNYkryZ446MHh2PX5pJmZusJ1mSuJ9gMnQwkfSqimT1HEOtJwFJJG4EZBGfYYGbvEPwIPUZw4O6ZBGd9L8FxhSfMbEOk/DsEP0LbCPqr2VOJzWwRcBnBD9Magv83unZ5E8Fa8b8kbSP48ZmcYHwtMrM1BGtqx0ZjNLOVBGvv3yc48L2S4Ae/ue9BDvBtgrXyjQQ/kk0nA1wLTCA40+nvBKdZNy1nN8EustMI1o5vAy5sZQXhjwTHCN4n2Aq+Iiy/i2ALcBXwJon9SH+NYFfFNoLPdLOnzIZbDT8CHiR4nw6l+c/kJoIzsJYQrED8CfifyIpRIn5KsEXzklq5FsrMVgCfIHgPNhIcGG/aYv86QXJeSvC5vofgB/ugmdkDBLsov0jw3q8N436openCFYiTCc5wej8yvAL8k+C3pTmbY/ZSfCuM5RmCEzw+Q/D+LAeOBo4zsyXhtFsJPtMrCI6P/Az4ajhtSig8sNIlSTKC3S7VceqOBOaF+2mdSxlJVcCfzOyOTMfiXCJ8iyOOcN/0OQQX6DnnnIvIy3QAHdQKgl0GnerWHc451x669K4q55xzbee7qpxzzrWJJw7nnHNtkrJjHJJmAqcD68zs8LDsPoJ7wUBwsd1mMxsf1l0FXEJwnvwVZvZIWD6V4NTNXOCOmIu14iotLbVhw4btHd++fTs9evRop/+sa/G+S573XfK875J3MH33yiuvbDCz/q02TNW9TAguHpsALGym/hfsu4fQWMJ71hDctO5dgkSRG74eQXCLgQVE7k/T3PChD33Iop588klzyfG+S573XfK875J3MH0HvGwJ/L6nbIvDzOZJGhavTpIITndtutDuLIIb4tUB70mqJriFN0C1hVeYhldEn0VwIZRzzrkMyNTpuMcDa23flY+D2f+K2Br23bBsZUx53KuLJU0nuB02ZWVlVFVV7a2rra3db9wlzvsued53yfO+S146+i5TieM8gttaNIl3AzYj/sH7uOcPm9kMgtttMHHiRKusrNxbV1VVRXTcJc77Lnned8nzvkteOvou7YlDwYNPPgN8KFJcw/53Bx1CcI8YWih3zjmXAZk4HXcK8LbtfzvlOYTP5pU0nOBpZy8CLwGjFDz4vYDg5mtz0h6xc865vVKWOCTdS3CH0jGSaiRdElady/67qZrupDqb4KD3P4HLzKzBzOqBywme0vYWMDts65xzLkNSeVbVec2UX9RM+fUEtxuPLZ9LeJvwVKurb+Cpxev3KwtOAIuMNzNttNl+r2OnUNyX+y1He8sOnE9TmWIaC+1XJ8WOh61i24R1+00fGV+5rZHF728Ly/bNV0BO+DonnDAnR+SE401tm8ZzJJTT9Hr/aXOb5hnT1865jslvchhRu6ue6X98JdNhdDzPzkvLYiTIDZNMTk74Okfk5mjfawXjebn7XjcNebk55DW9Dv/m5+aEf4PXeTk55OcG0+fn5lCQm0NBXg75uU2DKMjLoTAvJ/ybS0FuDoX5Qdtu+bnhEL7Oy6UwP2jvic91FZ44Inp1z+dvXz+u2frm7gdpkRO9om1im1uk0vYr339uzc2nqaxpPtHyvTFYUL63LRbWx0xnkbpI/b5pg7qFixYxduy4A+bTtMzGRmgMp8OC1422r6yxMbhgaG9Z+Leh6XXjvrqmoSGcZ0NjMDS9bvpb32g0Nv01o75hX3nwt5E9DY3s2B28rm8w9jQ0Ut9o+73eU9/I7oZgONh7fUrQPT+XooJcuhfkUpSfx55dO5mxZD49CvMoLsyjR2EuPQrz6FmYt7esV/d8enXLp1f3vPBvPj0L88jJ8STkOi5PHBH5uTkcPrh3psPoULp/sJjKIwdlOoyUa2gMEkpdfSO7mxJKfSN19Q3h30bq9gTju/Y0smtPA7uir/c0sGN3MOzcXc+O3Q2sXruDuvpGNm7fQW1dPdvr6tle18DuhsYWY5Ggd/d8SooK6FOUT9+iAvoUFVBSlE+/4kJKiwso7VlI/+JCSosL6VdcQH6u33bOpY8nDucg3N0V7IZqL8H59MceUF5X38D2ugZqd9Wzddcetu7cw9bo65172LxzDxu372bzjj28v3UXb7+/jY3bd7NzT0OcJUFJUT4De3dnUO9uDOzdjUG9wr+9uzO0b3fK+3T35OLajScO59KsMC+Xwrxc+vYoaPO0O3bXs2HbbtbX1rGhto7124K/67bVsXbLLtZs2cXrKzezcfvu/abLEQzq3Z2KvkVU9C1iaN/uHNKvByMHFDO8tEe7JkzX+XnicC6LFBXkUdEvj4p+RS2227WngXVb61i9ZScrN+4Ihk07WbFxB08sXsf6bXV720owtKSIkQOKg6F/MWPLezGqrJjCPE8o7kCeOJzrhLrl51LRr4iKfkUcM6LfAfU7dtezbMMOqtfX8u662r1/n6newO764BhMXo4YVdaTceW9OLy8F+MG92ZceS+KCvxno6vzT4BzXVBRQR5jy3sxtrzXfuUNjcbyD7bz5pqtLFodDE++vY4HXglu9JCbIw4v78WHDunLxGElTDykhAG9umXiX3AZ5InDObdXbo4Y0b+YEf2LOf3IciA4/Xrt1joWrtrCays38dKyTdz9wnJmPvseABV9i5g0vC8njO7P8SNLKUni2I3LLp44nHMtksTA8GytKWPLANhd38ii1Vt4edkmXl6+kcfeWssDr9QgwZFD+nDi6P6cOLo/44f2IdevSel0PHE459qsIC+HoytKOLqihC8zgoZG442azTz1znqeemc9v35iCTc/voTe3fOZclgZpx81iONGlvopwZ2EJw7n3EHLzdHeRPKNKaPZvGM3z1Rv4Im31/GvN9/nwVdr6FOUz6ljB3L6UYP4yIh+5HkSyVqeOJxz7a5PUQGnH1nO6UeWU1ffwNPvbOBvb6zmb2+s5r6XV9K3RwFnHlXO+ZMrGFXWM9PhujbyxOGcS6nCvFymjC1jytgydu1poGrxeh5+YzX3vLCCO59bxqThfbngmEOYOm4gBXm+FZINPHE459KmW34uUw8fyNTDB/JBbR33v1LD3S8s54p7X6O0uICzJw7l85MqMh2ma4UnDudcRvQrLuQrJx7K9ONHMG/Jev40fwX/99S7zJi3lGMG5lIxrpYR/YszHaaLwxOHcy6jcnJE5ZgBVI4ZwKrNO5n5zHv88bn3mPLLpzjzqHIuP2kUIwd4AulIPHE45zqMwX2686PTx3Jk/loW1pfxp/kreGjBas44spwrTh7JyAF+IL0j8CNRzrkOp3eh+MEnx/L09z7G9BNG8Nhba/n4r+Zx1Z/fOODOvy79UpY4JM2UtE7Swpjyr0taLGmRpJ9Fyq+SVB3WnRopnxqWVUu6MlXxOuc6ntLiQq467TCe+d5JXHzscGa/XMNJv6jiT/OX09B4kI9tdElL5RbHncDUaIGkjwFnAUea2Tjg52H5WOBcYFw4zW2SciXlArcCpwFjgfPCts65LqRvjwKuPmMsc684njFlPfnhXxfyqVuf5bUVmzIdWpeUssRhZvOAjTHFXwVuNLO6sM26sPwsYJaZ1ZnZe0A1MCkcqs1sqZntBmaFbZ1zXdCYgT2ZNf0Ybjp3PGu37uLTtz3H9x7w3Vfplu6D46OB4yVdD+wCvmNmLwGDgfmRdjVhGcDKmPLJ8WYsaTowHaCsrIyqqqq9dbW1tfuNu8R53yXP+y55rfVdb+Daybk8VJ3PA6+s5J//ruHLRxRyeKk/eCodn7t0J448oAQ4BvgwMFvSCCDe7TON+FtEcXdsmtkMYAbAxIkTrbKycm9d8OznyniTuVZ43yXP+y55ifbdaVPgrTVbueLe1/j5y7VceuIIvv3xMV36CvR0fO7S3bs1wJ8t8CLQCJSG5UMj7YYAq1sod845AA4b1Is5lx/H5ydX8H9PLeXs259j+QfbMx1Wp5buxPFX4CQASaOBAmADMAc4V1KhpOHAKOBF4CVglKThkgoIDqDPSXPMzrkOrntBLjd8+ghuv2AC723Yziduepq/vFaT6bA6rVSejnsv8DwwRlKNpEuAmcCI8BTdWcC0cOtjETAbeBP4J3CZmTWYWT1wOfAI8BYwO2zrnHMHmHr4IP7xjRMYV96bb963gG/PXkBdfUOmw+p0UnaMw8zOa6bqgmbaXw9cH6d8LjC3HUNzznVig/t0554vT+aWJ6q56fEl1GzawYwvTKR3UX6mQ+s0uu4RJOdcp5WXm8M3Pz6am84dz6srNvHZ25+jZtOOTIfVaXjicM51WmeNH8xdX5y895qPhau2ZDqkTsETh3OuU/vIof148KvHUpCbwzn/9zxPLl7X+kSuRZ44nHOd3uiynvz5a8cyvLQHX/rDy8x6cUWmQ8pqnjicc11CWa9u3HfpRzhuZClX/vnf/HH+8kyHlLU8cTjnuoziwjzumDaRk/9jAFc/tJC/vrYq0yFlJU8czrkuJT83h1vPn8Dk4X359v0LePTNtZkOKet44nDOdTnd8nO5Y9qHOby8F5fd8yrPVW/IdEhZxROHc65LKi7M486LJzG8Xw++dNfL/myPNvDE4Zzrskp6FPDHSyZRWlzIRb9/ibff35rpkLKCJw7nXJc2oFc37v7SZLrn53LBHS+ycqNfYd4aTxzOuS5vaN8i/vSlSdTVN/C1u19l1x6/MWJLPHE45xwwckBPfnXOeP69agvXPuw34W6JJw7nnAtNGVvGZR87lHtfXMnsl1e2PkEX5YnDOecivvXxMXx0ZD9+9NeFLFrtN0WMxxOHc85F5OaIm849mpKiAr76p1fZsmNPpkPqcDxxOOdcjNLiQm67YAJrtuzkW7Nfp7HRMh1Sh+KJwznn4phQUcIPPzmWx99ex2+eejfT4XQonjicc64ZF37kEM48qpxf/Gsx85d+kOlwOgxPHM451wxJ3PjZI6joW8SVD77h13eEUpY4JM2UtE7SwkjZjyWtkvR6OHwiUneVpGpJiyWdGimfGpZVS7oyVfE651w8RQV53PDpI1j2wQ5ufnxJpsPpEFK5xXEnMDVO+a/MbHw4zAWQNBY4FxgXTnObpFxJucCtwGnAWOC8sK1zzqXNsSNLOftDQ5gxbylvrfH7WaUscZjZPGBjgs3PAmaZWZ2ZvQdUA5PCodrMlprZbmBW2NY559LqB588jD5F+Vz54Bs0dPGzrPIysMzLJV0IvAx828w2AYOB+ZE2NWEZwMqY8snxZippOjAdoKysjKqqqr11tbW1+427xHnfJc/7Lnkdte8+d6i4fcEWfnTXY5wyLD/T4cSVjr5Ld+L4DfATwMK/vwC+CChOWyP+FlHcVG9mM4AZABMnTrTKysq9dVVVVUTHXeK875LnfZe8jtp3J5qxeNdL/HXpRr521kcZUlKU6ZAOkI6+S+tZVWa21swazKwR+C3BrigItiSGRpoOAVa3UO6cc2kniZ986nAAfvTXhZh1zV1WaU0ckgZFRj8NNJ1xNQc4V1KhpOHAKOBF4CVglKThkgoIDqDPSWfMzjkXNaSkiG+fMoYnF6/n4TfWZDqcjEjZripJ9wKVQKmkGuAaoFLSeILdTcuASwHMbJGk2cCbQD1wmZk1hPO5HHgEyAVmmpnf79g5l1EXHTuMOa+v4rqHF3HCqFL6FBVkOqS0SlniMLPz4hT/roX21wPXxymfC8xtx9Ccc+6g5OaIGz97JKff8gz/+9gSfnzmuEyHlFZ+5bhzziXhsEG9+NyEIdzzwgpWbd6Z6XDSyhOHc84l6YopowC4pYtdUe6JwznnkjS4T3fOP6aC+1+pYen62kyHkzatJg5J/SSdIelSSRdKmiAp3nUXzjnX5XytciQFuTn86rGus9XRbOKQdLykucCjBKfODgcmAD8FFkr6kaTi9ITpnHMdU/+ehXzxuGE8vGA1b67uGvexammL4zPA5WY2wcy+aGZXmtk3zOwTwNHAW8S/iaFzznUp048/lF7d8vjlo4szHUpatJQ4rjazpc3UHW5mD5jZA6kIyjnnsknvonwuPfFQHntrHa+u2JTpcFKupcTxmKTesYWSTgYeSl1IzjmXfS46dhilxQX8/JHOv9XRUuL4PfCkpH5NBZLOIbiI78xUB+acc9mkR2Eel31sJM+9+wHPVm/IdDgp1WziMLPbgVsIksfA8NYfPwVOMrPX0hWgc85li89PrqC8dzf+55HFnfoGiC2ejmtmvwf+G3gVuAj4aAvHPZxzrksrzMvl/5syitdXbubJxesyHU7KtHQ67muSXgW+C/QAegGPRMqdc87F+MyEIQzs1Y3fP7ss06GkTEs3Ofxc2qJwzrlOIj83hwuOqeDn/3qH6nW1jBzQ+S53a+kYx7stDekM0jnnssm5kyooyM3hrueXZTqUlGhpV9WTkr4qqTymPE/SCZJ+J+ni1IfonHPZpbS4kDOOKufBV2rYumtPpsNpdy0dHP8kkA/8RVKNpDckvQMsBS4GfhMePHfOORfjomOHsX13Aw+8XJPpUNpds8c4zGwHcDNws6RCYACw08w69wnKzjnXDo4Y0psJFX246/llXHTsMHJyOs+9YRO6rbqZ1ZnZSk8azjmXuGnHDmPZBzt4asn6TIfSrvx5HM45lyKnHT6IAT0LubOTnZqbssQhaaakdZIWxqn7jiSTVBqOS9LNkqrDYykTIm2nSVoSDtNSFa9zzrW3grwczp98CE+9s75TPegplVscdxLntuuShgIfB1ZEik8DRoXDdOA3Ydu+wDXAZGAScI2kkhTG7Jxz7eq8yUPJzxV3Pb8806G0m0SeAPhhSfMlbZG0S1KdpFafVmJm84CNcap+RXA1evRGLmcBd1lgPtBH0iDgVOBRM9toZpsIHirlzwBxzmWNAT27cfqR5TzwSg21dfWZDqddtHTleJPbgAuAWQRr/RcBQ5NZmKQzgVVmtiDm6bODgZWR8ZqwrLnyePOeTrC1QllZGVVVVXvramtr9xt3ifO+S573XfI6W98dUdjAX+rquXHWk0w5JD+ly0pH3yWSOHLMbLGkPDPbA/xW0nPA1W1ZkKQi4AfAKfGq45RZC+UHFprNAGYATJw40SorK/fWVVVVER13ifO+S573XfI6W99VAg+tepbnNuzhui+cmNJTc9PRd4kc49guqQBYIOkGSV8Hkrn5yqEEzy1fIGkZMAR4VdJAgi2J6FbMEGB1C+XOOZdVpn3kEJau387Ly7P/CYGJJI6LwnaXAw0EB7DbfANEM/u3mQ0ws2FmNowgKUwws/eBOcCF4dlVxwBbzGwN8AhwiqSS8KD4KWGZc85llVPHDaRbfg5zFqzKdCgHrcXEISkXuMbMdpnZZjP7kZldYWbvtDZjSfcCzwNjwluWXNJC87kEtzKpBn4LfA3AzDYCPwFeCofrwjLnnMsqPQrzOPmwMub++33qGxozHc5BafEYh5k1SBokKT88vpEwMzuvlfphkdcGXNZMu5nAzLYs2znnOqJi4oZ6AAAWeElEQVQzjyrn72+s4dl3P+DE0f0zHU7SEjk4vhR4WtJDwPamQjO7OWVROedcJ1Q5pj89u+Ux5/XVWZ04EjnGsZ7g+okioH9kcM451waFeblMHTeQfy16n117GjIdTtJa3eIwsx+lIxDnnOsKzjiqnPtfqaFq8TqmHj4o0+EkpdXEIelR4lw7YWbxrsdwzjnXgmMP7UdpcQEPL1jTeRMH8MPI627AZ4G61ITjnHOdW15uDp84YhD3vbSS2rp6igsT+RnuWFo9xmFmL0SGp8zsCoJbjzjnnEvCmUeVU1ffyKNvvp/pUJKSyE0Oe0WGPpJOBrJz+8o55zqACRUllPfuxpzXs/NGGIlsIy1i332j6oH3gC+nMijnnOvMcnLEGUeV87tn3mPT9t2U9CjIdEhtksjpuCPMrMLMhprZcDM7CXg21YE551xndsZR5dQ3GnMXrsl0KG2WSOJ4IU7Zi+0diHPOdSXjynsxon8PHl6Qfburmt1VJWkAwbGM7pKOYN8tznsRXAzonHMuSZI486hybnp8Ce9v2cXA3t0yHVLCWtri+CTwa4Jbmd8G3BoO3wf8okDnnDtIZxxVjhn87Y3s2upodovDzH4P/F7SOWY2O40xOedcl3Bo/2LGDurFPxa+z5eOH5HpcBKWyC1HZks6FRhHcAFgU/kNqQzMOee6gpMPG8CtT1azZcceehel9rGy7SWR6zhuA6YB3wK6Ezx/fGSK43LOuS7hxNH9aTR49t0NmQ4lYYmcVXWcmX0e+CC84eFkguMezjnnDtL4oX3o2S2Ppxavz3QoCUskcexq+hs+H3wXMCxlETnnXBeSl5vDcSNLmbdkPcEz7Tq+RBLHXEl9gJ8DrwPLgAdSGZRzznUlJ4zuz5otu1iyrjbToSSkxYPjknKAf5jZZuB+SX8Duvtzv51zrv2cED4NcN476xld1jPD0bSuxS0OM2sEboqM7/Sk4Zxz7Wtwn+6MHFDMU+9kx3GORHZVPSrprLbOWNJMSeskLYyU/UTSG5Jel/QvSeVhuSTdLKk6rJ8QmWaapCXhMK2tcTjnXDY4cXR/XnhvIzt3d/xHyiaSOC4H/iJpp6SNkjZJSmSr405gakzZ/5jZkWY2HvgbcHVYfhowKhymA78BkNQXuIbgTK5JwDWSShJYtnPOZZUTRvdnd30j89/7INOhtCqRxFEK5APFQP9wvH9rE5nZPGBjTNnWyGgP9j2S9izgLgvMB/pIGgScCjxqZhvNbBPwKAcmI+ecy3qTh/elMC+HeVmwuyqRK8cbJJ1LcHv1GyQNAcqAV5JZoKTrgQuBLcDHwuLBwMpIs5qwrLnyePOdTrC1QllZGVVVVXvramtr9xt3ifO+S573XfK6at+N6iP+8fpyTuyZfPJIR9+1mjgk/Zpgi+ME4AZgB3A78OFkFmhmPwB+IOkqgt1g17Dvzrv7NW2hPN58ZwAzACZOnGiVlZV766qqqoiOu8R53yXP+y55XbXv3s17j5/87U0OPXISQ/smdxPydPRdIruqjjWzSwkvBAzPqmqPx1XdA3w2fF0DDI3UDQFWt1DunHOdzomjSwGYt6Rj765KJHHsCa/nMABJ/YDGZBYmaVRk9Ezg7fD1HODC8OyqY4AtZrYGeAQ4RVJJeFD8lLDMOec6nUP7FzO4T/cOf5wjkWeO3wo8CPSXdC1wDnBtaxNJuheoBEol1RDskvqEpDEEiWc58JWw+VzgE0A1wa6wiyHYupH0E+ClsN11fh2Jc66zksQJo0t5eMEa9jQ0kp+byLp9+iVycPwuSa8AU8Kis81sYUvThNOdF6f4d820NeCyZupmAjNbW55zznUGJ47uz70vruS1FZuZNLxvpsOJK9F0lgvsAXa3YRrnnHNtdOzIUnJzxFPvrMt0KM1K5HkcPwDuBcoJDk7fE54R5Zxzrp316pbPhIo+zHun4z6fI5GthwuAD5vZD8NTaScRXIfhnHMuBU4Y1Z9/r9rCB7V1mQ4lrkQSx3L2PxaSByxNTTjOOeeOObQfAK+u2JzhSOJLJHHsABZJukPSb4F/A5sl/VLSL1MbnnPOdT1HDO5NXo54dcWmTIcSVyKn4/49HJrMT1EszjnngG75uYwt78Vr2Zo4zCzuKbTOOedSZ0JFCfe9tJL6hkbyOtj1HImcVTVV0kvhszXaclt155xzSTq6og879zTw9vvbMh3KARJJY78GLiW4K23Ct1V3zjmXvAkVwaOHOuLuqkQSRw3wupntMbOGpiHVgTnnXFc2pKQ7pcWFHfLMqkQOjn8XeFhSFbD3pGIzuzlVQTnnXFcniQkVfbJ2i+NaoAHoQ7CLqmlwzjmXQhMOKWHZBzs63IWAiWxxDDCzD6U8Euecc/vZd5xjM1PGlmU4mn0S2eJ4XNJJKY/EOefcfjrqhYCJJI4vA49JqvXTcZ1zLn26F+Ry2KBevNbBDpAnkjhKCZ453hs/Hdc559JqQkUfFtRspr4hqQevpkSriSM89fZs4Hvh60HA+FQH5pxzLjhAvmN3A4vXdpwLARO5cvzXwMeAL4RFO4DbUxmUc865wNFDgwPkHel6jkR2VR1rZpcCuyB4DjhQkNKonHPOATC0b3dKiws61PUciSSOPZJyAAOQ1A9odWebpJnh/a0WRsr+R9Lbkt6Q9BdJfSJ1V0mqlrRY0qmR8qlhWbWkK9v03znnXJaTxNEVJR3qAHmziUNS0zUetwIPAv0lXQs8A/z/Ccz7TmBqTNmjwOFmdiTwDnBVuKyxwLnAuHCa2yTlSsoNl38aMBY4L2zrnHNdxoSKEt7bsJ2N23dnOhSg5S2OFwHM7C7gh8DPgU3A2WY2q7UZm9k8YGNM2b/MrD4cnU/wDHOAs4BZZlZnZu8B1QSPqJ0EVJvZUjPbDcwK2zrnXJdxdEWwc6aj7K5q6cpxNb0ws0XAonZe9heB+8LXg9n/AVE1YRnAypjyyfFmJmk6MB2grKyMqqqqvXW1tbX7jbvEed8lz/sued53+6urN3IEf316AblrWz7EnI6+aylx9Jf0reYqzSzpx8ZK+gFQD9zdVBRvEcTfIrJm4pkBzACYOHGiVVZW7q2rqqoiOu4S532XPO+75HnfHWjsW0/zgfKprDymxXbp6LuWEkcuUEz8H/WkSZoGnA6cbGZNSaAGGBppNgRYHb5urtw557qMCRUlPPhKDQ2NRm5Ou/4st1lLiWONmV3XnguTNBX4HnCime2IVM0B7pH0S6AcGEVwjEXAKEnDgVUEB9A/354xOedcNji6og93Pb+cJeu28R8De2U0lpYOjh9USpN0L/A8MEZSjaRLCJ4m2BN4VNLrkm6HvcdQZgNvAv8ELgsfGFUPXA48ArwFzA7bOudclzJ2UG8AFneAR8m2tMVx8sHM2MzOi1P8uxbaXw9cH6d8LjD3YGJxzrlsN7y0B3k56hCJo9ktjvAKceeccx1AQV4OI/r34J0OcM+qRK4cd8451wGMLuvZIW526InDOeeyxJiynqzcuJPtdfWtN04hTxzOOZclRg/sCUD1utqMxuGJwznnssSYsiBxZHp3lScO55zLEkP7FtEtP4d3MnxmlScO55zLErk5YuSAYt/icM45l7jRZT0zfkquJw7nnMsiY8p6snZrHZt3ZO7ZHJ44nHMuizSdWfXO2sydWeWJwznnskhHOLPKE4dzzmWRQb270bMwL6NnVnnicM65LCKJ0QMze+sRTxzOOZdlms6s2vcsvPTyxOGcc1lmTFkxm3fsYf22uows3xOHc85lmUyfWeWJwznnskymz6zyxOGcc1mmX3EhpcUFGTuzyhOHc85loUw+1MkTh3POZaHRZT1ZsnYbjY3pP7MqZYlD0kxJ6yQtjJSdLWmRpEZJE2PaXyWpWtJiSadGyqeGZdWSrkxVvM45l01Gl/Vk++4GVm3emfZlp3KL405gakzZQuAzwLxooaSxwLnAuHCa2yTlSsoFbgVOA8YC54VtnXOuSxszsBggI3fKTVniMLN5wMaYsrfMbHGc5mcBs8yszszeA6qBSeFQbWZLzWw3MCts65xzXdqoDJ5ZlZf2JcY3GJgfGa8JywBWxpRPjjcDSdOB6QBlZWVUVVXtrautrd1v3CXO+y553nfJ875LTN9u4ukF1YylZm9ZOvquoyQOxSkz4m8RxT0SZGYzgBkAEydOtMrKyr11VVVVRMdd4rzvkud9lzzvu8QcsfRF1m+ro7Ly+L1l6ei7jnJWVQ0wNDI+BFjdQrlzznV5h/QrombTjrQvt6MkjjnAuZIKJQ0HRgEvAi8BoyQNl1RAcAB9TgbjdM65DmNISXe27qpny849aV1uynZVSboXqARKJdUA1xAcLL8F6A/8XdLrZnaqmS2SNBt4E6gHLjOzhnA+lwOPALnATDNblKqYnXMumwwtKQKgZtMOenfvnbblpixxmNl5zVT9pZn21wPXxymfC8xtx9Ccc65TGLI3cexkXHn6EkdH2VXlnHOujYaUdAeCxJFOnjiccy5L9SnKp0dBLis3pvcAuScO55zLUpIYUlLkWxzOOecSN7Rv97SfkuuJwznnstiQkiJWbdqZ1uePe+JwzrksNqSkO9vq6tm6sz5ty/TE4ZxzWazpzKqVadxd5YnDOeey2JDIRYDp4onDOeey2NDIRYDp4onDOeeyWK/uefQszPPE4ZxzLjGSGFyS3lNyPXE451yWG1JSxMqNvsXhnHMuQU0XAabrWg5PHM45l+WGlBSxfXcDm3ek57kcnjiccy7LpfsuuZ44nHMuy+1LHOk5QO6JwznnslzTRYDpunrcE4dzzmW53t3z6dUtfddyeOJwzrlOIJ3P5UhZ4pA0U9I6SQsjZX0lPSppSfi3JCyXpJslVUt6Q9KEyDTTwvZLJE1LVbzOOZfNhqTxIsBUbnHcCUyNKbsSeNzMRgGPh+MApwGjwmE68BsIEg1wDTAZmARc05RsnHPO7dN0EWA6ruVIWeIws3nAxpjis4A/hK//AHwqUn6XBeYDfSQNAk4FHjWzjWa2CXiUA5ORc851eUP7dmfnnga2peFSjrzUL2I/ZWa2BsDM1kgaEJYPBlZG2tWEZc2VH0DSdIKtFcrKyqiqqtpbV1tbu9+4S5z3XfK875Lnfdd2m9cFD3Ja+cH2lPdduhNHcxSnzFooP7DQbAYwA2DixIlWWVm5t66qqorouEuc913yvO+S533XdmVrtnLTq0+zXd1S3nfpPqtqbbgLivDvurC8BhgaaTcEWN1CuXPOuYjB4UWAG3Y2pnxZ6U4cc4CmM6OmAQ9Fyi8Mz646BtgS7tJ6BDhFUkl4UPyUsMw551xEr2759O6ez4YdqT84nrJdVZLuBSqBUkk1BGdH3QjMlnQJsAI4O2w+F/gEUA3sAC4GMLONkn4CvBS2u87MYg+4O+ecIzhAvmHn9pQvJ2WJw8zOa6bq5DhtDbismfnMBGa2Y2jOOdcpjSnrxbJVWZw4nHPOpdcvzjmKqqpNKV+O33LEOedcm3jicM451yaeOJxzzrWJJw7nnHNt4onDOedcm3jicM451yaeOJxzzrWJJw7nnHNtonQ89CPdJK0HNgNbwqJSYEMKF9k7sqxUTNdSu7bWxZa1Nt4V+y7R8s7cd621yea+y+RnrqX6RL6v8cras+8OMbP+rbYys045ADMir19O17JSMV1L7dpaF1uWwHiX67tEyztz37XWJpv7LpOfubb2XZLf4ZR+7sysU++qejgLlpXodC21a2tdbFlr46nWEfsu0fLO3HettcnmvsvkZ66l+kS+r/HK0v2565y7qmJJetnMJmY6jmzkfZc877vked8lLx1915m3OKJmZDqALOZ9lzzvu+R53yUv5X3XJbY4nHPOtZ+ussXhnHOunXjicM451yaeOJxzzrVJl00cknpIekXS6ZmOJZtIqpT0tKTbJVVmOp5sIilH0vWSbpE0LdPxZAtJx4eftzskPZfpeLKJpApJcyTNlHRle8230ySOsGPWSVoYUz5V0mJJ1TEd9z1gdnqj7Jja2HcG1ALdgJp0x9rRtLHvzgIGA3vo4n3Xln4zs6fN7CvA34A/ZCLejqSNn7nRwN/N7IvA2HYLItVXGKZrAE4AJgALI2W5wLvACKAAWBB23hTgXOAi4PRMx57poY19lxPWlwF3Zzr2TA9t7LsrgUvDNg9kOvZs6bdI/WygV6Zjz/TQxs9cP+BJ4Ang4vaKIS/JfNPhmNk8ScNiiicB1Wa2FEDSLIK1vmKgB0HH7pQ018wa0xhuh9KWvjOzN8P6TUBh2oLsoNr4uVsJ7A7bNKQrxo6ojf32pqQKYIuZbU1roB1QG/tuD3BNOM0DwO/bI4ZOkziaMZjgy9qkBphsZpcDSLoI2NCVk0YL4vadpM8ApwJ9gF9nIrAsELfvgJuAWyQdD8zLRGAdXHP9BnAJ7fSj10k113e3Az+W9HlgWXstrLMnDsUp23vFo5ndmb5Qsk7cvjOzPwN/TncwWaa5vttB8APo4mv2+2pm16Q5lmzT3GduIfC59l5Ypzk43owaYGhkfAiwOkOxZBvvu+R53yXH+y15ae27zp44XgJGSRouqYDggPicDMeULbzvkud9lxzvt+Slte86TeKQdC/wPDBGUo2kS8ysHrgceAR4C5htZosyGWdH5H2XPO+75Hi/Ja8j9J3f5NA551ybdJotDuecc+nhicM551ybeOJwzjnXJp44nHPOtYknDuecc23iicM551ybeOJwzjnXJp44nAMkNUh6XdJCSfdLKuoAMR0t6Y7w9UWS1ocxvinpy61Me4SkO9MSqOtyPHE4F9hpZuPN7HCCW59/JVqpQLq/L98HbomM32dm44FK4AZJZc1NaGb/BoaEtyN3rl154nDuQE8DIyUNk/SWpNuAV4Ghkk6R9LykV8Mtk2IASTeGWwJvSPp5WHaGpBckvSbpsaYfekl9Jf01bDtf0pGxAUjqCRxpZgti68xsHcFDew6RNEnSc+EynpM0JtL0YYJ7FjnXrjxxOBchKQ84Dfh3WDQGuMvMjga2Az8EppjZBOBl4FuS+gKfBsaZ2ZHAT8NpnwGOCaedBXw3LL8WeC1s+33grjihTAQWxilH0giCJ71VA28DJ4TLuBq4IdL0ZeD4tvWAc63r7M/jcC5R3SW9Hr5+GvgdUA4sN7P5YfkxBE+NfFYSBI/ofB7YCuwC7pD0d4JnY0Nwa+v7JA0K274Xlh8HfBbAzJ6Q1E9SbzPbEolnELA+Jsb/lHQcUEfwCNqNkoYCf5A0iuDZFfmR9uvC/8G5duWJw7nAzvD4wV5hctgeLQIeNbPzYieWNAk4mWDX0OXASQTHJ35pZnMkVQI/jswnVuzdRncC3WLK7mt6emXET4AnzezT4eNEqyJ13cL5ONeufFeVc4mbD3xU0kgASUWSRofHOXqb2VzgG0BTAuoNrApfT4vMZx5wfjiPSoLHF8c+S/stYGQCMUWXcVFM3Wia2d3l3MHwLQ7nEmRm68Pn1N8rqTAs/iGwDXhIUjeCrYlvhnU/Bu6XtIog6QyPlP9e0hvADvZPKk3LeltSb0k9zWxbC2H9jGBX1beAJ2LqPgb8vW3/pXOt8+dxONdBSfomsM3M7khi2kLgKeC48CE/zrUb31XlXMf1G4ID4cmoAK70pOFSwbc4nHPOtYlvcTjnnGsTTxzOOefaxBOHc865NvHE4Zxzrk08cTjnnGuT/wdJwNC0M6mHQAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "P_grid = np.zeros(len(V_grid))\n",
    "\n",
    "for j in range (len(V_grid)):\n",
    "    T_grid[j]=T_calc(V_grid[j])\n",
    "    P_grid[j] = pressure_total(T_grid[j],V_grid[j])\n",
    "\n",
    "plt.semilogx(P_grid,T_grid)\n",
    "plt.title('μJT=0 - Curva de inversão para o SH2 com CPA EOS ')\n",
    "plt.ylabel('Temperatura (K)') \n",
    "plt.xlabel('Presao (Pa)')\n",
    "plt.grid(True)\n",
    "plt.show()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}