I n d . E n g . C h e m . Res. 1987, 26, 1774-1781
1774
Simultaneous Optimization of Binary Phase Equilibrium and Thermodynamic Data for Organic Systems Paul K. Talley,* Christopher W. Bale, and Arthur D. Pelton Centre de Recherche en Calcul Thermochimique, Ecole Polytechnique de Montr&al, C.P. 6079, Montreal, QuQbec,Canada H3C 3A7
A technique is proposed which permits the simultaneous optimization of all available VLE, LLE, excess enthalpy, and other liquid excess property data for binary systems. T h e method permits the analysis of the complete range of experimental temperatures, pressures, and compositions. T h e result is a single set of self-consistent parameters for a power-series expression for the excess Gibbs energy which will reproduce experimental data with precision and can be used for extrapolations with confidence. T h e method is illustrated by examples including the following systems: hexane acetone; 2,3-dimethylbutane acetone; ethylbenzene + 2-ethoxyethanol; ethanol ethyl acetate; hexane hexadecane; water 1,4-dioxane; benzene + toluene; and water 2-butanone.
+
+
+ +
In the thermodynamic study of solutions of organic nonelectrolytes, there are many diverse types of experimental data available to the investigator. Vapor-liquid equilibria are a valuable source of excess Gibbs energy data. Calorimetric measurements of the liquid phase provide enthalpy of mixing data or excess heat capacities. Excess Gibbs energies can also be deduced from the boundaries of liquid immiscibility gaps, light-scattering measurements, freezing point depressions, activity coefficients at infinite dilution, and so on. A wide variety of equations have been proposed to express the thermodynamic properties of the liquid phase. Most of these equations have only two or three adjustable parameters. Such equations are often not flexible enough to represent very precise estimates of properties, especially for highly nonideal systems. In addition, extrapolation to temperatures and pressures outside the range of the fit cannot always be done with confidence. In the present article, an optimization procedure is proposed that enables one to treat all the different types of thermodynamic data simultaneously (phase equilibria data, heats of mixing, excess heat capacities, limiting activity coefficients, etc.) that have been measured under a variety of conditions. By means of a regression analysis, one analytical expression for the excess Gibbs energy of mixing for the liquid phase is obtained that represents all the diverse thermodynamic data a t the different temperatures, pressures, and compositions. The number of coefficients (ajustable parameters) in the expression is unconstrained, and thus, the data can be represented with precision. The proposed procedure is similar to a method used elsewhere (Bale and Pelton, 1983d) for the treatment of metallurgical systems. The strategy has been applied with considerable success to alloys (Ashtakala et al., 1981a-c; Bale and Pelton, 1983a-c), oxides (Blander and Pelton, 1983; Pelton and Blander, 1984), molten salts (Bale and Pelton, 1982; Gabriel et al., 1985; Gabriel and Pelton, 1985; Hatem et al., 1982; Lin et al., 1982; Pelton et al., 1982; Sangster and Pelton, 1984), and other inorganic systems where there may be several phases. The objective is t o provide one set of coefficients which enables the calculation of the various excess properties of mixing over the entire range of experimental conditions and permits reasonable extrapolations outside this region. In the proposed method, the integral molar excess Gibbs energy, gE,is expressed as the following power series in the mole fractions of the components, x1 and xg: g" = Za,x,xh + T C b , x l x $ + T In T C c k x l x $ ( I ) 0888-5885/87 / 2626- 1774$01.50/ 0
+
+
where Tis the temperature in kelvin and a,, b,, and ck are adjustable parameters. The temperature dependence is such as to permit enthalpies of mixing to be expressed as a linear function of temperature. The corresponding expressions for the enthalpy, hE, excess entropy, sE,and excess heat capacity. cpE, are h E = Calx& - TxckXIX; (2) sE = - x b , x l x
l2 -
(1 + In T ) Z c h x 1 x $
(3)
= -cckx1x$
(4) It is important to note that it is the excess properties which are fitted and not the activity coefficients or the logarithms of their ratios. The benefit of this has been discussed by Van Ness et al. (1973). An example illustrating the effective reproduction of experimental results is given in Figures 1-5 for the acetone + hexane system. The phase equilibria for this system have been measured in the ranges 20-900 mmHg and -20 to 55 "C. The heats of mixing have been measured from -30 to 40 "C. Each figure contains the experimental data points and calculated curves obtained from the proposed optimization procedure. The details of the analysis for this system are discussed later. All the curves in Figures 1-5 were calculated from the one single set of coefficients listed in Table I. The acetone + hexane system is one of over 4000 binary systems being evaluated by our group. The objective is to provide an on-line computing system with which one can perform phase equilibrium calculations for multicomponent systems as well as binaries. Extension of the power series expression to multicomponent systems will be discussed in a subsequent article. This on-line computing system will be coupled with an extensive bibliographic collection. The present work will supplement the Facility for the Analysis of Chemical Thermodynamics (F*A*C*T) (Bale et al., 1979) on-line database.
Theory Least-SquaresRegression Analysis. The basis of the optimization procedure is least-squares regression analysis. Equation 1 can be generalized into the form gE = c @ [ x I x ; (5) where the coefficients may be temperature-dependent. The partial excess Gibbs energies can be derived by differentiating eq 5. g? = X@[(l- iXl)Xl, (6) g: =
C@',(ZS~ - zX;)X1
0 1987 American Chemical Society
(7)
Ind. Eng. Chem. Res., Vol. 26, No. 9, 1987 1775 Table I. Power Series Coefficients of Equation 1 for Various Binary Liquids (cal/mol) system (1-2) 2 a,, cal/mol j b,, cali(mo1.K) hexane + acetone 1 1376.5 1 10.97 2 -5210.2 2 80.56 3 -70.93 4715.6 3 ethylbenzene + 2-ethoxyethanoP 1 1 2607.1 -33.37 2 -9356.0 2 174.96 3 6940.8 3 -136.96 1 989.1 2,3-dimethylbutane + acetone n -110.1 86.4 ethanol + ethyl acetate -41.96 3229.4 1 -4375.2 2 77.72 -3.62 1340.6 3 water + 2-butanone" 1 -2708.0 20.19 -23461.5 2 769.81 20645.1 3 -895.95 12492.4 a
k
cb. cali(mo1.K)
1 2 3 1 2 3
-2.115 -11.176 9.813 4.859 -25.444 20.039
1 2
5.870 -11.171
1
-0.909 -121.873 145.315 -7.245
2 3 4
Tsonopoulos (1974) binary interaction parameter ( k J = 0.10.
70
Hexane t Acetone e t 1 A t m
- Cgorodnikov
40
Hexane t Acetone a t 35 and 2o°C 500
et el
~
1
-0
Mole F r a c t i o n o f Hexane
0.0
1
Figure 1. Temperature-composition diagram for hexane a t 760.0 mmHg pressure. ioop"xane
1
119611
t Acetone
a t 55 and 45*c
,
,
.
,
00
+ acetone ,
~
Hole F r a i t l o n o f Hexane
Figure 3. Pressure-composition diagram for hexane 35 and 20 "C. BO
1
+ acetone a t
Hexane t Acetone a t -5 and -20.C
0 - Pnll end Schneffer 119591
i
3001 0.0
0
- Schaeffar
v
- Kudryavtseva and Susarav
and Rall 119581 Ii9631
1 Mole F r a c t i o n o f Hexane
Figure 2. Pressure-composition diagram for hexane 55 and 45 "C.
1
+ acetone a t
Note that the 4, coefficients are the same in eq 5-7. In the least-squares regression analysis, the objective function is written as N
M
5 = C ( y J- C@&2&J)2 ]=1
1=1
(8)
where [ is the variable to be minimized, N is the number of experimental points, M is the number of 9,coefficients, and YJcan have the value of an integral or a partial excess property at a selected point. The value of 2 , depends on the type of property in YJ(integral, partial of component 1, or partial of component 2). For example, if Y, is an integral property such as hE, it follows from eq 1 that 2 , is equal to xlxb. If Yl is a partial excess Gibbs energy for
10
J
1 00
Mole F r a c t i o n o f Hexane
Figure 4. Pressure-composition diagram for hexane -5 and -20 "C.
1
+ acetone a t
component 1, then 2, is (1 - ixl)x\. Thus, the partial properties of both components in solution can be treated simultaneously along with the integral excess properties in one simultaneous least squares regression. The temperature dependence of the coefficients in eq 5 is as shown in eq 1. In certain cases there are insufficient data to permit all three sets of coefficients to be calculated. For example, consider the case in which excess enthalpies have not been measured and the excess Gibbs energies (as derived from phase equilibrium data or otherwise) are sparse. In this instance, the bl and ck coefficients are generally set to zero and only the a,'s are calculated. In other systems, excess enthalpies as well as phase equilibrium data may be known but the excess heat
1776 Ind. Eng. Chem. Res., Vol. 26, No. 9, 1987 Table 11. Pure Component Data
cool
critical constants V,, cm3/mol 46.4 209 48.3 259 41.0 267 30.9 358 51.4 238 63.8 167 37.8 286 36.2 374 42.0 298 14.0 920 29.3 370 40.6 316 217.6 56
Pc, atm
T,, K 508.1 562.2 535.6 500.0 587.0 516.2 523.2 617.2 569.0 717.0 507.4 591.8 645.3
compound acetone benzene 2-butanone 2,3-dimethylbutane 1,4-dioxane ethanol ethyl acetate ethylbenzene 2-ethoxyethanol hexadecane hexane toluene water
0.309 0.210 0.329 0.247 0.288 0.635 0.363 0.301 0.75 0.742 0.296 0.264 0.344
>%
'9 '\
/ 3
svstem benzene + toluene water + 1,4-dioxane hexane + hexadecane
0.0
mole Frnctlon o f H+IIIxmn
Figure 5. Heats of mixing for hexane
+ acetone a t -30
Antoine coefficients" B 1210.595 1196.760 1261.339 1131.833 1554.679 1592.864 1244.951 1429.550 1801.900 1830.510 1189.640 1342.310 1730.630
C 229.664 219.161 221.969 229.462 240.337 226.189 217.881 213.767 230.000 154.450 226.280 219.187 233.420
Table 111. Power Series Coefficients for Heats of Mixing of Some Systems ( c a l h o l )
_ e . 30'C
A 7.117 14 6.879 87 7.063 56 6.818 32 7.431 55 8.11220 7.101 79 6.96580 7.819 10 7.028 70 6.910 58 6.950 87 8.071 31
W
1.0
i 1 1 2
3 4 1 2 3 4
a,, cal/mol 150.6 -7 405.0 27 736.6 -47 778.7 29 022.7 1312.8 -1 031.3 382.1 115.4
ckr
k
cal/(mol.K) 0.296 -18.016 64.655 -121.311 74.574 4.022 -3.535 1.713
where 11 and 12 refer to liquids 1 and 2, respectively. Equation 10 can be rearranged to give
and 30 "C.
capacities have not been determined. Here it is often necessary to set the ck coefficients to zero. In other systems for which a temperature dependence can be determined for the excess enthalpy, all three sets of coefficients can be optimized. Phase Equilibrium Data. Phase equilibrium data are divided into two classes: vapor-liquid equilibria, VLE, and liquid-liquid equilibria, LLE. The first class is further divided into two subclasses according to the variables reported, PTxy and PTx. These data are formulated in such a way as to allow all three types, or any combination thereof, to be used in the same regression analysis. Excess partial Gibbs energies from VLE data are determined from the relation c
Substituting eq 6 and simplifying yields
Expressions for Yj and 2, of eq 5 in this case are given as
Yj = RT In zij
bi,ii/Xi,d
(13)
= (1 - i ~ 1 , 1 2 ) ~ $ , 1-2 (1 - ~ ~ i , l J x i ? , ~ ~ (14)
Analogous expressions exist for the second component. where ,yLand x i are the vapor and liquid molar compositions, P and are the total and saturation pressures, and and are vapor and saturation fugacities. The Tsonopoulos correlation (Tsonopoulos, 1974) for the second virial coefficient is used for the calculation of the fugacities. The binary interaction parameters for the Tsonopoulos correlation have been set to zero here, except where otherwise indicated in Table I. The critical properties of the pure substances and Antoine coefficients used here are given in Table 11. PTxy data can be reduced explicitly by means of eq 9. LLE data are often available in the form of tie lines. The activities of component 1 in the two liquids can be equated as RT In xl,ll + (10) = RT 1n x1,12 + g f 1 2
Sample Calculations Seven systems are presented as sample calculations in Figures 1-22 in which diverse experimental data are compared to curves calculated from the coefficients listed in Tables I and 111. The ability of the present technique to represent many different types of data over wide ranges of T , P, and composition with one set of coefficients can be appreciated. In the present section, we discuss the procedures which were used to determine these coefficients with maximum utilization of the available data. The availability of data determines the procedure employed, but the method remains essentially the same. 2,3-Dimethylbutane + Acetone System. The phase diagram for the 2,3-dimethylbutane + acetone system has been measured at a constant pressure of 760 mmHg only
Ind. Eng. Chem. Res., Vol. 26, No. 9, 1987 1777 140
- Wrtl
and Van Winkla IlSnI
Liquid Liquid 40
H o l e F r a c t i o n o f 2.3-Dimethylbutane
0 0
1
Figure 6. Temperature-composition diagram for 2,3-dimethylbutane + acetone at 760 mmHg pressure.
(McConnell and Van Winkle, 1967; Willock and Van Winkle, 1970) (Figure 6). No information is available on the enthalpies for this system. The optimization is carried out for the a,'s alone, beginning with one coefficient and repeating with additional coefficients as deemed necessary. Generally, a satisfactory result is achieved when the rms deviation does not change significantly with the addition of one more coefficient. Care is taken to ensure that the resulting curve does not contain "wiggles". The results for the 2,3-dimethylbutane + acetone system indicated that three coefficients reproduced the data adequately. The phase diagram showing the experimental points along with the calculated phase boundaries is given in Figure 6. The calculated coefficients are given in Table I. Since no temperature dependence has been ascribed in this particular case, the coefficients should be used with caution a t other temperatures. Ethylbenzene 2-Ethoxyethanol System. The introduction of excess enthalpies into the optimization procedure increases the temperature range over which eq 1may be applied. In the ethylbenzene 2-ethoxyethanol binary system, enthalpies of mixing have been measured at 25 " C (Murti and Van Winkle, 1957). The most recent PTxy data have been measured at 50 mmHg (Fried et al., 1956) and 760 mmHg (Murti and Van Winkle, 1957). There are several possible approaches for the treatment of such data. The enthalpies could be fitted in a first step, thus determining the a, coefficients (the ck's being set to zero). Then the VLE data would be used to find the b 's in a subsequent step. Another approach is to optimize the a, and b, coefficients simultaneously using both the enthalpy and VLE data together. A third approach is to optimize all three sets of coefficients, a,, b,, c k , simultaneously. However, in this case, caution must be exercised to ensure that a reliable estimate of cpE can be obtained from the available data. In the present example, it was felt that the overall range of temperature of the data was sufficient to allow this type of analysis. An optimization using three a,, three b,, and three ck coefficients (see Table I) yielded average absolute deviations of 0.004 in vapor mole fractions, 0.41 "C in a bubble temperature calculation, and 2.7 cal/mol in the enthalpies of mixing. Phase diagrams for this system are given in Figures 7 and 8, showing experimental points and calculated phase boundaries. Figure 9 shows the calculated and experimental enthalpies. Hexane Acetone System. The binary systems which are of most interest are those for which a great deal of (possibly conflicting) data are available. A case in point is the hexane and acetone binary system mentioned earlier (Figures 1-5). As before, various possibilities exist in terms
+
120 0.0
YOIS
F r a c t i o n o f Ethylbenzene
Figure 7. Temperature-composition diagram for ethylbenzene 2-ethoxyethanol at 760 mmHg pressure. hylbenzene
t
2-Ethoxyethanol a t 50
0
-
nm
+
Hg
F r i e d e t sl. (19561
Liquid
Mole F r B c t i o n o f Ethylbenzene
u.0
Figure 8. Temperature-composition diagram for ethylbenzene 2-ethoxyethanol at 50 mmHg pressure. 200
+
Miti of Mlxlng for E t h y l b m z m i t 2-Ethoxyathinol
+
+
.
o
-b
t 1 and Vin Wlnkli IlSm
0
0.0 0.0
Woli F r i c t i o n o f E t h y l b i n z e n i
Figure 9. Heats of mixing for ethylbenzene 25 "C.
1.0
+ 2-ethoxyethanol at
of how the data may be fitted best. One of two procedures which we have found to yield consistantly good results is detailed in this section. In this procedure, the optimization is executed in two steps. First the excess enthalpies are treated in isolation to determine the a, and ck coefficients. This is discussed in detail later. Once this has been done, the phase equilibrium data are then used to calculate the b, coefficients. In the other procedure (also discussed in a later section), all three sets of coefficients are determined in one single simultaneous optimization. Fourteen phase equilibrium and nine enthalpies of mixing references for the hexane acetone system are recorded in our bibliographic collection. Enthalpies of
+
1778 Ind. Eng. Chem. Res., Vol. 26, No. 9, 1987 ~
P'
-
.:bar~:
-
~
E:?ilace:are
-
_ 1
at 1 A~mosprye
Ethanol
+
E t h y l a c e t a t e a t 40 ana 55'C
45C I
'
et al
J
3ai::h
1
%-t.
':96E1
a n d .an H:n