Relatório ELLB

Prévia do material em texto

UNIVERSIDADE DO ESTADO DO RIO DE JANEIRO
INSTITUTO DE QUÍMICA 
Departamento Físico Química 
Termodinâmica Experimental para Engenharia Química 
Prof. Pedro Alijó
Solubilidade Mútua de Dois Líquidos
Grupo 1: Gustavo Anciens
 Thaiane Nolasco
 Thainá Tussini
 Viviane Borges
Realização do experimento: 24/09/2015
Entrega: 
Rio de Janeiro
2015
Sumário
INTRODUÇÃO --------------------------------------------------------------------------2
OBJETIVO ------------------------------------------------------------------------------- 6
METODOLOGIA ------------------------------------------------------------------------6
MODELAGEM -------------------------------------------------------------------------- 7
RESULTADOS E DISCUSSÃO ----------------------------------------------------9
CONCLUSÃO --------------------------------------------------------------------------14
REFERÊNCIAS BIBLIOGRÁFICAS ---------------------------------------------14
APÊNDICE I – ROTINA DE CÁLCULOS ---------------------------------------14
Introdução teórica
Estudos sobre o equilíbrio termodinâmico entre duas fases líquidas em contato e seus equilíbrios são necessários para compreender diversos processos utilizados num processo químico, como destilação e extração.
Muitos pares de espécies químicas, quando misturadas em uma certa faixa de composição, a certas temperaturas e pressões, não satisfazem o critério de estabilidade (dGT,P ≤ 0), se dividindo em duas fases líquidas com composições diferentes, conforme a sua miscibilidade. Isto ocorre pois nesse caso, o estado bifásico se mostra mais estável que o estado monofásico.
A energia livre de mistura (ΔGmist) pode ser definida como a mudança na energia livre do sistema para levar a cabo o processo de mistura dos dois líquidos:
Δmist = (A+B) – [xA A + xB B]
O ΔGmist a T e P constantes, pode variar com a composição do sistema, ou seja, com a fração de mols dos componentes, como indicado na figura 1.
Figura 1: Variação do ΔGmist
No gráfico (a), ΔGmist é negativo em toda a faixa de composição, de modo que os dois líquidos são completamente miscíveis em uma determinada pressão e temperatura. Em (b), ΔGmist é positivo, de modo que os dois líquidos são imiscíveis a uma determinada P e T. Em (c), observamos uma situação mais complexa: ΔGmist <0, logo os dois líquidos são miscíveis. No entanto, se a mistura tem uma composição entre X1 e X2, o ΔGmist será menor. Assim, naquela P e T os líquidos são parcialmente miscíveis, ou seja, em composições X < X1 e X > X2, os líquidos são miscíveis, o que não se observa em composições intermediárias.
Pela regra das fases, temos dois graus de liberdade para uma mistura binária. Especificando a pressão, restará a temperatura ou a composição, havendo somente uma temperatura na qual teremos determinada composição. O ELL binário é representado em um diagrama de solubilidade, um gráfico de T vs X1.
A temperatura tem grande influência no processo de separação de fases do equilíbrio líquido-líquido, e para cada sistema de substâncias existem temperaturas limitantes, sendo uma superior e outra inferior, denominadas temperaturas consolutas. Acima ou abaixo dessas temperaturas, respectivamente, não é possível a formação das fases. Deste modo, a separação das fases ocorre dentro de determinada faixa de temperatura e para cada temperatura existe uma faixa de composição na qual a separação é possível, como ilustrado na figura 2.
Figura 2: Influência da temperatura no ELL
	
No diagrama A, há a presença de uma temperatura consoluta inferior, ou temperatura crítica inferior da solução (TCIS), que ocorre quando o diagrama intercepta a curva dos pontos de bolha do ELV.
Em B, o diagrama intercepta a curva de congelamento, e assim apenas a temperatura consoluta superior, ou temperatura crítica superior da solução (TCSS), ocorre.
	No diagrama C, há a presença de uma TCIS e uma TCSS, que consistem em estados limites do equilíbrio bifásico, nos quais todas as propriedades das duas fases em equilíbrio são idênticas.
As curvas binodais do ELL quando interceptam a curva dos pontos bolha do ELV, resulta no equilíbrio líquido-líquido-vapor (ELLV), que pela regra das fases, possui apenas um grau de liberdade. Com isso, para uma determinada pressão, a temperatura e as composições das três fases são especificadas. A figura 3 representa um típico diagrama temperatura versus composição do equilíbrio ELLV.
Figura 3: Diagrama T x composição típico do ELLV binário
Os pontos que representam os estados das três fases em equilíbrio estão em uma linha horizontal a T* (temperatura trifásica). Os pontos C e D representam as duas fases líquidas e o ponto E representa a fase vapor. Em temperaturas acima de T*, o sistema pode ser composto por apenas uma única fase líquida, duas fases (líquida e vapor) ou uma única fase vapor, o que irá depender da composição global do sistema. Na região α, o sistema é composto de um único líquido rico na espécie 2; na região β, ele é um único líquido rico na espécie 1. Na região α-V, líquido e vapor estão em equilíbrio. Os estados das fases individuais estão sobre a linha AC e AE. Na região β-V também existem fases líquida e vapor em equilíbrio, descritas pelas linhas BD e BE. Na região identificada por V, o sistema é descrito por uma única fase vapor e abaixo da T*, o sistema é totalmente liquido (região do ELL).
Objetivo
O experimento em questão visa realizar a modelagem do equilíbrio líquido-líquido utilizando a modelagem NRTL, em um sistema binário entre água e etileno glicol monobutil éter (EGME), observando como a miscibilidade destes dois líquidos varia em função da temperatura. Ainda, realizar a avaliação do equilíbrio líquido-líquido-vapor em questão.
Metodologia
3.1 Materiais
Bureta 
Termômetro digital
Pipetas
Pró-pipete
Erlenmeyer
Placa de aquecimento
Agitador magnético
3.2 Reagentes
Etileno glicol monobutil éter (EGME)
Água destilada
3.3 Procedimento Experimental
	Em um erlenmeyer de 125 mL, foi adicionada uma alíquota de 40 mL de água destilada. Após, foi adicionada uma quantidade determinada de etileno glicol (EGME). 
A solução formada, inicialmente límpida, foi tampada com uma rolha levemente apoiada no bocal do frasco, e foi aquecida com agitação constante. A temperatura foi controlada através de um termômetro encaixado. No momento em que foi observado o início do turvamento da solução, registrou-se a primeira temperatura. Após o registro, a mistura foi retirada da placa de aquecimento e mantida sob agitação manual, até que se tornasse límpida novamente. Neste ponto, registrou-se a segunda temperatura.
Este procedimento foi repetido para sete soluções, com variação da quantidade de EGME adicionada, como ilustrado na tabela 1.
	Volume de EGME (ml)
	Volume de água (ml)
	5
	40
	6
	40
	11
	40
	19
	40
	29
	40
	44
	40
	55
	40
Tabela1: Quantidades usadas dos reagentes.
Para realizar a verificação do ELLV, foi utilizada a última composição preparada. Esta foi aquecida até o aparecimento da primeira bolha, registrando a temperatura. Após, retirou-se a mistura da placa de aquecimento e observou-se a temperatura na qual a ebulição cessou. Esta temperatura também foi registrada.
Modelagem
O modelo de energia de gibbs em excesso, de composição local, a ser utilizado para a modelagem do equilíbrio líquido-líquido desse experimento será o NRTL (NonRandom Two Liquid). O modelo NRTL, de Renon e Prausnitz, é uma melhoria da equação de Wilson, introduzindo um terceiro parâmetro (α), para considerar o fato das misturas líquidas serem não-randômicas (efeito entrópico). 
A equação tem a seguinte forma: 
 (1) 
Derivando a equação, obtemos a expressão geral para o coeficiente de atividade: 
 (2) 
Onde:Para misturas binárias as seguintes equações são utilizadas:
 (3)
 (4)
Os coeficientes de interação binária (b12 e b21) e o o parâmetro de não-aleatoriedade  serão encontrados a partir de chutes iniciais para as concentrações nas fases, e a temperatura média obtida será usada como um valor de entrada. 
De início, usamos valores tabelados para b12, b21 e  encontrados na literatura (DECHEMA). Os valores podem ser observados na tabela abaixo:
	Temperatura (°C)
	b12
	b21
	60
	-1040,2
	2577,5
	110
	-1131,6
	2855,0
	120
	1220,5
	2960,5
Tabela 2: Valores de b12 e b21 encontrados na literatura.
O valor de foi mantido em 0,2 (de acordo com a literatura o valor pode variar entre 0,2 e 0,47).
Seguindo recomendações da literatura, a relação de dependência dos parâmetros bij com a temperatura pôde ser realizada através de um ajuste polinomial, no caso quadrático, como mostrado abaixo:
Gráficos 1 e 2 – Ajuste de segunda ordem dos dados da literatura (Dechema Liquid-Liquid Equilibrium)
Resultados e Discussão
Equilíbrio Líquido-Líquido
Os dados experimentais obtidos durantes a prática estão dispostos na tabela 3. Através desses valores, foi possível o cálculo da fração volumétrica (ɸ) de água e EGME, presentes na tabela 5.
	Volume EGME (mL)
	Volume de Água (mL)
	T1 (°C)
	T2 (°C)
	Tmédia (°C)
	5
	40
	62,4
	62
	62,2
	6
	40
	55,8
	55,5
	55,7
	11
	40
	49
	48,8
	48,9
	19
	40
	46,7
	46,3
	46,5
	29
	40
	49
	49
	49,0
	44
	40
	53,7
	53,5
	53,6
	55
	40
	63
	63
	63,0
Tabela 3: Dados experimentais.
	
	Água
	EGME
	Massa molar (g/mol)
	18,03
	118,18
	Massa Específica (g/ml)
	0,993
	0,9015
Tabela 4: Dados dos reagents.
	Fração volumétrica
EGME (ɸ)
	Fração volumétrica
Água(ɸ)
	Massa EGME (g)
	Massa Água (g)
	Número de mols EGME
	Número de mols Água
	Número de mols total
	Fração molar 
EGME 
	Fração molar
Água
	0,1111
	0,8889
	4,508
	39,72
	0,038141
	2,203
	2,24114
	0,0170
	0,982981
	0,1304
	0,8696
	5,409
	39,72
	0,045769
	2,203
	2,24876
	0,0204
	0,979647
	0,2157
	0,7843
	9,917
	39,72
	0,08391
	2,203
	2,28691
	0,0367
	0,963308
	0,3220
	0,6780
	17,13
	39,72
	0,144936
	2,203
	2,34793
	0,0617
	0,938271
	0,4203
	0,5797
	26,14
	39,72
	0,221218
	2,203
	2,42421
	0,0913
	0,908747
	0,5238
	0,4762
	39,67
	39,72
	0,335641
	2,203
	2,53864
	0,1322
	0,867787
	0,5789
	0,4211
	49,58
	39,72
	0,419551
	2,203
	2,62255
	0,1600
	0,840022
Tabela 5: Fração volumétrica e fração molar dos reagentes.
Através dos resultados de composição molar, foram realizados ajustes polinomiais de diversas ordens para observar o melhor ajuste dos dados, após os ajuste foi confirmado que o melhor é o ajuste de 4ª ordem, apresentado a seguir:
Gráfico 3: Composição molar de EGME vs. Temperatura do sistema, para ajuste polinomial de quarta ordem
Conforme se observa no gráfico 3, há uma temperatura consoluta inferior (TL) a, aproximadamente, 46,5 ºC. Para temperaturas menores que 46,5 ºC, obtém-se uma única fase líquida para qualquer composição; para temperaturas maiores que 46,5 ºC é possível existir o equilíbrio líquido-líquido dependendo da composição. Para esse sistema não há a formação da “ilha” com temperaturas consolutas inferior e superior ao mesmo tempo.
O polinômio de quarta ordem gerado possui quatro raízes, duas reais e duas imaginárias, sendo que apenas as raízes reais possuem significado físico. Elas representam as frações do componente 1 nas fases  e . A fase  é a fase rica em água e a  é a fase rica em EGME. 
Através da função “roots” no scilab, foram obtidas as raízes deste polinômio. Na tabela a seguir, são apresentadas as raízes x1α e x1β. 
	Temperatura (°C)
	x1
	x1
	62,2
	0,0155055
	0,1586215
	55,7
	0,0227688
	0,1430637
	48,9
	0,0355503
	0,0874192
	46,5
	0,0462508
	0,0656411
	49,0
	0,0352627
	0,0882965
	53,6
	0,0258087
	0,1332625
	63,0
	0,0147580
	0,1599303
Tabela 6: composições calculadas a partir do ajuste de quarta ordem.
Utilizando o modelo NRTL e os parâmetros ajustados da literatura (DECHEMA) foram calculadas as composições de cada fase a fim de comparar com os resultados experimentais.
	Temperatura (°C)
	x1
	x1
	62,2
	0.0172288 
	0.1464773 
	55,7
	0.0197113 
	0.1256138 
	48,9
	 0.0249302 
	0.0976440 
	46,5
	0.0283258 
	 0.0855036 
	49,0
	0.0248160 
	 0.0981153 
	53,6
	0.0209119 
	0.1177385 
	63,0
	0.0170137 
	0.1486987 
Tabela 7: composições calculadas a partir do modelo NRTL. 
Gráfico 4: Comparação entre modelado e experimental para as mesmas temperaturas.
Podemos notar que os dados experimentais em x1 estão coerentes com os apresentados pelo modelo utilizando os parâmetros DECHEMA, contudo, para os valores de x1 ocorreram alguns desvios, distorcendo a binodal experimental. Isso pode ser explicado pelo fato de que ocorreram possíveis erros durante a realização do experimento (medição da temperatura e momento exato de observação das duas fases formadas). Além disso, a obtenção de mais pontos para o parâmetro “b” poderia fornecer um ajuste melhor e eficaz que garantisse maior qualidade ao programa computacional e aos resultados.
5.2. Sensibilidade dos parâmetros 
No modelo NRTL os parâmetros que possuem dependência com a temperatura são os parâmetros  e b e a sensibilidade dos mesmos pode ser testada através de pequenas variações previamente estabelecidas, analisando-se posteriormente as alterações nos valores de x1α e x1β. Para o caso “α”, seu valor poderia ser alterado para o limite máximo recomendado, ou seja, de 0,45. Já o parâmetro “b” poderia ser variado em uma proporção de +10% e -10% em seu valor final.
Entretanto, não foi possível realizar esses cálculos uma vez que o procedimento computacional programado no SCILAB caminhou para uma solução trivial.
5.3. Equilíbrio Líquido-Líquido-Vapor
O sistema EGME/Água ao ser aquecido atinge um estágio onde há o aparecimento das primeiras bolhas, confirmando a entrada em uma região onde há a coexistência de 2 fases líquidas em equilíbrio com sua fase vapor, isto é, o Equilíbrio Líquido – Líquido – Vapor (ELLV). 
No experimento realizado, foi-se avaliado a temperatura na qual o sistema em estudo entrou na região de ELLV. A temperatura aferida experimentalmente, ao se verificar o aparecimento das primeiras bolhas, se encontra na tabela abaixo.
Consideramos que, segundo a regra das fases, um problema de ELLV apresenta um grau de liberdade. No entanto, o mesmo foi realizado a pressão constante (atmosférica), de forma a tornar o ELLV um problema invariante. 
Para a realização dessa previsão, a mistura binária foi estudada até o ponto no qual o mesmo entrava na região do ELLV, ou seja, o estágio a partir do qual há o aparecimento das primeiras bolhas. Assim, medimos duas temperaturas distintas durante o experimento: uma para o aquecimento (93ºC) e outra para o resfriamento (92,5ºC); ambas associadas ao aparecimento das primeiras bolhas. Os dados coletados se encontram presentem na tabela abaixo.
	Volume EGME (mL)
	Volume de Água (mL)
	Tmédia (°C)
	Fração volumétrica
EGME (ɸ)
	Fração volumétrica
Água (ɸ)
	Fração molar 
EGME
	Fração molar
Água
	55
	40
	92,75
	0,57894737
	0,42105263
	0,15998
	0,840022
Tabela 8: Dados do experimento do ELLV
Os valores da literatura encontrados para “b” foram novamente ajustados, desta vez utilizando-se uma regressão linear, uma vez que se verificou que o ajuste anteriormente feito não era adequado para previsão. Dessa forma, os novos gráficos obtidos podem ser observados abaixo. 
Gráficos 5 e 6 – Ajuste de primeira ordem dos dados da literatura (Dechema Liquid-Liquid Equilibrium.
	
	Temperatura ELLV
	x1 alfa
	x1 beta
	Composição do vapor (y1)
	Modelo
	94,09
	0.1776151
	0.1776151
	0.0052531
	Experimental
	92,75
	0.0044294
	0.1756074Não determinado
	DECHEMA
	90,89
	0.01545
	0.18030
	0.02689
Tabela 9: Comparação entre as frações molares na fase líquida (x1) e vapor (y1) experimentais do ELL, do modelo do ELLV e DECHEMA.
Gráfico 7: Comparação dos resultados do ELL e do ELLV.
Ao implementar os parâmetros na rotina de cálculo com a finalidade de se encontrar a composição em cada fase não foram obtidos bons resultados uma vez que os valores calculados para x1α e x1β foram iguais, mostrando que a modelagem caminhou para uma solução trivial. Com isso, conclui-se que os parâmetros utilizados não são indicados para prever este equilíbrio.
5.4. Equilíbrio Líquido-Vapor
O Equilíbrio Líquido – Vapor (ELV) está relacionado com o movimento relativo das moléculas em relação à interface entre a fase líquida e a fase vapor. As moléculas possuem um movimento aleatório (browniano), de modo que na interface líquido – vapor há moléculas que atravessam a mesma, indo da fase vapor para a fase líquida e vice-versa. Quando a taxa de moléculas que passam da fase líquida para a fase vapor se iguala a taxa das moléculas que passam da fase vapor para a fase líquida, entramos num estado de equilíbrio termodinâmico – o ELV. Para que o comportamento do ELV possa ser observado, é necessário um cálculo de BOL T.
	Temperatura (°C)
	Composição x1
	Composição y1
	166.67
	0.96
	0.9216516
Tabela 10: Resultados do cálculo de BOL T.
Conclusão
O experimento, quando realizado de forma correta, permite a determinação precisa das curvas de equilíbrio líquido-líquido entre misturas parcialmente miscíveis.
O modelo NRTL se mostrou uma boa ferramenta de modelagem termodinâmica do ELL e ELV, proporcionando a previsão desses equilíbrios de certa forma satisfatória aos dados obtidos da literatura (Dechema: liquid - liquid equilibrium), porém não foi obtido sucesso para a modelagem do ELLV, pois a esta caminhou para a obtenção da solução trivial. Este erro também ocorreu quando tentamos fazer o teste de sensibilidade dos parâmetros.
Alguns valores, gráficos e resultados poderiam ser melhorados e se tornarem mais coerentes mediante uma minimização dos erros durante o experimento (medição de temperatura, visualização das duas fases) e com valores tabelados mais acurados que permitissem uma regressão polinomial mais eficaz, com mais pontos a serem analisados
Referências Bibliográficas
W. ARLT, M.E.A. MACEDO, P. RASMUSSEN, J.M. SORENSEN. DECHEMA Chemistry Data Series, Volume V - Liquid-Liquid Equilibrium Data Collection. Parts 2 and 3 together ,1.290 pages. 1979.
SMITH, J. M.; VAN NESS, H. C.; ABBOTT, M. M. Introdução a termodinâmica da engenharia química ed.5 Rio de Janeiro: LTC - Livros Técnicos e Científicos Editora S. A., 2000;
Kackbart, L.M.: EQUILÍBRIO LÍQUIDO-LÍQUIDO DE SISTEMAS CONTENDO FENOL-ÁGUA-SOLVENTE: OBTENÇÃO E MODELAGEM TERMODINÂMICA DE LUCIANA MEIRELES. UNIVERSIDADE FEDERAL DO PARANÁ; 2007.
Apêndice I – Rotina de Cálculos
Rotina ELL:
//Componente 1 – EGME
//Componente 2 – Água
clc, clear
exec('functions_NRTL_2.sci');
//Dados de entrada.
alfa = 0.2;
R = 1.987;
coefb12= [-0.1177; 18.181; -1707.3]; 
coefb21= [0.0831; -8.5807; 2793.1]; 
//Chutes iniciais para as composições das fases alfa e beta
x1alfa = 0.018;
x1beta =0.18;
//Determinação das composições x1 na fase alfa e x1 na fase beta pelo modelo NRTL
for T = [335.35 328.85 322.05 319.65 322.15 326.75 336.15];
[x1alfabeta] = fsolve([x1alfa;x1beta],list(FobjELLB, T, alfa, coefb12, coefb21))
disp (T,"Temperatura= ")
disp(x1alfabeta(1),"x1 na fase alfa= ");
disp(x1alfabeta(2),"x1 na fase beta= ");
disp("*****************************************************")
end
Functions:
//Determinação dos coeficientes de atividade (gama).
function [g]=gama(x, T, alfa, coefb12, coefb21)
T1 = T-273.15 //Temperatura em °C
b12 = coefb12(1)*T1^2 + coefb12(2)*T1 + coefb12(3);
b21 = coefb21(1)*T1^2 + coefb21(2)*T1 + coefb21(3);
tal12 = b12/T;
tal21 = b21/T;
G12 = exp(-1*alfa*tal12);
G21 = exp(-1*alfa*tal21);
g =exp([(x(2)^2*(tal21*(G21/(x(1)+x(2)*G21))^2+((G12*tal12)/(x(2)+x(1)*G12)^2)));
(x(1)^2*(tal12*(G12/(x(2)+x(1)*G12))^2+((G21*tal21)/(x(1)+x(2)*G21)^2)))]);
endfunction
//Função Objetivo - Calcula as compoosições no ELL.
function F=FobjELLB(v, T, alfa, coefb12, coefb21)
xalfa=[v(1);1-v(1)];
xbeta=[v(2);1-v(2)];
galfa = gama(xalfa,T, alfa, coefb12, coefb21);
gbeta = gama(xbeta, T, alfa, coefb12, coefb21);
F = [xalfa(1)*galfa(1)-xbeta(1)*gbeta(1);xalfa(2)*galfa(2)-xbeta(2)*gbeta(2)];
endfunction
Rotina ELLV:
clc, clear 
exec('NRTL_2.sci'); 
//Dados de entrada // 
T = 365.9; 
T1 = T-273.15; 
alfa = 0.2; 
R = 1.987; 
cb12= [0;-5.0358;-610.69]; 
cb21= [0;16.755;1326.1]; 
P = 101.325; 
// Antoine EGME (Psat=atm e T=Kelvin) // 
A1 = 4.26; 
B1 = 1511.4; 
C1 = -88.2; 
// Antoine Água (Psat=KPa e T=Celsius) // 
A2 = 16.8; 
B2 = 3885.7; 
C2 = 231; 
//Chutes iniciais // 
x1alfa = 0.01; 
x1beta =0.1; 
y1 = 0.2; 
//Composições das fases - modelo NRTL 
[x1alfabeta]= fsolve([x1alfa;x1beta;y1;T],list(FOELLV, alfa, cb12, cb21)) 
disp(x1alfabeta(4),"Temp ELLV"); 
disp(x1alfabeta(1),"x1 alfa"); 
disp(x1alfabeta(2),"x1 beta"); 
disp(x1alfabeta(3),"y1");
Functions:
//Coeficientes de atividade . 
function [g]=gama(x, T, alfa, cb12, cb21) 
 T1 = T-273.15 //Temperatura em °C 
 b12 = cb12(1)*T1^2 + cb12(2)*T1 + cb12(3); 
 b21 = cb21(1)*T1^2 + cb21(2)*T1 + cb21(3); 
 tal12 = b12/T; 
 tal21 = b21/T; 
 G12 = exp(-1*alfa*tal12); 
 G21 = exp(-1*alfa*tal21); 
 g = exp([(x(2)^2*(tal21*(G21/(x(1)+x(2)*G21))^2+((G12*tal12)/(x(2)+x(1)*G12)^2))); (x(1)^2*(tal12*(G12/(x(2)+x(1)*G12))^2+((G21*tal21)/(x(1)+x(2)*G21)^2)))]); 
 endfunction 
//Equação de Antoine. 
function P=antoine(As, B, C, T) 
 P = exp(As - (B/(T+C))); 
 endfunction 
//Rotina do ELLV 
function F=FOELLV(v, alfa, cb12, cb21) 
 xalfa=[v(1);1-v(1)]; 
 xbeta=[v(2);1-v(2)]; 
 y1 = v(3); 
 y2 = 1-v(3); 
 T = v(4); 
 galfa = gama(xalfa, T, alfa, cb12, cb21); 
 gbeta = gama(xbeta, T, alfa, cb12, cb21); 
 P1sat = 101.325 * antoine(A1, B1, C1, T); 
 P2sat = antoine(A2, B2, C2, T1); 
 F = [xalfa(1)*galfa(1)-xbeta(1)*gbeta(1);xalfa(2)*galfa(2)-xbeta(2)*gbeta(2);xalfa(1)*galfa(1)*P1sat-y1*P;xalfa(2)*galfa(2)*P2sat-y2*P]; 
 endfunction
Rotina ELV:
clc, clear 
exec('BOLT.sci'); 
//Dados de entrada 
P= 101.325; 
alfa = 0.2; 
cb12= [0;-5.0358;-610.69]; 
cb21= [0;16.755;1326.1]; 
tol = 0.0001; 
// Antoine EGME (Psat=atm e T=Kelvin) // 
A1 = 4.26; 
B1 = 1511.4; 
C1 = -88.2; 
// Antoine Água (Psat=KPa e T=Celsius) // 
A2 = 16.8; 
B2 = 3885.7; 
C2 = 231; 
//Processo iterativo 
for x1 = 0.01:0.05:0.99; 
 x2 = 1-x1; 
 T1sat = ((B1/(A1-log(1)))- C1)-273.15; 
 T2sat =((B2/(A2-log(P)))- C2); 
 Tmi = (T1sat+T2sat)/2; 
 Tm = 1000; 
 N=1000; 
 
while (N>0) 
 P1sat = 101.325*(exp(A1-(B1/(C1+(Tmi+273.15))))); 
 P2sat = (exp(A2-(B2/(C2+Tmi)))); 
 alfa1 = P1sat/P2sat; 
 
g = gama([x1;x2], (Tmi+273.15), alfa, cb12, cb21) 
P2sat = P/(x1*g(1)*alfa1+x2*g(2)); 
Tm = ((B2/(A2-log(P2sat)))- C2); 
if abs(Tm-Tmi)<tol then 
 disp(10000-N,"Número de iterações = "); 
 y1 = 1-((x2*g(2)*P2sat)/P); 
 disp([Tm, x1, y1],"[T bolha; x1; y1] = "); 
 N=0 
else 
 N=N-1; 
 Tmi=Tm 
end; 
end; 
end
Functions:
// Coeficiente de atividade - modelo NRTL 
function [g]=gama(x, T, alfa, cb12, cb21) 
 T1 = T-273.15 //Temperatura em °C 
 b12 = cb12(1)*T1^2 + cb12(2)*T1 + cb12(3); 
 b21 = cb21(1)*T1^2 + cb21(2)*T1 + cb21(3); 
 tal12 = b12/T; 
 tal21 = b21/T; 
 G12 = exp(-1*alfa*tal12); 
 G21 = exp(-1*alfa*tal21); 
 g = exp([(x(2)^2*(tal21*(G21/(x(1)+x(2)*G21))^2+((G12*tal12)/(x(2)+x(1)*G12)^2))); (x(1)^2*(tal12*(G12/(x(2)+x(1)*G12))^2+((G21*tal21)/(x(1)+x(2)*G21)^2)))]); 
endfunction 
//Equação de Antoine. 
function P=antoine(As, B, C, T) 
 P = exp(As- (B/(T+C))); 
endfunction