*****************************************************************
*** Chapter 6
*** Test Problem 1
***
*** n-Butyl-Acetate and Water -- Gibbs energy minimization(NRTL)
*****************************************************************

SETS
    i     components                /1*2/
    k     phases                    /1*2/
    alias(i,j);

*****************************************************************
* P = pressure (atm)
* T = temperature (K)
* liqphase(k) = is phase k a liquid phase?
* ntot(i) = number of moles of component i in feed
* alpha(i,j) = binary interaction parameter (symmetric)
* tau(i,j) = binary interaction parameter (nonsymmetric)
* g(i,j) = binary interaction parameter (calculated)

PARAMETERS
P, T, liqphase(k), ntot(i), alpha(i,j), tau(i,j), g(i,j);

P = 1.0;
T = 298;

liqphase('1') = 1; liqphase('2') = 1;

ntot('1') = 0.5;
ntot('2') = 0.5;

alpha('1','1') = 0.0; alpha('1','2') = 0.391965;
alpha('2','1') = 0.391965; alpha('2','2') = 0.0;

tau('1','1') = 0.0; tau('1','2') = 3.00498;
tau('2','1') = 4.69071; tau('2','2') = 0.0;

g(i,j) = EXP(-alpha(i,j)*tau(i,j));

***************************************************************************
* gfe = Gibbs free energy
* n(i,k) = number of moles of component i in phase k
* psi(i,k) = substitution variable

VARIABLES gfe, n(i,k), psi(i,k);

***************************************************************************
* obj = objective function
* trans(i,k) = transformation constraint
* molesum(i) = material balance on component i

EQUATIONS
    obj
    trans(i,k)
    molesum(i);


obj.. gfe =e= SUM(k, SUM(i, n(i,k)*(LOG(n(i,k)) - LOG(SUM(j,n(j,k))))))
	     +SUM(k$liqphase(k), SUM(i, n(i,k)*SUM(j, g(i,j)*tau(i,j)*psi(j,k))));

trans(i,k)$liqphase(k).. psi(i,k)*SUM(j, g(j,i)*n(j,k)) - n(i,k) =e= 0.0;

molesum(i).. SUM(k, n(i,k)) =e= ntot(i);

MODEL gmin / all /;

n.lo(i,k) = 0.0000001; n.up(i,k) = ntot(i);
n.l('1','1') = 0.4993; n.l('1','2') = 0.0007;
n.l('2','1') = 0.3441; n.l('2','2') = 0.1559;
psi.lo(i,k) = 0.0;
psi.l(i,k) = n.l(i,k)/SUM(j, g(j,i)*n.l(j,k));

SOLVE gmin USING nlp MINIMIZING gfe;