- Email: [email protected]
SRK and PR alpha functions for ice
Fluid Phase Equilibria 356 (2013) 371–373
Contents lists available at ScienceDirect
Fluid Phase Equilibria journal homepage: www.elsevier.com/locate/fluid
Short communication
SRK and PR alpha functions for ice Claude F. Leibovici a,∗ , Dan Vladimir Nichita b a b
CFL Consultant Hélioparc, 2 Avenue Pierre Angot, 64053 Pau Cedex, France CNRS UMR 5150, Laboratoire des Fluides Complexes, Université de Pau et des Pays de l’Adour, B.P. 1155, 64013 Pau Cedex, France
a r t i c l e
i n f o
Article history: Received 14 June 2013 Received in revised form 29 July 2013 Accepted 31 July 2013 Available online 11 August 2013
a b s t r a c t In order to avoid possible phase inconsistencies in water saturation point calculations, this paper proposes to treat ice as a regular component. For SRK and PR cubic equations of state, the goal of this paper is to provide the parameters of an alpha function able to very accurately reproduce the vapor pressure over ice between 173 K and 273.16 K. © 2013 Elsevier B.V. All rights reserved.
Keywords: Vapor pressure Ice SRK and PR equations of state Alpha function
1. Introduction In natural gas processing, a major information is the water dew point. For many hydrocarbon supercritical systems, it happens that the water dew point temperature is predicted below the water ice point; when this is the case, it then implies that water droplets cannot form but instead, water directly precipitates as frost. Vapor pressure of ice being lower than the extrapolated vapor pressure of liquid water, the errors could be quite significant. In this paper, we try to provide at least a partial solution to this problem considering ice as a pure component for which an appropriate alpha function is proposed. 2. Proposed approach In the present work, between 173 K at 273 K, the vapor pressure over ice has been computed according to the equation proposed by Murphy and Koop [1]: log (P) = 9.550426 −
5723.265 T
+ 3.53068 log (T ) − 0.00728332 T
(1)
with temperature T in Kelvin, and pressure P in Pascal. Using Eq. (1), 101 equally spaced values were tabulated (step size of 1 K); for each temperature, the exact numerical value of the alpha function (which exactly cancels the fugacity difference between vapor and
∗ Corresponding author. Tel.: +33 684 62 89 48; fax: +33 559 84 42 96. E-mail address: cfl[email protected] (C.F. Leibovici). 0378-3812/$ – see front matter © 2013 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.fluid.2013.07.058
solid phases) has been generated using Soave–Redlich–Kwong [2] and Peng–Robinson [3] equations of state. In the first step, we compared these exact values to those estimated using the very simple alpha function proposed in 1993 [4].
˛(Tr ) = Tr 1 + kEoS log
T r
Pr
(2)
obtained by integration of the Clapeyron equation in the immediate vicinity of the critical point (kSRK = 0.327480, kPR = 0.298036); for ice, the selected “critical” constants are the coordinates of the triple point of water. Then, in this work, Tc = 273.16 K and Pc = 611.657 Pa; combined with Eq. (1), this would lead for ice to an apparent acentric factor of 3.18 (nine times larger than for water) which, in turn, probably explains the very large values obtained for the alpha function (up to 3.5 at 173 K). For this first step, rewritten in reduced coordinates, the experimental vapor pressure given by Eq. (1) was plugged into the parameter free Eq. (2). At a first glance, the agreement seemed to be quite good since, looking at the parity plots, we obtained for SRK ˛estimate = 0.987˛exact (R2 = 0.9995) and for PR ˛estimate = 0.974˛exact (R2 = 0.9999) but the sum of squares of relative differences between experimental and computed saturation pressures is 1.84 for SRK and 8.84 for PR. Using kEoS as the single tuning parameter reduces these sums of squares to 0.63 for SRK (kopt = 0.332802) and 0.51 for PR (kopt = 0.309685); even if this is an improvement, at the lowest temperatures the relative errors are as large as 20%. For the purposes of the applications of this work, these errors were considered to be unacceptable and ways for further improvements were investigated. However, in order to clarify the apparent contradiction mentioned above, tempera-
372
C.F. Leibovici, D.V. Nichita / Fluid Phase Equilibria 356 (2013) 371–373
Table 1 Influence of the number of terms on fitting performances.
a1 a1 , a2 a1 , a2 , a3 a1 , a2 , a3 , a4 a1 , a2 , a3 , a4 , a5
Soave–Redlich–Kwong
Peng–Robinson
SSQ
Max error
SSQ
Max error
5.012E-01 2.657E-02 9.205E-03 6.756E-04 8.281E-06
18.25% 3.18% 1.95% 0.66% 0.09%
4.010E-01 4.223E-02 1.283E-02 1.006E-03 1.902E-05
16.79% 3.95% 2.28% 0.79% 0.13%
0.00
ture was eliminated and we looked at the direct dependency of the experimental vapor pressure to the corresponding exact value of the alpha function; the following polynomial correla tions were obtained (R2 = 0.9999): for SRK Log Prexp = 10.513 −
Deviation on pressure, Pa
Used parameters
0.05
1 + a1 ×
1
Tr
− 1 + a2 log (Tr ) + a3 (Tr − 1)
+ a4 (Tr − 1)2 + a5 (Tr − 1)3
-0.15
SRK PR
-0.20
(3)
Coefficients a1 , to a5 were adjusted in order to minimize the sum of squares of relative differences between experimental and computed saturation pressures (their final values are reported in Table 2); the corresponding sum of squares of relative errors are 8.281E-06 for SRK and 1.902E-05 for PR. Just for comparison purposes, adjusting in the same way the single parameter appearing in Soave’s alpha function leads to sums of squares equal to 20.39 for SRK and 18.37 for PR, these very large numbers being mainly due to extremely large errors on the predicted vapor pressures below 220 K. One could notice that, without any tuning, the parameter free Eq. (2) led to significantly lower sums of squares. Because of the targeted and obtained accuracy levels, when plotted on the same graph and on any possible scale, it is impossible to notice any difference between calculated and experimental vapor pressure curves. So, as a function of temperature, the deviations on pressure (in Pa) are reported in Fig. 1 and the corresponding relative errors in Fig. 2. The proposed alpha function works very well in the temperature range over which its parameters have been adjusted; however, Table 2 Coefficients in Eq. (3) for SRK and PR EoS. Coefficient
SRK
PR
a1 a2 a3 a4 a5
−7.910222E+01 −4.623497E+02 3.761995E+02 −1.395614E+02 9.278779E+01
−9.214883E+01 −5.284656E+02 4.298903E+02 −1.591867E+02 1.049790E+02
-0.25
-0.30 170
180
190
200
210
220
230
240
250
260
270
280
Temperature, K Fig. 1. Deviation on vapor pressure vs. temperature, SRK and PR EoS.
0.15
SRK 0.10
Relative error,%
11.042 − 5.601˛exact + 1.479˛2exact − 0.501˛3exact This explains the extremely large sensitivity of predicted pressure to any even minor change of the alpha function value; taking the average value over the covered range (1 < ˛ < 3), we have more r or less P ≈ −6˛; so, to a change of 0.01 applied to the alpha Pr function corresponds a 6% change of the calculated pressure. None of the most classical alpha function formulations was able to give satisfactory results at the level of required accuracy. So, in the second step, the reduced pressure term which appears in Eq. (2) was replaced by a classical expansion derived from the assumption that heat of vaporization can be represented as a polynomial of temperature. After trying different numbers of terms in order to reach, as for any other component, a maximum relative error lower than 1% in the overall range of temperature (Table 1 summarizes the results as a function of the number of used terms), the final form of the proposed alpha function is:
-0.10
4.796˛exact + 1.089˛2exact − 0.383˛3exact and for PR Log Prexp =
˛(Tr ) = Tr
-0.05
PR
0.05
0.00
-0.05
-0.10 170
180
190
200
210
220
230
240
250
260
270
280
Temperature, K Fig. 2. Relative errors for vapor pressure vs. temperature, SRK and PR EoS.
outside this range, it poorly extrapolates; it goes through a minimum around 310 K and through a maximum around 120 K. Then, if during any iterative calculation temperature changes, it is strongly recommended to use, for any temperature above the freezing point of water, one of the following extrapolating functions: ˛(Tr ) = exp [b1 (1 − Tr b2 )]
(4a)
Table 3 Coefficients in Eqs. (4a) and (4b) for SRK and PR EoS. Coefficient
SRK
PR
b1 b2 c1 c2
3.087835E+01 6.041621E+00 1.166279E+00 5.180256E+00
2.261373E+01 5.416345E+00 1.376861E+00 3.933836E+00
C.F. Leibovici, D.V. Nichita / Fluid Phase Equilibria 356 (2013) 371–373
or ˛(Tr ) =
1 1 + c1 (Tr − 1) + c2 (Tr − 1)2
(4b)
The coefficients in Eqs. (4a) and (4b), obtained by matching the first and second derivatives of Eq. (3) at the critical point, are reported in Table 3. Finally, we note that what is proposed here does not change or improve in any manner the performances of cubic equations of state in water dew point calculations for natural gases;
373
it just extends the capabilities of such model to frost point calculations. References [1] D.M. Murphy, T. Koop, Quarterly Journal of the Royal Meteorological Society 131 (2005) 1539–1565. [2] G. Soave, Chemical Engineering Science 27 (1972) 1197–1203. [3] D.Y. Peng, D.B. Robinson, Industrial & Engineering Chemistry Fundamentals 15 (1976) 59–64. [4] C.F. Leibovici, Fluid Phase Equilibria 84 (1993) 1–8.
Contents lists available at ScienceDirect
Fluid Phase Equilibria journal homepage: www.elsevier.com/locate/fluid
Short communication
SRK and PR alpha functions for ice Claude F. Leibovici a,∗ , Dan Vladimir Nichita b a b
CFL Consultant Hélioparc, 2 Avenue Pierre Angot, 64053 Pau Cedex, France CNRS UMR 5150, Laboratoire des Fluides Complexes, Université de Pau et des Pays de l’Adour, B.P. 1155, 64013 Pau Cedex, France
a r t i c l e
i n f o
Article history: Received 14 June 2013 Received in revised form 29 July 2013 Accepted 31 July 2013 Available online 11 August 2013
a b s t r a c t In order to avoid possible phase inconsistencies in water saturation point calculations, this paper proposes to treat ice as a regular component. For SRK and PR cubic equations of state, the goal of this paper is to provide the parameters of an alpha function able to very accurately reproduce the vapor pressure over ice between 173 K and 273.16 K. © 2013 Elsevier B.V. All rights reserved.
Keywords: Vapor pressure Ice SRK and PR equations of state Alpha function
1. Introduction In natural gas processing, a major information is the water dew point. For many hydrocarbon supercritical systems, it happens that the water dew point temperature is predicted below the water ice point; when this is the case, it then implies that water droplets cannot form but instead, water directly precipitates as frost. Vapor pressure of ice being lower than the extrapolated vapor pressure of liquid water, the errors could be quite significant. In this paper, we try to provide at least a partial solution to this problem considering ice as a pure component for which an appropriate alpha function is proposed. 2. Proposed approach In the present work, between 173 K at 273 K, the vapor pressure over ice has been computed according to the equation proposed by Murphy and Koop [1]: log (P) = 9.550426 −
5723.265 T
+ 3.53068 log (T ) − 0.00728332 T
(1)
with temperature T in Kelvin, and pressure P in Pascal. Using Eq. (1), 101 equally spaced values were tabulated (step size of 1 K); for each temperature, the exact numerical value of the alpha function (which exactly cancels the fugacity difference between vapor and
∗ Corresponding author. Tel.: +33 684 62 89 48; fax: +33 559 84 42 96. E-mail address: cfl[email protected] (C.F. Leibovici). 0378-3812/$ – see front matter © 2013 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.fluid.2013.07.058
solid phases) has been generated using Soave–Redlich–Kwong [2] and Peng–Robinson [3] equations of state. In the first step, we compared these exact values to those estimated using the very simple alpha function proposed in 1993 [4].
˛(Tr ) = Tr 1 + kEoS log
T r
Pr
(2)
obtained by integration of the Clapeyron equation in the immediate vicinity of the critical point (kSRK = 0.327480, kPR = 0.298036); for ice, the selected “critical” constants are the coordinates of the triple point of water. Then, in this work, Tc = 273.16 K and Pc = 611.657 Pa; combined with Eq. (1), this would lead for ice to an apparent acentric factor of 3.18 (nine times larger than for water) which, in turn, probably explains the very large values obtained for the alpha function (up to 3.5 at 173 K). For this first step, rewritten in reduced coordinates, the experimental vapor pressure given by Eq. (1) was plugged into the parameter free Eq. (2). At a first glance, the agreement seemed to be quite good since, looking at the parity plots, we obtained for SRK ˛estimate = 0.987˛exact (R2 = 0.9995) and for PR ˛estimate = 0.974˛exact (R2 = 0.9999) but the sum of squares of relative differences between experimental and computed saturation pressures is 1.84 for SRK and 8.84 for PR. Using kEoS as the single tuning parameter reduces these sums of squares to 0.63 for SRK (kopt = 0.332802) and 0.51 for PR (kopt = 0.309685); even if this is an improvement, at the lowest temperatures the relative errors are as large as 20%. For the purposes of the applications of this work, these errors were considered to be unacceptable and ways for further improvements were investigated. However, in order to clarify the apparent contradiction mentioned above, tempera-
372
C.F. Leibovici, D.V. Nichita / Fluid Phase Equilibria 356 (2013) 371–373
Table 1 Influence of the number of terms on fitting performances.
a1 a1 , a2 a1 , a2 , a3 a1 , a2 , a3 , a4 a1 , a2 , a3 , a4 , a5
Soave–Redlich–Kwong
Peng–Robinson
SSQ
Max error
SSQ
Max error
5.012E-01 2.657E-02 9.205E-03 6.756E-04 8.281E-06
18.25% 3.18% 1.95% 0.66% 0.09%
4.010E-01 4.223E-02 1.283E-02 1.006E-03 1.902E-05
16.79% 3.95% 2.28% 0.79% 0.13%
0.00
ture was eliminated and we looked at the direct dependency of the experimental vapor pressure to the corresponding exact value of the alpha function; the following polynomial correla tions were obtained (R2 = 0.9999): for SRK Log Prexp = 10.513 −
Deviation on pressure, Pa
Used parameters
0.05
1 + a1 ×
1
Tr
− 1 + a2 log (Tr ) + a3 (Tr − 1)
+ a4 (Tr − 1)2 + a5 (Tr − 1)3
-0.15
SRK PR
-0.20
(3)
Coefficients a1 , to a5 were adjusted in order to minimize the sum of squares of relative differences between experimental and computed saturation pressures (their final values are reported in Table 2); the corresponding sum of squares of relative errors are 8.281E-06 for SRK and 1.902E-05 for PR. Just for comparison purposes, adjusting in the same way the single parameter appearing in Soave’s alpha function leads to sums of squares equal to 20.39 for SRK and 18.37 for PR, these very large numbers being mainly due to extremely large errors on the predicted vapor pressures below 220 K. One could notice that, without any tuning, the parameter free Eq. (2) led to significantly lower sums of squares. Because of the targeted and obtained accuracy levels, when plotted on the same graph and on any possible scale, it is impossible to notice any difference between calculated and experimental vapor pressure curves. So, as a function of temperature, the deviations on pressure (in Pa) are reported in Fig. 1 and the corresponding relative errors in Fig. 2. The proposed alpha function works very well in the temperature range over which its parameters have been adjusted; however, Table 2 Coefficients in Eq. (3) for SRK and PR EoS. Coefficient
SRK
PR
a1 a2 a3 a4 a5
−7.910222E+01 −4.623497E+02 3.761995E+02 −1.395614E+02 9.278779E+01
−9.214883E+01 −5.284656E+02 4.298903E+02 −1.591867E+02 1.049790E+02
-0.25
-0.30 170
180
190
200
210
220
230
240
250
260
270
280
Temperature, K Fig. 1. Deviation on vapor pressure vs. temperature, SRK and PR EoS.
0.15
SRK 0.10
Relative error,%
11.042 − 5.601˛exact + 1.479˛2exact − 0.501˛3exact This explains the extremely large sensitivity of predicted pressure to any even minor change of the alpha function value; taking the average value over the covered range (1 < ˛ < 3), we have more r or less P ≈ −6˛; so, to a change of 0.01 applied to the alpha Pr function corresponds a 6% change of the calculated pressure. None of the most classical alpha function formulations was able to give satisfactory results at the level of required accuracy. So, in the second step, the reduced pressure term which appears in Eq. (2) was replaced by a classical expansion derived from the assumption that heat of vaporization can be represented as a polynomial of temperature. After trying different numbers of terms in order to reach, as for any other component, a maximum relative error lower than 1% in the overall range of temperature (Table 1 summarizes the results as a function of the number of used terms), the final form of the proposed alpha function is:
-0.10
4.796˛exact + 1.089˛2exact − 0.383˛3exact and for PR Log Prexp =
˛(Tr ) = Tr
-0.05
PR
0.05
0.00
-0.05
-0.10 170
180
190
200
210
220
230
240
250
260
270
280
Temperature, K Fig. 2. Relative errors for vapor pressure vs. temperature, SRK and PR EoS.
outside this range, it poorly extrapolates; it goes through a minimum around 310 K and through a maximum around 120 K. Then, if during any iterative calculation temperature changes, it is strongly recommended to use, for any temperature above the freezing point of water, one of the following extrapolating functions: ˛(Tr ) = exp [b1 (1 − Tr b2 )]
(4a)
Table 3 Coefficients in Eqs. (4a) and (4b) for SRK and PR EoS. Coefficient
SRK
PR
b1 b2 c1 c2
3.087835E+01 6.041621E+00 1.166279E+00 5.180256E+00
2.261373E+01 5.416345E+00 1.376861E+00 3.933836E+00
C.F. Leibovici, D.V. Nichita / Fluid Phase Equilibria 356 (2013) 371–373
or ˛(Tr ) =
1 1 + c1 (Tr − 1) + c2 (Tr − 1)2
(4b)
The coefficients in Eqs. (4a) and (4b), obtained by matching the first and second derivatives of Eq. (3) at the critical point, are reported in Table 3. Finally, we note that what is proposed here does not change or improve in any manner the performances of cubic equations of state in water dew point calculations for natural gases;
373
it just extends the capabilities of such model to frost point calculations. References [1] D.M. Murphy, T. Koop, Quarterly Journal of the Royal Meteorological Society 131 (2005) 1539–1565. [2] G. Soave, Chemical Engineering Science 27 (1972) 1197–1203. [3] D.Y. Peng, D.B. Robinson, Industrial & Engineering Chemistry Fundamentals 15 (1976) 59–64. [4] C.F. Leibovici, Fluid Phase Equilibria 84 (1993) 1–8.