Distillation Curves for Ternary Systems
Case of Acetone - Chloroform - Benzene System
Authors' Data: Housam BINOUS*, Aïda WAKAD and Saïda Ben Achour
Department of Chemical Engineering
National Institute of Applied Sciences and Technology
Tunis, TUNISIA
Email: binoushousam@yahoo.com

Reference: L. T. Biegler, I. E. Grossmann and A. W. Westerberg, Systematic Methods of Chemical Process Design, Prentice Hall, 1999.

Acetone = component 1
Chloroform = component 2
Benzene = component 3

Off[General::"spell1"]

Gas Constant and Total Pressure

R = 1.987; P = 760;

Vapor Pressure Using Antoine Equation

A1 = 7.11714; B1 = 1210.595; C1 = 229.664;
A2 = 6.95465; B2 = 1170.966; C2 = 226.232;
A3 = 6.87987;B3 = 1196.760;C3 = 
      219.161;

PS1[i_] := 
        10^(A1 - B1/(C1 + T));
PS2[i_] := 
        10^(A2 - B2/(C2 + T));
PS3[i_] := 10^(A3 - B3/(C3 + T));

Binary inetraction Parameters and Activities for Wilson Model

l12 = 116.1171; l21 = -506.8519; 
    l13 = 682.40610;
l31 = -243.96510; l23 = -71.81089; 
    l32 = -11.82310;
V1 = 74.05;V2 = 80.67;V3 = 
      89.41;

A12 = 
        V2/V1 Exp[-l12/(R*(273.15 + 
                      T))];
A13 = 
        V3/V1 Exp[-l13/(R*(273.15 + 
                      T))];
A32 = 
        V2/V3 Exp[-l32/(R*(273.15 + 
                      T))];
A21 = 
        V1/V2 Exp[-l21/(R*(273.15 + 
                      T))];
A31 = 
        V1/V3 Exp[-l31/(R*(273.15 + 
                      T))];
A23 = V3/V2 Exp[-l23/(R*(273.15 + T))];

GAM1[i_] := 
        Exp[-Log[x1[i] + x2[i]*A12 + x3[i]*A13] + 
            1 - (x1[i]/(x1[i] + x2[i]*A12 + x3[i]*A13) + 
                x2[i]*A21/(x1[i]*A21 + x2[i] + x3[i]*A23) + 
                x3[i]*A31/(x1[i]*A31 + x2[i]*A32 + 
                        x3[i]))];
GAM2[i_] := 
        Exp[-Log[x1[i]*A21 + x2[i] + x3[i]*A23] + 
            1 - (x1[i]*A12/(x1[i] + x2[i]*A12 + x3[i]*A13) + 
                x2[i]/(x1[i]*A21 + x2[i] + x3[i]*A23) + 
                x3[i]*A32/(x1[i]*A31 + x2[i]*A32 + 
                        x3[i]))];
GAM3[i_] := 
        Exp[-Log[x1[i]*A31 + x2[i]*A32 + x3[i]] + 
            1 - (x1[i]*A13/(x1[i] + x2[i]*A12 + x3[i]*A13) + 
                x2[i]*A23/(x1[i]*A21 + x2[i] + x3[i]*A23) + 
                x3[i]/(x1[i]*A31 + x2[i]*A32 + x3[i]))];

GAMMA1[
          i_] := GAM1[i] /. x1[i] ->  x1 /. x2[i] ->  x2 /. 
          x3[i] ->  x3;

GAMMA2[
          i_] := GAM2[i] /. x1[i] ->  x1 /. x2[i] ->  x2 /. 
          x3[i] ->  x3;

GAMMA3[
        i_] := GAM3[i] /. x1[i] ->  x1 /. x2[i] ->  x2 /. 
        x3[i] ->  x3

Calculation of Distillation Curves for different initial conditions. We get the same figure as figure 14.8 page 466 of the above reference.

x1[0] = 0.2; x2[0] = 0.2; x3[0] = 0.6; i = 0;

While[
        i < 200, {T = FindRoot[
                P == 
                  Evaluate[
                    x1[i] PS1[i] GAM1[i] + x2[i] PS2[i] GAM2[i] + 
                      x3[i] PS3[i] GAM3[i]], {T, 100}][[1, 2]], 
          y1[i] = x1[i] PS1[i] GAM1[i]/P, 
y2[i] = 
            x2[i] PS2[i] GAM2[i]/P, 
y3[i] = 
            x3[i] PS3[i] GAM3[i]/P, x1[i + 1] = y1[i], x2[i + 1] = y2[i], 
          x3[i + 1] = y3[i], i++}];

t1 = Table[{x2[i], x3[i]}, {i, 0, 200}];

plt1 = 
        ListPlot[t1, PlotRange -> {{0, 1}, {0, 1}}, 
          PlotJoined ->  True, 
          PlotStyle ->  RGBColor[1, 0, 0]];

T =.

i = 0; x1[0] = 0.2; x2[0] = 0.2; x3[0] = 0.6; y1[0] = 0.2; y2[0] = 0.2; 
    y3[0] = 0.6;

While[
        i < 200, {sol = 
            FindRoot[{x1 + x2 + x3 ==  1, 
                y1[i] ==  
                  x1 PS1[i] GAMMA1[i]/P, 
y2[
                    i] ==  
                  x2 PS2[i] GAMMA2[i]/P, 
y3[
                    i] ==  x3 PS3[i] GAMMA3[i]/P}, {T, 60, 90}, {x1, 
                0.01, 0.2}, {x2, 0.01, 0.3}, {x3, 0.5, 0.9}, 
              MaxIterations ->  1000], 
          x1[i + 1] = Chop[sol[[2, 2]]], 
          x2[i + 1] = Chop[sol[[3, 2]]], 
          x3[i + 1] = Chop[sol[[4, 2]]], y1[i + 1] = x1[i], 
          y2[i + 1] = x2[i], y3[i + 1] = x3[i], i++}];

t2 = Table[{x2[i], x3[i]}, {i, 0, 200}];

plt2 = 
        ListPlot[t2, PlotRange -> {{0, 1}, {0, 1}}, 
          PlotJoined ->  True, 
          PlotStyle ->  RGBColor[1, 0, 0]];

x1[0] = 0.8; x2[0] = 0.1; x3[0] = 0.1; i = 0;

While[
        i < 200, {T = FindRoot[
                P == 
                  Evaluate[
                    x1[i] PS1[i] GAM1[i] + x2[i] PS2[i] GAM2[i] + 
                      x3[i] PS3[i] GAM3[i]], {T, 100}][[1, 2]], 
          y1[i] = x1[i] PS1[i] GAM1[i]/P, 
y2[i] = 
            x2[i] PS2[i] GAM2[i]/P, 
y3[i] = 
            x3[i] PS3[i] GAM3[i]/P, x1[i + 1] = y1[i], x2[i + 1] = y2[i], 
          x3[i + 1] = y3[i], i++}];

t3 = Table[{x2[i], x3[i]}, {i, 0, 200}];

plt3 = 
        ListPlot[t3, PlotRange -> {{0, 1}, {0, 1}}, 
          PlotJoined ->  True, 
          PlotStyle ->  RGBColor[1, 0, 0]];

T =.

i = 0; x1[0] = 0.8; x2[0] = 0.1; x3[0] = 0.1; y1[0] = 0.8; y2[0] = 0.1; 
    y3[0] = 0.1;

While[
        i < 200, {sol = 
            FindRoot[{x1 + x2 + x3 ==  1, 
                y1[i] ==  
                  x1 PS1[i] GAMMA1[i]/P, 
y2[
                    i] ==  
                  x2 PS2[i] GAMMA2[i]/P, 
y3[
                    i] ==  x3 PS3[i] GAMMA3[i]/P}, {T, 60, 90}, {x1, 
                0.01, 0.2}, {x2, 0.01, 0.3}, {x3, 0.5, 0.9}, 
              MaxIterations ->  1000], 
          x1[i + 1] = Chop[sol[[2, 2]]], 
          x2[i + 1] = Chop[sol[[3, 2]]], 
          x3[i + 1] = Chop[sol[[4, 2]]], y1[i + 1] = x1[i], 
          y2[i + 1] = x2[i], y3[i + 1] = x3[i], i++}];

t4 = Table[{x2[i], x3[i]}, {i, 0, 200}];

plt4 = 
        ListPlot[t4, PlotRange -> {{0, 1}, {0, 1}}, 
          PlotJoined ->  True, 
          PlotStyle ->  RGBColor[1, 0, 0]];

x1[0] = 0.1; x2[0] = 0.8; x3[0] = 0.1; i = 0;

While[
        i < 200, {T = FindRoot[
                P == 
                  Evaluate[
                    x1[i] PS1[i] GAM1[i] + x2[i] PS2[i] GAM2[i] + 
                      x3[i] PS3[i] GAM3[i]], {T, 100}][[1, 2]], 
          y1[i] = x1[i] PS1[i] GAM1[i]/P, 
y2[i] = 
            x2[i] PS2[i] GAM2[i]/P, 
y3[i] = 
            x3[i] PS3[i] GAM3[i]/P, x1[i + 1] = y1[i], x2[i + 1] = y2[i], 
          x3[i + 1] = y3[i], i++}];

t5 = Table[{x2[i], x3[i]}, {i, 0, 200}];

plt5 = 
        ListPlot[t5, PlotRange -> {{0, 1}, {0, 1}}, 
          PlotJoined ->  True, 
          PlotStyle ->  RGBColor[1, 0, 0]];

T =.

i = 0; x1[0] = 0.1; x2[0] = 0.8; x3[0] = 0.1; y1[0] = 0.1; y2[0] = 0.8; 
    y3[0] = 0.1;

While[
        i < 200, {sol = 
            FindRoot[{x1 + x2 + x3 ==  1, 
                y1[i] ==  
                  x1 PS1[i] GAMMA1[i]/P, 
y2[
                    i] ==  
                  x2 PS2[i] GAMMA2[i]/P, 
y3[
                    i] ==  x3 PS3[i] GAMMA3[i]/P}, {T, 60, 90}, {x1, 
                0.01, 0.2}, {x2, 0.01, 0.3}, {x3, 0.5, 0.9}, 
              MaxIterations ->  1000], 
          x1[i + 1] = Chop[sol[[2, 2]]], 
          x2[i + 1] = Chop[sol[[3, 2]]], 
          x3[i + 1] = Chop[sol[[4, 2]]], y1[i + 1] = x1[i], 
          y2[i + 1] = x2[i], y3[i + 1] = x3[i], i++}];

t6 = Table[{x2[i], x3[i]}, {i, 0, 200}];

plt6 = 
        ListPlot[t6, PlotRange -> {{0, 1}, {0, 1}}, 
          PlotJoined ->  True, 
          PlotStyle ->  RGBColor[1, 0, 0]];

x1[0] = 0.33; x2[0] = 0.65; x3[0] = 0.02; i = 0;

While[
        i < 200, {T = FindRoot[
                P == 
                  Evaluate[
                    x1[i] PS1[i] GAM1[i] + x2[i] PS2[i] GAM2[i] + 
                      x3[i] PS3[i] GAM3[i]], {T, 100}][[1, 2]], 
          y1[i] = x1[i] PS1[i] GAM1[i]/P, 
y2[i] = 
            x2[i] PS2[i] GAM2[i]/P, 
y3[i] = 
            x3[i] PS3[i] GAM3[i]/P, x1[i + 1] = y1[i], x2[i + 1] = y2[i], 
          x3[i + 1] = y3[i], i++}];

t7 = Table[{x2[i], x3[i]}, {i, 0, 200}];

plt7 = 
        ListPlot[t7, PlotRange -> {{0, 1}, {0, 1}}, 
          PlotJoined ->  True, 
          PlotStyle ->  RGBColor[0, 0, 1]];

T =.

i = 0; x1[0] = 0.33; x2[0] = 0.65; x3[0] = 0.02; y1[0] = 0.33; 
    y2[0] = 0.65; y3[0] = 0.02;

While[
        i < 200, {sol = 
            FindRoot[{x1 + x2 + x3 ==  1, 
                y1[i] ==  
                  x1 PS1[i] GAMMA1[i]/P, 
y2[
                    i] ==  
                  x2 PS2[i] GAMMA2[i]/P, 
y3[
                    i] ==  x3 PS3[i] GAMMA3[i]/P}, {T, 60, 90}, {x1, 
                0.01, 0.2}, {x2, 0.01, 0.3}, {x3, 0.5, 0.9}, 
              MaxIterations ->  1000], 
          x1[i + 1] = Chop[sol[[2, 2]]], 
          x2[i + 1] = Chop[sol[[3, 2]]], 
          x3[i + 1] = Chop[sol[[4, 2]]], y1[i + 1] = x1[i], 
          y2[i + 1] = x2[i], y3[i + 1] = x3[i], i++}];

t8 = Table[{x2[i], x3[i]}, {i, 0, 200}];

plt8 = 
        ListPlot[t8, PlotRange -> {{0, 1}, {0, 1}}, 
          PlotJoined ->  True, 
          PlotStyle ->  RGBColor[0, 0, 1]];

x1[0] = 0.31; x2[0] = 0.67; x3[0] = 0.02; i = 0;

While[
        i < 200, {T = FindRoot[
                P == 
                  Evaluate[
                    x1[i] PS1[i] GAM1[i] + x2[i] PS2[i] GAM2[i] + 
                      x3[i] PS3[i] GAM3[i]], {T, 100}][[1, 2]], 
          y1[i] = x1[i] PS1[i] GAM1[i]/P, 
y2[i] = 
            x2[i] PS2[i] GAM2[i]/P, 
y3[i] = 
            x3[i] PS3[i] GAM3[i]/P, x1[i + 1] = y1[i], x2[i + 1] = y2[i], 
          x3[i + 1] = y3[i], i++}];

t9 = Table[{x2[i], x3[i]}, {i, 0, 200}];

plt9 = 
        ListPlot[t9, PlotRange -> {{0, 1}, {0, 1}}, 
          PlotJoined ->  True, 
          PlotStyle ->  RGBColor[0, 0, 1]];

T =.

i = 0; x1[0] = 0.31; x2[0] = 0.67; x3[0] = 0.02; y1[0] = 0.31; 
    y2[0] = 0.67; y3[0] = 0.02;

While[
        i < 200, {sol = 
            FindRoot[{x1 + x2 + x3 ==  1, 
                y1[i] ==  
                  x1 PS1[i] GAMMA1[i]/P, 
y2[
                    i] ==  
                  x2 PS2[i] GAMMA2[i]/P, 
y3[
                    i] ==  x3 PS3[i] GAMMA3[i]/P}, {T, 60, 90}, {x1, 
                0.01, 0.2}, {x2, 0.01, 0.3}, {x3, 0.5, 0.9}, 
              MaxIterations ->  1000], 
          x1[i + 1] = Chop[sol[[2, 2]]], 
          x2[i + 1] = Chop[sol[[3, 2]]], 
          x3[i + 1] = Chop[sol[[4, 2]]], y1[i + 1] = x1[i], 
          y2[i + 1] = x2[i], y3[i + 1] = x3[i], i++}];

t10 = Table[{x2[i], x3[i]}, {i, 0, 200}];

plt10 = 
        ListPlot[t10, PlotRange -> {{0, 1}, {0, 1}}, 
          PlotJoined ->  True, 
          PlotStyle ->  RGBColor[0, 0, 1]];

Show[plt1, plt2, plt3, plt4, plt5, plt6, plt7, plt8, plt9, plt10, 
      AspectRatio -> 1, 
      Epilog -> {Line[{{1, 0}, {0, 1}}], RGBColor[0, 1, 0], 
          PointSize[0.03], Point[{0.6627290948414463`, 0}]}]

P = {{1, 0.5}, {0, 0.866025}};

t1eq = 
        Table[P . t1[[i]], {i, 1, 201}];
t2eq = 
        Table[P . t2[[i]], {i, 1, 201}];
t3eq = 
        Table[P . t3[[i]], {i, 1, 201}];
t4eq = 
        Table[P . t4[[i]], {i, 1, 201}];
t5eq = 
        Table[P . t5[[i]], {i, 1, 201}];
t6eq = 
        Table[P . t6[[i]], {i, 1, 201}];
t7eq = 
        Table[P . t7[[i]], {i, 1, 201}];
t8eq = 
        Table[P . t8[[i]], {i, 1, 201}];
t9eq = 
        Table[P . t9[[i]], {i, 1, 201}];
t10eq = 
        Table[P . t10[[i]], {i, 1, 201}];

plt1eq = 
        ListPlot[t1eq, PlotRange -> All, PlotJoined ->  True, 
          PlotStyle ->  RGBColor[0, 0, 1], Axes -> False, 
          DisplayFunction ->  Identity];
plt2eq = 
        ListPlot[t2eq, PlotRange -> All, PlotJoined ->  True, 
          PlotStyle ->  RGBColor[0, 0, 1], Axes -> False, 
          DisplayFunction ->  Identity];
plt3eq = 
        ListPlot[t3eq, PlotRange -> All, PlotJoined ->  True, 
          PlotStyle ->  RGBColor[0, 0, 1], Axes -> False, 
          DisplayFunction ->  Identity];
plt4eq = 
        ListPlot[t4eq, PlotRange -> All, PlotJoined ->  True, 
          PlotStyle ->  RGBColor[0, 0, 1], Axes -> False, 
          DisplayFunction ->  Identity];
plt5eq = 
        ListPlot[t5eq, PlotRange -> All, PlotJoined ->  True, 
          PlotStyle ->  RGBColor[0, 0, 1], Axes -> False, 
          DisplayFunction ->  Identity];
plt6eq = 
        ListPlot[t6eq, PlotRange -> All, PlotJoined ->  True, 
          PlotStyle ->  RGBColor[0, 0, 1], Axes -> False, 
          DisplayFunction ->  Identity];
plt7eq = 
        ListPlot[t7eq, PlotRange -> All, PlotJoined ->  True, 
          PlotStyle ->  RGBColor[0, 1, 0], Axes -> False, 
          DisplayFunction ->  Identity];
plt8eq = 
        ListPlot[t8eq, PlotRange -> All, PlotJoined ->  True, 
          PlotStyle ->  RGBColor[0, 1, 0], Axes -> False, 
          DisplayFunction ->  Identity];
plt9eq = 
        ListPlot[t9eq, PlotRange -> All, PlotJoined ->  True, 
          PlotStyle ->  RGBColor[0, 1, 0], Axes -> False, 
          DisplayFunction ->  Identity];
plt10eq = 
        ListPlot[t10eq, PlotRange -> All, PlotJoined ->  True, 
          PlotStyle ->  RGBColor[0, 1, 0], Axes -> False, 
          DisplayFunction ->  Identity];

Show[plt1eq, plt2eq, plt3eq, plt4eq, plt5eq, plt6eq, plt7eq, plt8eq, 
      plt9eq, plt10eq, 
      Epilog -> {RGBColor[1, 0, 0], Line[{{0, 0}, {1, 0}}], 
          Line[{{0, 0}, {0.5, 0.866025}}], Line[{{1, 0}, {0.5, 0.866025}}], 
          PointSize[0.02], RGBColor[0, 0, 1], 
          Text["Acetone", {-0.1, -0.1}], 
          Text["Chloroform", {1.1, -0.1}], 
          Text["Benzene", {0.5, 0.96}], RGBColor[0, 1, 0], 
          PointSize[0.015], Point[{0.6627290948414463`, 0}]}, 
      PlotRange -> {{-0.5, 1.5}, {-0.5, 1.5}}, 
      TextStyle ->  {FontSize -> 14, 
          FontWeight -> "Bold"}, 
      DisplayFunction -> $DisplayFunction, 
      AspectRatio -> 1]