Isobar Vapor-Liquid Diagram of 
Acetone-Water Mixture at 760 mmHg 
Author's Data: Housam BINOUS
Department of Chemical Engineering
National Institute of Applied Sciences and Technology
Tunis, TUNISIA
Email: binoushousam@yahoo.com

Activity Coefficients based on NRTL model
Component 1: Acetone
Component 2: Water
Model parameters are taken from Distil by Hyprotech Ltd.

Gam[i_]:=Exp[(((∑_(j = 1)^2)( tau[j, i] G[j, i] x[j]))/((∑_(l = 1)^2) (G[l, i] x[l])))+(∑_(j = 1)^2) (((G[i, j] x[j])/((∑_(l = 1)^2)( G[l, j] x[l])))(tau[i, j]-(((∑_(r = 1)^2)( tau[r, j] G[r, j] x[r]))/((∑_(l = 1)^2) (G[l, j] x[l])))))]/.x[1] -> X1/.x [2] -> X2

G[i_, j_] := Exp[-a[i, j] tau[i, j]]

tau[i_, j_] := g[i, j]/(R  (T + 273.15))

g[1, 1] = 0; g[2, 1] = 1299.3950;
g[1, 2] = 750.31810; g[2, 2] = 0;

a[1, 1] = 0; a[2, 1] = 0.5856;
a[1, 2] = 0.5856; a[2, 2] = 0;

Gas Constant and Total Pressure

R = 1.987; P = 760;

Vapor Pressure based on Antoine Equation and Normal Boiling Points

Psat[i_] := 10^((Aa[i] - Ba[i]/(T + Ca[i])))

Aa[1] = 7.63130; Aa[2] = 8.07131;
Ba[1] = 1566.69; Ba[2] = 1730.63;
Ca[1] = 273.419; Ca[2] = 233.426;

FindRoot[Psat[1] ==  760, {T, 100}]

FindRoot[Psat[2] ==  760, {T, 100}]

While loop to get bubble temperature and vapor mole fraction

i = 0; 
    While[i < 101, {xe[i] = i 0.01, X1 = xe[i], X2 = 1 - xe[i], 
        T[i] = Chop[NMinimize[(P -  Psat[1] Gam[1] X1 - 
                      Psat[2] Gam[2] (1 - X1))^2, {T}][[2, 1, 
                2]]], ye[i] = Psat[1] Gam[1] xe[i]/P /. T ->  T[i], 
        Print[i, " ", T[i], " ", xe[i], " ", ye[i], " ", 
          Gam[1] /. T ->  T[i], " ", 
          Gam[2] /. T ->  T[i]], i++}]

Ploting the VLE data

tb = Table[{xe[i], ye[i]}, {i, 0, 100}];

plt1 = 
      ListPlot[tb, PlotStyle -> RGBColor[1, 0, 0], 
        PlotJoined ->  True, PlotRange -> All, 
        Epilog ->  {RGBColor[0, 1, 0], 
            Line[{{0, 0}, {1, 1}}]}]

tb2 = Table[{xe[i], T[i]}, {i, 0, 100}];

tb3 = Table[{ye[i], T[i]}, {i, 0, 100}];

plt1 = 
      ListPlot[tb2, PlotStyle -> RGBColor[1, 0, 0], 
        PlotJoined ->  True, PlotRange -> All]

plt2 = 
      ListPlot[tb3, PlotStyle -> RGBColor[0, 0, 1], 
        PlotJoined ->  True, PlotRange -> All]

Show[plt1, plt2]

Plotting excess enthalpy and excess Gibbs free energy at 298.15K

X1 =. ; X2 = 1 - X1;

GEX[X1_] := 
      X1 R (273.15 + T) Log[
            Gam[1]] + (1 - X1) R (273.15 + T) Log[
            Gam[2]]

Plot[GEX[X1] /. T ->  25, {X1, 0, 1}, 
      PlotRange ->  All, 
PlotStyle ->  
        RGBColor[0, 0, 1]]

HEX[
            X1_] := -X1 R (273.15 +  T)^2 D[
                Log[X1 Gam[1]], {T, 1}] - (1 - 
                  X1) R (273.15 +  T)^2 D[
                Log[(1 - X1) Gam[2]], {T, 1}];;

Plot[HEX[X1] /. T ->  25, {X1, 0, 1}, 
      PlotRange ->  All, 
PlotStyle ->  
        RGBColor[1, 0, 0]]