THERMODYNAMICS OF SOLUTIONS Equation of State and Vapor Pressure OTTO R E D L I C H , F. J. A C K E R M A N , R. D. G U N N , M A X JACOBSON,AND SILVANUS LAU
Inorganic Materials Research Division, Lawrence Radiation Laboratory, and Department of Chemical Engineering, University of California, Berkeley, Calif.
A new improvement of the equation of state of Redlich and Kwong extends the application to the liquid state and vapor pressures. It leads to good results down to a reduced temperature of 0.65. The reasonable limiting relations for extreme temperatures or pressures of the original equation are maintained except for the case where both T and P at the same time assume extremely high values. Mean and individual fugacity coefficients are consistently derived from the equation of state without further assumptions or parameters. A flexible computer program provides convenient use for a variety of problems. HE practical interest in a n equation of state is mainly tied T t o the calculation of individual fugacity coefficients of mixtures since they are required for the investigation of chemical equilibria and phase equilibria. They can hardly be obtained by any other way but by derivation from a n equation of state. T h e results obtained by a n improvement (72) of the twoparameter equation of Redlich and Kwong (73) invited further development. T h e introduction of Pitzer's acentric factor (3, 70) as a third parameter in a "deviation function," added to the old equation for the compressibility factor, proved to be useful. But a considerable expansion of the scope and various further improvements appeared desirable. T h e use of automatic computation has removed almost all holds regarding the complexity of the deviation function, provided that no fourth individual parameter is introduced. I t is obvious that each new individual parameter improves the representation of a given set of d a t a and a t the same time tends to reduce the ability of a relation to predict the behavior of a substance outside the experimental range. Comparisons (7, 4, 5, 76) of the old two-parameter equation with the eightparameter relation of Benedict, Webb, and Rubin confirmed this reasoning much more strongly than was expected. A fourth parameter will be useful therefore, if a t all, only if it has a good theoretical meaning. Practically the only computational restriction to be imposed o n the deviation function is the algebraic integrability with respect to log P, desirable for the calculation of fugacity coefficients. T h e more important restrictions, those of a general nature, have been discussed (72). I n the preceding papers the scope was restricted to the gaseous state. The original equation, however, deliberately chosen to be of third degree, offered solutions for the liquid state that were reasonable though of course not accurate. Encouraged by these results as well as by the success of Pitzer's tables, we have now extended the program to liquids. T h e key to a n equation of state encompassing both phases is the vapor pressure, particularly if the attention is concentrated on fugacities. I t must be represented reasonably well. Conversely, vapor pressures are a n important help in the search for a n additional deviation function for the fugacity coefficients of the liquid phase. I n particular, the deviation function should closely satisfy the definition of the acentric factor
TV = -1 - logp,
(1)
based on the reduced vapor pressure, p7, a t the reduced temperature, T , = 0.7. The deviation functions for the gas and the liquid must of course become identical a t the critical temperature. I n addition, one wishes to refer the fugacity of the liquid to the same standard state as the gas. T h e combination of parameters chosen for mixtures differs somewhat from the usual linear and square-root combinations whenever the acentric factors of the components are different. The following empirical combination appears to represent the observed data more closely. I t entails some inconvenience in the differentiation of deviation functions with respect to the mole fractions, required in the calculation of individual fugacity coefficients. T h e combination rule will be expressed as a definition of the pseudocritical values, T3{ and PAr, and the acentric factor, T.V, , of the mixture. They are given as functions of the mole fractions, 2%)the critical constants, T,i and Pci, and the acentric factors, Wi, of the components with the aid of four auxiliary quantities by i
(3) (4)
Pitzer's tables provided considerable guidance in the development of suitable deviation functions. T h e present program goes beyond them by the extension to lower temperatures, by provision of reasonable limiting relations for high temperatures or high pressures, and by the consistent computation of individual fugacity coefficients without any further assumption or empirical parameter. T h e deviation functions fail, however, if both temperature and pressure are simultaneously increased to extreme values. T h e convenience of automatic computation overbalances the preparatory work of the user even in moderately frequent application. VOL. 4
NO. 4
NOVEMBER 1 9 6 5
369
Computer Program
The program, written in Fortran I1 language, can be used advantageously only in a high speed machine. A considerable effort has been made to achieve convenience in input a n d output as well as highest flexibility, so that one can solve a variety of problems without prescribing computations unnecessary for the specific problem. T h e input comprises the three parameters for each substance, a schedule of temperatures and pressures that may be prescribed in a variety of ways, and for mixtures the compositions in mole fractions. As an optional input, P-V-T data can be introduced for comparison with calculated values. T h e computation furnishes the compressibility factor and if, desired, the fugacity coefficient and the pressure-dependent terms of the entropy and enthalpy. For mixtures, in addition, individual fugacity coefficients and fugacities may be obtained. Vapor pressures of pure liquids can be prescribed by a special option; they are computed by a n efficient step-by-step approximation. In complicated integrations such as those of the deviation functions with respect to log P, the computer offers the advantage that the whole operation can be resolved in a series of intermediate steps, so that no really involved and lengthy expressions need be written down. Since the arc tan function during integration may change from one range to another, a special subroutine Angle has been added which eliminates its ambiguity. All integrations and differentiations have been carried out algebraically with the exception of the differentiation of deviation functions with respect to temperature. Such a differentiation was practically impossible. A numerical differentiation is sufficient if the temperature increment is chosen so low that the function is practically linear, and a t the same time so high that the division by the increment does not lead to a n intolerable loss of the number of significant figures. This condition can be easily checked by means of the program. T h e built-in standard increment of 0.1 will probably always be satisfactory for K. or R. T h e option provided for other increments should be used if reduced temperatures are introduced. T h e solution of the original equation of Redlich and Kwong leads sometimes to unsuitable results-Le., instable states or negative values. These results are automatically disregarded. T h e program does not allow computations for negative pressure (tension). If the computation for a pure substance leads to results for two phases, the comparison of the fugacity coefficients (mean fugacity coefficients for mixtures) indicates which is relatively stable and which is metastable. For mixtures, the relatively stable phase may still be metastable with respect to a twophase system. A special program is under consideration that furnishes the compositions and properties of phases in equilibrium a t a given temperature and pressure. Correctness and Accuracy
3. The correctness of the differentiation with respect to T is checked by comparing results for entropy and enthalpy with fugacity coefficients computed for a closely spaced set of temperatures. 4. A necessary and sufficient check of the computation of the individual fugacity coefficient, 41 and 42, of a binary mixture is given by their relations with the mean fugacity coefficient, $, and the mole fractions, ~1 and y~ log 9 = yi log
370
l&EC FUNDAMENTALS
+ yz log
d log $ldYl = log
(9)
$2
(10) These checks have been found efficient in disclosing program errors. They are believed to be sufficient. (91/$2)
Results
The vapor pressure comparisons in Table I are given for practically the whole range of the acentric factor W. The comparisons of Z in Figures 1 to 4 and the data in Table I1 Table 1.
Vapor Pressures (6)
p , Atm. T , ' K.
Obsd.
Calcd.
0 . 0 1 3)
113.15 125.15 137.15 149.15 161.15 173.15 185.15
1.129 2.686 5.474 9.906 16.501 25.744 38,382
1.001 2.597 5.488 10,104 16.849 26.121 38,385
Nitrogen ( W = 0.034)
81.15 89.15 97.15 105.15 113.15 121.15
1.526 3.317 6.306 10.844 17.342 26.295
1.493 3.314 6.353 10,960 17.466 26.133
83.15 95.15 107.15 119.15 127.15
1.183 3.686 8.797 17.692 26.413 1.01 5.38 17.57 42.94 62.55 1.180 4.238 J1.512 25.817 50.523 61.780 1,001 3.492 9.272 20.308 39.095
1.155 3.690 8.896 17.881 26.363 1.04 5.38 17.49 43.10 62.16 1.415 4.355 11.206 25.305 50.644 61.828 1 ,186 3.547 8.952 19,889 39.010
Methane ( W
=
Carbon monoxide ( W = 0.049)
Hydrogensulfide(W = 0.100)
213.15 253.15 293.15 333.15 353.15
Ammonia ( W
243.15 273.15 303.15 333.15 363.15 373.15
=
0.256)
Sulfur dioxide ( W = 0.256)
263.15 295.15 327.15 359.15 391.15
~~~~~
Table II. Vapor-liquid Equilibria (14) (v = vapor; 1 = liquid; c = critical) n
L
Hydrogen sulfide ( W = 0.100)
T h e detection of errors in a complicated computer program is by no means a trivial matter. The following checks have been developed.
1. T h e goodness of fit is of course checked by direct comparison with experimental data. For this purpose the program accepts P-V-T data of any conceivable kind expressed in any dimension and converts them to compressibility factors, Z. A few of the results are discussed below. 2. The computation of (mean) fugacity coefficients has been checked by the computation of isothermal sets of Z and numerical integration.
$1
Carbon dioxide ( W = 0.225)
T T
f'r
0.760
0.153
Obsd. v 0.845
1 0.0246 v 0.689 1 0.0749 v 0.529 1 0.139 v 0.420
0.886
0.459
0.959
0.766
0.988
0.919
1.ooo 0.931
1.000
0,957
0.748
0.981
0.888
v 0,4575
1.000
1.000
1 0.1641 c 0.2746
$9
Calcd.
Calcd.
0.872 0.0251 0.713 0.0769 0.561 0.153 0.458
0,886 0,884 0,778 0,782 0.707 0,705 0.676 0,673 0.658 0.743 0,740 0.710 0.700 0.680 0.669 0.657
0.608
0.476 0.193 0.288
are presented for substances with high values of W since good agreement with experimental data has never been a question for low W . For mixtures, examples of components with large differences in W were selected because their representation is in general more difficult. A satisfactory representation of the lower hydrocarbons can be taken for granted; but even for n-nonane (Figure 1) the
deviations are appreciable only for a very low reduced temperature or for very high reduced pressures. Even for a strongly polar substance such as ammonia (Figure 2) the discrepancies are large only for low temperatures. T h e remaining deviations for sulfur dioxide (Figure 3) are below 0.03 in 2. Although the deviation curve for the critical isotherm looks strange, the curve for 2 (Figure 4) appears to be reasonable.
0.02 0.0 I
0
D t -0.01
-0.02 -0.03
I I
0
IO
0
20
.
1
3
Pr
30 "
Deviation DZ (observed minus calculated)
Figure 3.
pr Figure 1
1
2
of the compressibility factor of sulfur dioxide (7)
Compressibility factor Z of n-nonane (1 4),
0 to 10,000 p.s.i.
v r, = 0.523 (559.680 R.)
0
r,
= 1.059 (799.680 R.) 0.6
1.0 0.5
0.9
z 0.4
z 0.8 0.3
0.7 0.2
I .o
0.5
0
I
1
0.3
1.0
Pr
Figure 2.
Pr
Compressibility factor Z of ammonia (2)
v rt = 0.796
rr = 0.91 9 A rI =
0
Table 111.
Figure
4.
P,P.S.I. 15.69 1000 2000 2455
679.68
94.91 1000
Compressibility Factor and Fugacities of Methane-+Pentane ( I 5)
779.68
2081 329.16 1025
x1
v 0.0 1 0.0 v 0.946 1 0.303 v 0.919 1 0.573 c 0.821 v 0.0
1 v 1 c v 1 c
Critical isotherm of sulfur dioxide
1.473
z T , OR. 559.68
1 1.1
0.0
0.795 0.247 0.665 0.0 0.0
0.290
Obsd. 0.957
fl
0.487 0.645
Calcd. 0.964 0.0052 0.864 0.265 0.740 0.458 0.609
0.840 0.0301 0.858 0.272 0.654 0.607 0.110 0.454
0.858 0.0296 0.851 0.265 0.608 0.556 0.103 0.363
Obsd. 0.0
861 1548 1773 0.0
808 1620 0.0 522
f2
Calcd. 0.0 0.0 861 864 1575 880 1890 0.0 0.0 804 821 1748 0.0 0.0 660
VOL. 4
NO. 4
Obsd.
15.10 17.9 18.1 20.2 82.02 85.0 96.4 238.68 236.9
Calcd. 15.14 15.59 19.7 13.6 16.6 27.1 13.9 83.0 80.7 92.5 82.5 79.1 232.6 210.8 208
NOVEMBER 1 9 6 5
371
P (PSI)
'
0.5; -
0
0
559.680 R.
919.680 R.
A V A
0
0
'
'
1
10000
(psi)
Temp.
559.68' R.
71 = 0.4055 y 1 = 0.8468
0
679.680R. 0
P
(atrn)
Figure 8. Compressibility nitrogen-ethylene (7)
740' R.
323.25' K.
V
A __-
0
0
y1 = 0.198 y1 = 0.401 yi = 0.796
yi 0.3 yi = 0.7
I
0
Figure 6. Compressibility factor of methanehydrogen sulfide (14 ) 620' R.
'
Figure 7. Compressibility factor of methane-carbon dioxide ( 7 7 )
P (psi)
Temp.
5000
P
Figure 5. Compressibility factor of methanedecane ( 7 4) Temp.
'
'
10000
5000
P (psi1
y1 = 0.1 yi = 0.3 yi = 0.7
"
Calcd.
__ - - - ----
factor
of
V 0
A
Table IV. Compressibility Factor and Fugacities of n-Butane-Carbon Dioxide ( 9 ) ( W l = 0.201; W2 = 0.225; v = vapor; 1 = liquid; c = critical)
Z T , R. 559.68
P , P.S.I.
x1
51.5
v 1.0 1 1.0 v 0.202 1 0.804 v 0.100 1 0.416 c 0.060 v 1.0 1 1.0 v 0.839 1 0.927 c 0.812
300 700 739.68
1095 436 600 758
372
I&EC F U N D A M E N T A L S
Obsd. 0.899 0,014 0.851 0.077 0.675 0.151 0.279 0.527 0.137 0.501 0.195 0,328
fz
f l
Calcd. 0.915 0.015 0.840 0.078 0.643 0.144 0.198 0.551 0.158 0.446 0,201 0.274
Obsd. 47.0
40.8 29.1 12.5 303 294 293
Calcd. 47.5 46.7 40.2 40.1 24.9 23.5 4.6 306 303 294 296 283
Obsd. 0.0
222.0 505 722 0.0
139 258
Calcd. 0.0 0.0 220 139 505 444 725 0.0 0.0
129 108 237
Results for the Critical Compressibility Factor
Table V.
T,, ’ K. 126.3 191 . O 373.6
Substance N2
CHI HqS
C 3-H8
370.0 304.2 430.66 64:. 3
con
SO2 HzO
W 0.040 0.013
Pc, Atm.
33.54 45.79 88.87 42.01 72.80 77.808 218.4
Obsd. 0.292 0.290 0.2833 0,2766 0.2746 0,2697 0.2276
0 100
0.152 0,225 0.2325 0.348
I t is especially the combination of the parameters that is being tested in Figures 4 to 8 and Tables I11 and IV. T h e phase equilibria in these tables represent a particularly severe test, since a t the same time the values of Z and the fugacity coefficients of both components are tested for both phases. These tests will furnish a general impression of the degree of accuracy to be expected in applications. Hydrogen, helium, and water are not satisfactorily represented. Tables I11 and I V show that the critical locus of mixtures is not well described by the program. A systematic study of critical points, critical condensation points, and vapor pressures of mixtures may well be a useful step in the future development of an equation of state.
Ai Az An A; A; Bi B2 B3 B1 B; B6 B7
Bs B9 Bio Bii Biz c 1
Appendix
C2
T h e proposed equation of state represents the cornpressibility factor as the sum of a number of terms
C4
c3
z = zo +
21
+ T.V x Z’ + L(Z3 + bV x Z,)
(11)
T h e quantity Zois a root of the Equation 2 proposed in 1949. T h e convenient representation for computation from the reduced temperature, T r ,and the reduced pressure, P,, is
+ [0.42?481 PrT,-’.5
Z03 - Zoz
(1
- 0.0866404 PrT,-l X
- 0.0866404 PrTr-’)]Z0 - 0.370371 PT2Tr-3.5 =
0
(12)
If this equation yields three real solutions, the computer presents the highest and lowest values. The acentric factor is denoted in Equation 11 by W . T h e quantity L discriminates between liquids and gases :
L
=
1 for’liquids; L
=
0 for gases
(13)
Z1
+ A2(Tr - 1.0)’ +
= -AiP,3/{1.0 A3[P7
- A4 - A5(7‘7 -
l.O)I4]
+
+ BaPrT,2)(1.0 - BJ’r BsTrPr)/ + B7(Tr - BS - B Q P -~ B I o P ~ T +~ ) ~ ] B1lTr3Pr3/(Tr4+ B12Pr4) (14)
BIP,(T~- BZ - B3Pr [1.0
22 = TrP,(Tr - 1.0 - 0.049 P,) X
+
(CI C2Pr - CaTrPr
+ CITr)/
[T: f C S ( T ~ CS - c7Pr Z3 = Di(1 2 4
= El(1
-
-
T,)T?Pr3/[(D2
TT)Pr/(l
- T,l3(1
C S T T P , ) ~ ](15)
+ D3Pr3)]
(16)
+ EZr’) +
+ E4Trz)(1 - Tr)PTZ/(1 + E5Pr4) E6(1 - Tr)Tr6Pr3/(1 + E7Pr6) + - Tr)Tr3(E9- Tr)Pr[I/(l + E10P72)EllPr/’(l + E12Pr4)I/(E13 - Tr’)
Coeficients f o r Equations 14- 17 0.035 C5 0.101212 14,137.6 c 6 2.46596 1.397.124 C, 0.220411 1.030 C* 0.0161963 D; 2.779 X 10-6 13.440 DZ 1.1 0.00260913 3.19325 D3 0.002 1 ,77486 E1 0.1068 0 434418 E4 0 001 0,144392 E; 0 02408 0 00704658 E4 0 . 0 1 124 616.830 E5 1 . 0 x 10-6 1.00122 E6 0.001686 0,0112141 E7 1.0 x 10-7 0.0495574 Es 0.698 0.000442593 E9 0.84 0.0602768 Eio 1. o 0.825714 Eii 0.090 0.00736587 Ei 2 0.1 0.00255204 E1 3 1.1 0.00115729
T h e authors express their appreciation of the advice extended to them in a part of their computational work by members of the Computation Center of the University of California in Berkeley. literature Cited (1) Ackerman, F. J., Redlich, O., J . Chem. Phys. 38, 2740 (1963). (2) Beattie, J. A., Lawrence, C. K., J . A m . Chem. SOC.52, 10 11930). (3)’ Edmister, W. C., Petrol. ReJner 37, No. 4, 173 (1958). (4) Edmister, W. C., Thompson, R. E., Yarborough, L., A.Z. Ch.E.J. 9, 116 (1963). , , (5) Ibid., p. 240. (6) Edwards. D. G.. “Collected and Smoothed Data.” Lawrence Radiation ‘Laboratory, University of California, Rept. UCRL7167, Rev. I (1964). (7) Hagenbach, W. P., Comings, E. W., Znd. Eng. Chem. 45, 606 (1953); (8) Kang, T. L., Hirth, L. J., Kobe, K. A., McKetta, J., J . Chem. Ene. Data 6. 220 (1961). (9) Olds, R. ‘H., Reamer, H. H., Sage, B. H., Lacey, W. N., Ind. Eng. Chem. 41, 475 (1949). (10) Pitzer, K. S., in “Thermodynamics,” by G. N. Lewis, M. Randall, K. S. Pitzer, and Leo Brewer, Appendix 1, McGrawHill, New York, 1061. (11) Reamer, H. H., Olds, R. H., Sage, B. H., Lacey, I V . N., Ind. Eng. Chem. 36, 88 (1944). (12) Redlich, O., Dunlop, A. K., Chem. Eng. Progr. Symp. Ser. 50. 95 (1963). (13) ’Rediich, O., Kwong, J. N. S., Chem. Reus. 44, 233 (1949). (14) Sage, B. H., Lacey, W.N., “Some Properties of the Lighter Hydrocarbons, Hydrogen Sulfide, and Carbon Dioxide,” American Petroleum Institute, New York, 1955. (15) Sage, B. H., Reamer, H. H., Olds, R. H., Lacey, W.N., Znd. Eng. Chem. 34,1108 (1942). (16) Shah, K. K., Thodos, George, Zbid., 57, 30 (1965). ,
RECEIVED for review March 29, 1965 ACCEPTED August 21, 1965
(-E3
E8(l
Present work 0.290 0,290 0.288 0,287 0,285 0,285 0.283
Acknowledgment
\
T h e additional ternis are given by
Calculated Redlich-Kwong Redlich-Dunlop 0.333 0.339 0.333 0.334 0.333 0.315 0.333 0.327 0.333 0.323 0.333 0.314 0.333 0.236
Pitzer 0.288 0.290 0.283 0.279 0.273 0.272 0.263
(17)
Work done under the auspices of the U. S. Atomic Energy Commission. The program has been submitted to SHARE as ES 64, Equation of State 1964. (Until the program is available from SHARE, the senior author wlll be glad to send a copy on request.) VOL. 4
NO. 4
NOVEMBER 1 9 6 5
373
