10080
Ind. Eng. Chem. Res. 2008, 47, 10080–10085
Perspective for Extended Group-Contribution Methods for the Prediction of Activity Coefficients Accounting for Sterical Effects. Linking UNIFAC and the Group Vector Space Method Stephan Machefer,*,† Martin Grasemann,‡ and Klaus Schnitzlein† Department of Chemical Reaction Engineering, Brandenburg UniVersity of Technology, Cottbus, Burger Chaussee 2, 03046 Cottbus, Germany, and Institute of Chemical Sciences and Engineering, Ecole polytechnique fe´de´rale de Lausanne, CH H3 594 Station 6, CH-1015 Lausanne, Switzerland
Group contribution methods for the calculation of an activity coefficient such as UNIFAC (universal factorial activity coefficient) currently do not take into account sterical effects. Consequently one of the major deficiencies of such methods is their inability to resolve isomeric systems and their often insufficient accuracy for systems including molecules of very different size. In this paper an approach is presented by which sterical information can be implemented into those methods by means of a modified group vector space (GVS) procedure. (Wen, X.; Qiang, Y. Group Vector Space method for estimating melting and boiling points of organic compounds. Ind. Eng. Chem. Res. 2002, 41, 5534-5537.) Exemplary calculations with UNIFAC/GVS underline the potential of this new approach resulting in a perspective for substantially improved group contributions methods with a wider range of applicability and even more reliable predictive performance. 1. Introduction During the last three decades, the UNIFAC model and its modifications have advanced to some kind of standard tool for the prediction of activity coefficients of nonideal, nonelectrolyte mixtures. As a group contribution methods, UNIFAC describes a system as mixture of functional groups. It is assumed that the nonideality of the mixture arises from combinatorial and residual contributions. The entropic part is calculated from Van der Waals (VdW) properties of the functional groups. For the calculation of the residual part, a semiempirical model is used that is based on binary interaction parameters of those functional groups. Consequently, UNIFAC has to be fitted to a database covering different sources of thermodynamic data in order to get representative group interaction parameters.1 The most powerful modification of UNIFAC (mUNIFAC), introduced an extended temperature relationship and modified the combinatorial part empirically.2 In this way, different kinds of phase equilibrium data (VLE, LLE, SLE,. . .) could be calculated with a single parameter set instead of case-specific sets (UNIFAC). Altogether the mUNIFAC has shown improved performance especially regarding temperature dependency and the calculation of asymmetric mixtures. However, all improvements were achieved at the expense of introducing more binary parameters as well as adjusting the combinatorial part and some of the VdW properties empirically. Consequently the mUNIFAC database has been successively extended, currently covering several tens of thousands of data sets of very different sources (VLE, LLE, and SLE data as well as heats of mixing and azeotropic data) for parameter adjustment and model validation.3,4 This huge and heterogeneous data background is the main advantage of mUNIFAC when compared to similar semiempirical group contribution methods like ASOG,5 DISQUAC,6 or other UNIFAC clones.7,8 However, the more empirical parameters that are introduced and the more the database is extended, the more one departs * To whom correspondence should be addressed. E-mail: s.machefer@ rt.tu-cottbus.de. Phone: +49 355 691178. Fax: +49 355 691110. † Brandenburg University of Technology. ‡ Ecole polytechnique fe´de´rale de Lausanne.
from the primary intention of group contribution methods. Instead, more substantial improvements should be the target of future developments in order to improve the performance of these methods in the ”knowledge vacuum”,9 i.e., for classes of mixtures not included in the parameter estimation procedure. For example, asymmetric mixtures containing molecules which are very different in size sometimes lead to poor predictions.3 Although this drawback has been improved empirically with mUNIFAC as mentioned earlier, the core of those problems is not approached: the lack of sterical information in UNIFAC, mUNIFAC, or any other of the existent group contribution methods. Consequently, isomeric mixtures cannot be resolved either. This is why a priori predictive but computationally extensive methods like COSMO-RS gain attractiveness.10 While mUNIFAC has recently been shown to still be superior to COSMO-RS in an extensive model-comparison,11 COSMO-RS is in rapid development and those results may not be representative for the current version of this method.12,13 It would be therefore highly desirable to introduce sterical information in group contribution methods. Appropriate approaches are very rare and mostly relatively complex. For example Monte Carlo simulation is applied for the modification of interaction partameters.14 Another approach focuses on isomers in terms of conformations of chainlike molecules for which the percentage of each conformer has to be assumed.15 Obviously this approach gets highly complicated for larger molecules. However, in this contribution we present a simple, staightforward, and universal concept: A group contribution method (i.e., UNIFAC) is extended and linked with the group vector space (GVS) method16 which weights a functional group’s sterical position inside a molecule. 2. Theory 2.1. UNIFAC. For the theoretical analysis of this work we chose the original UNIFAC1 as an exemplary group contribution method since its VdW properties have not been empirically adjusted (see above). The UNIFAC method is derived from the UNIQUAC local composition activity coefficient model.17 In
10.1021/ie800980m CCC: $40.75 2008 American Chemical Society Published on Web 11/19/2008
Ind. Eng. Chem. Res., Vol. 47, No. 24, 2008 10081
Figure 1. Sterical weighting factors ξ˜ l for some example molecules (nξ ) 1): top, p- and o-cresol; middle, 1-, 2-, 3-, and 4-octanol; bottom, 1- and 2-butanol.
UNIQUAC the activity coefficient of a mixture is derived from considerations of the local composition of a microstructure by means of statistical thermodynamics and certain assumptions. Combinatorial (superscript C) and residual aspects (superscript R) contribute to the overall activity coefficient of species i:
(
Vi Vi z ln γiC ) 1 - Vi + ln Vi - qi 1 - + ln 2 Fi Fi
)
(2)
z is usually set to 10. Volume fraction Vi and surface fraction Fi are related to the (normalized) VdW volumes and surface, respectively (ri and qi):
∑
Fi )
rx j j j
qi
∑
(3)
qx j j j
UNIFAC approximates the VdW properties from group increments (index k): qi )
∑
k
ri )
νk(i)Qk
∑
k
νk(i)Rk
(4)
The residual contribution is calculated from group residual activity coefficients Γk with the pure species (superscript (i)) as reference: ln γi ) R
∑
ν (ln Γk - ln Γk ) k k (i)
∑
m
Θm*Ψmk) -
∑ m
{
(1)
The combinatorial contribution is calculated from the number of configurations of a local composition in a lattice of coordination number z giving
ri
[
ln Γk ) Qk* 1 - ln(
Θm*Ψkm
∑
Θn*Ψnm
n
]
(6)
with
ln γi ) ln γiC + ln γiR
Vi )
where Γk is derived from the lattice energy of the mixture
(i)
(5)
Ψnm ) exp -
∆unm T
}
(7)
The interaction energy ∆unm between unlike functional groups has been used as empirical parameter and fitted to experimental data. Thus among interactive forces, different other aspects, such as sterical influences, are lumped together in this parameter. This is why a method is needed that explicitely accounts for sterical aspects. 2.2. Sterically Effective Surface of Interaction. In this contribution we propose a new approach to uncouple ∆unm from sterical impacts. We assume that the position of a functional group inside the molecule and its degree of shielding by surrounding groups affects its interaction with other groups. In other words, sterically undisturbed groups provide more surface for interaction than hidden groups. Consequently and in contrast to the combinatorial contribution, even identical functional groups are treated individually for the residual contribution of UNIFAC. Each functional group provides an effective surface fraction (asterisk): Θm* )
Qm*Xm
∑
Q *Xn n n
Xm
∑ ) ∑
(j)
j
νm
xj
xN j j j
UF
(8)
10082 Ind. Eng. Chem. Res., Vol. 47, No. 24, 2008
Figure 2. Activity coefficients for three binary isomeric mixtures: lines, UNIFAC/GVS calculations for the activity coefficient with different values for nξ and nR as well as two different group-matching procedures (av, max); circles, experimental data for γ1; triangles, experimental data for γ2.21
where NjUF is the overall number of functional groups in species j according to the UNIFAC definition of functional groups. Note that, for simplicity, all groups, even those with same Qk* value are treated individually according to eq 8. νm( j ) is therefore either 1 or 0. However, a computationally more effective definition would combine in kind and stericity identical groups. The effective VdW-surface parameter Qk* is calculated by means of the Group Vector Space procedure as follows. 2.3. Group Vector Space Method (GVS). The GVS has been shown to improve group contribution methods for the prediction of pure species properties, such as melting and boiling points,16 enthalpies of vaporization,18 or critical properties.19 In this work we apply a slightly modified version of the published GVS to the calculation of the sterically corrected VdW surface parameter Qk* for the residual part of the activity coefficient (eqs 5-8). A sterical weighting factor ξk is assigned to each functional group giving Qk* ) Qkξ˜ k
(9)
Herein ξ˜ is defined as a normalized value and is calculated from the GVS group position parameters λl: ξ˜ l )
λlnξNiGVS
∑
λ l l
nξ
(
1 Nl
∑
l
)
ξ˜ l ≡ 1
(10)
Note that the GVS related properties (index l) should be, strictly speaking, species specific (index i, e.g., ξ˜ l, i instead of ξ˜ l). However, in order to keep the equations simple, the species index i is omitted in the following. The power nξ g 0 introduced here should be treated as a global parameter. For nξ > 1 the sterical impact is enhanced whereas for nξ < 1 it is reduced.
Consequently for nξ ) 0 we return to the original UNIFAC. λl is calculated from the GVS module indices Rl: λR,l )
RlnRRl
√∑
l
nR ) 1 ⁄ 3, 1
(11)
(Rl Rl) nR
2
Note that in contrast to the original definition of λl,16 the VdW volume Rl has been introduced in order to take into account the groups different sizes. The global parameter nR expresses which dimension of the group is used (e.g., nR ) 1/3 for length or nR ) 1 for volume). Note that nR ) 0 corresponds to the original GVS procedure. The GVS module index Rl is calculated from the topological space matrix entries bly (l rows, y columns). Rl2 )
∑
y
bly2
(12)
The systematic procedure for setting up the topological space matrix from the chemical structure of a molecule is described in detail elsewhere.16 2.4. Linking UNIFAC and GVS. When applying eq 9, one faces the problem that the definition of functional groups is not the same for UNIFAC and GVS, respectively. GVS considers mainly basic groups whereas UNIFAC has a more differentiated definition. For example, functional groups which are attached to an aromatic ring are merged with the corresponding ring-atom in UNIFAC (AC-CH3, AC-OH, etc.) because of the special impact of delocalized electrons on energetic interactions. On the contrary, GVS would consider those groups separately (C, CH3, and C, OH, respectively). Consequently group-matching conditions have to be defined for those cases.
Ind. Eng. Chem. Res., Vol. 47, No. 24, 2008 10083
Figure 3. Activity coefficients for the o-cresol/m-divinylbenzene (left) and p-cresol/m-divinylbenzene (right). UNIFAC/GVS calculations for the activity coefficient with different values of nR and different group-matching procedures (av, max) for nξ ) 4: circles, experimental data for γ1; triangles, experimental data for γ2.21
If groups are identically defined for UNIFAC (index k) and GVS (index l), their weighting factors are also identical in first instance ξ˜ k ) ξ˜ l
∀
k≡l
(13)
In cases where UNIFAC and GVS groups are not identical, we discuss two possible matching definitions. First, the weighting factor for the UNIFAC group is calculated from an averaged value for the merged GVS groups (number Nw) and subsequent renormalization of all groups (subscript av): ξk,av )
(ξ˜ l1 + ξ˜ l2 + ... + ξ˜ lw) Nw
(14)
and ξ˜ k,av )
ξk,avNk
∑
k
ξk,av
for k ) l1, l2, ..., lw A second possibility consists of adopting the maximum weighting factor with subsequent renormalization (index max): ξk,max ) max(ξ˜ l1, ξ˜ l2, ..., ξ˜ lw)
(15)
and ξ˜ k,max )
ξk,maxNk
∑
k
ξk,max
for k ) l1, l2, ..., lw Note that in consequence of the renormalization of the matching procedure, the ξ˜ -value of identically defined groups (eq 13) changes afterward. 3. Examples Figure 1 shows the sterical weighting factors for different exemplary molecules as calculated with eqs 9-12 for nξ ) 1.
An analysis of the sterical weighting factors shows that different kinds of stericity are considered. The shielding of nearest neighbors can be exemplary demonstrated for the two cresols: The -OH and -CH3 groups in p-cresol have a higher sterical weighting factor than those in o-cresol. The impact of molecular size or chain length can be shown when comparing 2-octanol with 2-butanol. The sterical weight of the -OH group in relation to the furthermost -CH3 group decreases with increasing chain length: ξ˜ (1)/ξ˜ (5) ) 0.88 for 2-butanol whereas ξ˜ (1)/ξ(9) ) 0.76 for 2-octanol. The impact factors in Figure 1 are given for nR ) 1/3 (modified GVS procedure) and nR ) 0 (original GVS procedure). A slightly more realistic picture of the stericity is given for nR ) 1/3 (see group 3 in relation to group 4 in p-cresol for example). We would like to illustrate the impact of the UNIFAC/GVS combination on the calculation of the activity coefficient with some examples. The quantitative application of UNIFAC/GVS requires a readjustment of the interaction parameters (∆unm) giving values uncoupled from stericity. However, for the scope of this work we use the original UNIFAC interaction parameters.20 UNIFAC does not resolve isomers. This means that γ ) 1 is calculated for any kind and composition of isomeric mixtures. With UNIFAC/GVS it is now possible to calculate activity coefficients of isomeric mixtures as shown in Figure 2. The activity coefficients of three binary isomeric mixtures are correctly predicted with UNIFAC/GVS from a qualitative point of view when compared to experimental values which have been calculated from VLE data.21 The quantitative effect of UNIFAC/ GVS, however, strongly depends on the global parameters nξ and nR as well as on the choice of matching definition (av, eq 14; max, eq 15). The activity coefficients are mostly influenced by the parameter nξ. The power nR and the choice of matching definition have only minor impact. Since the qualitative effect on the activity coefficient seems to be similar, the latter parameters may be chosen arbitrarily, leaving nξ as single global parameter. Figure 3 shows the very different behavior of o- and p-cresol on the experimental activity coefficient in a binary mixture with
10084 Ind. Eng. Chem. Res., Vol. 47, No. 24, 2008
match UNIFAC and GVS functional groups (eqs 14 and 15). This way a new set of binary group interaction parameters could be generated which is uncoupled from sterical influences. Consequently we expect a better predictive performance of UNIFAC/GVS when compared to UNIFAC. This means that systems which have not been included into the parameter estimation procedure would be calculated with more precision. Furthermore, we believe that the UNIFAC/GVS would contribute to substantially improve two typical UNIFAC deficiencies: the inability or inaccuracy to calculate isomeric and asymmetric systems, respectively. Although the original UNIFAC was exemplary regarded in this work, the UNIFAC/GVS principles are transferable to any modification of UNIFAC or to other group contribution methods. Nomenclature
Figure 4. Activity coefficients at infinite dilution for different heptane/ octanol mixtures. UNIFAC/GVS calculations for the activity coefficient with different values of nR and nξ: circles, experimental data for γ1; triangles, experimental data for γ2.21
m-divinylbenzene.21 Again UNIFAC would virtually calculate the same activity coefficients for both systems. Small differences solely arise from the different boiling temperatures of p- and o-cresol (isobaric measurements). UNIFAC/GVS on the other hand clearly resolves the higher activity coefficients for the p-cresol/m-divinylbenzene system when compared to the ocresol/m-divinylbenzene system for all parameter variants. Finally the activity coefficient at infinite dilution for different octanol/n-heptane mixtures is shown in Figure 4. Since interactions with the -OH group mainly contribute to the activity coefficient, the impact of the -OH group decreases with increasing distance from the chain end (see Figure 1). Consequently UNIFAC/GVS calculates the highest activity coefficient for 1-octanol/n-heptane and systematically decreasing values for mixtures with 2-, 3-, and 4-octanol. Note that UNIFAC and GVS groups are all identical for this example; group matching is therefore not required. 4. Conclusions The UNIFAC method is combined with the principles of the group vector space (GVS) procedure16 providing a new approach (UNIFAC/GVS) which considers stericity. A sterically corrected surface of energetic interaction is calculated for each individual functional group of each molecule which takes into account the groups’ position inside the molecule. By means of exemplary binary systems, including isomers, the potential of the UNIFAC/ GVS method is demonstrated. Although the original UNIFAC binary interaction parameters were used, the UNIFAC/GVS method shows qualitative agreement with experimental observations. However, readjustment of the UNIFAC binary interaction parameters (eq 8) will be inevitable in order to evaluate the actual quantitative improvement of UNIFAC/GVS compared to UNIFAC. For this purpose two global parameters (nξ, eq 10 and nR, eq 11) were introduced which increase the flexibility of UNIFAC/GVS. Furthermore two methods were proposed to
F ) surface fraction (m2/m2) N ) overall number of functional groups Q ) normalized VdW surface (group) R ) normalized VdW volume (group) T ) temperature (K) V ) volume fraction (m3/m3) X ) molar fraction (group) (mol/mol) b ) topological GVS space matrix nξ ) sterical impact power (GVS) nR ) group volume impact power (GVS) q ) normalized VdW surface (species) r ) normalized VdW volume (species) ∆u ) interaction energy (T independent) (K) x ) mole fraction (species) (mol/mol) z ) lattice coordination number R ) GVS module γ ) activity coefficient λ ) GVS group position ξ ) sterical weighting factor ξ˜ ) sterical weighting factor (normalized) Θ ) surface fraction (groups) (m2/m2) Ψ ) interaction parameter (T dependent) Subscripts i, j ) species indices k, m, n ) group indice (UNIFAC) l ) group index (GVS) w ) matched group index y ) column index (matrix b) Superscripts * ) sterically corrected property 0 ) infinite dilution C ) combinatorial R ) residual (i) ) pure species property UF ) UNIFAC specific property GVS ) GVS specific property AbbreViations ASOG ) analytical solution of groups COSMO-RS ) conductor-like screening models for real solvents DISQUAC ) dispersive quasi-chemical GVS ) group vector space LLE ) liquid-liquid equilibrium SLE ) solid-liquid equilibrium UNIFAC ) universal factorial activity coefficient
Ind. Eng. Chem. Res., Vol. 47, No. 24, 2008 10085 mUNIFAC ) modified UNIFAC (Dortmund) UNIQUAC ) universal quasi-chemical VdW ) Van der Waals VLE ) vapor-liquid equilibrium av ) group matching by averaging max ) group matching by maximal impact
Literature Cited (1) Fredenslund, A.; Jones, L. J.; Prausnitz, J. M. Group-contribution estimation of activity coefficient in nonideal liquid mixtures. AIChE J. 1975, 21, 1086–1099. (2) Weidlich, J.; Gmehling, J. A modified UNIFAC model. 1. Prediction of VLE, hE and γ∞. Ind. Eng. Chem. Res. 1987, 26, 1372–1381. (3) Lohmann, J.; Joh, R.; Gmehling, J. From UNIFAC to modified UNIFAC (Dortmund). Ind. Eng. Chem. Res. 2001, 40, 957–964. (4) Voutsas, E. C.; Tassios, D. P. Prediction of infinite-dilution activity coefficients in binary mixtures with UNIFAC. A critical evaluation. Ind. Eng. Chem. Res. 1996, 35, 1438. (5) Wilson, G. M.; Deal, C. H. From activity coefficients and molecular structure. Activity coefficients in changing environments-solutions of groups. Ind. Eng. Chem. Fundam. 1962, 1, 20–23. (6) Kehiaian, H. V. Thermodynamics of binary organic liquid mixtures. Pure Appl. Chem. 1985, 57, 15–30. (7) Larsen, B. L.; Rasmussen, P. A modified UNIFAC group-contribution model for prediction of phase equilibria and heats of mixing. Ind. Eng. Chem. Res. 1987, 26, 2274–2286. (8) Torres-Marchal, C.; Cantalino, A. L. Industrial applications of UNIFAC. Fluid Phase Equilib. 1986, 29, 69–76. (9) Gupta, S.; Olson, J. D. Industrial Needs in Physical Properties. Ind. Eng. Chem. Res. 2003, 42, 6359–6374. (10) Arlt, W.; Spuhl, O.; Klamt, A. Challenges in thermodynamics. Chem. Eng. Process. 2004, 43, 221–238. (11) Grensemann, H.; Gmehling, J. Performance of a conductor-like screening model for real solvents model in comparison to classical group contribution methods. Ind. Eng. Chem. Res. 2005, 44, 1610–1624.
(12) Klamt, A. Comments on ”Performance of a conductor-like screening model for real solvents model in comparison to classical group contribution methods”. Ind. Eng. Chem. Res. 2005, 44, 7042–7042. (13) Grensemann, H.; Gmehling, J. Rebuttal to the comments of Andreas Klamt on “Performance of a COSMO-RS model in comparison to classical group contribution methods”. Ind. Eng. Chem. Res. 2005, 44, 7043–7044. (14) Iwai, Y.; Mori, Y.; Uchida, H.; Arai, Y. Modification of parameters in group contribution method to distinguish isomers based on information by Monte Carlo simulation. Fluid Phase Equilib. 1999, 158-160, 357– 366. (15) Chen, J.; Li, Z.; Li, Y. A method for prediction VLE and HE of chain like isomers using a modified UNIFAC model. Fluid Phase Equilib. 1992, 74, 47–65. (16) Wen, X.; Qiang, Y. Group Vector Space method for estimating melting and boiling points of organic compounds. Ind. Eng. Chem. Res. 2002, 41, 5534–5537. (17) Abrams, D. S.; Prausnitz, J. M. Statistical thermodynamics of liquid mixtures: A new expression for the excess Gibbs energy of partly or completely miscible systems. AIChE J. 1975, 21, 116–128. (18) Wenying, W.; Jinyu, Y.; Wen, X. Group Vector Space method for estimating the normal boiling temperature and enthalpy of vaporization of organic compounds. J. Chem. Eng. Data 2004, 49, 1249–1253. (19) Wen, X.; Wenying, Y. Group Vector Space method for estimating critical properties of organic compounds. Ind. Eng. Chem. Res. 2003, 42, 6258–6262. (20) Poling, B. E.; Prausnitz, J. M.; O’Connel, J. P. The Properties of Gases and Liquids; McGraw-Hill: Boston, 2001. (21) Gmehling, J.; Onken, U.; Arlt, W.; Grenzhauser, P.; Weidlich, U.; Kolbe, B.; Rarey, J. Vapor-Liquid Equilibrium Data Collection; Chemistry Data Series, Volume I; Dechema: Frankfurt (Main), 1978-2006.
ReceiVed for reView June 24, 2008 ReVised manuscript receiVed October 6, 2008 Accepted October 7, 2008 IE800980M
