I n d . Eng. Chem. Res. 1992, 31, 1203-1211
Auwers, K.; Jacob, A. Ueber stereoisomere Butantetracarbonsauren. Ber. Dtsch. Chem. Ces. 1894,27,1114-1132. Boistelle, R. Fundamentals of Nucleation and Crystal Growth. In Crystallization and Polymorphism of Fats and Fatty Acids; Garti, N., Sato, K., Eds.; Marcel Dekker: New York, 1988; pp 189-226. Brodmann, G. L. Performance of Nonformaldehyde Cellulose Reactants. Text. Chem. Color. 1990,22(ll),13-16. Brotherton, D. L.; Fung, K. W.; Addison, A. L. An Alternative Zero Formaldehyde Durable Press System. Book Pap.-Znt. Conf. Exhib., AATCC 1989,170-175. Moms, C. E.; Harper, R. J., Jr. Effects of Wet Pickup Variation and Swelling Pretreatment on Properties of Cotton Fabrics Crosslinked with Butanetetracarboxylic Acid. Am. Dyest. Rep. 1991, 80 (9),58, 60,62,64,66. Rowland, S.P.; Welch, C. M.; Brannan, M. A. F.; Gallagher, D. M. Introduction of Ester Cross Links into Cotton Cellulose by a
1203
Rapid Curing Process. Text. Res. J. 1967,37,933-941. Rowland, S. P.; Welch, C. M.; Brannan, M. A. F. US. Patent 3,526,048,1970. Welch, C. M. Tetracarboxylic Acids as Formaldehyde-Free Durable Press Finishing Agents, Part I Catalyst, Additive, and Durability Studies. Text. Res. J. 1988,58,480-486. Welch, C. M. Durable Press Finishing without Formaldehyde. Text. Chem. Color. 1990,22(5),13-16. Welch, C. M.; Andrews, B. A. K. Ester Crosslinks: A Route to High Performance Nonformaldehyde Finishing of Cotton. Text. Chem. Color. 19898,21 (2),13-17. Welch, Clark M.; Andrews, B. K. U.S. Patent 4,820,307,1989b. Received for review September 10,1991 Revised manuscript received December 23, 1991 Accepted January 9, 1992
Thermodynamic Inconsistencies in and Accuracy of Chemical Equations of State for Associating Fluidst Ioannis G.Economou and Marc D. Donohue* Department of Chemical Engineering, The Johns Hopkins University, Baltimore, Maryland 21218
Hydrogen bonding is usually taken into account through the use of chemical theory. Thermodynamic properties of associating systems are calculated by solving the equations for chemical and physical equilibria simultaneously. Recently, many equations of state have been proposed for associating components with a built-in chemical theory based on the approach of Heidemann and Prausnitz. Many of these theories assume that the chemical and physical interactions can be separated in the equation of state. In this work, we review equations of state for associating fluids using chemical theory and show that an arbitrary separation of chemical and physical interactions can be thermodynamically inconsistent. Comparisons are made for the various equations of state with Monte Carlo simulation data.
Introduction Hydrogen bonding has a significant effect on the thermodynamic properties of pure components as well as mixtures. In pure compounds, hydrogen bonding results in the formation of dimers, as, for example, in carboxylic acids, linear chains or closed rings as in alcohols, or three-dimensional structures, as in water. The association decreases the vapor pressure and increases the boiling point of the associating component compared with a nonassociating component of equal size and mass. In mixtures, hydrogen bonding can occur between molecules of the same compound (self-association) and/or between molecules of different compounds (solvation). Strongly selfassociating components show limited mutual miscibility (liquid-liquid equilibria) with nonassociating components whereas solvation between species that do not self-associate results in negative deviations from Raoult’s law (for example in the chloroform-acetone mixture). Therefore, it is not surprising that numerous approaches have been tried to describe the thermodynamic properties and phase behavior of associating and solvating systems. Traditional cubic equations of state, such as the Redlich-Kwong-Soave and the Peng-Robinson equations of state, do not account explicitly for hydrogen bonding, and strong specific interactions are taken into account through the use of temperature- and density-dependent binary parameters. However, this requires a large amount of
* Author to whom correspondence should be addressed. ‘This paper is dedicated to the memory of Professor Basil Iatridis, National Technical University of Athens, Greece.
experimental data and phase predictions based on these models often are unreliable outside the range of experimental data. In addition, the use of many adjustable parameters introduces a thermodynamic inconsistency. In phase equilibrium calculations the fugacity of the component in the liquid phase is calculated by an integral from zero density to the liquid density of the system and the use of different temperature- and density-dependent parameters cause the inconsistency (Robinson et d., 1985). A more powerful approach is to assume that hydrogen bonding leads to the formation of new species and to treat the pure component as a mixture of monomers and oligomers. This method, known as “chemical theory”, was originally introduced near the turn of the century (Dolezalek, 1908), and initially it was used for phase equilibrium calculations through activity coefficient models (Renon and Prausnitz, 1967; Wiehe and Bagley, 1967). Recent efforts using chemical theory in activity coefficient models have been made by Nath and Bender (1981) and Nagata and co-workers (Nagata, 1985; Cibulka and Nagata, 1987a,b; Hofman and Nagata, 1987). This method however has several disadvantages. An arbitrary assumption must be made on the type of species resulting from the associations in the system. In addition, the number of parameters needed increases considerably since the association species are treated as real species. This usually requires a trialand-error procedure in order to solve the chemical equilibrium and phase equilibrium equations simultaneously. Chemical theory also can be used with an equation of state. To eliminate the need of many parameters for the equation of state, assumptions usually are made relating the parameters for the oligomers to the parameters for the
o s a ~ - s a a 5 / 9 2 ~ 2 ~ 3 i - i ~ ~ ~ ~0o 1992 3 . o oAmerican ~o Chemical Society
1204 Ind. Eng. Chem. Res., Vol. 31, No. 4,1992
monomer. Gmehling et al. (1979)used the perturbedhard-chain theory (PHCT) with a monomer-dimer equilibrium scheme to calculate vapor-liquid equilibria of several associating systems. The equation of state parameters for the dimer species were given as a function of the parameters for the monomer. However, this approach is not readily generalized to other chemical equilibrium models. Wenzel et al. (1982)combined a linear association model leading to the formation of species up to dodecamer with a three-parameter cubic equation of state. For pure associating components, the equation of state needs 3(N + 1)nonassociation parameters, where N is the number of different oligomers, and 3N association parameters. It is apparent that this model is impractical for phase equilibrium calculations of more than one associating component because the solution of the nonlinear algebraic equations becomes intractable unless several unrealistic assumptions are made. A more elegant approach is to solve for the chemical equilibria and incorporate the result into the equation of state. Initially demonstrated by Heidemann and Prausnitz (1976)for a simple van der Waals type equation of state, this approach has been used for various equilibrium schemes and for several equations of state. Following Heidemann and Prausnitz (1976),Hu et al. (1984)and Junyang and Ying (1990)used chemical theory combined with the Carnahan-Starling-van der Waals equation of state to calculate Henry’s constants for nonpolar solutes in water and excess properties of mixtures of alcohols, respectively. Hong and Hu (1989)combined the monomer-dimer model with the Redlich-Kwong (RK)equation of state for their phase equilibrium calculations. In these models, chemical and physical interactions are coupled in a complicated way. Ikonomou and Donohue (1986)incorporated the infinite equilibrium model and the monomer-dimer model into the perturbed-anisotropic-chain theory (PACT) and derived a new equation of state, the associated PACT (APACT),that is capable of predicting thermodynamic properties of pure associating componenh as well as mixtures of more than one associating component (Ikonomou and Donohue, 1988). In the APACT, the expressions used to calculate the parameters of the associating species from the parameters of the monomer and the mixing rules in the equation of state were chosen so that the chemical and physical terms in the equation of state can be separated and that the physical terms, namely, the repulsive and the attractive terms, are association independent. However, this is not a general result and it cannot be assumed a priori that the contributions from the chemical and the physical interactions are separable. For other equations of state and for different mixing rules, the repulsive and the attractive terms become association dependent. Elliott et al. (1990),in their simplified version of APACT, relaxed some of the assumptions made in the original derivation of APACT (assumptions that were made specifically to keep the equation simple) and obtained the result that there is a contribution from the repulsive term in the association term of the equation of state. Recently, Economou and Donohue (1991)showed also that if a more realistic assumption is made about the variation of c parameter with the size of an i-mer, both the repulsive and the attractive terms of the equation of state become association dependent and there is a contribution in the association term from both of the physical terms. Unfortunately, many recent attempts to incorporate chemical theory into cubic equations of state have made the a priori assumption that the chemical and physical
terms in the equation of state can be separated (Anderko, 1989a-c;Wenzel and Krop, 1990). In these equations, the chemical terms are evaluated independently and then incorporated into the physical equation of state. Chemical theory also has been combined with lattice theory. Panayiotou (1988,1990)incorporated chemical theory into the Sanchez-Lacombe equation of state (Lacombe and Sanchez, 1976;Sanchez and Lacombe, 1976), a lattice fluid theory, to describe thermodynamic properties and phase equilibria of pure alcohols and mixtures of alcohols with alkanes. The hydrogen bonding expressions are essentially the same as the expressions derived by Heidemann and Prausnitz (1976)and by Ikonomou and Donohue (1986). In this work, we attempt to show how chemical theory should and should not be combined into the equation of state calculations. We review the derivation of Heidemann and Prausnitz (1976)for a simple equation of state and the derivation of Ikonomou and Donohue (1986)for a more complicated equation of state. We also present Anderko’s derivation based on the assumption that the physical and chemical terms in the equation of state are separable (Anderko, 1989a),and we show that this approach leads to a thermodynamic inconsistency. Wenzel and Krop’s equation (Wenzel and Krop, 1990) is similar to that of Anderko, and so it is not presented. Comparisons are made between the different theories and simulation data for hard spheres with one and two bonding sites.
Equations of State Using Chemical Theory As discussed in the Introduction, a number of different researchers have used chemical theory with an equation of state for hydrogen bonding systems. This usually is accomplished by solving the equations for chemical equilibria and physical equilibria simultaneously by a trial and error procedure. Heidemann and Prausnitz (1976) were the first to show that it is possible to solve for the chemical equilibria analytically and combine the result into the equation of state to solve for phase equilibria. Infinite Equilibrium Model. The Heidemann and Prausnitz approach has been applied to a number of different equations of state and for a number of different reaction schemes. We present here a brief review of the derivation for a component with one and two hydrogen bonding sites per molecule. For a pure associating component M with two bonding sites per molecule, it is assumed that monomers “react” to form new species according to the chemical equilibria Mi + Mi + Mz Mz + Mi + M3
and as a result the pure component is treated as a mixture of the associated species. In addition to linear species, cyclic species can be formed as well. The true number of moles, +,and the number of moles that would exist in the system in the absence of any association, no,are defined as m
nT = nl
+ n2 + n3 + ... + ni = Eni
(2)
i=l m
no = n, + 2n2 + 3n3 + ... + ini = Cini i= 1
(3)
where ni is the number of moles of the i-mer. The mole fractions of the associated species are defined as
Ind. Eng. Chem. Res., Vol. 31, No. 4,1992 1205
(4) Conservation of mass requires that the total number of moles, %, and the superficial number of moles, no,are related through the expression
(5) The quantity nT/nOis a measure of the extent of the association. In the absence of any association in the system, %/no is unity whereas for strong associating species %/no is always less than unity. Kempter and Mecke (1940)expressed the equilibrium constant for this chemical equilibria in terms of the volume fractions of the reacting species:
where 8 denotes the volume fraction of the species. Alternatively, Flory (1944)expressed this chemical equilibria in terms of volume fractions as
(7) or in terms of molar concentrations as
The differences between the two expressions for K and K J in terms of volume fractions are due to the difference in the reference states (Painter et al., 1989). According to Flory (1944),an equilibrium constant expressed in terms of mole fractions is not acceptable based on statistical mechanical arguments. This equilibrium constant is a function of the extent of the association and as a result is not a constant. In all the above expressions, it is assumed that the equilibrium constant is independent of the length of the linear species. However, it has been suggested that the equilibrium constant for the alcohols depends on the length of the chain (Barker, 1952). The attenuated equilibrium constant (AEC) model (Garland and Christian, 1975)assumes that the equilibrium constant for the formation of species Mitl depends on the length of the chain according to the expression
of the fugacities of the species according to the expression
where 4 is the fugacity coefficient. Heidemann and Prausnitz (1976)made the simplifying assumption that the equilibrium constant Ki is independent of the degree of association (K2= K3 = ... = K ) . The fugacity coefficient is calculated through classical thermodynamics (Prausnitz et al., 1986) according to the equation In di =
s,T -( -) 1 ap RT ani
T,V,nj,,
$1
d V - In Z
(11)
In order to evaluate the partial derivative in eq 11, an equation of state is needed. Heidemann and Prausnitz used a Carnahan-Starling-van der Waals type equation of state to calculate the fugacity coefficient in eq 11. The following simplifying assumptions were used to calculate the energy parameter ai and the size parameter bi of the associating species Mi:
ai = i2a
(12)
bi = ib (13) where a and b refer to the monomer. The usual Berthelot mixing rules were used to calculate the mixture parameters ( a ) and ( b ) for the pseudomixture (Prausnitz et al., 1986). The ratio P414i/$i+lfinally is obtained
where p is the density, p = l/uo. In eq 14,g is a function of density that depends on the equation of state used to evaluate the fugacity coefficients and on the expressions used to calculate the equation of state parameters for the associating species from the parameters for the monomer species and is given from the expression
(g)T,V,nj+i
-
(
E)T,V,nj+i+l]
-
+)
dV (15)
Equations 12 and 13 and the Berthelot mixing rules used by Heidemann and Prausnitz (1976)give the result that only repulsive interactions contribute to eq 15 and the following expression is obtained:
49 - 392 where KM is a constant. A detailed discussion of the similarities and differences among these theories is given by Acree (1984)and by Painter et al. (1989). In all these models, as well as in the Heidemann and Prausnitz analysis (HP), in APACT, and in more recent models for associating fluids, it is assumed that only linear associating species are formed according to the chemical equilibria described by eq 1. Nonetheless, there is spectroscopic evidence that the alcohols form cyclic species in addition to the linear species (Liddel and Becker, 1957;Becker et al., 1958). However, the addition of another chemical equilibrium and subsequently another equilibrium constant complicates the problem and thus usually is avoided. For real systems where nonidealities are due to the physical interactions in addition to the chemical interactions, the equilibrium constant is written in terms of the activities of the species (Prausnitz et al., 1986)or in terms
(16) (1 - d2 where 9 = nT(b ) / V is the reduced density. On the basis of mixing rules for ( b ) , 9 = b/uo. The expression from eq 14 is substituted back to eq 10. From eq 10 an expression for zi in terms of z1 can be obtained by applying this equation to the formation of dimer, trimer, etc. The result is g=-
( 2 )'
zi = K-RTpeg
zli
The expression for zi is substituted back to the material balance expression (eq 5):
1206 Ind. Eng. Chem. Res., Vol. 31, No. 4, 1992
Another equation can be written based on the constraint C'LlZi = 1: (19)
i=l
Equations 18 and 19 are solved simultaneously for z1and nT/nO,and the result is nT 21=--
no
-
2 1 + (1 + 4KRTpefl)1/2
(21)
Finally the equation of state is given in terms of the compressibility factor, 2,from the expression (Heidemann and Prausnitz, 1976)
z = - '-V
- 1+
p s o c
+ Zrep + Zattr
noRT
(22)
where
For compounds with two hydrogen bonding sites per molecule +/no = (nl/no)1/2 and so 2- is obtained from eqs 21 and 23:
In eq 16, only repulsive forces contribute to g. This is a result of the mixing rules used for ( a ) and the assumptions used for ai (eq 12). However, Hu et al. (1984) showed that eq 12 is not accurate and it overestimates the effect of the chain length on ai. They proposed that a more reasonable assumption is the following: ai= il+wa
(27)
where 0 Iw I1. For w # 1 there is an attractive contribution to the g function in addition to the repulsive one. Junyang and Ying (1990) showed that for the CarnahanStarling-van der Waals equation of state, the best value for w is zero, and so g for this equation is given by the expression 47 - 37= a g=---(1 - ,# bRT'
(28)
Hong and Hu (1989) incorporated the Heidemann and Prausnitz approach into the RK equation of state. In their analysis, Hong and Hu calculate ai and bi from eqs 27 and 13with w = 0. Using Berthelot mixing rules, the following expression results for g: U
+
(30)
(20)
The mole fraction of monomers based on the superficial number of moles is
nl nT 2 =z1= no no 1 + 2KRTpeg + (1 + 4KRTpeg)'l2
equation of state. If one uses the expression for ai proposed by Heidemann and Prausnitz (w= l),this term vanishes. The RK equation of state for an associating component with two bonding sites per molecule is given by eq 22 with 2- calculated from eq 26. +/no and nl/n,, are calculated from eqs 20 and 21, and the repulsive and attractive terms are given by the expressions
g = -In (1 - 1) - -In (1 1) (29) bRT In eq 29, the second term in the right-hand side of the equation is a contribution from the attractive term of the
In eq 31, Zattris association dependent because of the assumption that w = 0. By setting w = 1,Zatt'becomes association independent as in Heidemann and Prausnitz analysis. In the equations of state described above, the physical and the chemical terms in the equation are coupled in a complicated way. We present next an equation of state where the contributions from the association and the physical interactions are separated because of the particular assumptions made to keep the equation simple. In APACT the derivation is similar to the derivation described for the HP equation of state. The number of moles are defined according to eqs 2 and 3, and the equilibrium constant is defined as in eq 10. The perturbed anisotropic chain theory (PACT) of Vimalchand and Donohue (1985) is used to calculate the fugacity coefficient in eq 11. The molecular parameters for the association species are calculated on the basis of the parameters for the monomers using simple linear relations. One fluid van der Waals type mixing rules are used for mixtures. These particular expressions used in APACT result in a zero value for g in eq 14 over the whole range of density. Finally the expression for z1 and nT/nOis given from eq 20 and the expression for nl/no is given from eq 21 where g = 0. However, if one relaxes any of the assumptions made in the original derivation of the APACT, a nonzero function for g is obtained. The equation of state in terms of the compressibility factor, 2, is given from eq 22 where 2"" is calculated from eqs 23 and 26 and
zattr
= ZU + zani
(33) Lennard-Jones and anisotropic interactions are calculated using a perturbation expansion (Vimalchand, 1985). The parameter c takes into account the nonspherical shape of the molecules. The introduction of this parameter was based initially on phenomenological arguments, and following Prigogine, the parameter c was defined as onethird the number of external degrees of freedom of the molecules. However, the generalized Flow equation of state (Dickman and Hall, 1986),which is identical in form to the repulsive term in PHCT and PACT, shows that the parameter c can be related to the molecular geometry and the ratio of the excluded volume for a chain molecule over the excluded volume for a sphere. A linear mixing rule is used to calculate ( c ) : (c)=
&Ci
i=l
(34)
In order to keep the equation simple, Ikonomou and Donohue (1986) made the assumption that the c parameter of an i-mer is simply ci = ic and so the repulsive compressibility factor becomes association independent:
Ind. Eng. Chem. Res., Vol. 31, No. 4, 1992 1207 (35) Panayiotou (1988) derived eq 20 for b / n oand eq 21 for nl/noindependently using a completely different thermodynamic formalism. However, his derivation was based on the same principles as the derivation presented above. The equilibrium constant was independent of the chain length, the associated species were treated as real species, and as a result the equation of state can be used to calculate thermodynamic properties of these species. The resulting expressions were incorporated into the Sanchez-Lacombe equation of state (Lacombe and Sanchez, 1976; Sanchez and Lacombe, 1976). Monomer-Dimer Model. For compounds with one hydrogen bonding site per molecule, chemical equilibria between the reacting sites result in the formation only of dimer species. The equilibrium between monomers and dimers is described by the expression (36) The ratio of fugacity coefficients again is evaluated using eq 11,and z2 can be expressed in terms of zl.The material balance is of the form no = 21 + 222 (37) nT and is solved simultaneously with the equation 21 + 22 = 1 (38) to obtain the result 2(KRTpeS - 1) nT - 4KRTpeg - 1- (1 + 8 K R T ~ e g ) l / ~ (39) no
--
From eqs 39 and 40, an expression for nl/nois obtained nl 2 _ -- 2 - nT- 1 = (41) no no 1 (1 8KRTpeB)1/2
+ +
The association compressibility factor, Z"", for one hydrogen bonding site per molecule is obtained through eqs 23, 39, and 41:
Equations 36-42 were derived on the basis of the APACT (Ikonomou and Donohue, 1986). However, if one applies the monomer-dimer equilibrium to the H P and to the RK equations of state, the same expressions for +/no (eq 39), for z1(eq 40), for nl/no(eq 41), and for (eq 42) are obtained. The only difference between these equations of state is that in APACT g is zero whereas for HP and RK g has a nonzero value. In APACT, the repulsive and attractive compressibility factors are association independent (Ikonomou and Donohue, 1986) whereas in the H P and in RK equations of state P P is association dependent and Zattrcan be association dependent or independent based on the assumption used for w. Recent Equations of State for Associating Fluids. Elliott et al. (1990) relaxed the assumption about the variation of the parameter c with the extent of association made in APACT. In their equation of state, referred to as ESD hereafter, they did not make an a priori assumption about the variation of ci with i. Instead, ci was cal-
culated from the expressions used for the other molecular parameters for the associating chains. Finally, the ESD equation of state is given from eq 22 where ZBssoc is calculated from eq 23 and (43) Zettr
=-
9.49qqY 1 + 1.7745Yq
(44)
where q is a shape parameter and Y is given by Elliott et al. (1990). Due to the mixing rules used, both the repulsive and the attractive terms in the ESD equation of state are association independent. The ESD equation of state derivation for associating systems was based on the approach used in APACT, and as a result the expressions for %/no, nl/no,and z1are identical to the expressions derived for APACT and are not repeated. The only difference is the functional form of g which depends on the physical terms of the equation of state and on the mixing rules. In ESD, g is a nonzero function given by g = -In (1- 1.9s) (45) Recently Economou and Donohue (1991) relaxed the assumption about the ci parameter and rederived APACT. Valuea for the parameter c for chains of various length were determined from the generalized Flory-dimer (GF-D) equation of state of Honnel and Hall (1989). In the new equation of state, referred to as APACT-2, g is a nonzero function and P P and ZEk become association dependent. Although more complicated than APACT and for engineering calculations intractable, APACT-2 is in quantitative agreement with molecular simulation data and with the thermodynamic perturbation theory of Wertheim (Wertheim, 1986; Jackson, et al., 1988) for associating fluids. By setting the fugacity coefficients equal to unity in eq 10, the ideal gas equilibrium constant is obtained. Similarly, in ideal solution theory the equilibrium constant is written in terms of molar concentrations according to the expression
For both cases, the derivation for nT/nOand nl/nofor componenta with one and two bonding sites per molecule results in the same expressions as the other chemical theories. However,g has a zero value over the whole range of density. As a result, the expressions obtained from APACT, ideal gas law, and ideal solution theory are identical for these quantities. However, there is a fundamental difference between these approaches. The result g = 0 in APACT was not made as an a priori assumption, but it is a result of the mixing rules and of the particular assumptions made in the derivation of the equation. By changing any of these assumptions a nonzero vdue for g is obtained. In addition, the ratio of fugacity coefficienta is equal to unity but the fugacity coefficients themselves are not unity. This result has been misunderstood by Wenzel and Krop (1990). On the other hand, in the ideal gas calculations the fugacity coefficients are unity and this also results in a zero value for g. Finally, in the ideal gas law the repulsive and the attractive contributions to the compressibility factor, P P and Zaw,are zero and only 2"" has a nonzero value. Anderko (1989a) derived an equation of state based on the a priori assumption that chemical and physical terms can be separated in the equation for associating compo-
1208 Ind. Eng. Chem. Res., Vol. 31, No. 4, 1992
nents. Accordingly, he obtained the expressions presented above for %/noand nl/nowith eg = 1. To account for the physical interactions, the three-parameter cubic equation of state of Yu et al. (1987) (YLI), was used. His equation of state is given from eq 22 where FeP and Zatt'are association independent and given by the expressions p
p
=
9 1-9
(47)
1.50
1
-IC,
= -- U bRT 1
t
1.25
9
'
1
'
1
'
RK. Y L I , ESD
+
-
3
N
\ s_ 1.00
s,
-.
E?,
+ v ( d / b + 3) + v 2 d / b
1
- - - Anderko
