Use of transformed correlations to help screen and populate properties within databanks

Fluid Phase Equilibria 256 (2007) 54–61

Use of transformed correlations to help screen and populate properties within databanks David L. Morgan ∗ Dow Corning Corporation, P.O. Box 994, M/S C043D1, Midland, MI 48686, USA Received 25 September 2006; received in revised form 22 January 2007; accepted 23 January 2007 Available online 30 January 2007

Abstract In-house thermophysical property databanks are often sparsely populated with properties and contain legacy correlations which impede using database information in other simulation or database tools. As part of their flatfile origin, they often contain errors which eluded earlier detection. While many organizations, such as DIPPRs 801 database consortium, employ quality checks to help eliminate defects within their databanks, what is proposed here is a more systematic use of “transformed” correlations that automatically “fit” properties to desired equation forms used within a database. This approach is most applicable to properties that can be estimated by relatively simple group contribution (GC) or corresponding states principle (CSP) methods. Dow Corning (DC) is applying this approach to help find errors and to convert non-standard correlation forms to more widely accepted DIPPR equation forms. Three transforms of tabular correlations are presented in this paper: (1) Halm and Stiel’s [R.L. Halm, L.I. Stiel, A fourth parameter for the vapor pressure and entropy of vaporization of polar fluids, AIChE J. 13 (2) (1967) 351–355] four-parameter (Tc , Pc , w, x) CSP vapor pressure (VP) correlation to DIPPR equation 101 (Riedel); (2) Pitzer’s [K.S. Pitzer, D.Z. Lippman, R.F. Curl Jr., C.M. Huggins, D.E. Petersen, The volumetric and thermodynamic properties of fluids. I. Theoretical basis and virial coefficients. II. Compressibility factor, vapor pressure and entropy of vaporization, J. Am. Chem. Soc. 77 (1955) 3427–3440] saturated entropy of vaporization tables to DIPPR equation 106 for heat of vaporization (Hv ); (3) Harrison and Seaton’s [B.K. Harrison, W.H. Seaton, Solution to missing group problem for estimation of ideal gas heat capacities, Ind. Eng. Chem. Res. 27 (8) (1988) 1536–1540] ideal gas heat capacity (ICP) GC method with Si, B, and Al group updates from Myers and Danner [K.H. Myers, R.P. Danner, Prediction of properties of silicon, boron, and aluminum compounds, J. Chem. Eng. Data 38 (2) (1993) 175–200] to DIPPR equation 107 (Aly and Lee). Statistical measures of the VP re-correlation in the ranges, 0.44 ≤ Tr ≤ 1, 0 ≤ w ≤ 1.2, |x| ≤ 0.25 are S.D. = 4.2%, AAD = 2.6%, and maximum deviation = 27.5%. The VP transform represents Pitzer’s original VP correlation from which the Halm and Stiel method is derived with RMS = 0.66% and AAD = 0.37%. The Hv transform was fitted in the 0.56 ≤ Tr ≤ 1 and 0 ≤ w ≤ 1.2 ranges with S.D. = 0.26%. The ICP transform reproduces the original Harrison and Seaton method with a S.D. = 0.7%, and AAD = 0.3%. It has been used to identify problems within DC databanks, with a literature ICP evaluation, and with siloxane ICP data in the DIPPR database which were estimated by ab initio methods. © 2007 Elsevier B.V. All rights reserved. Keywords: Databases; Correlations; Vapor pressure; Heat of vaporization; Ideal gas heat capacity; Aly–Lee equation; Corresponding states principle; Groupcontribution method; Estimation methods; DIPPR; PPDS; ASPEN Plus; Siloxanes; Silanes

1. Introduction Dow Corning (DC) maintains an in-house physical properties database for pure components. Currently the DC databank has 389 materials which are predominately silicon-containing. To expand the range of pure substances, an additional 1864 components are obtained from a Sponsor version of the DIPPR database

∗

Tel.: +1 989 496 8399; fax: +1 989 496 5956. E-mail address: [email protected].

0378-3812/$ – see front matter © 2007 Elsevier B.V. All rights reserved. doi:10.1016/j.fluid.2007.01.016

[1]. Since 77 compounds overlap these two databanks, the total number of compounds is 2176. Originally the DC database was stored in a flatfile and accessed via an in-house FORTRAN program called PHYPROP. In recent years the database has been moved from flatfile into a relational, Excel database to improve its maintainability. The in-house program has been replaced by TUV NELs PPDS program. While NEL has been very accommodating of Dow Corning’s legacy equations, these forms continue to impose their limitations on users and to complicate database maintenance; e.g., a parallel databank, INHSPCD, which is based on an older version of DCs databank, is used in

D.L. Morgan / Fluid Phase Equilibria 256 (2007) 54–61

55

ASPEN Plus. This database uses different equation forms and contains only a small subset of the properties which are available in PPDS. We are currently using an NEL utility called PINST3 that can export PPDS databanks into an intermediate DFMS file which ASPEN can import.1 ASPEN converts this file into a second pure-component user databank called “PPDS”. Users of ASPEN can then obtain new or updated component properties by first searching PPDS and then the legacy INHSPCD. Unfortunately, NELs utility will only export properties to ASPEN which are represented by equation forms which ASPEN understands. DCs choice to replace its legacy equations with DIPPR equation forms follows from our sponsorship of the DIPPR 801 Project. The transforms discussed in this paper could also be adapted to other equation forms such as NELs. Use of transformed correlations is very common. Perhaps one of the better known examples is Lee and Kesler’s [2] reformation of Pitzer’s et al. [3] tabular CSP method into one which uses simple and n-octane reference fluids which are fitted to an analytical BWR equation of state. Their auxiliary vapor pressure correlation: ln(Pr ) = f0 (Tr ) + w f1 (Tr )

(1)

where Tr = T/Tc and Pr = P/Pc are the reduced temperature and pressure, w the Pitzer’s acentric factor, and f0 and f1 are functions of Tr : f0 (Tr ) =

5.92714−6.09648/Tr −1.28862 ln(Tr ) + 0.169347Tr6 ,

f1 (Tr ) =

15.2518−15.6875/Tr −13.4721 ln(Tr ) + 0.43577Tr6

can be made into a transform which takes three input parameters (Tc , Pc , w) and generates a Riedel vapor pressure equation fit (DIPPRs equation 101 with a T6 term): ln(P) = c1 +

c2 + c3 ln(T ) + c4 T 6 T

(2)

where c1 = ln(Pc ) + (5.92714 + 15.2518w) + (1.28862 + 13.4721w) ln(Tc ), c2 = (−6.09648 − 15.6875w)Tc , c3 = (−1.28862 − 13.4721w), c4 = (0.169347 + 0.43577w)/Tc6 .

2. Transformation of Halm and Stiel vapor pressure method (H&S-1) Most vapor pressure data in Dow Corning’s databank are represented by using four parameters (Tc , Pc , w, x). These were fitted from vapor pressure data to minimize deviations from the

Fig. 1. Halm and Stiel’s 1967 extension [4] of the Pitzer et al. 1955 CSP method [3] for vapor pressure.

Halm and Stiel (1967) CSP method [4]: log(Pr ) = log Pr(0) + w log Pr(1) + x log Pr(2)

(3)

(0)

where the term log Pr represents the reduced vapor pressure (1) of a simple spherical fluid having w = 0, log Pr a perturbation term which represents a contribution caused by a molecule’s (2) non-spherical shape or “acentricity”, and the log Pr term represents contributions caused by a molecule’s polarity. The parameter x is called the Stiel polarity factor. For nonpolar molecules, x = 0, and H&S-1 reduces to the original 1955 Pitzer et al. Three-parameter (Tc , Pc , w) method from which it is derived. Each term is tabulated as a function of reduced temperature; these are reproduced in Table 1 and in Fig. 1. Halm and Stiel extended the lower temperature range of the original Pitzer correlation from Tr = 0.56 to 0.44. As can be seen in Fig. 1, (2) the log Pr term goes to zero for Tr ≥ 0.7. This has the mathematical consequence of leading to an undefined vapor pressure derivative at Tr = 0.7. The abrupt change in vapor pressure slope at Tr = 0.7, especially when the absolute magnitude of x is large, can severely limit the accuracy of re-fits of the H&S-1 method to analytical equations such as DIPPRs equation 101. To re-correlate H&S-1, Eq. (2) was first written in terms of reduced variables: ln(Pr ) = c0 +

c1 + c2 ln(Tr ) + c3 Tr6 Tr

(4)

and each coefficient ck was assumed to be of the form: ck = bk,0 + bk,1 w + bk,2 x

(5)

1

NEL also has a Cape-Open exchange DLL which allows PPDS information to be directly exchanged with ASPEN and avoid the use of intermediate DFMS files. DC does not currently use it.

The bk,j were fitted using non-linear regression analysis of 5408 H&S-1 grid points over the variable ranges which occur

56

D.L. Morgan / Fluid Phase Equilibria 256 (2007) 54–61

Table 1 Vapor pressure and entropy of vaporization correlations of Pitzer et al. [3] and of Halm and Stiel [4] (0)

(1)

(2)

Tr

Pitzer log Pr

Pitzer log Pr

H&S log Pr

Pitzera,b S(0)

Pitzera,b S(1)

Pitzera,b S(2)

H&S-1a,c S(2)*

0.44 0.46 0.48 0.50 0.52 0.54 0.56 0.58 0.60 0.62 0.64 0.66 0.68 0.70 0.72 0.74 0.76 0.78 0.80 0.82 0.84 0.86 0.88 0.90 0.92 0.94 0.95 0.96 0.97 0.98 0.99 1.00

−3.011 −2.745 −2.516 −2.321 −2.141 −1.974 −1.834 −1.688 −1.552 −1.426 −1.308 −1.198 −1.096 −1.000 −0.909 −0.823 −0.742 −0.665 −0.592 −0.522 −0.455 −0.391 −0.330 −0.270 −0.212 −0.156 −0.129 −0.102 −0.076 −0.050 −0.025 0.000

−3.810 −3.540 −3.250 −2.900 −2.610 −2.340 −2.080 −1.880 −1.700 −1.540 −1.390 −1.250 −1.120 −1.000 −0.895 −0.800 −0.705 −0.620 −0.545 −0.475 −0.405 −0.345 −0.285 −0.230 −0.180 −0.133 −0.109 −0.086 −0.064 −0.042 −0.021 0.000

4.600 3.780 3.440 2.870 2.410 2.110 1.670 1.310 1.000 0.745 0.520 0.324 0.156 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000

18.64 17.74 16.92 16.12 15.36 14.62 13.89 13.19 12.49 11.84 11.20 10.57 9.97 9.37 8.79 8.19 7.58 6.95 6.23 5.44 5.00 4.52 4.00 3.38 2.57 0.00

27.80 26.20 24.60 23.20 21.80 20.50 19.30 18.10 17.00 16.00 14.90 13.90 13.00 12.10 11.20 10.30 9.39 8.53 7.54 6.51 5.96 5.39 4.72 3.91 2.83 0.00

2.10 2.20 2.30 2.40 2.50 2.60 2.70 2.80 2.80 2.90 2.90 2.80 2.70 2.60 2.50 2.40 2.20 2.00 1.80 1.50 1.40 1.30 1.10 0.90 0.60 0.00

−44.70 −42.90 −41.10 −39.10 −36.90 −34.50 −31.80 −28.60 −24.80 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00

a b c

Units are calorie/mol K, where 1 calorie = 4.186 J. Pitzer et al.: S = S (0) + wS (1) + w2S (2) . Halm and Stiel: S = S (0) + wS (1) + xS (2)∗ .

within DCs database: 0.44 ≤ Tr ≤ 1, 0 ≤ w ≤ 1.2, |x| ≤ 0.25. The fitting was done using the Marquardt–Levenberg algorithm in Systat Software’s SigmaPlot program. After re-arranging the bk,j and Eq. (4), the resulting coefficients for Eq. (2) are given by: c1 = ln(Pc ) + (6.1526 + 13.296w − 37.857x) + (1.787 + 10.278w − 49.396x) ln(Tc ), c2 = Tc (−6.3825 − 13.51w + 39.265x),

comparison statistics are RMS = 3.0% and AAD = 1.4%. When transform fits are restricted to the approximate range of the Pitzer’s original VP correlation, Tr ≥ 0.56, w ≤ 0.8, and x = 0, RMS = 0.66% and AAD = 0.37%. Fig. 2 shows the effect of (2) discontinuity in the derivative of log Pr at Tr = 0.7 for three hypothetical compounds. As the absolute value of x becomes larger, the deviation at Tr = 0.7 sharply increases. The solid circles represent the worst case compound which had the largest deviation in the transform’s regression. The overall worst case comparison is AAD = 5.8%.

c3 = −1.787 − 10.278w + 49.396x, c4 =

0.23735 + 0.20225w − 1.2942x Tc6

(6)

Fit statistics expressed in terms of relative pressure deviations (%) are summarized in Table 2. Overall the fit had a S.D. of 4.2% and an AAD of 2.7%. The maximum deviation of 27.5% occurred at (Tr , w, x) = (0.44, 1.2, −0.25). It is apparent that the largest deviations between original correlation and transform fit occur at low temperatures for compounds which have large w and x because vapor pressure has a very low magnitude in this region. In the range 0.56 ≤ Tr ≤ 1, w ≤ 0.8, and |x| ≤ 0.05,

3. Transformation of Pitzer’s heat of vaporization method Although Halm and Stiel [4] presented a CSP method to calculate heat of vaporization for polar fluids which was based on Pitzer’s work [3], this method was designed primarily as an alternative way to evaluate parameter x. Instead, Dow Corning adapted Pitzer’s entropy of vaporization correlation to evaluate heats of vaporization: S = S (0) + wS (1) + w2 S (2) ,

(7)

D.L. Morgan / Fluid Phase Equilibria 256 (2007) 54–61

57

Table 2 Statistics for fitting the transform of Halm and Stiel vapor pressure method [4] to Riedel’s Eq. (2)

Tr -range w-range x-range N AAD (%) Maximum deviation (%) RMS (%) S.D. (%) Maximum deviation Tr Maximum deviation w Maximum deviation x

All points used in regression

H&S CSP restricted x

Pitzer CSP

Pitzer CSP w-extrapa

Pitzer CSP Tr -extended w-extrapa

0.44–1 0–1.2 ≤|0.25| 5408 2.6 27.5 4.2 4.2 0.44 1.2 −0.25

0.44–1 0–1.2 ≤|0.05| 2080 1.4 21.8 3.0 N/A 0.44 1.2 −0.05

0.56–1 0–0.8 0 234 0.37 3.8 0.66 N/A 0.56 0.8 0

0.56–1 0–1.2 0 338 0.58 6.9 1.4 N/A 0.56 1.2 0

0.44–1 0–1.2 0 416 1.2 20.4 2.9 N/A 0.44 1.2 0

Pitzer’s 1955 method [3] used n-heptane as the compound with the largest acentric factor, w = 0.352. n-Octadecane (C18) and n-octacosane (C28) have acentric factors of approximately 0.8 and 1.2, respectively. a

Perturbation terms for Eq. (7), which were originally tabulated as a function of Tr in the range 0.56 ≤ Tr ≤ 1.0, are reproduced in Table 1. The Halm and Stiel version of S(2) is also included in this table. To transform Eq. (7) into DIPPRs Hv equation 107:

The bi,j were fitted using non-linear regression analysis of 338Hv∗ grid points in the ranges 0.56 ≤ Tr ≤ 1.0 and 0 ≤ w ≤ 1.20. The fit’s overall S.D. and AAD were 0.26% and 0.25%, respectively. The maximum deviation, 0.89%, occurred at (Tr , w) = (0.9, 0). The resulting coefficients for Eq. (8) are:

Hv = d1 (1 − Tr )(d2 +d3 Tr +d4 Tr +d5 Tr )

d1 = (RTc )(7.8149 + 11.409w + 2.1674w2 − 0.65342w3 ),

2

3

(8)

it was first re-arranged into dimensionless heat vaporization: Hv∗ =

Hv TS Tr (S (0) + wS (1) + w2 S (2) ) = = RTc RTc R

(9)

d3 = −0.84408 + 1.8297w − 3.2435w2 + 1.1449w3 , d4 = 0.41923 − 1.0892w + 1.9138w2 − 0.65758w3 , d5 = 0.0

where R is the gas constant, R = 8314.41 J/(K kmol) = 1987.19 cal/(K kmol). Each coefficient of Eq. (8) was then made a cubic function of w: di = bi,0 + bi,1 w + bi,2 w2 + bi,3 w3

d2 = 0.81892 − 0.67637w + 1.2798w2 − 0.47594w3 ,

(10)

(11)

When the Pitzer Hv transform was compared with Hv departure function results for the simple and n-octane reference fluids of Lee and Kesler [2], the AAD and maximum deviation were found to 1.3% and 4.0%, respectively. The Lee and Kesler results are given in Morgan [5], who accurately fitted Lee and Kesler BWR-EOS results in the 0.36 ≤ Tr ≤ 1.0 range to the scaled equation of Torquato and Stell [6]. In doing the comparisons, the Pitzer transform and the fits of Lee and Kesler data were extrapolated down to Tr = 0.3. 4. Transformation of Harrison and Seaton ideal gas heat capacity method (H&S–2) Harrison and Seaton’s 1988 group-contribution method [7] predicts ICP over the temperature range from 300 to 1500 K. It was developed for use in the ASTMs chemical thermodynamic and energy release computer program, CHETAH [8]. This method requires as inputs only types of atoms in a molecule and their numbers; i.e., an empirical formula. Although much simpler to use than 1st order GC methods such as Benson’s method [9], it often gives predictions which are only slightly inferior. When more elaborate methods cannot be used because of missing groups (such as for most silicon chemicals), it may be the best method to use. The form of this correlation is:

Fig. 2. Deviation plot comparing re-correlated Halm and Stiel vapor pressure fits with original correlation for three sets of input parameters (w and x).

ICP = φConstant +

14  k=1

n k φk

(12a)

58

D.L. Morgan / Fluid Phase Equilibria 256 (2007) 54–61

Table 3 Harrison and Seaton [7] ideal gas heat capacity group parametersa , φk (J/(kmol K)), for use in Eqs. (12a) and (12b) Temperature (K)

Constant

C

H

O

Si

N

F

Cl

Br

I

S

Al

B

P

Other

300 400 500 600 800 1000 1500

4860 864 −1850 −4610 −7490 −8530 −7370

9,040 12,600 15,500 17,500 20,100 21,600 23,900

5,690 7,370 8,890 10,500 13,100 15,200 17,900

11,400 13,900 15,700 17,500 19,400 20,400 20,600

17,500 19,800 22,300 22,100 24,500 24,600 24,000

11,900 14,000 16,000 17,300 19,400 20,400 21,100

12,700 16,200 17,900 20,100 21,500 22,400 22,100

16,800 18,900 20,200 21,400 22,400 22,800 22,600

17,800 19,900 21,200 22,400 23,400 23,800 23,000

18,700 20,500 22,100 23,300 25,000 25,400 24,600

15,300 17,000 19,400 20,300 22,300 22,900 22,500

16,600 19,100 21,000 22,500 24,000 24,500 25,000

10,500 13,700 16,100 17,800 20,500 21,600 22,700

18,000 20,900 21,600 22,800 23,000 23,400 24,200

19,500 20,800 21,700 22,100 23,000 23,300 23,300

a

The parameters for Si, B, and Al groups are updates from Myers and Danner [10]. Units in the original publications, J/(mol K), were changed to J/(kmol K).

ICP = φConstant + nC φC + nH φH + nO φO + nSi φSi + nN φN + nF φF + nCl φCl + nBr φBr + nI φI + nS φS + nAl φAl + nB φB + nP φP + nOther φOther

(12b)

where nk is the number of the kth atom type and the function φ is tabulated over the 300–1500 K range for each atom group. The 14 atom groups are explicitly identified in Eq. (12b). “Other” is used for elements other than C, H, O, Si, N, F, Cl, Br, I, S, Al, B, or P. The constant term, φConstant , is only a function of temperature. Table 3 gives the values of each φ over the temperature range. The groups for Si, Al, and B in this table are updates of the H&S-2 method given in Myers and Danner [10]. Units were changed from J/(mol K) in the original papers to J/(kmol K) to be consistent with DIPPRs units. Most components in DCs databank have an ICP fit. While some of the fits are based on evaluations of published spectroscopic or speed-of-sound data, most were obtained by fitting H&S-2 data to a polynomial equation. Polynomials accurately represent H&S-2 data, but when extrapolated to low temperatures in a simulator such as ASPEN Plus, they can give bad results. To keep PHYPROPs polynomial ICP fits from exhibiting weird low-temperature behavior in ASPEN Plus, an auxiliary extrapolation function, C9 + C10 T1.5 for T < Tmin , is used, where coefficients C9 and C10 are evaluated from the value and derivative of ICP at Tmin . This ASPEN add-on function, along with the tabular Halm and Stiel and Pitzer correlations used in PHYPROP, are some of the reasons why DC has parallel databanks. A much better solution would be to replace the legacy equations with forms that can be used in both platforms. For ICP, DIPPRs equation 107 is a good choice. This equation was developed from statistical mechanical considerations by Aly and Lee [11] and has much better extrapolation behavior2 outside of the range of H&S-2:  ICP = a0 + a1

a2 /T sinh(a2 T )



2 + a3

a4 /T cosh(a4 /T )

2 (13a)

To convert the H&S-2 atomic groups and constants in Table 3 into an Aly and Lee equation transform, each coefficient ai of Eq. (13a) was represented as a summation that followed the form

2 DIPPR has noted that the Aly–Lee equation sometimes levels off to an unrealistic constant value when extrapolated below 200 K.

of Eq. (12a): ai = ai,0 +

14 

(13b)

nj ai,j

j=1

where ai,0 are fitted constants, and the ai,j are fitted coefficients which are multiplied by nj , the number of the jth atom in a molecule. The net effect of doing a regression analysis was to reduce the number of functions and constants in Table 3 from 105 to 75 coefficients in Table 4. As the footnote in Table 4 explains, another 28 of the small a2,j and a4,j coefficients in this table can be set to zero with little affect on the fit. During the initial regression, coefficients for the C, H, O, and Si atoms were determined by fitting a table of H&S-2 data for 163 compounds. The other groups were subsequently fitted by using the H&S-2 data of other compounds. The overall fit, which used 3115 H&S2 data points generated for 445 compounds, had a S.D. = 0.66%, AAD = 0.31%, and a maximum deviation of 13% (P at 1000 K). These statistics are summarized in Table 5. In Table 6, the 10 deviations larger than 4% are listed; these occurred for the monatomic species of Other, P, B, and Al. These relatively large Table 4 Coefficientsa , ai,j , for the transform of the Harrison and Seaton ICP method [7] to the Aly and Lee Eqs. (13a) and (13b) Group ID

j

a0,j

a1,j

a2,j

a3,j

a4,j

Constant C H O Si N F Cl Br I S Al B P Other

0 1 2 3 4 5 6 7 8 9 10 11 12 13 14

11,980 3,561 3,898 7,705 14,030 8,342 7,899 13,060 13,960 14,770 11,410 12,420 5,745 13,880 17,490

−20,890 22,470 17,040 14,100 10,160 14,250 14,760 10,060 9,440 10,740 12,080 13,350 18,690 10,560 6,204

1358 0.0 0.0 −0.3 −0.2 −0.4 −0.4 0.1 1.5 1.9 0.0 0.3 0.4 1.3 0.6

−27,180 20,500 6,379 14,870 14,880 13,290 19,760 13,970 14,580 14,390 14,540 16,050 18,410 15,310 7,814

633.9 0.0 −0.2 1.2 0.9 0.3 1.0 0.1 −1.9 −5.0 −0.6 −2.0 −1.3 −0.7 −1.3

a ICP units are J/(kmol K). To accurately retain original fit statistics, four significant digits were retained for each coefficient. Coefficients with absolute magnitudes <1000 were then rounded to nearest 0.1. If all a2,j and a4,j coefficients in this table with magnitudes of five and less are set to zero, the calculated S.D. of the fit degrades only slightly from 0.66% to 0.69%. The AAD statistic degrades from 0.31% to 0.35%.

D.L. Morgan / Fluid Phase Equilibria 256 (2007) 54–61

59

Table 5 Fit statistics for the re-correlated Harrison and Seaton ICP method Without correction AAD (%) Maximum deviation (%) S.D. (%) No. Pts. No. compounds

Note

0.3 13.1 0.7 3115 445

Temperature (K)

Without error correction absolute deviation (%)

0.3 4.2 0.5 3115 445

P at 1000 K

Table 6 ICP transform fit deviations from Harrision and Seaton method greater than 4% Compound

With Eq. (14) error correction

With Eq. (14) error correction absolute deviation (%)

Other

800 1000

6.6 7.9

1.3 1.1

P

400 800 1000

4.9 12.1 13.1

1.8 4.2 4.2

B

500 800 1000

4.6 5.4 7.7

0.8 4.1 2.4

Al

800 1000

5.0 7.2

2.4 1.1

% deviations appear to be the result of small errors in the fit in representing the H&S-2 constants in Table 3, and of the relatively small magnitude of ICP for mono-atomics (single atoms have no rotational or vibrational modes and their high-temperature limiting ICP approaches 2.5R). The calculated constants, which were obtained by setting all nj in Eqs. (13a) and (13b) to zero, are compared with the H&S-2 constants in Table 7. The errors in this table were fitted to the quartic polynomial: error = −13, 020 + 83.54T − 0.1725T 2 + 1.381E − 4T 3 − 3.742E − 8T 4

(14)

and applied as corrections to calculated H&S-2 constants. The correction reduced the error in representing the constants from AAD = 22% to 4.7%. When applied to all the calculated ICP data used in the transform’s regression, the maximum deviation

Note P at 800 K

Table 8 Compounds with negative Aly and Lee equation coefficients produced by the ICP transform, Eqs. (13a) and (13b) Compound

Negative coefficients

Al B Br2 C HBr HCl HF HI Other Other2 P Si

a1 , a3 a1 , a3 a1 a3 a3 a3 a3 a3 a1 , a3 a1 , a3 a1 , a3 a1 , a3

from the original H&S-2 method was reduced from 13% to 4.2%. For ordinary applications, use of the quartic polynomial error fit should be unnecessary. Eq. (14) should never be extrapolated outside of the range 298.15 < T < 1500 K. Coefficients a2 and a4 represent limiting temperatures in the Aly and Lee Eq. (13a); ideally all five coefficients fitted to this equation should be positive. On applying the H&S-2 transform to generate coefficients for the 445 regression compounds, 12 compounds were found to have negative coefficients, a1 or a3 . All of these compounds were found to be monatomic or diatomic molecules which have relatively small and constant ICP compared to the ICP of more complicated polyatomic molecules. The “abnormal” fits are identified in Table 8. Fig. 3 illustrates ICP of a polyatomic molecule, hexane, and ICP of a simple diatomic, HCl, which was identified as “abnormal” in Table 8. As can be seen in the plot, both the original H&S-2 method and

Table 7 Error of ICP transform in representing constants of the Harrison and Seaton method Temperature (K)

H&S-2 constant (J/(kmol K))

Calc constant, Eqs. (13a) and (13b) (J/(kmol K))

Error (J/(kmol K))

Absolute deviation (%)

Error fit, Eq. (14) (J/(kmol K))

Calc constant + error fit (J/(kmol K))

Absolute deviation (%)

300 400 500 600 800 1000 1500

4860 864 −1850 −4610 −7490 −8530 −7370

4890 325 −2627 −4430 −6306 −7198 −8114

−30 539 777 −180 −1184 −1332 744

0.6 62.3 42.0 3.9 15.8 15.6 10.1

−59 674 544 −24 −1224 −1327 740

4831 999 −2083 −4453 −7530 −8525 −7374

0.6 15.6 12.6 3.4 0.5 0.1 0.1

AAD (%)

21.5

4.7

60

D.L. Morgan / Fluid Phase Equilibria 256 (2007) 54–61

Fig. 3. ICP of a polyatomic molecule (hexane) and of a diatomic molecule (HCl) compared with predictions of the Harrison and Seaton (H&S-2) method.

its transform fits appear to agree well with the accepted data used by DIPPR. DIPPRs accepted HCl data are from TRC [12], McBride et al. [13], and Stull [14], while its hexane data are from TRC [15]. 5. Discussion Converting legacy equations and correlations within databases to alternative equation forms is time-consuming because it is important to maintain traceability within the database and not to significantly degrade the quality of old fits without review of original data. It is an opportunity to systematically examine the database for defects and to improve some of most questionable fits. Unlike the Thermodynamics Research Center (TRC) SOURCE Data System, which stores experimental data from the literature, Dong et al. [16], DCs database only contains recommended point and temperature-dependent property fits along with temperature ranges, and several codes that represent quality and references. The original property reviews and raw data for each compound are recorded in company reports or in Excel workbooks that must be reviewed and updated whenever significant changes are made. There are tentative plans to populate an Access database with recommended, experimental, and literature review results; such a database could possibly be accessed via BYUs (Brigham Young University) DIADEM program.3 Recently the transform of Harrison and Seaton’s method was applied to check the integrity of ICP correlations in DCs INHSPCD databank and to several siloxanes which were being reviewed in the DIPPR database. It was easy to automate this procedure in Excel since the only input parameter used is an empirical formula which is available within the database. When fits in INHSPCD were compared with Aly–Lee transform fits in 3

http://dippr.byu.edu/.

the 300–1500 K range to generate AAD statistics, it was found that 32 of 371 components were missing ICP fits, and that 15 had an AAD which exceeded 10%. Of these largest discrepancies, most were attributed to typing errors made while entering coefficients into the database, errors in applying H&S-2, use of estimation methods less reliable than H&S-2, or to a unit conversion factor error of 1000. Possibly the most significant discovery was that Mosin’s [17] evaluation of ICP for hexamethylcyclotrisiloxane (D3) averaged 46% lower than H&S-2. This was significant because ICP of only one other siloxane, hexamethylsiloxane (MM), had been rigorously evaluated before by Scott et al. [18] and by Mosin and Mikhailov [19]. Was the Mosin review of D3 bad or did Mosin’s evaluation show that the H&S-2 method exhibits large deviations for ring-strained cyclic siloxanes like it does for ring-strained, carbon compounds? Use of the ICP transform to screen siloxanes within the DIPPR revealed that fits of three siloxanes (MD7M, MD10M, D8) averaged about 20% lower than H&S-2. These fits were based on ab initio calculations that were done at BYU. To test the questionable ICP results for D3, D8, MD7M, and MD10M, a non-rigorous, thermodynamic self-consistency test was applied to estimate ICP data for each siloxane. This test uses liquid heat capacity (LCP) data and the temperature derivative of heat of vaporization (dHv /dT) to estimate ICP; viz., ICP ≈ LCP + dHv /dT. The resulting estimated data agreed better with H&S-2 than with the suspect ICP results. The whole puzzle was sent to the DIPPR team at BYU to further investigate. Professor Rowley concluded that the spectroscopic frequency assignments made by Mosin in evaluating D3 thermodynamic properties were unrealistic. The group at BYU worked with a team led by Professor Colonna at Delft University, Colonna et al. [20], who eventually had speed of sound measurements done on two cyclic siloxanes, D4 & D5. The new ICP results, which were derived from these results, agree very well with H&S-2 and will be published in a future edition of Fluid Phase Equilibria. The BYU team adjusted the scaling parameters and basis sets used in their ab initio calculations for siloxanes and will also publish their results in an upcoming issue of Fluid Phase Equilibria. Generally the uncertainty code assigned to H&S-2 method by DIPPR is <25%. Our experience with reviewing a couple dozen ICP evaluations for silanes and siloxanes suggests that the uncertainty code could possibly be set to <10% for siliconcontaining materials. The uncertainty code which DC assigns to Pitzer’s Hv method is 10%. This is the same code as what DIPPR assigns to Clapeyron-derived Hv . When results of the Pitzer and Clapyron methods for predicting Hv are compared, usually the differences fall within this level of uncertainty. Uncertainty codes assigned to the H&S-1 VP method are much more variable and will depend on reviewer assessments of the quality of vapor pressure used to fit parameters, Tc , Pc , w, and x. 6. Conclusions Transforms, which can convert legacy physical property equations or correlation methods in databases to more acceptable equation forms, can sometimes be developed by fitting tabular

D.L. Morgan / Fluid Phase Equilibria 256 (2007) 54–61

data or by re-arranging existing equations. These are powerful quality-control tools for locating database entry errors, unit conversion mistakes, thermodynamic inconsistencies, or problems with published literature results. Transforms which are sufficiently faithful to original methods can be used to help populate a database with missing properties. Statistics obtained from fitting the transformations over the entire range of databank input parameters can help database administrators rank conversion methods to be used on particular components. Sometimes a more time-consuming re-analysis of original data is preferable over use of automated fitting methods. Acknowledgements We wish to thank AnnaMaria Morgan for her patience and understanding while this manuscript was prepared, and the Dow Corning Corporation for supporting this work. References [1] R.L. Rowley, W.V. Wilding, J.L. Oscarson, Y. Yang, T.E. Daubert, R.P. Danner, DIPPR® Data Compilation of Pure Chemical Properties, Design Institute for Physical Properties, AIChE, New York, NY, 2006. [2] B.I. Lee, M.G. Kesler, A generalized thermodynamic correlation based on three-parameter corresponding states, AIChE J. 21 (3) (1975) 510– 527. [3] K.S. Pitzer, D.Z. Lippman, R.F. Curl Jr., C.M. Huggins, D.E. Petersen, The volumetric and thermodynamic properties of fluids. I. Theoretical basis and virial coefficients. II. Compressibility factor, vapor pressure and entropy of vaporization, J. Am. Chem. Soc. 77 (1955) 3427–3440. [4] R.L. Halm, L.I. Stiel, A fourth parameter for the vapor pressure and entropy of vaporization of polar fluids, AIChE J. 13 (2) (1967) 351–355. [5] D.L. Morgan, Extension of Pitzer Corresponding States Correlations using New Vapor Pressure Measurements of the n-Alkanes C10 to C28, Ph.D. Dissertation, Rice University, Houston, TX, 1990. [6] S. Torquato, G.R. Stell, An equation for the latent heat of vaporization, Ind. Eng. Chem. Fundam. 21 (1982) 202–205.

61

[7] B.K. Harrison, W.H. Seaton, Solution to missing group problem for estimation of ideal gas heat capacities, Ind. Eng. Chem. Res. 27 (8) (1988) 1536–1540. [8] B.K. Harrison, A. Madas, A. Sharma, The ASTM Computer Program for Chemical Thermodynamic and Energy Release Evaluation CHETAH 8.0 User Guide, ASTM International, West Conshohocken, Pennsylvania, 2005. [9] S.W. Benson, Thermochemical Kinetics, 2nd ed., John Wiley & Sons, New York, 1976. [10] K.H. Myers, R.P. Danner, Prediction of properties of silicon, boron, and aluminum compounds, J. Chem. Eng. Data 38 (2) (1993) 175–200. [11] F.A. Aly, L.L. Lee, Self-consistent equations for calculating the ideal gas heat capacity, enthalpy, and entropy, Fluid Phase Equilib. 6 (3/4) (1981) 169–179. [12] Thermodynamics Research Center (TRC), Selected Values of Properties of Chemical Compounds, Data Project, Texas A&M University, College Station, TX, 1980. [13] B.J. McBride, S. Heimel, J.G. Ehlers, S. Gordon, Thermodynamic Properties to 6000 K for 210 Substances Involving the First 18 Elements, Office of Scientific and Technical Information, National Aeronautics and Space Administration, Washington, DC, 1963. [14] D.R. Stull, JANAF Thermochemical Tables, 2nd ed., U.S. Department of Commerce, Washington, DC, 1971, NSRDS-NBS 37. [15] Thermodynamics Research Center (TRC), Selected Values of Hydrocarbons and Related Compounds, Hydrocarbon Project, Texas A&M University, College Station, Texas, 1984. [16] Q. Dong, X. Yan, R.C. Wilhoit, X. Hong, R.D. Chirico, V.V. Diky, M. Frenkel, Data quality assurance for thermophysical property databases— applications to the TRC SOURCE data system, J. Chem. Inf. Comput. Sci. 42 (2002) 473–480. [17] A.M. Mosin, Thermodynamic functions of dimethylsiloxane and hexamethylcyclotrisiloxane, Russ. J. Phys. Chem. 49 (3) (1975) 437. [18] D.N. Scott, J.F. Messerly, S.S. Todd, G.B. Guthrie, I.A. Hossenlopp, R.T. Moore, A. Osborn, W.T. Berg, J.P. McCullough, Hexamethyldisiloxane: chemical thermodynamic properties and internal rotation about the siloxane linkage, J. Phys. Chem. 65 (1961) 1320. [19] A.M. Mosin, A.M. Mikhailov, Thermodynamic functions of hexamethyldisiloxane, Zh. Fiz. Khim. 46 (1972) 537. [20] P. Colonna, N.R. Nannan, A. Guardone, E.W. Lemmon, Multiparameter equations of state for selected siloxanes, Fluid Phase Equilib. 244 (2006) 193–211.