Calculation of Chemical Potentials and Occupancies in Clathrate

Aug 12, 2013 - Department of Chemical Engineering, Indian Institute of Science, Bangalore 560012, Kamataka, India ... Copyright © 2013 American Chemi...
0 downloads 0 Views 491KB Size
Article pubs.acs.org/JPCC

Calculation of Chemical Potentials and Occupancies in Clathrate Hydrates through Monte Carlo Molecular Simulations Srikanth Ravipati and Sudeep N. Punnathanam* Department of Chemical Engineering, Indian Institute of Science, Bangalore 560012, Kamataka, India ABSTRACT: The flexibility of the water lattice in clathrate hydrates and guest−guest interactions has been shown in previous studies to significantly affect the values of the thermodynamic properties, such as chemical potentials and free energies. Here we describe methods for computing occupancies, chemical potentials, and free energies that account for the flexibility of water lattice and guest−guest interactions in the hydrate phase. The methods are validated for a wide variety of guest molecules, such as methane, ethane, carbon dioxide, and tetrahydrodfuran by comparing the predicted occupancy values of guest molecules with those obtained from isothermal isobaric semigrand Monte Carlo simulations. The proposed methods extend the van der Waals and Platteuw theory for clathrate hydrates, and the Langmuir constant is calculated based on the structure of the empty hydrate lattice. These methods in combination with development of advanced molecular models for water and guest molecules should lead to a more thermodynamically consistent theory for clathrate hydrates.

1. INTRODUCTION Clathrate hydrates are inclusion compounds in which the host lattice is made up of water molecules connected to each other in a tetrahedral manner via hydrogen bonds. The water lattice contains cavities which are occupied by guest molecules such as methane, ethane, etc. The presence of these guest molecules stabilizes the hydrate by lowering the chemical potential of water in the hydrate phase. The study of clathrate hydrates has been historically important due to their role in blockage of natural gas pipelines.1 In recent years, their importance has increased, as they play an important role in many scientific and technological areas, such as gas storage and transportation,2,3 climate change,4 and desalination of water.5−7 Natural gas hydrates containing methane are also considered as a potential fuel source for the future.8−11 The current theoretical understanding of clathrate hydrates is based on the statistical thermodynamic model proposed by van der Waals and Platteeuw known as vdWP theory.8,12 The vdWP theory uses the Langmuir adsorption model to describe the thermodynamic properties of clathrate hydrates. Accordingly, the water lattice is assumed to be rigid; each cavity is assumed to contain at most one guest molecule; interaction of the guest molecules with the surrounding lattice is limited to the immediate neighboring water molecules; and interactions with other guest molecules are ignored. Understandably, this model has its drawbacks,13−19 and a lot of efforts have been made to improve the theory, such as (i) inclusion of long-range guest−host interactions,20−22 (ii) rigorous Monte Carlo integration of the Langmuir constant,23,24 (iii) inclusion of guest−guest interactions using the mean-field approximation,14,15,25 and (iv) inclusion of multiple cage occupancy26 for small molecules. In recent years, molecular simulations14−16,18,27−45 have been widely used to compute various properties of clathrate hydrates. Molecular simulations provide exact results for a given © 2013 American Chemical Society

molecular model and are able to overcome the limitations of various theories for clathrate hydrates. However, molecular simulations of clathrate hydrates are still computationally intensive, and developments leading to improvement of current theories for clatharate hydrates are highly desirable. In this study, our focus is on improving the vdWP theory for predicting thermodynamic properties of clathrate hydrates. Within the framework of the vdWP theory, the Langmuir constant is an important quantity that describes the guest− water interactions in the hydrate. Among the widely used versions of vdWP theory,8,12,46−49 the Langmuir constant is computed after assuming the water lattice to be rigid. However, recent studies13,14,25 have shown that the assumption of a rigid water lattice is incorrect and can lead to significant errors in predictions of phase equilibria for clathrate hydrates. This suggests that the flexibility of the water lattice has to be accounted for during the calculation of the Langmuir constant. To date, none of the suggested improvements for the vdWP theory accounts for the flexibility of the water lattice. In this paper, we propose a method for computing the Langmuir constant that accounts for the flexibility of the water lattice in the hydrate phase. In addition to the calculation of the Langmuir constant, we also propose a method for computing guest−guest interactions. The efficacy of these methods is tested by comparing the predictions of occupancy of clathrate cavities by various guest molecules with values obtained from Monte Carlo simulations of clathrate hydrates in the isothermal−isobaric semigrand ensemble.14,15 The clathrate guest molecules chosen for this study are methane, ethane, carbon dioxide (CO2), and tetrahydrofuran (THF). Of these, methane and ethane are hydrophobic molecules with very low Received: June 11, 2013 Revised: July 29, 2013 Published: August 12, 2013 18549

dx.doi.org/10.1021/jp405771u | J. Phys. Chem. C 2013, 117, 18549−18555

The Journal of Physical Chemistry C

Article

solubility in water, whereas carbon dioxide and THF have significant solubility in water. Also, THF molecules can form hydrogen bonds with water molecules even in the hydrate phase and thus affect the structure of the water lattice upon inclusion.

Z(N , Nw , P , T ) = Z 0(N , Nw , P , T ) Z(Nw , P , T ) 1 = N!

∫ ... ∫ e−β[U



2. THEORY The thermodynamic stability of clathrate hydrates is due to the presence of guest molecules in the cavities. The inclusion of guest molecules inside the cavities of clathrate lowers the fugacity (or chemical potential) of water in the hydrate phase as follows:8,12 ln fwH = ln f w0 −

∫0

fg

⎛ xg ⎞ dfg ⎜ ⎟ ⎝ xw ⎠ fg

gg

e − βU ⟨ N!

∫ ... ∫ e−βU

Z 0(N , Nw , P , T ) =

dr1...drN ⟩(Nw , P , T )

where the RHS is now the average configurational integral for the guest molecules due to interaction only with the surrounding water lattice. The second assumption is that this average can be replaced by a product of average configurational integrals for each cavity whose values are independent of the number of guest molecules, N, inside the hydrate phase. The clathrate hydrate contains cavities of different types. Assuming an occupancy of at most one molecule per cavity and rewriting the guest−guest interaction as Ugg = ∑iNiUgg i , eq 5 can be simplified to the following equation 1 N!

Z 0(N , Nw , P , T ) =

∏ i

(2)

[e

M1

Mi

∑*

...

N1= 0

∑*

···

Ni = 0

Mi ! Ni! (Mi − Ni)!

−βUigg

gw i

∫V e−βU



dr⟩]Ni

(6)

i

where

where Mi is the number of cavities of type i, Ni is the number of guest molecules in cavity of type i, Vi is the volume of cavity of type i, and Ugg i is the interaction between a guest molecule in cavity of type i and other guest molecules. The asterisk on the summations denotes the restriction ∑iNi = N. The subscript (Nw,P,T) has been removed for clarity. All the averages denoted henceforth are computed in the isothermal−isobaric ensemble. Substitution of eq 6 into eq 3 and summing over N gives

Δ(N , Nw , P , T ) = qwNwq N Z(N , Nw , P , T )

is the isothermal−isobaric partition function for the hydrate; qw and q are the contributions to the hydrate molecular partition function from the internal degrees of freedom of water and guest molecules; λ = exp(μ/kBT) is the activity of the guest molecule; Z(N,Nw,P,T) is the configurational integral of the hydrate; Nw and N are the number of molecules of water and guest in the hydrate, respectively; kB is Boltzmann’s constant; and μ is the chemical potential of guest molecule. P and T represent the pressure and absolute temperature of the system, respectively. Equation 2 can be rewritten as follows: ∞

Γ(Nw , μ , P , T ) = Δ(Nw , P , T )

gw

(5)



N =0

dr1...drN ⟩(Nw , P , T )

where Ugg and Ugw denote the guest−guest and guest−water interactions and β = 1/kBT. The RHS is the average configurational integral for the guest molecules due to interaction among themselves and with the surrounding water lattice. To proceed further, we make two key assumptions. The first assumption is that the guest−guest interaction, Ugg, is constant and can be taken outside the average. As a result eq 4 becomes

(1)

∑ λ N Δ(N , Nw , P , T )

+ U gg ]

(4)

where f Hw is the fugacity of water in the hydrate phase, f 0w is the fugacity of water in the hydrate phase in the absence of any guest molecules, fg is the fugacity of the guest molecule, xw is the mole fraction of water, and xg is the mole fraction of guest molecule. When the fugacity of water in the hydrate phase becomes lower than that in the aqueous phase, the clathrate hydrate becomes thermodynamically stable. The extent of lowering of the fugacity of water depends on the variation of the occupancy of the guest molecules in the clathrate cavities with fugacity of the guest molecule. The lattice structure of clathrate hydrates contains cavities of different types. We define occupancy as the ratio of the number of clathrate cavities occupied by a guest molecule to the total number of cavities. An expression for the occupancy of cavities of the clathrate hydrates by guest molecules is derived as follows. The isothermal−isobaric semigrand ensemble partition function for clathrate hydrate with one type of guest molecule can be written as follows: Γ(μ , Nw , P , T ) =

gw

∑ N=0

zN

Γ(Nw , μ , P , T ) = Δ(Nw , P , T ) gg i

∏ [1 + ze βU i

gw i

∫V e−βU



i

dr⟩]Mi (7)

Rewriting z in terms of fugacity, fg, we have z = fg/kBT. The Langmuir constant, Ci, for a cavity of type i is defined as

Z(N , Nw , P , T ) Z(Nw , P , T ) (3)

Ci =

where Δ(Nw,P,T) is the isothermal−isobaric partition function and Z(Nw,P,T) is the configurational integral of the empty hydrate lattice; and z = λq. The ratio of the integrals on the right hand side (RHS) of eq 3 can be written as

1 ⟨ kBT

∫V e−u(r)/k T dr⟩ B

i

(8)

The guest−guest interaction contribution for the guest molecule in cavity of type i, Ugg i , is obtained using the meanfield approximation as follows: 18550

dx.doi.org/10.1021/jp405771u | J. Phys. Chem. C 2013, 117, 18549−18555

The Journal of Physical Chemistry C Uigg =

Article

∑ aijθj

Table 1. Force Field Parameters for Various Species Used in Molecular Simulations

j

molecule

where aij is the guest−guest interaction obtained by a guest molecule in cavity of type i from guest molecules in cavity of type j, when all cavities of type j are filled with the guest molecules and θj is the occupancy of guest molecules in cavity of type i. Substitution of these expression into eq 7 then gives

water (H2O)

methane (CH4) ethane (C2H6) carbon dioxide (CO2)

Γ(Nw , μ , P , T ) = Δ(Nw , P , T )

∏ [1 + Cifg

exp(β ∑ aijθj)]Mi

i

tetrahydrofuran ((CH2)4O)

j

The average number of guest molecules in the hydrate can be obtained from the above equation as follows: ⎞ ⎛ ∂ ln Γ ⎟ ⟨N ⟩ = ⎜⎜ ⎟ = ⎝ ∂ ln fg ⎠ P , T

∑ Mi i

Cifg exp( −β ∑j aijθj)

⟨N ⟩ = M

∑ νθi i i

(9)

where M is the total number of cavities, νi = Mi/M is the fraction of cavities of type i, and θi =

Cifg exp( −β ∑j aijθj) 1 + Cifg exp( −β ∑j aijθj)

ε/kB (K)

σ (Å)

q (e)

106.1

3.1668

148.0 98.0 28.129 80.507 190.0 56.3 56.3

3.73 3.75 2.757 3.033 2.20 3.88 3.88

0.0 0.5897 −1.1794 0.0 0.0 0.6512 −0.3256 −0.410 0.160 0.045

cutoff radius 11.5 Å. For sII hydrate a system comprising one unit cell was taken with the cutoff radius 8.4 Å. The initial structures of the water lattice in both sI and sII structures were obtained through simulated annealing simulations as described in refs 13 and 15. Each simulation run consisted of 105 equilibrium cycles and 4 × 105 production cycles. A cycle consists of an average of one particle move (Translation move, Rotation move, Insertion and Deletion with equal probability) per molecule and one volume move. Methane hydrate simulations are done at 30 bar and 273 K, ethane hydrate simulations are done at 5 bar and 273 K, CO2 hydrate simulations are done at 20 bar and 273 K, and THF hydrate simulations are done at 5 kPa and 273 K.54 These conditions have been chosen such that they lie close to the triple point lines for these clathrate hydrates. 3.2. Calculation of the Langmuir Constant. To calculate Ci, simulations of the empty hydrate lattice were performed in the isothermal−isobaric ensemble. The values of the temperature and pressure were the same as those used to compute occupancy values. These simulations were run for a total of 105 equilibrium cycles and 105 production cycles. During the simulations, 106 test insertions of the guest molecules were performed for every 25th cycle, and the average on the RHS of eq 8 was computed. This method is similar in spirit to the Widom’s test particle insertion method55 for computing excess chemical potentials. In addition to water molecules, the simulations also contained hard spheres placed in all cavities that are not of type i. The size of the hard sphere is chosen to be large enough such that it blocks the insertion of guest molecules into these cavities but not too large to distort the lattice structure of the clathrate hydrate. In our study, we have chosen 3.73 Å as contact distance for the hard sphere interaction. 3.3. Guest−Guest Interactions. The mean-field values for the guest−guest interactions were computed by performing simulations in a canonical ensemble consisting of only guest molecules. The volume of the system was kept equal to the standard volume of the empty clathrate hydrate, i.e., a unit cell length of 12 Å for sI hydrate and 17.3 Å for sII hydrate. The temperature was kept equal to that used to compute the occupancy curves. The guest molecules were placed such that their centers of mass coincided with the centers of the clathrate cavities. During simulations, only rotation moves around the center of mass are allowed. To compute aij, every cavity center of type j and one cavity center of type i contains the guest molecules. The simulations consisted of 5 × 105 equilibrium cycles and 5 × 105 production cycles.

1 + Cifg exp( −β ∑j aijθj)

The total occupancy is then obtained as θ=

group O H pair CH4 CH3 C O O (CH2)α (CH2)β

(10)

3. METHODOLOGY In this study we compute the value of θ through simulations of the clathrate hydrate in the isothermal−isobaric semigrand ensemble. We also compute the value of Ci through simulations of the empty hydrate in the isothermal−isobaric ensemble. The occupancy data obtained from simulations is then compared with the predictions from eqs 9 and 10. 3.1. Simulation Details. Methane, ethane, and THF are modeled using the united atom (TraPPE-UA)50,51 force field. According to this force field, methane is modeled as a united atom with a single site; ethane is modeled as a rigid two-site dumbell with each site representing a CH3 united atom; and THF is represented with a five-site model. Although the original model treats THF as a flexible molecule, in this study THF is kept as rigid for ease of computation. CO2 is modeled using the EPM252 force field. The water molecules are modeled using the TIP4P/Ice53 model. The interactions between the sites consist of van der Waals interactions modeled using the Lennard-Jones equation and electrostatic interactions modeled using the Coulomb’s Law. Table 1 lists the various force field parameters for all the molecules used in this study. The Lorrentz−Berthelot rules are used to obtain Lennard-Jones parameters for interactions between unlike sites. The electrostatic interactions are computed using the Ewald summation method. Pure methane, ethane, and CO2 form sI clathrate hydrate, whereas pure THF forms sII clathrate hydrate. Also at certain compositions, a mixture of methane and ethane forms sII clathrate hydrate. Hence occupancy data were obtained for methane and ethane in both sI and sII clathrate hydrate, for CO2 in sI hydrate, and for THF in sII hydrate. For sI hydrate, a system comprising 2 × 2 × 2 unit cells was taken with the 18551

dx.doi.org/10.1021/jp405771u | J. Phys. Chem. C 2013, 117, 18549−18555

The Journal of Physical Chemistry C

Article

4. RESULTS AND DISCUSSION As mentioned earlier, in this study, the method for calculating the occupancy is tested on clathrate hydrates containing methane, ethane, carbon dioxide, and THF as guest molecules. Earlier studies40−42,56−61 have shown that THF can form hydrogen bonds with the water molecules of the hydrate lattice. To check whether our molecular simulations show the formation of a hydrogen bond between THF and water, we compute the radial distribution function (RDF) for the oxygen atom of the THF molecule with the hydrogen atom of the water molecules. Figure 1 shows the RDF calculated for the

of methane and ethane forms sII hydrate where methane occupies both small and large cavities of sII hydrate and ethane only occupies large cavities of sII hydrate. Table 2 lists the values of Ci and aij for all the above combinations of guest molecules and cavity types. In Figure 2, we compare the predictions of θ using eqs 9 and 10 with those obtained from semigrand isothermal−isobaric simulations of clathrate hydrates. The comparison shows an excellent agreement between the simulations and predictions of theory. The accuracy of the predictions for a wide variety of guest molecules attests to the robustness of our method. The close match between the predictions of θ from eqs 9 and 10 and simulations also allows us to make another important observation. The expression for the Langmuir constant as given in eq 8 was derived under the assumption that the net average interaction between the guest and water molecules, that is, the cavity potential u(r), is independent of the occupancy of neighboring cavities. Strictly speaking, this is not true, as changes is the structure of the water lattice can affect the cavity potential. The two possible ways by which this happens are (i) expansion of the hydrate phase due to inclusion of guest molecules and (ii) disruption of coordination of the water lattice due to hydrogen bond formation with the guest molecules. Inclusion of any of the guest molecules always causes an expansion of volume of the hydrate phase, and inclusion of THF molecules results in the formation of a hydrogen bond with the surrounding water molecules (Figure 1). The estimates of Ci given in Table 2 are based on the structure of the empty hydrate lattice. The close agreement between the predictions of θ and values obtained from simulations indicates that the effects of changes in structure of the water lattice on the thermodynamic properties are minimal and that the assumptions used to derive eq 8 are justified for almost all kinds of molecules. This is a fortunate development, and our proposed method for computing θ in conjunction with recent developments in molecular modeling of both water44,53 and guest molecules50,51,62−66 should lead to a more thermodynamically consistent theory for clathrate hydrates.

Figure 1. Radial distribution function of hydrogen atoms surrounding the oxygen atom of the THF molecule inside the sII hydrate at 5 KPa and 273 K.

THF oxygen atom with the water molecule hydrogen atom, in which there is a peak around 1.5 Å. This peak indicates the formation of a hydrogen bond between the THF oxygen atom and the water hydrogen atom.56 In pure clathrate hydrates, that is, hydrates with only one type of guest molecules, the methane molecule occupies both small and large cavities of sI hydrate; ethane occupies only the large cavities of sI; CO2 occupies both small and large cavities of sI hydrate and THF occupies only the large cavities of sII hydrate. In addition, for certain compositions, a binary mixture

Table 2. Langmuir Constants and Guest−Guest Interaction Contributions for the Guest Molecules Used in the Study in Different Cavities of sI and sII Structures guest molecule

structure

methane (T = 273 K, P = 30 × 10 Pa) 5

methane (T = 273 K, P = 30 × 105 Pa)

ethane (T = 273 K, P = 5 × 105 Pa) ethane (T = 273K, P = 5 × 105 Pa) CO2 (T = 273 K, P = 20 × 105 Pa)

THF (T = 273 K, P = 5 × 103Pa)

sI

sII

sI sII sI

sII

cavity 12

small(5 )

Langmuir constant (Pa−1) 8.5625 × 10

−7

large(51262)

4.4855 × 10−6

small(512)

8.4854 × 10−7

large(51264)

7.1806 × 10−6

large(51262) large(51264) small(512)

2.870 × 10−5 2.657 × 10−4 1.402 × 10−7

large(51262)

1.189 × 10−5

large(51264)

0.097

18552

guest−guest interaction (1/K) aSS = −15.16 aSL = −230.92 aLS = −76.97 aLL = −171.68 aSS = −208.04 aSL = −86.25 aLS = −172.51 aLL = −52.90 aLL = −513.57 aLL = −150.90 aSS = −16.98 aSL = −284.30 aLS = −94.78 aLL = −212.78 aLL = −649.63

dx.doi.org/10.1021/jp405771u | J. Phys. Chem. C 2013, 117, 18549−18555

The Journal of Physical Chemistry C

Article

Figure 2. Comparison of predicted occupancy values for methane in sI, ethane in sI, methane in sII, ethane in sII, CO2 in sI, and THF in sII hydrate structures obtained using eq 9 with simulations. The symbols denote occupancy values obtained from simulations whereas the lines denote occupancy values obtained from eq 9.

5. CONCLUSIONS In this study, we have successfully extended the van der Waals and Platteeuw theory for clathrate hydrates to predict the variation of occupancy of clathrate cavities with fugacities for a wide variety of guest molecules. The only exceptions are flexible molecules with internal degrees of freedom, such as propane, nbutane, etc. Recently,13−15 it has been shown that the assumption of the rigid water lattice leads to significant errors in the predictions of the hydrate phase diagram. The proposed extensions allow one to account for the flexibility of the water lattice in the hydrate phase and also guest−guest interactions in the prediction of the occupancy. In addition, the methods for obtaining estimates of the Langmuir constants are based on the structure of the water lattice in an empty hydrate. This indicates that the effect of changes in the structure of the water lattice upon inclusion of guest molecules is minimal as far as the thermodynamic properties of the clathrate hydrates are concerned. This is a fortunate development which in combination with advanced molecular modeling can lead to

the development of thermodynamically consistent theories for clathrate hydrates.



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.



REFERENCES

(1) Hammerschmidt, E. G. Formation of Gas Hydrates in Natural Gas Transmission Lines. Ind. Eng. Chem. 1934, 26, 851−855. (2) Shi, G.; Jing, Y.; Zhang, X.; Zheng, G.4th IEEE Conference on Prospects of natural gas storage and transportation using hydrate technology in China. Industrial Electronics and Applications, 2009. ICIEA 2009; 2009; pp 530−534. (3) Gudmundsson, J.; Borrehaug, A. Frozen Hydrate for Transport of Natural Gas. Proc. 2nd Int. Conf. Natural Gas Hydrates 1996, 415−422.

18553

dx.doi.org/10.1021/jp405771u | J. Phys. Chem. C 2013, 117, 18549−18555

The Journal of Physical Chemistry C

Article

(4) Kennett, J.; Cannariato, K.; Hendy, I.; Behl, R. Methane Hydrates in Quaternary Climate Change: The Clathrate Gun Hypothesis; AGU: Washington, DC, 2003. (5) Barduhn, A. J.; Towlson, H. E.; Hu, Y. C. The properties of some new gas hydrates and their use in demineralizing sea water. AIChE J. 1962, 8, 176−183. (6) Barduhn, A. J.; Roux, G. M.; Richard, H. A.; Giuliano, G. B.; Stern, S. A. The rate of formation of the hydrates of F-31(CH2ClF) and F-142b(CH3CClF2) in a stirred tank reactor. Desalination 1976, 18, 59−69. (7) Barduhn, A. J.; Wen-Ching, L. Thermodynamic properties of F22(CHClF2) hydrate for use in desalting. Desalination 1978, 25, 151− 162. (8) Sloan, E. D.; Koh, C. A. Clathrate Hydrates of Natural Gases, 3rd ed.; CRC Press: Boca Raton, FL, 2007. (9) Lee, S.-Y.; Holder, G. D. Methane Hydrates Potential as a Future Energy Source. Fuel Process. Technol. 2001, 71, 181−186. (10) Demirbas, A. Methane hydrates as potential energy resource: Part 1 -Importance, resource and recovery facilities. Energy Convers. Manage. 2010, 51, 1547−1561. (11) Demirbas, A. Methane Hydrates as Potential Energy Resource: Part 2Methane Production Processes from Gas Hydrates. Energy Convers. Manage. 2010, 51, 1562−1571. (12) van der Waals, J. H.; Platteeuw, J. C. Clathrate Solutions. Adv. Chem. Phys. 1959, 2, 1−57. (13) Ravipati, S.; Punnathanam, S. N. Analysis of Parameter Values in the van der Waals and Platteeuw Theory for Methane Hydrates Using Monte Carlo Molecular Simulations. Ind. Eng. Chem. Res. 2012, 51, 9419−9426. (14) Wierzchowski, S. J.; Monson, P. A. Calculation of Free Energies and Chemical Potentials for Gas Hydrates Using Monte Carlo Simulations. J. Phys. Chem. B 2007, 111, 7274−7282. (15) Pimpalgaonkar, H.; Veesam, S. K.; Punnathanam, S. N. Theory of Gas Hydrates: Effect of the Approximation of Rigid Water Lattice. J. Phys. Chem. B 2011, 115, 10018−10026. (16) Chialvo, A. A.; Houssa, M.; Cummings, P. T. Molecular Dynamics Study of the Structure and Thermophysical Properties of Model sI Clathrate Hydrates. J. Phys. Chem. B 2002, 106, 442−451. (17) Cao, Z.; Tester, J. W.; Sparks, K. A.; Trout, B. L. Molecular Computations Using Robust Hydrocarbon Water Potentials for Predicting Gas Hydrate Phase Equilibria. J. Phys. Chem. B 2001, 105, 10950−10960. (18) Rodger, P. M. Stability of Gas Hydrates. J. Phys. Chem. 1990, 94, 6080−6089. (19) Kvamme, B.; Lund, A.; Hertzberg, T. The Influence of Gas−Gas Interactions on the Langmuir Constants for Some Natural Gas Hydrates. Fluid Phase Equilib. 1993, 90, 15−44. (20) John, V. T.; Holder, G. D. Choice of Cell Size in the Cell Theory of Hydrate Phase Gas−Water Interactions. J. Phys. Chem. 1981, 85, 1811−1814. (21) John, V. T.; Holder, G. D. Contribution of the Second and Subsequent Water Shells to the Potential Energy of Guest−Host Interactions in Clathrate Hydrates. J. Phys. Chem. 1982, 86, 455−459. (22) John, V. T.; Holder, G. D. Langmuir Constants for Spherical and Linear Molecules in Clathrate Hydrates. Validity of the Cell Theory. J. Phys. Chem. 1985, 89, 3279−3285. (23) Sparks, K. A.; Tester, J. W. Intermolecular Potential Energy of Water Clathrates: The Inadequacy of the Nearest-Neighbor Approximation. J. Phys. Chem. 1992, 96, 11022−11029. (24) Sparks, K. A.; Tester, J. W.; Cao, Z.; Trout, B. L. Configurational Properties of Water Clathrates: Monte Carlo and Multidimensional Integration versus the Lennard-Jones and Devonshire Approximation. J. Phys. Chem. B 1999, 103, 6300−6308. (25) Klauda, J. B.; Sandler, S. I. Phase Behavior of Clathrate Hydrates: A Model for Single and Multiple Gas Component Hydrates. Chem. Eng. Sci. 2003, 58, 27−41. (26) Martin, A. A Simplified van der Waals-Platteeuw Model of Clathrate Hydrates with Multiple Occupancy of Cavities. J. Phys. Chem. B 2010, 114, 9602−9607.

(27) Rodger, P. M. Lattice Relaxation in Type I Gas Hydrates. AIChE J. 1991, 37, 1511−1516. (28) Sizov, V. V.; Piotrovskaya, E. M. Computer Simulation of Methane Hydrate Cage Occupancy. J. Phys. Chem. B 2007, 111, 2886− 2890. (29) Zhang, J.; Hawtin, R. W.; Yang, Y.; Nakagava, E.; Rivero, M.; Choi, S. K.; Rodger, P. M. Molecular Dynamics Study of Methane Hydrate Formation at a Water/Methane Interface. J. Phys. Chem. B 2008, 112, 10608−10618 (PMID: 18671369). (30) Alavi, S.; Ripmeester, J. A. Nonequilibrium Adiabatic Molecular Dynamics Simulations of Methane Clathrate Hydrate Decomposition. J. Chem. Phys. 2010, 132, 144703. (31) Bai, D.; Chen, G.; Zhang, X.; Wang, W. Microsecond Molecular Dynamics Simulations of the Kinetic Pathways of Gas Hydrate Formation from Solid Surfaces. Langmuir 2011, 27, 5961−5967. (32) Sarupria, S.; Debenedetti, P. G. Molecular Dynamics Study of Carbon Dioxide Hydrate Dissociation. J. Phys. Chem. A 2011, 115, 6102−6111. (33) Tse, J. S.; Klein, M. L.; McDonald, I. R. Molecular Dynamics Studies of Ice Ic and the Structure I Clathrate Hydrate of Methane. J. Phys. Chem. 1983, 87, 4198−4203. (34) Tse, J. S.; Klein, M. L.; McDonald, I. R. Computer Simulation Studies of the Structure I Clathrate Hydrates of Methane, Tetrafluoromethane, Cyclopropane, and Ethylene Oxide. J. Chem. Phys. 1984, 81, 6146−6153. (35) Moon, C.; Taylor, P. C.; Rodger, P. M. Molecular Dynamics Study of Gas Hydrate Formation. J. Am. Chem. Soc. 2003, 125, 4706− 4707. (36) Moon, C.; Taylor, P. C.; Rodger, P. M. Clathrate Nucleation and Inhibition from a Molecular Perspective. Can. J. Phys. 2003, 81, 451−457. (37) Moon, C.; Hawtin, R. W.; Rodger, P. M. Nucleation and Control of Clathrate Hydrates: Insights from Simulation. Faraday Discuss. 2007, 136, 367−382. (38) Walsh, M. R.; Koh, C. A.; Sloan, E. D.; Sum, A. K.; Wu, D. T. Microsecond Simulations of Spontaneous Methane Hydrate Nucleation and Growth. Science 2009, 326, 1095−1098. (39) Wierzchowski, S. J.; Monson, P. A. Calculating the Phase Behavior of Gas-Hydrate-Forming Systems from Molecular Models. Ind. Eng. Chem. Res. 2006, 45, 424−431. (40) Alavi, S.; Udachin, K.; Ripmeester, J. Effect of Guest-Host Hydrogen Bonding on the Structures and Properties of Clathrate Hydrates. Chem.Eur. J. 2010, 16, 1017−1025. (41) Alavi, S.; Susilo, R.; Ripmeester, J. A. Linking Microscopic Guest Properties to Macroscopic Observables in Clathrate Hydrates: Guest− Host Hydrogen Bonding. J. Chem. Phys. 2009, 130, 174501. (42) Susilo, R.; Alavi, S.; Moudrakovski, I. L.; Englezos, P.; Ripmeester, J. A. Guest−Host Hydrogen Bonding in Structure H Clathrate Hydrates. ChemPhysChem 2009, 10, 824−829. (43) Jensen, L.; Thomsen, K.; van Solms, N.; Wierzchowski, S.; Walsh, M. R.; Koh, C. A.; Sloan, E. D.; Wu, D. T.; Sum, A. K. Calculation of Liquid Water-Hydrate-Methane Vapor Phase Equilibria from Molecular Simulations. J. Phys. Chem. B 2010, 114, 5775−5782. (44) Conde, M. M.; Vega, C. Determining the Three-Phase Coexistence Line in Methane Hydrates Using Computer Simulations. J. Chem. Phys. 2010, 133, 064507. (45) Conde, M. M.; Vega, C. Note: A Simple Correlation to Locate the Three Phase Coexistence Line in Methane-Hydrate Simulations. J. Chem. Phys. 2013, 138, 056101. (46) Ballard, L.; Sloan, E. D., Jr. The Next Generation of Hydrate Prediction: I. Hydrate Standard States and Incorporation of Spectroscopy. Fluid Phase Equilib. 2002, 194−197, 371−383. (47) Jager, M.; Ballard, L.; Sloan, E. D., Jr. The Next Generation of Hydrate Prediction: II. Dedicated Aqueous Phase Fugacity Model for Hydrate Prediction. Fluid Phase Equilib. 2003, 211, 85−107. (48) Ballard, L.; Sloan, E. D., Jr. The Next Generation of Hydrate Prediction: Part III. Gibbs Energy Minimization Formalism. Fluid Phase Equilib. 2004, 218, 15−31. 18554

dx.doi.org/10.1021/jp405771u | J. Phys. Chem. C 2013, 117, 18549−18555

The Journal of Physical Chemistry C

Article

(49) Ballard, L.; Jr, E. S. The next generation of hydrate prediction IV: A comparison of available hydrate prediction programs. Fluid Phase Equilib. 2004, 216, 257−270. (50) Stubbs, J. M.; Potoff, J. J.; Siepmann, J. I. Transferable Potentials for Phase Equilibria. 6. United-Atom Description for Ethers, Glycols, Ketones, and Aldehydes. J. Phys. Chem. B 2004, 108, 17596−17605. (51) Keasler, S. J.; Charan, S. M.; Wick, C. D.; Economou, I. G.; Siepmann, J. I. Transferable Potentials for Phase Equilibria−United Atom Description of Five- and Six-Membered Cyclic Alkanes and Ethers. J. Phys. Chem. B 2012, 116, 11234−11246. (52) Harris, J. G.; Yung, K. H. Carbon Dioxide’s Liquid-Vapor Coexistence Curve And Critical Properties as Predicted by a Simple Molecular Model. J. Phys. Chem. 1995, 99, 12021−12024. (53) Abascal, J. L. F.; Sanz, E.; Fernandez, R. G.; Vega, C. A Potential Model for the Study of Ices and Amorphous Water: TIP4P/Ice. J. Chem. Phys. 2005, 122, 234511. (54) Makino, T.; Sugahara, T.; Ohgaki, K. Stability Boundaries of Tetrahydrofuran + Water System. J. Chem. Eng. Data 2005, 50, 2058− 2060. (55) Widom, B. Some Topics in the Theory of Fluids. J. Chem. Phys. 1963, 39, 2808−2812. (56) Alavi, S.; Ripmeester, J. A. Effect of Small Cage Guests on Hydrogen Bonding of Tetrahydrofuran in Binary Structure II Clathrate Hydrates. J. Chem. Phys. 2012, 137, 054712. (57) Udachin, K.; Alavi, S.; Ripmeester, J. A. Communication: Single Crystal X-ray Diffraction Observation of Hydrogen Bonding between 1-Propanol and Water in a Structure II Clathrate Hydrate. J. Chem. Phys. 2011, 134, 121104. (58) Buch, V.; Devlin, J. P.; Monreal, I. A.; Jagoda-Cwiklik, B.; UrasAytemiz, N.; Cwiklik, L. Clathrate Hydrates with Hydrogen-Bonding Guests. Phys. Chem. Chem. Phys. 2009, 11, 10245−10265. (59) Monreal, I. A.; Cwiklik, L.; Jagoda-Cwiklik, B.; Devlin, J. P. Classical to Nonclassical Transition of Ether−HCN Clathrate Hydrates at Low Temperature. J. Phys. Chem. Lett. 2010, 1, 290−294. (60) Hernandez, J.; Uras, N.; Devlin, J. P. Coated Ice Nanocrystals from Water−Adsorbate Vapor Mixtures: Formation of Ether−CO2 Clathrate Hydrate Nanocrystals at 120 K. J. Phys. Chem. B 1998, 102, 4526−4535. (61) Alavi, S.; Ohmura, R.; Ripmeester, J. A. A Molecular Dynamics Study of Ethanol-Water Hydrogen Bonding in Binary Structure I Clathrate Hydrate with CO2. J. Chem. Phys. 2011, 134, 054702. (62) Zhu, X.; Lopes, P. E. M.; MacKerell, A. D. Recent Developments and Applications of the CHARMM Force Fields. Wiley Interdiscip. Rev.: Comput. Mol. Sci. 2012, 2, 167−185. (63) Cornell, W. D.; Cieplak, P.; Bayly, C. I.; Gould, I. R.; Merz, K. M.; Ferguson, D. M.; Spellmeyer, D. C.; Fox, T.; Caldwell, J. W.; Kollman, P. A. A Second Generation Force Field for the Simulation of Proteins, Nucleic Acids, and Organic Molecules. J. Am. Chem. Soc. 1995, 117, 5179−5197. (64) Jorgensen, W. L.; Maxwell, D. S.; Tirado-Rives, J. Development and Testing of the OPLS All-Atom Force Field on Conformational Energetics and Properties of Organic Liquids. J. Am. Chem. Soc. 1996, 118, 11225−11236. (65) Kaminski, G. A.; Friesner, R. A.; Tirado-Rives, J.; Jorgensen, W. L. Evaluation and Reparametrization of the OPLS-AA Force Field for Proteins via Comparison with Accurate Quantum Chemical Calculations on Peptides. J. Phys. Chem. B 2001, 105, 6474−6487. (66) Mayo, S. L.; Olafson, B. D.; Goddard, W. A. DREIDING: A Generic Force Field for Molecular Simulations. J. Phys. Chem. 1990, 94, 8897−8909.

18555

dx.doi.org/10.1021/jp405771u | J. Phys. Chem. C 2013, 117, 18549−18555