2332
J. Phys. Chem. B 2006, 110, 2332-2337
Prediction of Hydrogen Hydrate Equilibrium by Integrating ab Initio Calculations with Statistical Thermodynamics Jae W. Lee,† Prasad Yedlapalli,† and Sangyong Lee*,‡ Department of Chemical Engineering, The City College of the CUNY, New York, New York 10031, and Department of Chemical and Natural Gas Engineering, Texas A&M UniVersity-KingsVille, KingsVille, Texas 78363 ReceiVed: June 10, 2005; In Final Form: December 5, 2005
This paper addresses a new calculation approach for the prediction of hydrogen hydrate equilibrium by introducing the concept of a single hydrogen cluster in one cavity. By integrating ab initio calculations with classical statistical thermodynamics, this approach enables the van der Waals model to predict the dissociation pressure of hydrogen hydrates. Compared to hydrates formed by light hydrocarbon gases, structure II (sII) hydrogen hydrates stably encage two and four hydrogen molecules in the small and large cavities, respectively. By treating two hydrogen molecules or four hydrogen molecules as one rigid body cluster, we determine ab initio binding energies between water molecules and hydrogen clusters at the MP2 level with the 6-31++G(2d,2p) basis set. These binding energies will be used to determine the parameters of the Exp-6 potential function from which the smooth cell potential and the Langmuir constant of each cluster are calculated. Then, the dissociation pressure is determined using the Zele-Lee-Holder cell distortion model: 105, 625, and 2000 bar at 150, 200, and 250 K, respectively.
1. Introduction In the past few decades, the research and development of new hydrogen storage methods and applications have been driven by the need for a new energy source that will replace the limited amount of petroleum.1 Consequently, the technology of using hydrogen as an environmentally clean and efficient fuel became an active research area. To achieve a secure supply of hydrogen, a safe and economical storage technique needs to be developed. Hydrogen gas hydrates are one such promising medium for hydrogen storage.2 Hydrogen gas hydrates’ storage density is 50 g/L (5.3 wt %) around 200 MPa and 250 K.3 It was experimentally3 and theoretically4 confirmed that two hydrogen molecules and four hydrogen molecules are stably encaged in the small and large cavities (512 and 51264), respectively, of type II gas hydrates. Recently, Sloan and his colleauges5 synthesized low-pressure hydrogen hydrates with the addition of a second component (tetrahydrofuran, THF), and Lee et al.6 increased the hydrogen storage density up to 4.03 wt % by manipulating THF molar concentrations at 120 bar and 270 K. Despite this series of important discoveries about hydrogen hydrates, there is no unified thermodynamic model to predict their equilibrium behavior over a wide range of conditions. Understanding the thermodynamic stability of hydrogen hydrates is essential for designing suitable storage and transport systems. The main difficulty is that the conventional van der Waals and Platteuuw’s statistical mechanics model7 assumes the single occupancy of a gas molecule in one cavity and cannot calculate the chemical potentials of H2 hydrates due to the multiple * To whom correspondence should be addressed. MSC 193, 700 University Blvd, Department of Chemical and Natural Gas Engineering, Texas A&M University, TX 78363. Tel: 361-593-2629. Fax: 361-5934026. E-mail:
[email protected]. † The City College of CUNY. ‡ Texas A&M University-Kingsville.
occupancy of H2 in one cavity. For the group of two or four H2 molecules, no experimental or theoretical values of parameters in pair potential functions are available to calculate their Langmuir constants. In this work, we will overcome this difficulty by incorporating ab initio calculations into the van der Waals model modified by cavity distortion.8,9 Ab initio calculations will proceed by considering two or four hydrogen molecules as one rigid body. In other words, the 4H2 molecules in the large cavity will be treated as a single, rigid tetrahedron cluster and the 2H2 molecules in the small cavity will be considered as a single, rigid cluster.4 From the ab initio calculations, the intermolecular potential between a hydrogen cluster and a water molecule will be determined and will be fit to the Exp-6 potential form. Then, a smooth cell potential will be derived from the Exp-6 potential, and the Langumir constants for each H2 cluster will be determined. These Langmuir constants will be used in the ZeleLee-Holder model8,9 to determine the equilibrium conditions of the hydrogen hydrates. 2. Determination of Intermolecular Potential by ab Initio Calculations Hydrogen hydrates form type II gas hydrates with 136 water molecules to form 24 cavities per unit cell. The hydrogen hydrate contains 61 ( 7 hydrogen molecules in a unit cell and the general formula of a hydrogen hydrate is H2(H2O)23. Two hydrogen molecules are stably enclathrated in a small cavity of 12 pentagonal faces (512, dodecahedron) and 4 hydrogen molecules in a large cavity with 12 pentagonal and 4 hexagonal faces (51264, hexakaidecahedron) as shown in Figure 1. To determine the optimized geometries of H2 molecules in an H2 cluster, Patchkovskii and Tse4 employed the density functional theory (DFT). The optimal distance between the center of mass of H2 molecules was determined as 2.58 Å in
10.1021/jp0531311 CCC: $33.50 © 2006 American Chemical Society Published on Web 01/17/2006
Hydrogen Hydrate Equilibrium
J. Phys. Chem. B, Vol. 110, No. 5, 2006 2333
Figure 1. Two H2 molecules and four H2 molecules in the small and large cavities of sII hydrates.
the 2H2 cluster of the small cavity. The 4H2 molecules in the large cavity form a slightly distorted tetrahedron with the H2H2 distance varying between 2.90 and 3.13 Å. To calculate binding energies between the hydrogen clusters (bi-hydrogen or tetra-hydrogen cluster) and one water molecule, the ab initio second-order Mo¨ller-Plesset (MP2) calculations will be performed using the GAMESS program.10 The basis set of 6-31++G(2d,2p) is used for the calculations. The distances between hydrogen molecules are selected as 2.58 Å for the bi-hydrogen cluster and 2.90 Å for the tetra-hydrogen cluster to allow more free rotation and more geometrically valid configurations within the water cavity (that is, not lying outside the cavity) than the case of 3.13 Å. These two hydrogen clusters are assumed to be rigid bodies in the calculations. By varying the distance (r) between the water molecule and the hydrogen clusters, the binding energies can be calculated. Following Cao et al.’s pair potential calculation of methane hydrates,11 the binding energy between the H2 cluster and water is given in eq 1 as a function of the distance (r) between a water molecule (oxygen atom) and the center of mass of the H2 cluster.
∆Ecp ) Ehw - Eh,cp - Ew,cp
#ofconfi.
〈∆Ecp(r)〉 )
∆Ecp(r) exp(-∆Ecp(r)/kT) (2)
#ofconfi.
∑
TABLE 1: Number of Configurations for ab Initio Calculations small cavity parameter
range
'r, A0 ξ, deg φ, deg R, deg β, deg γ, deg no. of configurations
1.6-6.0 -40 to +40 -40 to +40 0-180 0-180 0-180
no. of points 12 5 5 3 3 3 8100 × 2 ) 16 200
large cavity range 1.6-6.0 -40 to +40 -40 to +40 0-120 0-360 0-120
no. of points 12 5 5 3 4 3 10 800 × 2 ) 21 600
(1)
where ∆Ecp is the binding energy (or pair potential), Ehw is the total energy of the hydrogen cluster-water system at a specific configuration, and Eh,cp and Ew,cp are the counterpoise-corrected potentials12 of the hydrogen cluster and water molecules. While Anderson et al.13 employed 50% of the counterpoise correction with a larger basis set (aug-cc-pVQZ) for the accurate prediction of the molecular potentials of methane and water, this study uses full counterpoise correction with the 6-31++G(2d,2p) basis set to achieve reasonable accuracy. With respect to all of the configurations possible at each H2 cluster-oxygen distance (r), the average potentials will be calculated as the following Boltzmann-weighted average11
∑
Figure 2. Geometry for the calculation of interatomic potentials between a bi-hydrogen cluster and water. For the tetra-hydrogen cluster, a tetrahedral hydrogen cluster will replace the linear bi-hydrogen cluster.
exp(-∆Ecp(r)/kT)
where k is the Boltzmann constant and T is the absolute temperature. Here, 250 K is used for the calculation of the Boltzmann-weighted average since this temperature is one of the experimental values3 that is relatively easy to obtain in the lab. Using eq 2 predicts a fractional occupancy for methane larger than the experimental values as the temperature increases above 273 K14 as reported by Anderson et al.13 However, this work employs this Boltzmann average because the temperature range of interest is 78-250 K (below 273 K) and the prediction error for the H2 fractional occupancy would be smaller in this temperature range.
The average value of the pair potentials will include all of the possible configurations by modifying the geometric variations11 of methane hydrates: (1) The distance between the oxygen atom of water and the center of mass of the hydrogen clusters will be chosen between 1.6 and 6 Å. Distances less than 1.6 Å are not taken to avoid very high repulsive binding energies. The range of polar angles (φ and ξ) as shown in Figure 2 will be chosen between -40° and +40° so that hydrogen clusters are not outside the cavity structure at the minimum distance (r ) 1.6 Å). (2) The rotation of hydrogen clusters will be considered in terms of Euler angles R, β, and γ with respect to x′, y′, and z′. For the bi-hydrogen cluster, these Euler angles of R, β, and γ will change from 0 to 180°. For the tetra-hydrogen cluster, the hydrogen molecules form a tetrahedron,4 R and γ will change between 0 and 120°, and β will change from 0 to 360°. (3) We will have 12 points in the distance (r), 3 points in each Euler angle, and 5 points in each polar angle (φ and ξ). So, for one water orientation (parallel to the y-z plane as shown in Figure 2), we calculate 8100 () 12 × 3 × 3 × 3 × 5 × 5) different configurations. When another water orientation (parallel to the x-z plane) is included, the total number of configurations will be 16 200 (8100 × 2). For a distance (r), we will take an average value of 1350 interaction energies using eq 2. For the large cavity, we sample 4 points from Euler angle β. Thus, its total number of configurations is 21 600 () 12 × 3 × 4 × 3 × 5 × 5 × 2). Table 1 shows the detailed configurations for ab initio calculations. Figure 3 shows the calculated ab initio binding energies (shown by circles) between a water molecule and two types of H2 clusters in the small and the large cavities, respectively. The
2334 J. Phys. Chem. B, Vol. 110, No. 5, 2006
Lee et al. TABLE 2: Exp-6 Curve Fitting Results interaction pair
2H2-H2O
4H2-H2O
, j/mol R rm, A°
694.70 4.0434 3.28
674.9456 4.2628 4.0418
A smooth cell potential is obtained by calculating the total interaction energy between the H2 cluster and the water cage using the Exp-6 potential and is given by
〈W〉 )
[ [
z 6eR rm ((a - R)e-R(a-R)/rm 2(R - 6) aR R
(a + R)e-R(a+R)/rm) +
r2m R2 R
]
[e-R(a-R)/rm - e-R(a+R)/rm] -
(
)]
r6m 1 1 4 4aR (a - R) (a + R)4
(4)
where z is the coordination number of the cavity, a is the cell (cavity) radius, and R is the distance between a hydrogen cluster and the center of the cavity. The detailed derivation from eq 3 to eq 4 is available in the Appendix. With the smooth cell potential, the Langmuir constant for gas molecule j in a cavity of type i is determined by taking a volume integral.
Cij )
4π kT
∫0a exp(-
)
〈W(R)〉 2 R dR kT
(5)
4. Water Chemical Potentials in the Water and Hydrate Phases At thermodynamic equilibrium, the chemical potentials of water in the hydrate and in the other coexisting phases are equal.7 Then Figure 3. Binding energy curves from ab initio calculations: (a) 2H2H2O, (b) 4H2-H2O.
attractive force goes to zero for distances larger than 6 Å in both cases. The minimum binding potentials occur around the distances of 3.2 and 4.0 Å for the small and large cavities. The repulsive part becomes significant when the distance is less than 2 Å. The real time required for the calculation of one configuration is around 0.5-1 min on a 1.2 GHz PC, and it takes approximately one week for each cavity. 3. Smooth Cell Potential from Exp-6 Pair Potential and Langmuir Constant We have first tried to fit the calculated potential curves in Figure 3 to the Kihara pair potential model. But one of the parameters (σ) in the Kihahra potential turns out to be negative. Even if we calculate very high potentials with the distance (r) of 0.5-1.5 Å and use those values for the fitting, the parameter is still negative. Thus, the following Exp-6 pair potential is used to fit the potential curves in Figure 3.
〈∆Ecp(r)〉 ) φ )
[ ((
)) ( ) ]
r 6 exp R 1 R-6 rm
-R
rm r
µwR (T,P) ) µwH (T,P)
where µwR (T,P) is the water chemical potential in the pure liquid or ice phase, and µwH (T,P) is the water chemical potential in the hydrate phase at temperature T and pressure P. By the use of the chemical potential of the theoretical empty hydrate lattice, µβ, as the reference state7 (273.15 K and zero pressure), the equilibrium equation becomes
∆µwR ) ∆µwH
where the negative value of represents the minimum potential energy at a separation distance of rm. R describes the steepness of the repulsive part. The curve fitting results are shown in Table 2.
(7)
where ∆µwR ) µβ - µwR and ∆µwH ) µβ - µwH. The statistical thermodynamic model for the hydrate phase as stated by van der Waals and Platteeuw7 is
∆µwH
2
)-
RT
∑ Vi ln(1 - ∑ θij) i)1
(8)
where Vi is the number of i-type cavities per water molecule and θij is the fractional occupancy of i-type cavities with j-type molecules. The fractional occupancy is expressed as
6
(3)
(6)
θij )
Cij fj (1 +
∑j
(9) Cijfj)
where f j is the fugacity of the gas. Peng-Robinson’s correlation15 is used for the calculation of the fugacity in this paper
Hydrogen Hydrate Equilibrium
J. Phys. Chem. B, Vol. 110, No. 5, 2006 2335 TABLE 3: Reference Parameters and Optimized Cavity Radii of H2 Hydrates ∆µw0 (J/mol) ∆hw0 (J/mol) small cavity radius (Å) large cavity radius (Å)
Figure 4. Flowchart for determining reference values.
by assuming that it can be extrapolated from 8016 or 10017 to 200 MPa. It is also assumed that the fugacity of gas is independent of the occupancy of H2 molecules. Cij is the Langmuir constant that can be estimated as shown in the previous section. The model developed by Holder et al.18 is used to predict the chemical potential of water in the aqueous phase that is in equilibrium with the hydrates. The model is given as
∆µwR ∆µw0 ) RT RT0
∆h
∫TT RTw2 dT + ∫0P 0
∆Vw dP - ln γwXw (10) RT
The first term on the right-hand side is the reference chemical potential difference, which is experimentally determined based on the distortion model.8,9 The second term gives the temperature dependence of enthalpy at constant pressure. The third term accounts for the change in the chemical potential difference due to pressure. The fourth term gives the change in chemical potential due to the activity change of water. 5. Calculation Procedure for Determining Equilibrium Pressure of H2 Hydrates The temperature range used for this study is 70 to 250 K, which is below the ice point and is the same range as experimental and theoretical studies.3-4 Optimal radii of distorted cavities and corresponding reference chemical potential and enthalpy differences (∆µw0 and ∆hw0) are determined using Zele et al.’s model8 and the Lee-Holder model.9
Figure 5. Dissociation pressures of H2 hydrates.
1486.49 693.66 3.89 4.73
The procedure for determining these reference values with the cavity radii is in the following as shown in Figure 4: (1) An initial guess of the small cavity radius is made with reported values between 3.83 and 3.96 Å.4 With the initial guess of the radius, the reference chemical potential difference is calculated in terms of a ) R + β × ∆µw0, where there are two sets of parameters for small (R ) 3.7547 and β ) 3.8441 × 10-4) and large (R ) 4.5628 and β ) 4.6878 × 10-4) water cavities.8 The initial guess of the large cavity radius is determined by this formula with its R and β. (2) With these two radii and a given temperature, eq 5 gives the Langmuir constants (Cij) for bi-hydrogen and tetra-hydrogen clusters in the small and large cavities. The fractional occupancy in eq 9 is determined by using the Langmuir constants and the calculated fugacity from the Peng-Robinson equation of state.13 Then, the water chemical potential in the hydrate phase (∆µH) is calculated in terms of eq 8. (3) This hydrate phase chemical potential is equal to the righthand side of eq 10, and a new reference chemical potential (∆µw0) is determined by using a reported equilibrium value4 of 2000 bar at 250 K. Here, the reference enthalpy difference (∆hw0) in ∆hw of eq 10 is updated in the inner loop with respect to the new chemical potential (∆µw0) and three pairs of experimental data:3,4 0.1 bar, 1 bar, and 2000 bar at 78, 80, and 250 K. To update ∆hw in the inner loop, we use a Y ) aX form in the Lee-Holder distortion model9 (refer to eq 11 in ref 9). An ideal solution is assumed with the diluted H2 in the aqueous phase (γw ) 1, Xw ) 1 in eq 10). (4) The ∆µw0 and ∆hw0 are to be continuously updated to match with the experimental dissociation pressures for the previously given cavity radii. With the converged ∆µw0, the new cavity radii will be obtained by a ) R + β × ∆µw0. Then, the Cij values are updated. With these updated radii, Cij values, and fugacity, the chemical potential of the hydrate phase (∆µwH) is updated. (5) If ∆µwH is equal to the right-hand side of eq 10 with the previous ∆µw0 and ∆hw0, then the final values of ∆µw0, ∆hw0, and the two radii are determined. Otherwise, repeat steps (3) through (5) until converged to ∆µw,new0 ) ∆µw0. Table 3 shows the calculated reference chemical potential difference, reference enthalpy difference, and two cavity radii by using the cavity distortion model.8-9 These values, along
2336 J. Phys. Chem. B, Vol. 110, No. 5, 2006 with the pair potential parameters in Table 2, allow us to predict dissociation pressures over a wide range of temperatures as shown in Figure 5. The experimental values of 0.1 and 1 bar at 78 and 80 K are collected from the sharp peaks of the Raman and IR spectroscopy analyses.3 Another reported value of 2000 bar at 250 K is borrowed from the theoretical stability calculation.4 In real experiments, this value is given as 18002200 bar at 249 K. So, the average value of 2000 bar at 250 K can be considered to be one of the experimental values. Some other experimental values19,20 are available but were not employed for this study because multiple equilibrium pressures exist at a single temperature and these values are too delicate to use. The calculated dissociation pressures are 0.29 and 0.39 bar at 78 and 80 K, respectively, and they show only a slight mismatch with the experimental values. The calculated dissociation pressures at 150, 200, and 250 K are 105, 625, and 2000 bar, respectively. While the Boltzmann average gives higher occupancy in the small cavity than in the large cavity for methane hydrates,11,21 employing the Boltzmann average11 in our calculation scheme for H2 hydrates does not cause such an aphysical situation. As shown in Figure 5, the occupancy of the tetra-hydrogen cluster in the large cavity is always a little bit higher than that of the bi-hydrogen cluster in the small cavity. The possible reasons why we do not have methane’s aphysical situation with the Boltzmann average would be (1) ab initio intermolecular potentials between water and the H2 cluster are different from those between methane and water;21 (2) Exp-6 is used in this study while Cao et al.21 employed Lennard-Jones and Kihara potential forms for the calculation of methane hydrate equilibrium. Figure 5 shows that the predicted dissociation pressure of H2 hydrate with respect to temperature follows the same trend as other light hydrocarbon gas hydrates. In particular, we have the capability to predict the dissociation pressure from 100 to 250 K where experimental data are absent. Conclusion A new calculation procedure for H2 hydrate equilibria has been proposed by combining ab initio calculations and the ZeleLee-Holder distortion model. The pair binding energies between water molecules and hydrogen clusters are determined by ab initio calculations. The Exp-6 model gives the best fit for these binding energies. The smooth cell potential is derived from the Exp-6 model and is used for calculating the Langmuir constants for the bi-hydrogen and tetra-hydrogen clusters in the small and large cavities. Using these Langmuir constants, we can determine the water chemical potentials in the hydrate and aqueous phases. The values of reference chemical potential and enthalpy differences are determined in terms of the Zele-LeeHolder model and are used for the prediction of H2 hydrate dissociation pressures.
Lee et al. fj ) fugacity of gas component j hwL, hwR, hwβ ) molar enthalpy of pure liquid water, ice, and empty hydrate lattice, respectively, J/mol k ) Boltzmann constant, 1.38 × 10-23 J/K T ) absolute temperature, K z ) coordination number of the cavity r ) distance between the water molecule (oxygen atom) and the center of mass of the hydrogen cluster, Å rm ) parameter in the Exp-6 form, Å R ) radial coordinate of the H2 cluster, Å Vw ) molar volume of water 〈W(R)〉 ) spherically symmetric smooth cell potential, J xw ) water mole fraction in liquid phase Greek Letters φ, ξ ) polar coordinates of the H2 cluster as shown in Figure 2 ) depth of intermolecular potential well, J R ) parameter in Exp-6 form γw ) water activity coefficient in liquid phase µ ) chemical potential, J/mol µwH, µwL, µwR, µwβ ) chemical potential of water in hydrate, liquid, ice, and empty hydrate phases, respectively, J/mol Vi ) number of i-type cavities per water molecule θij ) fractional occupancy of i-type cavities with j-type molecules Superscripts 0 ) property at 273.15 K and 0 atm R ) ice phase β ) empty hydrate phase i ) cavity type i j ) component w ) water Appendix: Derivation of Spherically Symmetric Smooth Cell Potential (〈W〉) Here, the cavity is treated as a sphere of radius a. It is assumed that all of the water molecules are evenly smeared over the surface of the sphere. Chart 1 refers to a slice of that sphere.22 With the above assumptions, the number of water molecules in dA is
dA z 4πa2
(A-1)
The Exp-6 intermolecular potential energy form is used to represent the interaction energy between a water molecule and an H2 cluster. Then, the interaction energy between this element CHART 1 a
Nomenclature a ) cavity radius, Å Cij ) Langmuir constant for component j in cavity I ∆Ecp ) binding energy of the water molecule and the H2 cluster, J Ehw ) total energy of the hydrogen cluster-water system at a specific configuration, J Eh,cp ) counterpoise-corrected potential of the hydrogen cluster, J Ew,cp ) counterpoise-corrected potential of the water molecule, J
a Key: a, radius of a spherical cavity; O, origin; P, position of the center of mass of the H2 cluster; R, distance between O and P; W, position of differential area element; dA, surface area of the differential element containing point W, dA ) (a sin θ dφ)(a dθ); r, distance between P and W; z, coordination number of the cell.
Hydrogen Hydrate Equilibrium
J. Phys. Chem. B, Vol. 110, No. 5, 2006 2337
and the H2 cluster is
φ(r) )
[ ( ( )) ( ) ]
r 6 exp R 1 R-6 rm
-R
rm r
6
z dA 4πa2
(A-2)
With the assumption of pair-wise additivity, the total intermolecular potential energies between the H2 cluster and the water molecules are
〈W〉 ) 〈W〉 )
[ ( ( )) ( ) ] [ ( ( )) ( )]
∫ R - 6
6 exp R 1 -
∫02π ∫0π R - 6
r rm
6 exp R 1 R
rm r
6
-R
r rm
rm r
6
z dA 4πa2
-
za2 sin θ dθ dφ (A-3) 4πa2
In the ∆OPW, we have r2 ) R2 + a2- 2aR cos θ. Differentiating gives 2r dr ) 2aR sin θ dθ. Inserting this into eq A-3 and integrating with the limits of θ ) 0 ∼ π, r ) a - R ∼ a + R, and φ ) 0 ∼ 2π gives
〈W〉 )
〈W〉 )
z 2(R - 6)
[ ( ( )) ( ) ]
∫aa-+RR
6 exp R 1 -
[ [
r rm
-R
6eR rm z ((a - R)e-R(a-R)/rm 2(R - 6) aR R
(a + R)e
-R(a+R)/rm
rm2
6
r dr aR
] )]
[e-R(a-R)/rm - e-R(a+R)/rm] R2 rm6 1 1 R (A-4) 4aR (a - R)4 (a + R)4
)+
(
References and Notes
rm r
(1) Conte, M.; Prosini, P. P.; Passerini, S. Mater. Sci. Eng. B 2004, 108, 2.
(2) Mao, W. L.; Mao, H. K. Proc. Natl. Acad. Sci. 2004, 101, 708. (3) Mao, W. L.; Mao, H. K.; Goncharov, A. F.; Struzhkin, V. V.; Guo, Q.; Hu, J.; Shu, J.; Hemley, R. J.; Somayazulu, M.; Zhao, Y. Science 2002, 297, 2247. (4) Patchkovskii, S.; Tse, J. S. Proc. Natl. Acad. Sci. 2003, 100, 14645. (5) Florusse, L. J.; Peters, C. J.; Schoonman, J.; Hester, K. C.; Koh, C. A.; Dec, S. F.; Marsh, K. N.; Sloan, E. D. Science 2004, 306, 469. (6) Lee, H.; Lee, J.-W.; Kim, D. Y.; Park, J.; Seo, Y.-T.; Zeng, H. I.; Moudrakovski, L.; Ratcliffe, C. I.; Ripmeester, J. A. Nature 2005, 434, 743. (7) van der Waals, J. H.; Platteeuw, J. C. AdV. Chem. Phys. 1959, 2, 1. (8) Zele, S. R.; Lee, S.-Y.; Holder, G. D. J. Phys. Chem. B 1999, 43, 10250. (9) Lee, S.-Y.; Holder, G. D. AIChE J. 2002, 48, 161. (10) Schmidt, M. W.; Baldridge, K. K.; Boatz, J. A.; Elbert, S. T.; Gordon, M. S.; Jensen, J. H.; Koseki, S.; Matsunaga, N.; Nguyen, K. A.; Su, S.; Windus, T. L.; Dupuis, M.; Montgomery, J. A., Jr. J. Comput. Chem. 1993, 14, 1347. (11) Cao, Z.; Tester, J. W.; Trout, B. L. J. Chem. Phys. 2001, 115, 2250. (12) Boys, S. F.; Bernardi, F. Mol. Phys. 1970, 19, 553. (13) Anderson, B. J.; Tester, J. W.; Trout, B. L. J. Phys. Chem. B 2004, 108, 18705. (14) Cao, Z.; Tester, J. W.; Sparks, K. A.; Trout, B. L. J. Phys. Chem. B 2001, 105, 10950. (15) Peng, D. Y.; Robinson, D. B. Ind. Eng. Chem. Fundam. 1976, 15, 59. (16) Sloan, E. D. Clathrate Hydrates of Natural Gases, 2nd ed.; Mercel Dekker: New York, 1998; p 265. (17) Rudoph, E. S. J.; Langeveld, J. H.; de Loos, Th. W.; de Swaan Arons, J. Fluid Phase Equilibr. 2000, 173, 81. (18) Holder, G. D.; Corbin, G.; Papadopoulos, K. D. Ind. Eng. Chem. Fundam. 1980, 19, 282. (19) Dyadin, Y. A.; Larinov, E. G.; Manakov, Y. A.; Zhurko, F. V.; Aladko, E. Y.; Mikina, T. V.; Komarov, V. Y. MendeleeV Commun. 1999, 9, 209. (20) Vos, W. L.; Finger, L. W.; Hemley, R. J.; Mao, H.-K. Phys. ReV. Lett. 1993, 71, 3150. (21) Cao, Z.; Tester, J. W.; Sparks, K. A.; Trout, B. L. J. Phys. Chem. B 2001, 105, 10950. (22) John, V. T.; Holder, G. D. J. Phys. Chem. 1985, 89, 3279.
