4 Diagrama de fases para sustancias puras

En esta sección se presentan los diagramas de fases comunes para sustancias puras. Como son:

  1. Envolvente de fases liquido-vapor
  2. Isoterma
  3. Isobara
  4. Sólido-líquido

In [ ]:


In [ ]:


In [ ]:


In [ ]:


In [ ]:


In [4]:
%load_ext fortranmagic  # activating magic


---------------------------------------------------------------------------
ImportError                               Traceback (most recent call last)
<ipython-input-4-430be5188958> in <module>()
----> 1 get_ipython().magic('load_ext fortranmagic  # activating magic')

/home/andres-python/anaconda3/lib/python3.5/site-packages/IPython/core/interactiveshell.py in magic(self, arg_s)
   2156         magic_name, _, magic_arg_s = arg_s.partition(' ')
   2157         magic_name = magic_name.lstrip(prefilter.ESC_MAGIC)
-> 2158         return self.run_line_magic(magic_name, magic_arg_s)
   2159 
   2160     #-------------------------------------------------------------------------

/home/andres-python/anaconda3/lib/python3.5/site-packages/IPython/core/interactiveshell.py in run_line_magic(self, magic_name, line)
   2077                 kwargs['local_ns'] = sys._getframe(stack_depth).f_locals
   2078             with self.builtin_trap:
-> 2079                 result = fn(*args,**kwargs)
   2080             return result
   2081 

<decorator-gen-62> in load_ext(self, module_str)

/home/andres-python/anaconda3/lib/python3.5/site-packages/IPython/core/magic.py in <lambda>(f, *a, **k)
    186     # but it's overkill for just that one bit of state.
    187     def magic_deco(arg):
--> 188         call = lambda f, *a, **k: f(*a, **k)
    189 
    190         if callable(arg):

/home/andres-python/anaconda3/lib/python3.5/site-packages/IPython/core/magics/extension.py in load_ext(self, module_str)
     35         if not module_str:
     36             raise UsageError('Missing module name.')
---> 37         res = self.shell.extension_manager.load_extension(module_str)
     38 
     39         if res == 'already loaded':

/home/andres-python/anaconda3/lib/python3.5/site-packages/IPython/core/extensions.py in load_extension(self, module_str)
     81             if module_str not in sys.modules:
     82                 with prepended_to_syspath(self.ipython_extension_dir):
---> 83                     __import__(module_str)
     84             mod = sys.modules[module_str]
     85             if self._call_load_ipython_extension(mod):

ImportError: No module named 'fortranmagic  # activating magic'

In [ ]:


In [ ]:

Envolvente de fases para sustancias puras

En esta sección se presenta el desarrollo y manipulación de las ecuaciones para establecer un algortimo que permita calcular los puntos del diagrama de fases para una sustancia pura utilizando ecuaciones de estado para modelar el comportamiento de fases.

$$ P^l(T, V^l) - P^v(T, V^v) = 0$$
$$ ln \phi(T, V^l) - ln \phi(T, V^v) = 0 $$
$$ g(T, V^l, V^v) = X_S - S = 0 $$
$$ F = \begin{bmatrix} P^l(T, V^l) - P^v(T, V^v)\\ ln \phi(T, V^l) - ln \phi(T, V^v)\\ X_S - S \end{bmatrix} = 0$$
$$ F = \begin{bmatrix} ln \left( \frac{P^l(T, V^l)} {P^v(T, V^v)} \right)\\ ln f_l(T, V^l) - ln f_v(T, V^v)\\ 0 \end{bmatrix} $$
$$ J \begin{bmatrix} \Delta ln T\\ \Delta ln V^l\\ \Delta ln V^v\\ \end{bmatrix} + F = 0$$

Se inicia con el cálculo de la envolvente de fases de uns sustancia pura.

$$ n \left(\frac{\partial ln \hat\phi_i}{\partial n_j}\right)_{T,P} = n \left(\frac{\partial ln \hat\phi_j}{\partial n_j}\right)_{T,P} = n \left(\frac{\partial^2 F} {\partial n_j \partial n_i} \right)_{T,V} + 1 + \frac{n}{RT} + \frac{ \left(\frac{\partial P} {\partial n_j}\right)_{T,V} \left(\frac{\partial P} {\partial n_i}\right)_{T,V} } {\left(\frac{\partial P} {\partial V} \right)_{T,n}}$$

Relación entre derivadas de la fugacidad y el coeficiente de fugacidad

$$ \left(\frac{\partial ln\hat f_i }{\partial n_i}\right)_{T,P} = \left(\frac{\partial \hat \phi_i}{\partial n_i}\right)_{T,P} + \left(\frac{\delta_{ij}}{n} - \frac{1}{n} \right)$$
$$RJAC(1,1)=DPDTx/Pl-DPDTy/Pv$$$$ RJAC(1,1)=T*RJAC(1,1) $$$$ RJAC(1,1) = T \left( \frac {\left(\frac{\partial P_{x} }{\partial T}\right)} {P_l} - \frac {\left(\frac{\partial P_{y} }{\partial T}\right)} {P_v} \right) $$$$RJAC(1,2)=Vl*DPDVx/Pl$$$$ RJAC(1,2) = -V_l \left( \frac {\left(\frac{\partial P_{x} }{\partial V}\right)} {P_l} \right) $$$$RJAC(1,3)=-Vv*DPDVy/Pv$$$$ RJAC(1,2) = -V_v \left( \frac {\left(\frac{\partial P_{y} }{\partial V}\right)} {P_y} \right) $$$$RJAC(2,1)=FUGTx(icomp)-FUGTy(icomp)$$$$ RJAC(2,1)= \left(\frac{\partial f_{ix} } {\partial T} \right) - \left(\frac{\partial f_{iy} } {\partial T} \right)$$

$$RJAC(2,1)=T*RJAC(2,1)$$$$RJAC(2,1)=T* \left(\left(\frac{\partial f_{ix} } {\partial T} \right) - \left(\frac{\partial f_{iy} } {\partial T} \right) \right) $$$$RJAC(2,2)=Vl*FUGVx(icomp)$$$$ RJAC(2,1)= V_l \left(\frac{\partial f_{ix} } {\partial V} \right) $$$$RJAC(2,3)=-Vv*FUGVy(icomp)$$$$ RJAC(2,3)= V_v \left(\frac{\partial f_{iy} } {\partial V_{y}} \right) $$$$DLFUGT(I)=(ArTn(I)-Arn(I)/T)/RT+1.D0/T ! term DPDT/P is cancelled out $$$$DLFUGT(I)=(ArTn(I)-Arn(I)/T)/RT+1.D0/T $$$$DLFUGT(I) = \frac{ ArTn(I) - Arn(I)} {R_GT^2} + \frac{1} {T} $$$$DLFUGT(I) = \frac{ \frac{\partial Ar}{\partial T \partial n} - \frac{\partial Ar}{\partial n }} {R_GT^2} + \frac{1} {T} $$$$ DPDT = -ArTV+TOTN*RGAS/V $$$$ \frac{\partial P}{\partial T} = -\frac{\partial Ar}{\partial T \partial V} + N_T \frac{R_{G}} {V} $$
$$RJAC(1,1)=DPDTx/Pl-DPDTy/Pv$$$$ RJAC(1,1)=T*RJAC(1,1) $$$$ RJAC(1,1) = T \left( \frac {\left(\frac{\partial P_{x} }{\partial T}\right)} {P_l} - \frac {\left(\frac{\partial P_{y} }{\partial T}\right)} {P_v} \right) $$$$RJAC(1,2)=Vl*DPDVx/Pl$$$$ RJAC(1,2) = -V_l \left( \frac {\left(\frac{\partial P_{x} }{\partial V}\right)} {P_l} \right) $$$$RJAC(1,3)=-Vv*DPDVy/Pv$$$$ RJAC(1,3) = -V_v \left( \frac {\left(\frac{\partial P_{y} }{\partial V}\right)} {P_y} \right) $$$$RJAC(2,1)=FUGTx(icomp)-FUGTy(icomp)$$$$ RJAC(2,1)= \left(\frac{\partial f_{ix} } {\partial T} \right) - \left(\frac{\partial f_{iy} } {\partial T} \right)$$

$$RJAC(2,1)=T*RJAC(2,1)$$$$RJAC(2,1)=T* \left(\left(\frac{\partial f_{ix} } {\partial T} \right) - \left(\frac{\partial f_{iy} } {\partial T} \right) \right) $$$$RJAC(2,2)=Vl*FUGVx(icomp)$$$$ RJAC(2,1)= V_l \left(\frac{\partial f_{ix} } {\partial V} \right) $$$$RJAC(2,3)=-Vv*FUGVy(icomp)$$$$ RJAC(2,3)= V_v \left(\frac{\partial f_{iy} } {\partial V_{y}} \right) $$$$DLFUGT(I)=(ArTn(I)-Arn(I)/T)/RT+1.D0/T ! term DPDT/P is cancelled out $$$$DLFUGT(I)=(ArTn(I)-Arn(I)/T)/RT+1.D0/T $$$$DLFUGT(I) = \frac{ \frac{ArTn(I) - Arn(I)} {T}} {R_GT} + \frac{1} {T} $$$$ DPDT = -ArTV+TOTN*RGAS/V $$$$ \frac{\partial P}{\partial T} = -\frac{\partial Ar}{\partial T \partial V} + N_T \frac{R_{G}} {V} $$

_Note:

$$\begin{bmatrix} T \left( \frac {\left(\frac{\partial P_{x} }{\partial T}\right)} {P_l} - \frac {\left(\frac{\partial P_{y} }{\partial T}\right)} {P_v} \right) & x_{12} & x_{13} & \dots & x_{1n} \\ x_{21} & x_{22} & x_{23} & \dots & x_{2n} \\ \hdotsfor {5} \\ x_{d1} & x_{d2} & x_{d3} & \dots & x_{dn} \end{bmatrix}$$
$$J_x = \begin{bmatrix} T \left( \frac {\left(\frac{\partial P_{x} }{\partial T}\right)} {P_l} - \frac {\left(\frac{\partial P_{y} }{\partial T}\right)} {P_v} \right) & -V_l \left( \frac {\left(\frac{\partial P_{x} }{\partial V}\right)} {P_l} \right) & -V_v \left( \frac {\left(\frac{\partial P_{y} }{\partial V}\right)} {P_y} \right) \\ \left(\frac{\partial f_{ix} } {\partial T} \right) - \left(\frac{\partial f_{iy} } {\partial T} \right) & V_l \left(\frac{\partial f_{ix} } {\partial V} \right) & V_v \left(\frac{\partial f_{iy} } {\partial V_{y}} \right) & \\ 0 & 0 & 0 & \end{bmatrix}$$

In [1]:
J_x =  \begin{bmatrix}
T \left( \frac {\left(\frac{\partial P_{x} }{\partial T}\right)} {P_l} - \frac {\left(\frac{\partial P_{y} }{\partial T}\right)} {P_v}  \right) & 
-V_l \left( \frac {\left(\frac{\partial P_{x} }{\partial V}\right)} {P_l} \right) & 
-V_v \left( \frac {\left(\frac{\partial P_{y} }{\partial V}\right)} {P_y} \right) \\
    \left(\frac{\partial f_{ix} } {\partial T} \right) - \left(\frac{\partial f_{iy} } {\partial T} \right)       & V_l \left(\frac{\partial f_{ix} } {\partial V} \right) & V_v \left(\frac{\partial f_{iy} } {\partial V_{y}} \right) &   \\
    0       & 0 & 0 &  
\end{bmatrix}$$


  File "<ipython-input-1-72792e2ef634>", line 1
    J_x =  \begin{bmatrix}
                          ^
SyntaxError: unexpected character after line continuation character

In [ ]:

Envolvente para mezclas

$$ f_i = ln K_i + ln \hat\phi_i^v(T,P,y) - ln \hat\phi_i^l(T,P,x) = 0 $$$$ i = 1,2,... C $$$$ f_{C+1} = \sum_{i=1}^C(y_i - x_i) = 0 $$$$ f_{C+2} = X -X_{spec} = 0 $$
$$ x_i = \frac{z_i} {1-\beta+ \beta K_i}$$$$ y_i = \frac{K_iz_i} {1-\beta+ \beta K_i}$$
$$ J_{ij} = \frac{\partial f_i}{\partial ln K_j} = \frac{\partial ln K_i}{\partial ln K_j} + \frac{\partial \hat \phi_i^v}{\partial ln K_j} - \frac{\partial \hat \phi_i^l}{\partial ln K_j} $$
$$ \frac{\partial ln K_i}{\partial ln K_j} = \left\{ \begin{array}{lcc} 1 & i = j \\ \\ 0 & i \neq j \\ \end{array} \right.$$
$$ \frac{\partial ln \hat \phi_i^v}{\partial ln K_j} = \sum_{k=1}^C\frac{\partial ln \hat \phi_i^v}{\partial y_k} \frac{\partial y_k}{\partial ln K_j} $$
$$ \frac{\partial y_k}{\partial ln K_j} = 0 $$$$ k \neq j $$
$$ \frac{\partial ln \hat \phi_i^v}{\partial ln K_j} = \frac{\partial ln \hat \phi_i^v}{\partial y_k} \frac{\partial y_k}{\partial ln K_j} $$
$$ \frac{\partial x_i}{\partial ln K_j} = K_j \frac{\partial x_i}{\partial K_j} = \frac{ \beta K_j z_j}{(1 - \beta + \beta K_j)^2} = -\beta \frac{x_i y_i}{z_i} $$
$$ \frac{\partial y_i}{\partial ln K_j} = (1 - \beta) \frac{x_i y_i}{z_i} $$

Finalmente, el termino $ \frac{\partial ln \hat \phi_i^v}{\partial y_k} $

4. Descripción del algoritmo

Donde $J_F$ es la matriz jacobiana de la función vectorial $F$, $Λ$ es el vector de variables del sistema $F=0$, $S_{Spec}$ es el valor asignado a una de las variables del vector $Λ$, $\frac{dΛ}{ dS_{Spec}}$ es la derivada, manteniendo la condición $F=0$, del vector de variables con respecto al parámetro $S_{spec}$. Observe que si $S_{spec}=Λ_i$, entonces $\frac{dΛi} {dS_{Spec}} =1$. El vector $\frac{dΛ}{ dS_{Spec}}$ es llamado “vector de sensitividades”.

$\frac{\partial F} {\partial S_{Spec}}$ es la derivada parcial del vector de funciones $F$ con respecto la variable $S_{spec}$.

En la ecuación A.3-1 la matriz jacobiana $J_F$ debe ser valuada en un punto ya convergido que es solución del sistema de ecuaciones $F=0$. Observe en los distintos sistemas de ecuaciones presentados en el capítulo 3, que sólo una componente del vector $F$ depende explícitamente de $S_{spec}$. Por tanto, las componentes del vector $\frac{\partial F} {\partial S_{Spec}}$ son todas iguales a cero, excepto la que depende de $S_{spec}$, en esta tesis el valor de dicha componente es siempre $“-1”$.

Conocidos $J_F$ y $\frac{\partial F} {\partial S_{Spec}}$ es posible calcular todas las componentes del vector $\frac{dΛ}{ dS_{Spec}}$ .

Con dΛ dSSpec conocido es posible predecir los valores de todas las variables del vector Λ para el siguiente punto de la “hiper-línea que se está calculando, aplicando la siguiente ecuación:

$$ A_{next point}^0 = A _{conve. pont} + \left(\frac{dA}{dS_{Spec}}\right) \Delta S_{Spec} $$

Aquí Λ0 next point corresponde al valor inicial del vector Λ para el próximo punto a ser calculado. Λconv. point es el valor del vector Λ en el punto ya convergido.

Por otra parte, el vector de sensitividades dΛ dSSpec provee información sobre la próxima variable que debe ser especificada en el próximo punto a ser calculado. La variable a especificar corresponderá a la componente del vector dΛ dSSpec de mayor valor absoluto. Supongamos que la variable especificada para el punto convergido fue la presión P, es decir en el punto convergido Sspec = P. Luego al calcular el vector de sensitividades para el punto convergido usando A.3-1, supongamos que se determina que la componente de mayor valor absoluto de dicho vector es la correspondiente a dT dP , esto implica que en el próximo punto a ser calculado la variable que se debe especificar ya no es P sino T. Esto es Sspec = T. Cuando existe un cambio de variables especificadas como el caso del ejemplo anterior, el vector de sensitividades se normaliza dividiendo todas las componentes del vector por la de mayor valor absoluto, en el caso anterior se deberían dividir todas las componentes por dT dP . Finalmente se aplica A.3-2 para encontrar los valores de Λ0 next point

La variable ∆SSpec se computa en esta tesis de la siguiente manera:


In [ ]:

Donde ∆SSpec _Old es el paso que se dio en la variable que se especifico para obtener el punto convergido. “N” es una constante impuesta por el usuario, e ITER es el número de iteraciones requeridas para el punto convergido. Note que para calcular ∆SSpec se utiliza la componente de dΛ dSSpec de mayor valor absoluto,

max   Λ dSSpec d , antes de que el vector dΛ dSSpec haya sido normalizado.

4. Referencias

[1] key

[#] E.L. Allgower, K. Georg, Introduction to Numerical Continuation Methods, SIAM. Classics in Applied Mathematics, Philadelphia, 2003.

[#] M. Cismondi, M.L. Michelsen, Global phase equilibrium calculations: Critical lines, critical end points and liquid-liquid-vapour equilibrium in binary mixtures, Journal of Supercritical Fluids, 39 (2007) 287-295.

[#] M. Cismondi, M.L. Michelsen, M.S. Zabaloy, Automated generation of phase diagrams for binary systems with azeotropic behavior, Industrial and Engineering Chemistry Research, 47 (2008) 9728-9743.