THE JOURNAL OF CHEMICAL PHYSICS 122, 074503 共2005兲

On the thermodynamic stability and structural transition of clathrate
hydrates
Yuji Koyama, Hideki Tanaka,a兲 and Kenichiro Koga
Department of Chemistry, Faculty of Science, Okayama University, 3-1-1 Tsushima, Okayama 700-8530,
Japan

共Received 15 October 2004; accepted 30 November 2004; published online 8 February 2005兲
Gas mixtures of methane and ethane form structure II clathrate hydrates despite the fact that eachof pure methane and pure ethane gases forms the structure I hydrate. Optimization of the interaction
potential parameters for methane and ethane is attempted so as to reproduce the dissociationpressures of each simple hydrate containing either methane or ethane alone. An account for the
structural transitions between type I and type II hydrates upon changing the mole fraction of the gas
mixture is given on the basis of the van der Waals and Platteeuw theory with these optimized
potentials. Cage occupancies of the two kinds of hydrates are also calculated as functions of themole fraction at the dissociation pressure and at a fixed pressure well above the dissociationpressure. © 2005 American Institute of Physics. 关DOI: 10.1063/1.1850904兴
I. INTRODUCTION

A gas hydrate consists of guest molecules and host watermolecules forming a hydrogen bonded network whose constituents are planar pentagonal and hexagonal rings only.
Two common forms of clathrate hydrates are called structure
I and structure II.1,2 Another crystalline form called structure
H, which was found more recently, falls into a category of
mixed hydrates, which contain more than one guest species.3
The clathrate hydrate is thermodynamically stable only when
guest molecules are encapsulated in the host cages. The
guest-host interactions, though much weaker than the hydrogen bonds of the host lattice, play a crucial role for clathrate
hydrates to be stabilized. The thermodynamic stability of theclathrate hydrates has long been explained by the van der
Waals and Platteeuw 共vdWP兲 theory.4 This has been the only
rigorous statistical mechanical theory to account for the stability of clathrate hydrates and has been widely used to predict phase behavior of many hydrates. However, the theory
requires some empirical parameters for all practical
purposes.1,2 When applied directly, with a simple free energy
evaluation method, to estimation of clathrate hydrate stability for large guest species, it leads to an incorrect phase
diagram due to its built-in assumptions that are invalid forhost lattices occupied with large guests.
We proposed a method to calculate the free energy ofclathrate hydrates encaging molecules comparable in size tothe vacant space of the larger cages5 and of hydrates inwhich both larger and smaller cages are occupied by a single
guest species.6 The fixed lattice approximation, the most serious shortcoming in the original vdWP theory, was removed
by introducing the modulation of the phonon frequencies of
the host molecules. Thus, more accurate evaluation of the
dissociation pressure was made possible, which requires onlyintermolecular interactions but no other empirical parameters. It was revealed that the modulation of the frequenciesa兲Electronic mail: htanakaa@cc.okayama-u.ac.jp

0021-9606/2005/122共7兲/074503/7/$22.50

lifts the free energy and thereby destabilizes the clathrate
hydrates encaging large guest molecules such as propane.5,6The free energy calculation was further sophisticated in order
to treat clathrate hydrates encaging nonspherical guest suchas propane and ethane.7 The free energy associated with therotational degrees of freedom of the guest molecules wasevaluated as a perturbation from the spherical guest. However, evaluation of the anharmonic free energy by the thermodynamic integration demands a large amount of CPU
time. It is therefore desirable for the free energy calculation,
whenever possible, to employ an effective spherical potential
for a nonspherical guest that incorporates all the contributions from the anisotropy of the guest.
It is well known that both pure methane and pure ethane
form the clathrate hydrate structure I, i.e., the each simple
hydrate, which contains only one guest species, is type I.
However it was found that temperature dependence of the
dissociation pressure of a mixed hydrate formed by methane
mixed with a small amount of ethane did not fit well to theprediction assuming structure I,8 which naturally suggested atransition from structure I to II.9 Convincing evidences havebeen reported recently that the gas mixture of methane andethane stabilizes structure II rather than structure I.10,11 The
transition has been accounted for by minimizing the relevant
Gibbs free energy of model systems with parameters optimized to reproduce thermodynamic properties of clathrate
hydrates.12,13 Here, we show that the structural transition canbe accounted for on the basis of the vdWP theory with thefree energy of cage occupancy evaluated from intermolecularinteractions only. Our method treats a clathrate hydrate of agas mixture as an open system with respect to the guest
species and therefore is a more advantageous way to examine its thermodynamic stability and cage occupancy.
We first obtain sets of the potential parameters for methane and those for ethane, each optimized for reproducing theobserved dissociation pressure of each simple hydrate. It isessential that these optimized potential parameters, when

122, 074503-1© 2005 American Institute of Physics

Downloaded 01 Nov 2005 to 150.46.100.222. Redistribution subject to AIP license or copyright, see http://jcp.aip.org/jcp/copyright.jsp

074503-2

J. Chem. Phys. 122, 074503 共2005兲

Koyama, Tanaka, and Koga

employed for the free energy calculation based on the vdWP
theory, always give rise to structure I being more stable thanstructure II for each simple hydrate. Second, we examine
whether or not structure II becomes more stable than structure I in some range of the mole fraction of ethane when they
are in equilibrium with a methane+ ethane gas mixture.
Third, the cage occupancies are calculated with the sets ofthe potential parameters at dissociation pressure and at
higher constant pressure.
II. MODEL AND METHOD
A. Intermolecular interaction

The water-water intermolecular interaction is describedby the TIP4P potential.14 This is one of the most frequently
used potentials for water because of its ability to reproduce
many important properties of pure water within the framework of pair potential. The TIP4P model consists of fourinteraction sites: a positive charge qH on each of two hydrogen atoms, a negative charge 共−2qH兲 on the bisector of two
OH bonds, and a soft core interaction between oxygen atoms. The methane-methane and the ethane-ethane interactions are both described by the Lennard–Jones 共LJ兲 potential,
which is spherical. For each species, we determine its LJ
energy and size parameters that give the observed dissociation pressure of the simple hydrate of that species at 263.15
K. The set of such two parameters is not unique, i.e., it is nota single point but forms a curve in the parameter space. The
range of the LJ size parameter thus determined is from 0.37
to 0.39 nm for methane and from 0.44 to 0.46 nm for ethane.
Adoption of the LJ potential for methane and ethane meansthat rotational motions of the two guest species are not considered. In general, this is a reasonable approximation formethane, for its anisotropy is very small. For ethane, however, the LJ potential with an optimized set of parameters
should be regarded as an effective potential for the simple
and mixed hydrates of ethane. Some remarks will be givenon this effective potential. The water-guest intermolecularinteraction is taken to be the LJ potential between the oxygensite of the TIP4P water and the center of each guest molecule, with the parameter set given by the Lorentz–Berthelot
mixing rule from that of the oxygen site 共⑀OO
= 0.6487 kJ mol−1 and ␴OO = 0.3154 nm兲 and that of theguest. The interaction potentials for all pairs of molecules aretruncated smoothly at 0.8655 nm as was done in the previous
calculation.5

and the structure are given, Nl and Ns are determined automatically for each structure.兲 The corresponding grand partition function ⌶ is written as冋 兺⌶ = exp共− ␤Aw0 兲 1 +

冋 兺

⫻ 1+

exp兵␤共␮ j − f sj 兲其

jexp兵␤共␮ j − f lj兲其

j

册

Nl

册

Ns

共1兲

,

where ␮ j stands for the chemical potential of guest species
j , f j is its free energy of cage occupancy with respect toeither a larger cage 共l兲 or a smaller cage 共s兲 at temperature
T , Aw0 is the free energy of the empty hydrate, and ␤ is 1 / kBTwith Boltzmann’s constant kB. The free energy of cage occupancy has mostly been treated as input parameters. One of us
showed that the free energy of cage occupancy f j can becalculated for hydrates encaging smaller guests such as argon
and methane with an assumption that all the host vibrationalmodes are independent of the presence of guest molecules.
This assumption leads to the following expression of f j:
exp共− ␤ f j兲 = 共2␲m jkBT/h2兲3/2

冕exp关− ␤w共r兲兴dr, 共2兲

vcage

where h is Planck’s constant, m j is molecular mass of guest
j , w共r兲 is the interaction potential between the guest molecule at r and all the water molecules with fixed configuration, and the integration spans the entire volume of the single
cage vcage.
Cage occupancy y i for guest i is given by
yi ⬅
=

具Ni典
= 共Nl + Ns兲−1 ⳵ ln⌶/⳵ 共␤␮i兲
共Nl + Ns兲
Ns
exp关␤共␮i − f si 兲兴
共Nl + Ns兲

冋 兺

⫻ 1+

j

exp兵␤共␮ j − f sj 兲其冋 兺⫻exp关␤共␮i − f li兲兴 1 +

册

−1

+

Nl
共Nl + Ns兲exp兵␤共␮ j − f lj兲其

j

册

−1

,

共3兲

where 具Ni典 is the mean number of cages occupied by guest
species i. The chemical potential ␮c of water in the hydrate isgiven by冠 再␮c = − kBT ⳵ ln⌶/⳵ Nw = ␮0c − kBT ␣sln 1
B. The vdWP theory

The vdWP theory is applicable to any sort of hydrate,either type I or II and either simple or mixed. It is assumed in
the vdWP theory that first the cage structure is not distortedby the incorporation of guest molecules, second a guest molecule inside a cage moves in the force field of water molecules fixed at lattice sites and there is no coupling between
host and guest molecular motions. Here, we describe only anessential part of it for convenience of the later discussion.
Let us consider a clathrate hydrate which contains Nw watermolecules, Nl larger cages, and Ns smaller cages. 共Once Nw+冎 再兺j exp关␤共␮ j − f sj 兲兴 + ␣lln 1 + 兺j⫻exp关␤共␮ j − f lj兲兴

冎冡

,

共4兲

where ␮0c is the chemical potential of water in the hypothetical empty hydrate and ␣l and ␣s are the ratios of the numbersof the larger and smaller cages to the number of water molecules, respectively. Equations 共1兲–共4兲 describe the essential
part of the vdWP theory.

Downloaded 01 Nov 2005 to 150.46.100.222. Redistribution subject to AIP license or copyright, see http://jcp.aip.org/jcp/copyright.jsp

074503-3

J. Chem. Phys. 122, 074503 共2005兲

Stability of clathrate hydrates

TABLE I. Chemical potential ␮c0 of empty hydrate and chemical potential
difference ⌬␮ between ice and empty hydrate for structures I and II at
263.15 K, where the anharmonic free energy and the residual entropy terms
are omitted.

TABLE II. Potential energy u at structures of minimum potential energy andfree energy f of cage occupancy at 273.15 K 共Ref. 6兲. Those free energy andthe potential energy are calculated for larger and smaller cage occupancies,
denoted by superscripts l and s, respectively.

kJ mol−1

Structure I

Structure II

kJ mol−1

Methane in hydrate I

Methane in hydrate II␮c0⌬␮

−47.93−0.834−48.03−0.732

us
ul
fs
fl

−22.13−20.26−27.83−30.30−27.72−16.67−28.13−31.19Argon in hydrate I

Argon in hydrate II−17.01−14.25−28.73−29.77−17.61−11.43−28.97−30.48

C. Free energies of ice and empty clathrate hydrates

The equilibrium condition between ice and clathrate hydrate is given by

␮i = ␮c ,共5兲

where ␮i and ␮c denote the chemical potentials of ice andclathrate hydrate. By subtracting the chemical potential ␮0c of
the corresponding hypothetical empty clathrate hydrate fromboth sides of Eq. 共5兲,

␮i − ␮0c = ␮c − ␮0c .

共6兲

The left-hand side of the above equation is calculated assuming that the free energy is expressed by the sum of the potential energy of a given system at mechanically equilibrium
positions of molecules, the harmonic vibrational free energy,
the anharmonic vibrational free energy, and the residual entropy term, the last of which is common to all the protondisordered ice forms. It is advantageous to calculate thechemical potentials relative to that of the empty hydrate 关Eq.
共6兲兴 instead of the absolute chemical potentials 关Eq. 共5兲兴 because, on both sides of Eq. 共6兲, the residual entropy terms
cancel out exactly and the difference in the anharmonic vibrational free energy terms is expected to be negligible, although each contribution may not be so small. We generate
100 proton-disordered structures for each of clathrate hydrate
I and II and for ice Ih. The chemical potentials of clathrate
hydrate I and II and that of ice Ih are calculated according to
the above mentioned way and the differences ⌬␮’s between
ice and empty hydrate are tabulated in Table I. The chemical
potential of empty hydrate II is always lower than that ofempty hydrate I, which is expected from the preferential formation of structure II for argon and other smaller molecules.
III. INFERENCES FROM EXPERIMENTAL
OBSERVATIONS

It would be worthwhile here to summarize the experimental facts that are relevant to the structural transition of
the methane+ ethane system and what we can learn fromthem.
共1兲 The structure II is more stable than the structure I
when the guests are very small molecules such as argon,
krypton, and nitrogen.15
共2兲 A methane molecule is encaged in both larger andsmaller cages of the structure I.2
共3兲 An ethane molecule is encapsulated in only larger
cages of the structure I except for a small amount in smaller
cages near the dissociation pressure.16
The individual observations lead naturally to the following inferences.

us
ul
fs
fl

共1兲 The free energy of cage occupancy for smaller guest
molecules is rather insensitive to the type of the cage. In theprevious work, we calculated the free energy for argon using
the most reliable intermolecular interaction, which is listed inTable II.6 Even if we choose a different intermolecular interaction, we obtain a similar tendency. This fact together withthe preferential formation of the structure II hydrate necessarily indicates that the chemical potential of the empty
structure II hydrate is lower than that of the empty structure
I hydrate at an ambient pressure and temperature. Although
there are some exceptions,2 this view is supported by experiments tabulated in Ref. 2.
共2兲 The free energy difference, f s − f l becomes larger
with increasing guest size. A methane molecule is smallenough to be accommodated in any type of cage but large
enough to have fairly different free energy values for different types of cages as shown in Table II for each structure.
The free energy of cage occupancy in the larger cage is lower
than that in the smaller cage for both structure I and structure
II. However, ␣l 共ratio of the number of larger cages to thenumber of water molecules兲 for structure I 共6 / 46兲 is muchlarger than ␣l for structure II 共8 / 136兲. Thus, structure I is
more stabilized by accommodating a large number of thelarger cages. The free energies of cage occupancy for structure II are lower than the corresponding those for structure I.
The preferential formation of structure I may occur under thesubtle balance between the larger ␣l and the less stable
empty hydrate of structure I as realized in methane hydrate.
共3兲 Only larger cages can accommodate very large guest
molecules such as ethane. Such guests prefer structure I to
structure II since ␣l for structure I is much larger than thatfor structure II.
The above inferences conduct to the following speculation on the mixed hydrate. If a small amount of ethane ismixed with a methane gas, ethane molecules occupy thelarger cages exclusively while methane molecules can occupy both the smaller and larger cages. That is, both types ofcages can be potentially occupied unlike in the case of thesimple hydrate of ethane. The free energy of such hydrates
would be lower than that of the simple hydrate of ethane.
Under the methane-rich condition, structure II is more favorable because the chemical potential of empty structure II

Downloaded 01 Nov 2005 to 150.46.100.222. Redistribution subject to AIP license or copyright, see http://jcp.aip.org/jcp/copyright.jsp

074503-4

J. Chem. Phys. 122, 074503 共2005兲

Koyama, Tanaka, and Koga

TABLE III. Free energy f of cage occupancy at 263.15 K calculated by
optimized LJ parameters ␴ and ⑀ which can reproduce dissociation pressure
pd 共␴ , ⑀, and pd are also listed here兲. Superscripts s and l stand for smaller
and larger cages, and subscripts, I and II, stand for hydrate I and hydrate II.
Linear ethane molecule means a rigid rotor composed two LJ interaction
site.

pd 共MPa兲␴ 共nm兲
⑀ 共kJ mol−1兲
f Is 共kJ mol−1兲
f Il 共kJ mol−1兲
f IIs 共kJ mol−1兲
f IIl 共kJ mol−1兲

MethaneEthaneEthane 共linear兲

1.84
0.38501.277−27.38−30.47−27.43−31.340.306
0.45201.736…
−37.48…
−41.76…
0.37751.720…
−37.62…
…

hydrate is lower than that of empty structure I hydrate. Thisis also supported by our direct calculation of the chemical
potential of the empty hydrates relative to that of ice asshown in Table I. Under ethane-rich condition, on the contrary, structure I becomes more favorable. This is because theproportion of larger cages, which are available for ethane
molecules, is much larger in structure I than in structure II.

IV. RESULTS AND DISCUSSION

We determine the LJ parameters for each guest species
that reproduce the dissociation pressure of simple clathrate
hydrate of that species. Ref. 2 accumulates a huge number ofdata on dissociation pressure. A simple interpolation provides a way to estimate the dissociation pressure at a given
temperature. Conversely, once the dissociation pressure pd isgiven as listed in Table III, the LJ parameters of a guest isdetermined. Since there are two parameters ⑀ and ␴ to be
fitted and a single datum pd to fit, the optimized LJ parameters are not determined uniquely but their relation is determined, e.g., ⑀ is determined as a function of ␴. For bothethane and methane, the LJ parameters are fitted to the dissociation pressure of structure I hydrate at T = 263.15 K. The
free energy of cage occupancy in structure II is calculatedusing the parameter set thus obtained. It is revealed that, for

FIG. 1. Appropriate energy parameter for methane against size parameter at
263.15 K.

FIG. 2. Appropriate energy parameter for ethane against size parameter at
263.15 K.

both methane and ethane, structure I simple hydrate is always more stable than structure II within the parameter space
we have examined.
The energy parameter thus obtained is plotted against thesize parameter in Fig. 1 for methane and Fig. 2 for ethane.
Note that the region above the curve of the optimized parameter set in Fig. 1 and that in Fig. 2 give rise to the free energy
lower than the value to fit. The energy parameter decreases
with increasing the size parameter for methane whereas thereverse holds for ethane. The figure contrasts the role ofvacant volume in the larger cage for the two guest species.
To account for the difference, the potential energy w共r兲 of aguest molecule in a larger cage 共interacting with surrounding
water molecules fixed to the lattice sites兲 is calculated as afunction of radial distance r from the center of the cage,
which is the average over all the directions. Figure 3 shows
the potential energy curve 共solid line兲 for each guest species
obtained from an optimized parameter set 共␴ = 0.38 nm and⑀ = 1.3 kJ mol−1 for methane and ␴ = 0.45 nm and ⑀

FIG. 3. Potential energy curves of methane and ethane molecules for the
optimized size parameter 共␴ = 0.38 nm for methane, ␴ = 0.45 nm for ethane兲
as a function of radial distance r from the center of the larger cage 共solidline兲. The dotted and dot-dashed lines are corresponding to those for larger
and smaller size parameters 共␴ ± 0.01 nm兲, respectively. The three potential
energy curves are calculated with the same energy parameter 共⑀= 1.3 kJ mol−1 for methane, ⑀ = 1.7 kJ mol−1 for ethane兲.

Downloaded 01 Nov 2005 to 150.46.100.222. Redistribution subject to AIP license or copyright, see http://jcp.aip.org/jcp/copyright.jsp

074503-5

Stability of clathrate hydrates= 1.7 kJ mol−1 for ethane兲. The others show the potential energy curves for the optimized ␴ ± 0.01 nm with the common
energy parameter. Both the potential energy curves aremonotonically increasing functions of r. The potential curve
for methane is fairly flat around the center, indicating that amethane molecule can move around a wide space. This relatively large region of lower energy seems to contribute tolowering the free energy of cage occupancy. The curves obtained from three size parameters 共0.37–0.39 nm兲 show thatthe potential energy becomes lower as the size parameterincreases in the region from the center to r = 0.09 nm. Because this region is wide for methane, the free energy islowered with increasing the size parameter if the energy parameter is fixed. This explains why the curve of the optimized parameter set has a negative slope as shown in Fig. 1.
On the other hand, the potential energy of ethane exhibits aquadratic character that confines an ethane molecule in a
relatively small space. It is shown by the curves obtained
from three size parameters 共0.44, 0.45, 0.46 nm兲 that thepotential energy near the center of the cage is lowered as thesize parameter increases as found for methane. However, thistendency is limited to a very small region from the center tor = 0.04 nm, where the energy is still low, and the reverse
holds for the outer region. The difference in the energy isintensified with increasing the distance from the crossover
point. This character of the potential energy curve for ethanegives rise to increase in free energy with increasing the size
parameter if the energy parameter is fixed, and therefore explains why the slope of the optimized parameter curve ispositive as shown in Fig. 2. This completes the explanationof the difference between Figs. 1 and 2.
Table III lists a set of the free energies of cage occupancy calculated from the optimized LJ parameters for structure I and II at 263.15 K 共these are optimized not only for thesimple hydrates but also for the mixed hydrates as shown
below兲. The size parameter for ethane is fairly large compared to that given in a classical textbook.17 This is justified
by considering oblate shape of the larger cage in clathrate
structure I and anisotropic nature of an ethane molecule. If a
two-site model of ethane is adopted,18 a reasonable value forthe free energy as that in Table III would be obtained without
invoking the fitting. Figure 4 shows to what extent structure
I hydrate encaging methane is stabilized compared to structure II by plotting both sides of Eq. 共6兲 against gas pressure.
The left-hand side of Eq. 共6兲 is constant with a good approximation while the right-hand side is heavily dependent on thegas pressure. The intersection gives a dissociation pressure.
Although the difference, ␮i − ␮0c is not preferable to formation of structure I compared to structure II, the slope of ␮c− ␮0c for structure I is steeper than that for structure II, which
gives rise to the lower dissociation pressure of structure I.
This is true for the methane hydrate with the parameters
examined here.
With those LJ parameters, we examine whether or not
structural transitions take place at 263.15 K. Now we consider a three-component system 共water+ ethane+ methane兲 atfour-phase equilibrium where an ice phase, a methane+ ethane gas mixture, and the two types of clathrate hydrates
coexist. As Gibbs’s phase rule tells us 共F = 2 + C − P, where

J. Chem. Phys. 122, 074503 共2005兲

FIG. 4. The chemical potential changes due to the incorporation of methane
molecules for hydrate I 共thick line兲 and II 共thin line兲 as a function of thecoexistence pressure at 263.15 K. The dotted line shows the chemical potential difference between ice and empty hydrate I 共below, thick兲 and II
共above, thin兲.

F , C, and P stand for the number of degree of freedom, thenumber of components, and the number of phases in thesystem兲, we are left with F = 2 + 3 − 4 = 1. Since the temperature is fixed, the composition of the gas mixture is determined uniquely. Thus, once the temperature is given, the fourphases can coexist only at a unique composition of the system, not in a finite range of composition. This equilibrium
condition is derived from the vdWP theory alone. The two
sets of the LJ parameters shown in Figs. 1 and 2 are subject
to further optimization in order to reproduce the structuraltransitions. Since the optimized LJ parameter set for each
guest species has one degree of freedom, the two sets of suchparameters have two degrees of freedom. So we first calculate the dissociation pressures of structure I and II mixed
hydrates over the two-dimensional optimized parameter
space, varying the methane mole fraction of the gas mixture.
The structure with lower dissociation pressure is more stable
at a given condition. For each given LJ parameter set, thestructural transition 共i.e., the four-phase equilibrium兲 is observed if the two dissociation pressures of structures I and II
become equal at some mole fraction of methane. Such mole
fraction is found for a wide range of the two-dimensional
optimized parameter space. In Fig. 5, we show the methanemole fraction at which the structural transition occurs at
263.15 K 共or the mole fraction at four-phase equilibrium兲 asa function of the size parameters of methane and ethane.
There are two surfaces in the three-dimensional space, whichindicates the structural transition occurs twice as the methanemole fraction is increased from 0 to 1. The structure II clathrate hydrate is the stable phase in the space between the two
surfaces and the structure I is more stable elsewhere. Now
we can further optimize the LJ parameter sets for the twoguest species such that they give the experimentally observed composition range, 78%–98%, in which structure II
forms at 263.2 K.13 Given in Table III are the free energies
obtained from those doubly optimized LJ parameters. The
free energy of larger cage occupancy for structure II is lower
than that for structure I for both guest species. But the difference is more prominent in the case of ethane, which is the

Downloaded 01 Nov 2005 to 150.46.100.222. Redistribution subject to AIP license or copyright, see http://jcp.aip.org/jcp/copyright.jsp

074503-6

Koyama, Tanaka, and Koga

J. Chem. Phys. 122, 074503 共2005兲

FIG. 5. Methane mole fraction for theoccurrence of structural transition
against both methane and ethane size
parameters at 263.15 K.

dominant factor for inducing the transition from structure I to
structure II together with the difference in relative number ofthe larger cages to the smaller cages. Figure 6 shows occupancy of methane and ethane along the dissociation pressure
at 263.15 K. Two discontinuities in each occupancy arecaused by the structural transition. The figure also shows thateach guest does not fully occupy the cages of clathrate hydrate at any mole fraction of methane. The mean cage occupancies averaged over the mole fraction for structures I and
II are about 76% and 71%, respectively. The transition to
structure II does not necessarily gives rise to larger cageoccupancy. Thus, the most dominant contribution to the stability of structure II is higher stability in its empty state
共Table I兲. We also show cage occupancy in Fig. 7 at a fixedpressure, 10 MPa, which is far higher than the dissociationpressure. Under this high pressure, the chemical potential ofwater ␮c given by Eq. 共4兲 determines which is the morestable structure of clathrate hydrate. As seen from Fig. 7, amethane+ ethane gas mixture almost fully occupies the cagesof either structure I or II; the mean cage occupancies ofstructure I and II are about 95%. The structural transition
from I to II occurs in higher methane mole fraction at 10
MPa compared with the case around the dissociation pressure although this feature is slightly different from the previous work.19 Comparing the results at 10 MPa 共Fig. 7兲 withthose at the dissociation pressure 共Fig. 6兲, we note that cageoccupancy of each guest increases with increasing gas pressure except for that of ethane in the structure II hydrate,
which hardly increases with increasing pressure because almost all the available cages are already occupied by ethane

FIG. 6. Occupancy of cages for methane 共solid line兲 and ethane 共dotted line兲
as a function of methane mole fraction at dissociation pressure and at
263.15 K.

FIG. 7. Occupancy of cages for methane 共solid line兲 and ethane 共dotted line兲
as a function of methane mole fraction at constant pressure, 10 MPa and at
263.15 K.

Downloaded 01 Nov 2005 to 150.46.100.222. Redistribution subject to AIP license or copyright, see http://jcp.aip.org/jcp/copyright.jsp

074503-7

J. Chem. Phys. 122, 074503 共2005兲

Stability of clathrate hydrates

molecules at the dissociation pressure. Hence it follows thatthe structure II hydrate cannot gain a significant amount of
the extra 共negative兲 free energy of larger cage occupancy
with increasing ethane gas pressure. This explains why thelower limit of the methane mole fraction for structure II hydrate is much higher at 10 MPa than at the dissociation pressure.
V. CONCLUSIONS

The appropriate LJ parameters are so determined as toreproduce the dissociation pressure for simple clathrates,
which contain only a single component of guest species, either methane or ethane. With those appropriately chosen LJ
parameters, we examine stability of structure I and II hydrates over the entire range of the mole fraction. Many parameter sets reproduce the expected structural transitions in acertain range of methane and ethane composition. We choose
among many the two parameter sets that reproduces the compositions of the phase transitions experimentally observed. In
this way, the LJ parameter sets for methane and ethane are
doubly optimized. The structural transition for the mixed
clathrate hydrate is accounted for by only the intermolecularinteractions, which are deduced from the simple clathrate
hydrates of methane and ethane.
There are three important factors to explain which form,
structure I or II, of the methane+ ethane mixed hydrate ismore stable. First empty structure II hydrate takes the lower
chemical potential value than empty structure I hydrate, second the ratio of the number of larger cages to the number ofwater molecules for structure I is much larger than that forstructure II, third the free energy of larger cage occupancy
for structure II is lower than that for structure I in the bothcases of methane and ethane. The first and third factors contribute to stabilizing structure II hydrate, while the second
one contributes to stabilizing structure I hydrate. The preferential structure at a given condition is determined by competition among these three factors.
The seemingly gradual transformation from structure I to
II with increasing the ethane mole fraction observed byexperiments12,13 needs to be examined carefully. Coexistenceof four phases in a certain range of the mole fraction of thegas mixture should be due to a kinetic factor, because it
violates the Gibbs phase rule and therefore cannot be a true
equilibrium state. Such observation, however, may be helpful
to provide some information on the role of guest species tostabilize cage structure.
ACKNOWLEDGMENTS

This work was supported by grant-in-aid from JSPS andMinistry of Education in Japan, by Okayama Foundation for
Science and Technology, and by NAREGI nanoscience
project. The authors are grateful to M. Horino for some help
in the early stage of the present study.
1

D. W. Davidson, in Water—A Comprehensive Treatise, edited by F. Franks
共Plenum, New York, 1973兲, Vol. 2.
2
E. D. Sloan, Clathrate Hydrates of Natural Gases, 2nd ed. 共Marcel Dekker, New York, 1998兲.
3
J. A. Ripmeester, J. S. Tse, C. I. Ratcliffe, and B. M. Powell, Nature
共London兲 325, 135 共1987兲.
4
J. H. van der Waals and J. C. Platteeuw, Adv. Chem. Phys., 2, 1 共1959兲.
5
H. Tanaka and K. Kiyohara, J. Chem. Phys. 98, 4098 共1993兲.
6
H. Tanaka and K. Kiyohara, J. Chem. Phys. 98, 8110 共1993兲.
7
H. Tanaka, J. Chem. Phys. 101, 10833 共1994兲.
8
G. D. Holder and J. H. Hand, AIChE J. 28, 440 共1982兲.
9E. M. Hendriks, B. Edmonds, R. A. S. Moorwood, and R. Szezepanks,
Fluid Phase Equilib. 117, 193 共1996兲.
10
S. Subramanian, R. A. Kini, S. F. Dec, and E. D. Sloan, Chem. Eng. Sci.
55, 1981 共2000兲.
11
S. Takeya, Y. Kamata, T. Uchida, J. Nagao, T. Ebinuma, H. Narita, A.
Hori, and T. Hondoh, Can. J. Phys. 81, 479 共2003兲.
12
A. L. Ballard and E. D. Sloan, Chem. Eng. Sci., 55, 5773 共2000兲.
13
T. Uchida, S. Takeya, Y. Kamata et al., J. Phys. Chem. B 106, 12426
共2002兲.
14
W. L. Jorgensen, J. Chandrasekhar, J. D. Madura, R. W. Impey, and M. L.
Klein, J. Chem. Phys. 79, 926 共1983兲.
15
D. W. Davidson, Y. P. Handa, C. I. Ratcliffe, J. S. Tse, and B. M. Powell,
Nature 共London兲 311, 142 共1984兲.
16
We calculate the free energy of cage occupancy by ethane in the smaller
cage of structure I. It is about −15 kJ mol−1, which is much higher thanthat in the larger cage, −37 kJ mol−1. Therefore, the probability for an
ethane molecule to be encaged in the smaller cage is extremely low when
the pressure is close to the dissociation pressure.
17
J. O. Hirschfelder, C. F. Curtiss, and R. B. Bird, Molecular Theory ofGases and Liquids 共Wiley, New York, 1954兲.
18
W. L. Jorgensen, J. D. Madura, and C. J. Swenson, J. Am. Chem. Soc.
106, 6638 共1984兲.
19
S. Subramanian, A. L. Ballard, R. A. Kini, S. F. Dec, and E. D. Sloan,
Chem. Eng. Sci. 55, 5763 共2000兲.

Downloaded 01 Nov 2005 to 150.46.100.222. Redistribution subject to AIP license or copyright, see http://jcp.aip.org/jcp/copyright.jsp