RCM Computation using UNIQUAC Model
Case of a Quaternary Reactive System : Acetic Acid - Ethanol - Water- Ethyl Acetate​
​Authors: Housam Binous and Ahmed Bellagi

The esterification reaction between ethanol and acetic acid produces water and ethyl acetate:
​
CH
3
COOH+
CH
3
CH
2
OH⇌
H
2
O+
CH
3
COOCH
2
CH
3
.
This reaction extent is limited by the equilibrium constant, taken to be 13.4. In the liquid phase esterification occurs. The residue curves are computed using the UNIQUAC model and gas-phase acetic acid dimerization.
We plot three residue curves. There are two nonreactive azeotropes in the system: (1) azeotrope between ethanol and ethyl acetate labeled A (an unstable node with a low boiling temperature show in orange); (2) azeotrope between water and ethanol labeled B (a saddle point with an intermediate boiling temperature). Finally, pure water and ethyl acetate are two saddle points and pure ethanol is a stable node.
The transformed compositions of the two azeotropes are the following:
A:
X
1
=0.5336 and
X
2
=1
B:
X
1
=0 and
X
2
=0.9070
​
​
Reference : H. Binous, Residue Curve Map for Homogeneous Reactive Quaternary Mixtures, Computer Applications in Engineering Education Journal, 15(1), pp. 73-77, 2007.
In[]:=
Off[General::"spell1"]
Vapor pressure and normal boiling points
In[]:=
Psat1=10^3Exp[A1+B1/(C1+T)+D1Log[T]+E1T^F1]
Out[]=
1000
A1+E1
F1
T
+
B1
C1+T
+D1Log[T]

In[]:=
Psat2=10^3Exp[A2+B2/(C2+T)+D2Log[T]+E2T^F2]
Out[]=
1000
A2+E2
F2
T
+
B2
C2+T
+D2Log[T]

In[]:=
Psat3=10^3Exp[A3+B3/(C3+T)+D3Log[T]+E3T^F3]
Out[]=
1000
A3+E3
F3
T
+
B3
C3+T
+D3Log[T]

In[]:=
Psat4=10^3Exp[A4+B4/(C4+T)+D4Log[T]+E4T^F4]
Out[]=
1000
A4+E4
F4
T
+
B4
C4+T
+D4Log[T]

In[]:=
A1=61.34090;B1=-6768.88;C1=0;D1=-6.72663;E1=4.8410^-6;F1=2;
In[]:=
A2=86.4860;B2=-7931.10;C2=0;D2=-10.24980;E2=6.3910^-6;F2=2;
In[]:=
A3=88.3760;B3=-7147.86;C3=0;D3=-10.99170;E3=8.5510^-6;F3=2;
In[]:=
A4=65.92780;B4=-7227.53;C4=0;D4=-7.17695;E4=4.0310^-6;F4=2;
In[]:=
FindRoot[Psat11.0132510^5,{T,380}]
Out[]=
{T390.951}
In[]:=
390.951-273.15
Out[]=
117.801
In[]:=
FindRoot[Psat21.0132510^5,{T,380}]
Out[]=
{T351.32}
In[]:=
351.32-273.15
Out[]=
78.17
In[]:=
FindRoot[Psat31.0132510^5,{T,380}]
Out[]=
{T350.184}
In[]:=
350.184-273.15
Out[]=
77.034
In[]:=
FindRoot[Psat41.0132510^5,{T,380}]
Out[]=
{T373.153}
In[]:=
373.153-273.15
Out[]=
100.003
Activity coefficients using UNIQUAC model
Component 1 = Acetic Acid
Component 2 = Ethanol
Component 3 = Ethyl Acetate
Component 4 = Water
In[]:=
tau[j_,i_]:=Exp[-(U[j,i]-U[i,i])/(RT)]
In[]:=
U[1,1]=0.0;U[1,2]=122.3620;U[1,3]=-373.4015;U[1,4]=427.7405;
In[]:=
U[2,1]=-162.7242;U[2,2]=0.0;U[2,3]=-62.7027;U[2,4]=50.8846;
In[]:=
U[3,1]=690.2536;U[3,2]=363.9664;U[3,3]=0.0;U[3,4]=757.5853;
In[]:=
U[4,1]=-305.4517;U[4,2]=232.0091;U[4,3]=151.6420;U[4,4]=0.0;
In[]:=
Gam1=Exp[Log[Gam1C]+Log[Gam1R]]
Out[]=
Gam1CGam1R
In[]:=
Gam2=Exp[Log[Gam2C]+Log[Gam2R]]
Out[]=
Gam2CGam2R
In[]:=
Gam3=Exp[Log[Gam3C]+Log[Gam3R]]
Out[]=
Gam3CGam3R
In[]:=
Gam4=Exp[Log[Gam4C]+Log[Gam4R]]
Out[]=
Gam4CGam4R
In[]:=
Gam1C=ExpLog[phi[1]/x[1]]+z2q[1]Log[v[1]/phi[1]]+l[1]-phi[1]x[1]
4
∑
j=1
x[j]l[j];
In[]:=
Gam2C=ExpLog[phi[2]/x[2]]+z2q[2]Log[v[2]/phi[2]]+l[2]-phi[2]x[2]
4
∑
j=1
x[j]l[j];
In[]:=
Gam3C=ExpLog[phi[3]/x[3]]+z2q[3]Log[v[3]/phi[3]]+l[3]-phi[3]x[3]
4
∑
j=1
x[j]l[j];
In[]:=
Gam4C=ExpLog[phi[4]/x[4]]+z2q[4]Log[v[4]/phi[4]]+l[4]-phi[4]x[4]
4
∑
j=1
x[j]l[j];
In[]:=
r[1]=2.20230;r[2]=2.10540;r[3]=3.47860;r[4]=0.920000;
In[]:=
q[1]=2.07200;q[2]=1.97200;q[3]=3.11590;q[4]=1.39970;
In[]:=
l[i_]:=z/2(r[i]-q[i])-(r[i]-1)
In[]:=
z=10;
In[]:=
Gam1R=Expq[1]1-Log
4
∑
j=1
v[j]tau[j,1]-
4
∑
j=1
v[j]tau[1,j]
4
∑
k=1
v[k]tau[k,j]
;
In[]:=
Gam2R=Expq[2]1-Log
4
∑
j=1
v[j]tau[j,2]-
4
∑
j=1
v[j]tau[2,j]
4
∑
k=1
v[k]tau[k,j]
;
In[]:=
Gam3R=Expq[3]1-Log
4
∑
j=1
v[j]tau[j,3]-
4
∑
j=1
v[j]tau[3,j]
4
∑
k=1
v[k]tau[k,j]
;
In[]:=
Gam4R=Expq[4]1-Log
4
∑
j=1
v[j]tau[j,4]-
4
∑
j=1
v[j]tau[4,j]
4
∑
k=1
v[k]tau[k,j]
;
In[]:=
phi[i_]:=r[i]x[i]
4
∑
j=1
r[j]x[j]
In[]:=
v[i_]:=q[i]x[i]
4
∑
j=1
q[j]x[j]
In[]:=
R=1.98721;
Dimerization constant for acetic acid
In[]:=
K=10^(-12.5454+3166/(T));
Total pressure
In[]:=
P=1.0132510^5;
Two azeotropes exits: Water-Ethanol and Ethanol-Ethyl Acetate
In[]:=
Quiet@FindRoot[{Psat1Gam1x[1]/((1+Sqrt[(1+4KPsat1)])/(1+Sqrt[(1+4KPy[1](2-y[1]))]))+Psat2Gam2x[2]/((2(1-y[1]+Sqrt[(1+4KPy[1](2-y[1]))]))/((2-y[1])(1+Sqrt[(1+4KPy[1](2-y[1]))])))+Psat3Gam3x[3]/((2(1-y[1]+Sqrt[(1+4KPy[1](2-y[1]))]))/((2-y[1])(1+Sqrt[(1+4KPy[1](2-y[1]))])))+Psat4Gam4x[4]/((2(1-y[1]+Sqrt[(1+4KPy[1](2-y[1]))]))/((2-y[1])(1+Sqrt[(1+4KPy[1](2-y[1]))])))P,​​y[1]Psat1Gam1x[1]/((1+Sqrt[(1+4KPsat1)])/(1+Sqrt[(1+4KPy[1](2-y[1]))]))1/P,​​x[1]0,x[4]0,x[2]==Psat2Gam2x[2]/((2(1-y[1]+Sqrt[(1+4KPy[1](2-y[1]))]))/((2-y[1])(1+Sqrt[(1+4KPy[1](2-y[1]))])))1/P,x[3]==Psat3Gam3x[3]/((2(1-y[1]+Sqrt[(1+4KPy[1](2-y[1]))]))/((2-y[1])(1+Sqrt[(1+4KPy[1](2-y[1]))])))1/P},{x[1],0},{x[2],0.5},{x[3],0.5},{x[4],0},{T,350},{y[1],0}]
Out[]=
{x[1]1.53942×
-24
10
,x[2]0.466059,x[3]0.533941,x[4]0.,T345.324,y[1]5.29563×
-25
10
}
In[]:=
345.324-273.15
Out[]=
72.174
In[]:=
FindRoot[{Psat1Gam1x[1]/((1+Sqrt[(1+4KPsat1)])/(1+Sqrt[(1+4KPy[1](2-y[1]))]))+Psat2Gam2x[2]/((2(1-y[1]+Sqrt[(1+4KPy[1](2-y[1]))]))/((2-y[1])(1+Sqrt[(1+4KPy[1](2-y[1]))])))+Psat3Gam3x[3]/((2(1-y[1]+Sqrt[(1+4KPy[1](2-y[1]))]))/((2-y[1])(1+Sqrt[(1+4KPy[1](2-y[1]))])))+Psat4Gam4x[4]/((2(1-y[1]+Sqrt[(1+4KPy[1](2-y[1]))]))/((2-y[1])(1+Sqrt[(1+4KPy[1](2-y[1]))])))P,​​y[1]Psat1Gam1x[1]/((1+Sqrt[(1+4KPsat1)])/(1+Sqrt[(1+4KPy[1](2-y[1]))]))1/P,​​x[1]0,x[3]0,x[2]==Psat2Gam2x[2]/((2(1-y[1]+Sqrt[(1+4KPy[1](2-y[1]))]))/((2-y[1])(1+Sqrt[(1+4KPy[1](2-y[1]))])))1/P,x[4]==Psat4Gam4x[4]/((2(1-y[1]+Sqrt[(1+4KPy[1](2-y[1]))]))/((2-y[1])(1+Sqrt[(1+4KPy[1](2-y[1]))])))1/P},{x[1],0},{x[2],0.5},{x[3],0},{x[4],0.5},{T,350},{y[1],0}]
Out[]=
{x[1]1.34226×
-26
10
,x[2]0.907561,x[3]-1.68043×
-26
10
,x[4]0.0924388,T351.197,y[1]-1.62406×
-27
10
}
Chemical reaction equilibrium constant