H-xy Diagram for Ethanol-Water 
at 760mmHg and 76 mmHg
Author's Data: Housam BINOUS
Department of Chemical Engineering
National Institute of Applied Sciences and Technology
Tunis, TUNISIA
Email: binoushousam@yahoo.com

Supression error and warning messages

Off[General::spell1]

Off[FindRoot::lstol]

Vapor Pressures, Normal boiling points and total pressure = 760 mmHg 
Component 1= ethanol
Component 2= Water

A2 = 8.07131; B2 = 1730.630; C2 = 233.426;
A1 = 8.11220; B1 = 1592.864; C1 = 226.184;

PS1 = 10^(A1 - B1/(C1 + T));
PS2 = 10^(A2 - B2/(C2 + T));

FindRoot[PS1 ==  760,  {T, 100}]

FindRoot[PS2 ==  760,  {T, 100}]

P = 760;

Activity coefficients using the Van Larr model

G1 = 
      Exp[A12 (A21 (1 - x)/(A12 x + A21 (1 - x)))^2] /. 
        x -> x[i]
G2 = 
      Exp[A21 (A12 x/(A12 x + A21 (1 - x)))^2] /. 
        x -> x[i]

A21 = 0.9227; A12 = 1.6798;

Vapor-liquid Equilibrium curve

h = 0.01; i = 0;

Do[{x[i] = h i, 
        T[i] = FindRoot[
              P ==  PS1 G1 x[i] + PS2 G2 (1 - x[i]), {T, 
                100}][[1, 2]], 
        y[i] =  PS1 G1 x[i]/P /. T -> T[i]}, {i, 0, 100}]

tbl1 = Table[{x[i], y[i]}, {i, 0, 100}];

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

Saturated liquid and saturated vapor Enthalpy. h-x and H-y curves.

M1 = 46.07; M2 = 18.016;

cl2[T_] := 
        8.314 (8.712  + 1.25 10^-3 T - 
              1.8 10^-7 T^2);
cl1[T_] := 
        8.314 (33.866 - 1.7260 10^-1 T + 
              3.4917 10^-4 T^2);

Do[{hl[
            i] = (x[i] 
                    Integrate[
                      cl1[T], {T, 273.15, T[i] + 273.15}] + (1 - x[i]) 
                    Integrate[cl2[T], {T, 273.15, T[i] + 273.15}] )/ (x[
                      i] M1 + (1 - x[i]) M2)  0.43028, 
        hlm[i] = (x[i] M1/(x[i] M1 + (1 - x[i]) M2) 
                    Integrate[
                      cl1[T], {T, 273.15, T[i] + 273.15}] + (1 - 
                        x[i]) M2/(x[i] M1 + (1 - x[i]) M2) 
                    Integrate[
                      cl2[T], {T, 273.15, T[i] + 273.15}] ) / (x[
                      i] M1 + (1 - x[i]) M2) 0.43028}, {i, 0, 
        100}]

tbl2 = Table[{x[i], hl[i]}, {i, 0, 100}];

tbl2bis = 
        Table[{x[i] M1/(x[i] M1 + (1 - x[i]) M2), hlm[i]}, {i, 0, 
            100}];

G1 = 
        Exp[A12 (A21 (1 - x)/(A12 x + A21 (1 - x)))^2] /. 
          x -> xh;
G2 = 
        Exp[A21 (A12 x/(A12 x + A21 (1 - x)))^2] /. 
          x -> xh;

Do[{yh[i] = h i, 
        Th[i] = FindRoot[{P ==  PS1 G1 xh + PS2 G2 (1 - xh), 
                yh[i] ==   PS1 G1 xh/P}, {T, 100}, {xh, 
                0.5}][[1, 2]]}, {i, 0, 100}]

tbl3 = Table[{yh[i], Th[i]}, {i, 0, 100}];

CG1[T_] := 
        9.014 + 2.141 10^-1 T - 8.39 10^-5 T^2 + 
          1.373 10^-9 T^3;
CG2[T_] := 
        32.24 + 1.924 10^-3 T + 1.055  10^-5 T^2 - 
          3.596 10^-9 T^3;

dh1 = 
      1000 38.58 ((516.3 - 273.15)/(516.3 - 78 - 
                  273.15))^0.38
dh2 = 
      1000 40.656 ((647.4 - 273.15)/(647.4 - 100 - 
                  273.15))^0.38

Do[{hv[
            i] = (yh[i] dh1 + (1 - yh[i]) 
                  dh2 + (yh[i] 
                      Integrate[
                        CG1[T], {T, 273.15, Th[i] + 273.15}] + (1 - 
                          yh[i]) 
                      Integrate[CG2[T], {T, 273.15, Th[i] + 273.15}]) ) 
            0.43028, 
        hvm[i] = (yh[i] 
                    M1/(yh[i] M1 + (1 - yh[i]) M2) dh1 + (1 - 
                        yh[i]) M2/(yh[i] M1 + (1 - yh[i]) M2) 
                    dh2 + yh[i] M1/(yh[i] M1 + (1 - yh[i]) M2) 
                    Integrate[
                      CG1[T], {T, 273.15, Th[i] + 273.15}] + (1 - yh[i]) 
                    M2/(yh[i] M1 + (1 - yh[i]) M2) 
                    Integrate[
                      CG2[T], {T, 273.15, Th[i] + 273.15}] )/ (x[
                      i] M1 + (1 - x[i]) M2)  0.43028}, {i, 0, 
        100}]

tbl4 = Table[{yh[i], hv[i]}, {i, 0, 100}];

tbl4bis = 
        Table[{yh[i] M1/(yh[i] M1 + (1 - yh[i]) M2), hvm[i]}, {i, 
            0, 100}];

plt2 = 
      ListPlot[tbl2, PlotJoined -> True, 
        PlotStyle -> RGBColor[1, 0, 0], 
        DisplayFunction -> Identity]

plt2bis = 
      ListPlot[tbl2bis, PlotJoined -> True, 
        PlotStyle -> RGBColor[1, 0, 0], 
        DisplayFunction -> Identity]

plt3 = 
      ListPlot[tbl4, PlotJoined -> True, 
        PlotStyle -> RGBColor[0, 0, 1], 
        DisplayFunction -> Identity]

plt3bis = 
      ListPlot[tbl4bis, PlotJoined -> True, 
        PlotStyle -> RGBColor[0, 0, 1], 
        DisplayFunction -> Identity]

vapor = Interpolation[tbl4, InterpolationOrder ->  1];

vaporbis = 
        Interpolation[tbl4bis, InterpolationOrder ->  1];

liquid = Interpolation[tbl2, InterpolationOrder ->  1];

liquidbis = 
        Interpolation[tbl2bis, InterpolationOrder ->  1];

tbl5 = 
        Table[Line[{{x[i], liquid[x[i]]}, {y[i], vapor[y[i]]}}], {i, 2, 90, 
            10}];

tbl5bis = 
        Table[Line[{{x[i] M1/(x[i] M1 + (1 - x[i]) M2), 
                liquidbis[
                  x[i] M1/(x[i] M1 + (1 - x[i]) M2)]}, {y[
                    i] M1/(y[i] M1 + (1 - y[i]) M2), 
                vaporbis[
                  y[i] M1/(y[i] M1 + (1 - y[i]) M2)]}}], {i, 2, 
            90, 10}];

Show[plt3, plt2, Graphics[{RGBColor[0, 1, 0], tbl5}], 
      DisplayFunction -> $DisplayFunction]

Show[plt3bis, plt2bis, Graphics[{RGBColor[0, 1, 0], tbl5bis}], 
      DisplayFunction -> $DisplayFunction]

tbl6 = Table[{x[i], vapor[y[i]]}, {i, 0, 100}];

tbl6bis = 
        Table[{x[i] M1/(x[i] M1 + (1 - x[i]) M2), 
            vaporbis[y[i] M1/(y[i] M1 + (1 - y[i]) M2)]}, {i, 0, 
            100}];

plt4 = 
      ListPlot[tbl6, PlotStyle -> RGBColor[1, 0, 1], 
        PlotJoined -> True, DisplayFunction -> Identity]

plt4bis = 
      ListPlot[tbl6bis, PlotStyle -> RGBColor[1, 0, 1], 
        PlotJoined -> True, DisplayFunction -> Identity]

plt5 = 
      Show[plt3, plt2, plt4, Graphics[{RGBColor[0, 1, 0], tbl5}], 
        PlotRange -> All, AspectRatio -> 1.5, 
        Epilog -> {Line[{{x[20], liquid[x[20]]}, {x[20], 
                  vapor[y[20]]}}], 
Line[{{x[20], 
                  vapor[y[20]]}, {y[20], vapor[y[20]]}}]}, 
        DisplayFunction -> $DisplayFunction]

plt5bis = 
      Show[plt3bis, plt2bis, plt4bis, Graphics[{RGBColor[0, 1, 0], tbl5bis}], 
        PlotRange -> {-100, 1200}, AspectRatio -> 2, 
        Epilog -> {Line[{{x[
                      17] M1/(x[17] M1 + (1 - x[17]) M2), 
                  liquidbis[
                    x[17] M1/(x[17] M1 + (1 - x[17]) M2)]}, {x[
                      17] M1/(x[17] M1 + (1 - x[17]) M2), 
                  vaporbis[
                    y[17] M1/(y[
                                17] M1 + (1 - 
                                  y[17]) M2)]}}], 
Line[{{x[17] M1/(x[17] M1 + (1 - x[17]) M2), 
                  vaporbis[
                    y[17] M1/(y[17] M1 + (1 - y[17]) M2)]}, {y[
                      17] M1/(y[17] M1 + (1 - y[17]) M2), 
                  vaporbis[
                    y[17] M1/(y[17] M1 + (1 - y[17]) M2)]}}]}, 
        DisplayFunction -> $DisplayFunction]

Total pressure 76 mmHg
Component 1= ethanol
Component 2= Water

P = 76;

i =. ;

Activity coefficients using the Van Larr model

G1 = 
        Exp[A12 (A21 (1 - x)/(A12 x + A21 (1 - x)))^2] /. 
          x -> x[i];
G2 = 
        Exp[A21 (A12 x/(A12 x + A21 (1 - x)))^2] /. 
          x -> x[i];

A21 = 0.9227; A12 = 1.6798;

Vapor-liquid Equilibrium curve

h = 0.01; i = 0;

Do[{x[i] = h i, 
        T[i] = FindRoot[
              P ==  PS1 G1 x[i] + PS2 G2 (1 - x[i]), {T, 
                90}][[1, 2]], 
        y[i] =  PS1 G1 x[i]/P /. T -> T[i]}, {i, 0, 100}]

Tbl1 = Table[{x[i], y[i]}, {i, 0, 100}];

Plt1 = 
      ListPlot[Tbl1, PlotJoined -> True, 
        PlotRange -> {{0, 1}, {0, 1}}, 
        PlotStyle -> RGBColor[1, 0, 0], 
        Epilog -> {RGBColor[0, 1, 0], 
            Line[{{0, 0}, {1, 1}}]}]

Saturated liquid and saturated vapor Enthalpy. h-x and H-y curves.

Do[{
hl[
            i] = (x[i] 
                    Integrate[
                      cl1[T], {T, 273.15, T[i] + 273.15}] + (1 - x[i]) 
                    Integrate[cl2[T], {T, 273.15, T[i] + 273.15}] )/ (x[
                      i] M1 + (1 - x[i]) M2)  0.43028, 
        hlm[i] = (x[i] M1/(x[i] M1 + (1 - x[i]) M2) 
                    Integrate[
                      cl1[T], {T, 273.15, T[i] + 273.15}] + (1 - 
                        x[i]) M2/(x[i] M1 + (1 - x[i]) M2) 
                    Integrate[
                      cl2[T], {T, 273.15, T[i] + 273.15}] ) / (x[
                      i] M1 + (1 - x[i]) M2) 0.43028}, {i, 0, 
        100}]

Table[{x[i], T[i]}, {i, 0, 100}];

Tbl2 = Table[{x[i], hl[i]}, {i, 0, 100}];

Tbl2bis = 
        Table[{x[i] M1/(x[i] M1 + (1 - x[i]) M2), hlm[i]}, {i, 0, 
            100}];

G1 = 
        Exp[A12 (A21 (1 - x)/(A12 x + A21 (1 - x)))^2] /. 
          x -> xh;
G2 = 
        Exp[A21 (A12 x/(A12 x + A21 (1 - x)))^2] /. 
          x -> xh;

Do[{yh[i] = h i, 
        Th[i] = FindRoot[{P ==  PS1 G1 xh + PS2 G2 (1 - xh), 
                yh[i] ==   PS1 G1 xh/P}, {T, 30}, {xh, 
                0.001}][[1, 2]]}, {i, 0, 100}]

Tbl3 = Table[{yh[i], Th[i]}, {i, 0, 100}];

Do[{hv[
            i] = (yh[i] dh1 + (1 - yh[i]) 
                  dh2 + (yh[i] 
                      Integrate[
                        CG1[T], {T, 273.15, Th[i] + 273.15}] + (1 - 
                          yh[i]) 
                      Integrate[CG2[T], {T, 273.15, Th[i] + 273.15}])) 
            0.43028, 
        hvm[i] = (yh[i] 
                    M1/(yh[i] M1 + (1 - yh[i]) M2) dh1 + (1 - 
                        yh[i]) M2/(yh[i] M1 + (1 - yh[i]) M2) 
                    dh2 + yh[i] M1/(yh[i] M1 + (1 - yh[i]) M2) 
                    Integrate[
                      CG1[T], {T, 273.15, Th[i] + 273.15}] + (1 - yh[i]) 
                    M2/(yh[i] M1 + (1 - yh[i]) M2) 
                    Integrate[CG2[T], {T, 273.15, Th[i] + 273.15}] )/(x[
                      i] M1 + (1 - x[i]) M2) 0.43028}, {i, 0, 
        100}]

Tbl4 = Table[{yh[i], hv[i]}, {i, 0, 100}];

Tbl4bis = 
        Table[{yh[i] M1/(yh[i] M1 + (1 - yh[i]) M2), hvm[i]}, {i, 
            0, 100}];

Plt2 = 
      ListPlot[Tbl2, PlotJoined -> True, 
        PlotStyle -> RGBColor[1, 0, 0], 
        DisplayFunction -> Identity]

Plt2bis = 
      ListPlot[Tbl2bis, PlotJoined -> True, 
        PlotStyle -> RGBColor[1, 0, 0], 
        DisplayFunction -> Identity]

Plt3 = 
      ListPlot[Tbl4, PlotJoined -> True, 
        PlotStyle -> RGBColor[0, 0, 1], 
        DisplayFunction -> Identity]

Plt3bis = 
      ListPlot[Tbl4bis, PlotJoined -> True, 
        PlotStyle -> RGBColor[0, 0, 1], 
        DisplayFunction -> Identity]

Vapor = Interpolation[Tbl4, InterpolationOrder ->  1];

Vaporbis = 
        Interpolation[Tbl4bis, InterpolationOrder ->  1];

Liquid = Interpolation[Tbl2, InterpolationOrder ->  1];

Liquidbis = 
        Interpolation[Tbl2bis, InterpolationOrder ->  1];

Tbl5 = 
        Table[Line[{{x[i], Liquid[x[i]]}, {y[i], Vapor[y[i]]}}], {i, 2, 90, 
            10}];

Tbl5bis = 
        Table[Line[{{x[i] M1/(x[i] M1 + (1 - x[i]) M2), 
                Liquidbis[
                  x[i] M1/(x[i] M1 + (1 - x[i]) M2)]}, {y[
                    i] M1/(y[i] M1 + (1 - y[i]) M2), 
                Vaporbis[
                  y[i] M1/(y[i] M1 + (1 - y[i]) M2)]}}], {i, 2, 
            90, 10}];

Show[Plt3, Plt2, Graphics[{RGBColor[0, 1, 0], Tbl5}], 
      DisplayFunction -> $DisplayFunction]

Show[Plt3bis, Plt2bis, Graphics[{RGBColor[0, 1, 0], Tbl5bis}], 
      DisplayFunction -> $DisplayFunction]

Tbl6 = Table[{x[i], Vapor[y[i]]}, {i, 0, 100}];

Tbl6bis = 
        Table[{x[i] M1/(x[i] M1 + (1 - x[i]) M2), 
            Vaporbis[y[i] M1/(y[i] M1 + (1 - y[i]) M2)]}, {i, 0, 
            100}];

Plt4 = 
      ListPlot[Tbl6, PlotStyle -> RGBColor[1, 0, 1], 
        PlotJoined -> True, DisplayFunction -> Identity]

Plt4bis = 
      ListPlot[Tbl6bis, PlotStyle -> RGBColor[1, 0, 1], 
        PlotJoined -> True, DisplayFunction -> Identity]

Plt5 = 
      Show[Plt3, Plt2, Plt4, Graphics[{RGBColor[0, 1, 0], Tbl5}], 
        PlotRange -> All, AspectRatio -> 1.5, 
        Epilog -> {Line[{{x[20], Liquid[x[20]]}, {x[20], 
                  Vapor[y[20]]}}], 
Line[{{x[20], 
                  Vapor[y[20]]}, {y[20], Vapor[y[20]]}}]}, 
        DisplayFunction -> $DisplayFunction]

Plt5bis = 
      Show[Plt3bis, Plt2bis, Plt4bis, Graphics[{RGBColor[0, 0, 1], Tbl5bis}], 
        PlotRange -> {-100, 1200}, AspectRatio -> 2, 
        Epilog -> {Line[{{x[
                      12] M1/(x[12] M1 + (1 - x[12]) M2), 
                  Liquidbis[
                    x[12] M1/(x[12] M1 + (1 - x[12]) M2)]}, {x[
                      12] M1/(x[12] M1 + (1 - x[12]) M2), 
                  Vaporbis[
                    y[12] M1/(y[
                                12] M1 + (1 - 
                                  y[12]) M2)]}}], 
Line[{{x[12] M1/(x[12] M1 + (1 - x[12]) M2), 
                  Vaporbis[
                    y[12] M1/(y[12] M1 + (1 - y[12]) M2)]}, {y[
                      12] M1/(y[12] M1 + (1 - y[12]) M2), 
                  Vaporbis[
                    y[12] M1/(y[12] M1 + (1 - y[12]) M2)]}}]}, 
        DisplayFunction -> $DisplayFunction]

Enthalpy composition diagram of ethanol water system at 760 mmHg and 76 mmHg

Show[Plt5bis, plt5bis, 
      DisplayFunction -> $DisplayFunction]