Liquid Liquid Equilibrium Curve and Tie Lines 
for ternary system: water - ethyl acetate - acetone at 15 °C
Author's Data: Housam BINOUS
Department of Chemical Engineering
National Institute of Applied Sciences and Technology
Tunis, TUNISIA
Email: binoushousam@yahoo.com

Off[General::"spell1"]

Activities using the NRTL model

GAMMA1 = 
      Exp[(t21/(R T) G21 x2 + t31/(R T) G31 x3)/(x1 + 
                G21 x2 + G31 x3) + 
          x1/(x1 + G21 x2 + 
                  G31 x3) (-(x2 t21/(R T) G21 + 
                      x3 t31/(R T) G31)/(x1 + G21 x2 + 
                    G31 x3)) + 
          x2 G12/(G12 x1 + x2 + 
                  G32 x3) (t12/(R T) - (x1 t12/(R T) G12 + x3 t32/(R T) G32)/(G12 x1 + x2 + 
                      G32 x3)) + 
x3 G13/(G13 x1 + 
                  G23 x2  + 
                  x3) (t13/(R T) - (x1 t13/(R T) G13 + 
                      x2 t23/(R T) G23)/(G13 x1 + G23 x2  + 
                      x3))]

GAMMA1p = 
      Exp[(t21/(R T) G21 xp2 + t31/(R T) G31 xp3)/(xp1 + 
                G21 xp2 + G31 xp3) + 
          xp1/(xp1 + G21 xp2 + 
                  G31 xp3) (-(xp2 t21/(R T) G21 + 
                      xp3 t31/(R T) G31)/(xp1 + G21 xp2 + 
                    G31 xp3)) + 
          xp2 G12/(G12 xp1 + xp2 + 
                  G32 xp3) (t12/(R T) - (xp1 t12/(R T) G12 + xp3 t32/(R T) G32)/(G12 xp1 + xp2 + 
                      G32 xp3)) + 
xp3 G13/(G13 xp1 + G23 xp2  + 
                  xp3) (t13/(R T) - (xp1 t13/(R T) G13 + 
                      xp2 t23/(R T) G23)/(G13 xp1 + G23 xp2  + 
                      xp3))]

GAMMA2 = 
      Exp[(t12/(R T) G12 x1 + 
                t32/(R T) G32 x3)/(G12 x1 + x2 + G32 x3) + 
          x1 G21/(x1 + G21 x2 + 
                  G31 x3) (t21/(R T) - (x2 t21/(R T) G21 + x3 t31/(R T) G31)/(x1 + G21 x2 + G31 x3)) + 
          x2/(G12 x1 + x2 + 
                  G32 x3) (-(x1 t12/(R T) G12 + 
                      x3 t32/(R T) G32)/(G12 x1 + x2 + 
                    G32 x3)) + 
          x3 G23/(G13 x1 + G23 x2  + 
                  x3) (t23/(R T) - (x1 t13/(R T) G13 + 
                      x2 t23/(R T) G23)/(G13 x1 + G23 x2  + 
                      x3))]

GAMMA2p = 
      Exp[(t12/(R T) G12 xp1 + 
                t32/(R T) G32 xp3)/(G12 xp1 + xp2 + G32 xp3) + 
          xp1 G21/(xp1 + G21 xp2 + 
                  G31 xp3) (t21/(R T) - (xp2 t21/(R T) G21 + xp3 t31/(R T) G31)/(xp1 + G21 xp2 + G31 xp3)) + 
          xp2/(G12 xp1 + xp2 + 
                  G32 xp3) (-(xp1 t12/(R T) G12 + 
                      xp3 t32/(R T) G32)/(G12 xp1 + xp2 + 
                    G32 xp3)) + 
          xp3 G23/(G13 xp1 + G23 xp2  + 
                  xp3) (t23/(R T) - (xp1 t13/(R T) G13 + 
                      xp2 t23/(R T) G23)/(G13 xp1 + G23 xp2  + 
                      xp3))]

GAMMA3 = 
      Exp[(t13/(R T) G13 x1 + 
                t23/(R T) G23 x2)/(G13 x1 + G23 x2 + x3) + 
          x1 G31/(x1 + G21 x2 + 
                  G31 x3) (t31/(R T) - (x2 t21/(R T) G21 + x3 t31/(R T) G31)/(x1 + G21 x2 + G31 x3)) + 
          x2 G32/(G12 x1 + x2 + 
                  G32 x3) (t32/(R T) - (x1 t12/(R T) G12 + x3 t32/(R T) G32)/(G12 x1 + x2 + G32 x3)) + 
          x3/(G13 x1 + G23 x2  + 
                  x3) (-(x1 t13/(R T) G13 + 
                      x2 t23/(R T) G23)/(G13 x1 + G23 x2  + 
                    x3))]

GAMMA3p = 
      Exp[(t13/(R T) G13 xp1 + 
                t23/(R T) G23 xp2)/(G13 xp1 + G23 xp2 + xp3) + 
          xp1 G31/(xp1 + G21 xp2 + 
                  G31 xp3) (t31/(R T) - (xp2 t21/(R T) G21 + xp3 t31/(R T) G31)/(xp1 + G21 xp2 + G31 xp3)) + 
          xp2 G32/(G12 xp1 + xp2 + 
                  G32 xp3) (t32/(R T) - (xp1 t12/(R T) G12 + xp3 t32/(R T) G32)/(G12 xp1 + xp2 + G32 xp3)) + 
          xp3/(G13 xp1 + G23 xp2  + 
                  xp3) (-(xp1 t13/(R T) G13 + 
                      xp2 t23/(R T) G23)/(G13 xp1 + 
                    G23 xp2  + xp3))]

G12 = Exp[-a12 t12/(R T)]; 
    G21 = Exp[-a12 t21/(R T)]; G13 = Exp[-a13 t13/(R T)]; 
    G31 = Exp[-a13 t31/(R T)]; G32 = Exp[-a32 t32/(R T)]; 
    G23 = Exp[-a32 t23/(R T)];

Gas constant, temperature and binary interaction parameters for the NRTL model

R = 1.987;
a13 = 0.585600018501282; 
    a12 = 0.410400003194809;
a32 = 0.291200011968613;
t31 = 750.318115234375; 
    t13 = 1299.39501953125;
t23 = 467.125701904297; 
    t32 = -304.988403320313;
t12 = 2316.36303710938;
t21 = 935.687927246094;

T = 273.15 + 15

Computing Tie lines and LLE curve 
water = component 1
ethylacetate = component 2
acetone = component 3

X1 = 0.7; X2 = 0.2;

\(sol[1] = 
      FindRoot[{\((GAMMA1\ x1 == \ 
              GAMMA1p\ xp1)\), \((GAMMA2\ x2 == \ 
              GAMMA2p\ xp2)\), 
\((GAMMA3\ x3 == \ 
              GAMMA3p\ xp3)\), 
\((x1\ L + 
                xp1\ Lp == \ X1)\), \((x2\ L + xp2\ Lp == \ 
              X2)\), \((1 == \ xp1 + xp2 + xp3)\), \((1 == \ 
              x1 + x2 + x3)\), \((L + Lp == 1)\)}, {x1, 0.95}, {x2, 
          0.01}, {xp1, 0.5}, {xp2, 0.4}, {L, 0.05}, {Lp, 0.04}, {x3, 
          0.01}, {xp3, 0.1}, 
	MaxIterations -> 10000]\)

X1 = 0.8; X2 = 0.1;

\(sol[2] = 
      FindRoot[{\((GAMMA1\ x1 == \ 
              GAMMA1p\ xp1)\), \((GAMMA2\ x2 == \ 
              GAMMA2p\ xp2)\), 
\((GAMMA3\ x3 == \ 
              GAMMA3p\ xp3)\), 
\((x1\ L + 
                xp1\ Lp == \ X1)\), \((x2\ L + xp2\ Lp == \ 
              X2)\), \((1 == \ xp1 + xp2 + xp3)\), \((1 == \ 
              x1 + x2 + x3)\), \((L + Lp == 1)\)}, {x1, 0.95}, {x2, 
          0.01}, {xp1, 0.5}, {xp2, 0.4}, {L, 0.05}, {Lp, 0.04}, {x3, 
          0.01}, {xp3, 0.1}, 
	MaxIterations -> 10000]\)

X1 = 0.85; X2 = 0.1;

sol[3] = 
      FindRoot[{(GAMMA1 x1 ==  
              GAMMA1p xp1), (GAMMA2 x2 ==  
              GAMMA2p xp2), 
(GAMMA3 x3 ==  
              GAMMA3p xp3), 
(x1 L + 
                xp1 Lp ==  X1), (x2 L + xp2 Lp ==  
              X2), (1 ==  xp1 + xp2 + xp3), (1 ==  
              x1 + x2 + x3), (L + Lp == 1)}, {x1, 0.95}, {x2, 
          0.01}, {xp1, 0.5}, {xp2, 0.4}, {L, 0.05}, {Lp, 0.04}, {x3, 
          0.01}, {xp3, 0.1}, MaxIterations -> 10000]

X1 = 0.9; X2 = 0.1;

sol[4] = 
      FindRoot[{(GAMMA1 x1 ==  
              GAMMA1p xp1), (GAMMA2 x2 ==  
              GAMMA2p xp2), 
(GAMMA3 x3 ==  
              GAMMA3p xp3), 
(x1 L + 
                xp1 Lp ==  X1), (x2 L + xp2 Lp ==  
              X2), (1 ==  xp1 + xp2 + xp3), (1 ==  
              x1 + x2 + x3), (L + Lp == 1)}, {x1, 0.95}, {x2, 
          0.01}, {xp1, 0.5}, {xp2, 0.4}, {L, 0.05}, {Lp, 0.04}, {x3, 
          0.01}, {xp3, 0.1}, MaxIterations -> 100000]

X1 = 0.65; X2 = 0.25;

\(sol[5] = 
      FindRoot[{\((GAMMA1\ x1 == \ 
              GAMMA1p\ xp1)\), \((GAMMA2\ x2 == \ 
              GAMMA2p\ xp2)\), 
\((GAMMA3\ x3 == \ 
              GAMMA3p\ xp3)\), 
\((x1\ L + 
                xp1\ Lp == \ X1)\), \((x2\ L + xp2\ Lp == \ 
              X2)\), \((1 == \ xp1 + xp2 + xp3)\), \((1 == \ 
              x1 + x2 + x3)\), \((L + Lp == 1)\)}, {x1, 0.95}, {x2, 
          0.01}, {xp1, 0.5}, {xp2, 0.4}, {L, 0.05}, {Lp, 0.04}, {x3, 
          0.01}, {xp3, 0.1}, AccuracyGoal -> 10, 
	
        MaxIterations -> 10000]\)

X1 = 0.65; X2 = 0.35;

\(sol[6] = 
      FindRoot[{\((GAMMA1\ x1 == \ 
              GAMMA1p\ xp1)\), \((GAMMA2\ x2 == \ 
              GAMMA2p\ xp2)\), 
\((GAMMA3\ x3 == \ 
              GAMMA3p\ xp3)\), 
\((x1\ L + 
                xp1\ Lp == \ X1)\), \((x2\ L + xp2\ Lp == \ 
              X2)\), \((1 == \ xp1 + xp2 + xp3)\), \((1 == \ 
              x1 + x2 + x3)\), \((L + Lp == 1)\)}, {x1, 0.95}, {x2, 
          0.01}, {xp1, 0.5}, {xp2, 0.4}, {L, 0.05}, {Lp, 0.04}, {x3, 
          0.01}, {xp3, 0.1}, AccuracyGoal -> 10, 
	
        MaxIterations -> 10000]\)

X1 = 0.95; X2 = 0.01;

\(sol[7] = 
      FindRoot[{\((GAMMA1\ x1 == \ 
              GAMMA1p\ xp1)\), \((GAMMA2\ x2 == \ 
              GAMMA2p\ xp2)\), 
\((GAMMA3\ x3 == \ 
              GAMMA3p\ xp3)\), 
\((x1\ L + 
                xp1\ Lp == \ X1)\), \((x2\ L + xp2\ Lp == \ 
              X2)\), \((1 == \ xp1 + xp2 + xp3)\), \((1 == \ 
              x1 + x2 + x3)\), \((L + Lp == 1)\)}, {x1, 0.95}, {x2, 
          0.01}, {xp1, 0.5}, {xp2, 0.4}, {L, 0.05}, {Lp, 0.04}, {x3, 
          0.01}, {xp3, 0.1}, AccuracyGoal -> 10, 
	
        MaxIterations -> 10000]\)

X1 = 0.75; X2 = 0.2;

\(sol[8] = 
      FindRoot[{\((GAMMA1\ x1 == \ 
              GAMMA1p\ xp1)\), \((GAMMA2\ x2 == \ 
              GAMMA2p\ xp2)\), 
\((GAMMA3\ x3 == \ 
              GAMMA3p\ xp3)\), 
\((x1\ L + 
                xp1\ Lp == \ X1)\), \((x2\ L + xp2\ Lp == \ 
              X2)\), \((1 == \ xp1 + xp2 + xp3)\), \((1 == \ 
              x1 + x2 + x3)\), \((L + Lp == 1)\)}, {x1, 0.95}, {x2, 
          0.01}, {xp1, 0.5}, {xp2, 0.4}, {L, 0.05}, {Lp, 0.04}, {x3, 
          0.01}, {xp3, 0.1}, AccuracyGoal -> 10, 
	
        MaxIterations -> 10000]\)

X1 = 0.75; X2 = 0.15;

\(sol[9] = 
      FindRoot[{\((GAMMA1\ x1 == \ 
              GAMMA1p\ xp1)\), \((GAMMA2\ x2 == \ 
              GAMMA2p\ xp2)\), 
\((GAMMA3\ x3 == \ 
              GAMMA3p\ xp3)\), 
\((x1\ L + 
                xp1\ Lp == \ X1)\), \((x2\ L + xp2\ Lp == \ 
              X2)\), \((1 == \ xp1 + xp2 + xp3)\), \((1 == \ 
              x1 + x2 + x3)\), \((L + Lp == 1)\)}, {x1, 0.95}, {x2, 
          0.01}, {xp1, 0.5}, {xp2, 0.4}, {L, 0.05}, {Lp, 0.04}, {x3, 
          0.01}, {xp3, 0.1}, AccuracyGoal -> 10, 
	
        MaxIterations -> 10000]\)

X1 = 0.65; X2 = 0.3;

\(sol[10] = 
      FindRoot[{\((GAMMA1\ x1 == \ 
              GAMMA1p\ xp1)\), \((GAMMA2\ x2 == \ 
              GAMMA2p\ xp2)\), 
\((GAMMA3\ x3 == \ 
              GAMMA3p\ xp3)\), 
\((x1\ L + 
                xp1\ Lp == \ X1)\), \((x2\ L + xp2\ Lp == \ 
              X2)\), \((1 == \ xp1 + xp2 + xp3)\), \((1 == \ 
              x1 + x2 + x3)\), \((L + Lp == 1)\)}, {x1, 0.95}, {x2, 
          0.01}, {xp1, 0.5}, {xp2, 0.4}, {L, 0.05}, {Lp, 0.04}, {x3, 
          0.01}, {xp3, 0.1}, AccuracyGoal -> 10, 
	
        MaxIterations -> 10000]\)

X1 = 0.55; X2 = 0.4;

\(sol[11] = 
      FindRoot[{\((GAMMA1\ x1 == \ 
              GAMMA1p\ xp1)\), \((GAMMA2\ x2 == \ 
              GAMMA2p\ xp2)\), 
\((GAMMA3\ x3 == \ 
              GAMMA3p\ xp3)\), 
\((x1\ L + 
                xp1\ Lp == \ X1)\), \((x2\ L + xp2\ Lp == \ 
              X2)\), \((1 == \ xp1 + xp2 + xp3)\), \((1 == \ 
              x1 + x2 + x3)\), \((L + Lp == 1)\)}, {x1, 0.95}, {x2, 
          0.01}, {xp1, 0.5}, {xp2, 0.4}, {L, 0.05}, {Lp, 0.04}, {x3, 
          0.01}, {xp3, 0.1}, AccuracyGoal -> 10, 
	
        MaxIterations -> 10000]\)

Plots of Liquid-Liquid Equilibrium Curve and Tie Lines for Water-Benzene-Ethanol Ternary System at 25°C. Plait Point is shown as a red dot.

tbl1 = Table[{x1 /. sol[i], x2 /. sol[i]}, {i, 1, 11}]

tbl2 = Table[{xp1 /. sol[i], xp2 /. sol[i]}, {i, 1, 11}]

tbl1final = Sort[tbl1]; tbl2final = Sort[tbl2];

plt1 = 
      ListPlot[tbl1final, AspectRatio -> 1, PlotJoined -> True, 
        DisplayFunction -> Identity]

plt2 = 
      ListPlot[tbl2final, AspectRatio -> 1, PlotJoined -> True, 
        DisplayFunction -> Identity]

\(Show[plt1, plt2, 
      Graphics[{RGBColor[1, 0, 0], PointSize[0.015], 
          Point[{0.848637294208835, 
0.0357681230623852}], 
          RGBColor[0, 0, 1], Text[acetone, {0, \(-0.1\)}], 
          Text[water, {1.0, \(-0.1\)}], Text[ethylacetate, {0, 1.1}], 
          RGBColor[0, 0, 1]}], 
      PlotRange -> {{\(-0.2\), 1.2}, {\(-0.2\), 1.2}}, 
      Epilog -> {{RGBColor[0, 0, 
              1], 
Line[{{x1 /. sol[1], 
                  x2 /. sol[1]}, {xp1 /. sol[1], 
                  xp2 /. sol[1]}}], 
Line[{{x1 /. sol[2], 
                  x2 /. sol[2]}, {xp1 /. sol[2], xp2 /. sol[2]}}], 
            Line[{{x1 /. sol[3], x2 /. sol[3]}, {xp1 /. sol[3], 
                  xp2 /. sol[3]}}], 
            Line[{{x1 /. sol[4], x2 /. sol[4]}, {xp1 /. sol[4], 
                  xp2 /. sol[4]}}], 
            Line[{{x1 /. sol[5], x2 /. sol[5]}, {xp1 /. sol[5], 
                  xp2 /. sol[5]}}], 
            Line[{{x1 /. sol[6], x2 /. sol[6]}, {xp1 /. sol[6], 
                  xp2 /. sol[6]}}], 
            Line[{{x1 /. sol[7], x2 /. sol[7]}, {xp1 /. sol[7], 
                  xp2 /. sol[7]}}], 
            Line[{{x1 /. sol[8], x2 /. sol[8]}, {xp1 /. sol[8], 
                  xp2 /. sol[8]}}], 
            Line[{{x1 /. sol[9], x2 /. sol[9]}, {xp1 /. sol[9], 
                  xp2 /. sol[9]}}], 
            Line[{{x1 /. sol[10], x2 /. sol[10]}, {xp1 /. sol[10], 
                  xp2 /. sol[10]}}], 
            Line[{{x1 /. sol[11], x2 /. sol[11]}, {xp1 /. sol[11], 
                  xp2 /. sol[11]}}], {RGBColor[0, 1, 0], 
              Line[{{1, 0}, {0, 1}}]}}}, 
      DisplayFunction -> $DisplayFunction]\)