Effect of Guest Size and Conformation on Crystal Structure and

May 2, 2013 - We first quantitatively analyze the dependence of the experimental lattice constants and formation pressure on the average maximum lengt...
1 downloads 11 Views 3MB Size
Article pubs.acs.org/JPCC

Effect of Guest Size and Conformation on Crystal Structure and Stability of Structure H Clathrate Hydrates: Experimental and Molecular Dynamics Simulation Studies Kyoichi Tezuka,† Kotaro Murayama,† Satoshi Takeya,‡ Saman Alavi,*,†,§ and Ryo Ohmura*,† †

Department of Mechanical Engineering, Keio University, 3-14-1 Hiyoshi, Kohoku-Ku,Yokohama 223-8522, Japan National Institute of Advanced Industrial Science and Technology (AIST), 1-1-1 Higashi, Tsukuba 305-8565, Japan § Department of Chemistry, University of Ottawa, Ottawa, Ontrario, Canada K1N 6N5 ‡

S Supporting Information *

ABSTRACT: To better understand the effect of size and flexibility of large molecule guest substances (LMGSs) on crystal lattice structure and thermodynamic stability of structure H (sH) clathrate hydrates, we performed powder X-ray diffraction (PXRD) measurements and Parrinello− Rahman molecular dynamics (MD) simulations on six alkane and cycloalkane LMGS in the presence of methane help gas. We first quantitatively analyze the dependence of the experimental lattice constants and formation pressure on the average maximum length of the LMGS guests as determined by MD simulation. The PXRD results show that within a family of LMGSs with similar molecular structures molecules of optimum size give better stability of sH hydrate phase. The most stable hydrate in each class has larger lattice constants along the a-axis and smaller lattice constant along the c-axis. The lattice expansion/shrinkage is well reproduced in the MD simulations. From the MD simulations, we determine the changes in the conformation as a result of encapsulation and the tilt angle with respect to the long axis of the sH large cages of the LMGSs. The results indicate that the molecular shapes inside the sH large cages can significantly differ from those of the most stable molecular structure in the gas phase. In the case of flexible molecules, such as 2methylbutane, the 1−4 dihedral angles and effective molecular sizes change upon encapsulation. Molecules with shorter length generally have larger tilt angles in the large cages; however, the effective width dimension of the LMGS also affects the tilt angle. Understanding the stability of sH hydrates of various LMGSs requires a consideration of guest molecule size, structural flexibility, and tilt angle in the cages. None of these quantities alone explain trends in the stability. During simulation trajectories, we observe changes in conformation in the LMGS molecules in the large cages. The effects of these different factors make a priori structural determinations of the large cages guests extremely complex.



INTRODUCTION

and development of heat pump/refrigeration system utilizing the heat for hydrate formation/dissociation.8 The formation conditions for sH hydrates vary widely, depending on the chemical species of LMGSs and help gases. For example, at 276 K, the methane gas pressure required for formation of sH clathrate hydrates containing a series of LMGSs ranges between 1.5 and 3.5 MPa.9 The details of how molecular factors affect the formation conditions are not clearly understood despite that the recent advances on understanding guest−host interactions on the lattice structure, molecular conformation of LMGSs in the lattice, and the formation conditions. The previous studies related to these points are briefly reviewed below.9−12 Udachin et al. revealed that lowering the temperature below 167 K induced the distortion of the sH hydrate lattice containing cyclooctane.10 They speculated that the distortion can due to the change of cyclooctane orientation at the low temperature. Lee et al. observed that 13C chemical shifts of the LMGSs in the sH

Clathrate hydrates, or simply hydrates, are ice-like crystalline inclusion compounds consisting of hydrogen-bonded host water molecules forming cages that contain guest molecules.1 Depending on the size and chemical structure of the guest substances, there are variations in the crystallographic structures of hydrates, leading to forms designated as structure I (sI), structure II (sII), and structure H (sH).1,2 The formation of sH hydrates requires two different guest substances: one is a relatively large molecule guest substance (LMGS), and the other is a small help gas, such as methane or xenon. Ripmeester and Ratcliffe have discussed a wide variety of LMGSs, including 2-methylbutane and 2,2-dimethylbutane.3 The LMGSs must generally have one molecular dimension of ∼8 Å (including the van der Waals radii) in order for the sH hydrate to be formed. Phase-equilibrium conditions of sH hydrates are milder (i.e., higher temperature/lower pressure) than those of sI or sII hydrates formed only with small molecule guest substances. The milder phase equilibrium conditions in the sH hydrate forming systems is advantageous for applications such as storage media for natural gas4 and hydrogen,5,6 CO2 capture,7 © XXXX American Chemical Society

Received: January 18, 2013 Revised: April 30, 2013

A

dx.doi.org/10.1021/jp4005899 | J. Phys. Chem. C XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry C

Article

Figure 1. Chosen structures of the LMGSs of the sH hydrates along with Newman projections along the 1−4 carbon chain: (a) 2-methylbutane with two conformers; (b) pinacolone; (c) 2,2-dimethylbutane with one conformer; (d) 2,3-dimethylbutane with two conformers; (e) methylcyclohexane; (f) methylcyclopentane.

constant along the c-axis decreases. Nevertheless, their data and analyses showed that there are exceptions to the general correlation. The exceptions could be partially caused by the development of the correlation without considering the presence of different chemical functional groups on the LMGSs. To understand this discrepancy and rationalize the relationship between the molecular volume and lattice constant systematically, further analysis concerning the conformations of the large molecules in the cage is needed. To better understand the effect of LMGSs on lattice constants and thermodynamic stability of sH hydrates, Tezuka et al. linked the equilibrium methane pressure of sH hydrate phases with the guest gas phase molecular volume, largest molecule dimension, and molar mass.9 These correlations are partially successful as the molecular flexibility was not considered when making these links.

hydrates move upfield in the spectrum compared to the same peaks of the molecule in gas phase.11 They suggest that such differences can be attributed to conformational changes required for fitting the LMGSs in the large hydrate cages and/or to guest−cage water interactions. On the basis of the results of powder X-ray diffraction (PXRD) measurements and semiempirical molecular orbital calculations, Takeya et al. describe the changes in lattice size and cohesive energies which reflect equilibrium conditions of sH hydrates are due to the changes in chemical structures of LMGSs.12 In their study the existence of some relations and some exceptions between the molecular volume and lattice constant were reported. They found that as a general correlation structure H hydrates with the greater molecular volume of LMGSs are more stable. The lattice constant along the a-axis also generally increases with the increase in the molecular volume of LMGSs, while the lattice B

dx.doi.org/10.1021/jp4005899 | J. Phys. Chem. C XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry C

Article

liquid nitrogen pool, and the hydrate sample was obtained from the vessel. The hydrate sample was subjected to PXRD measurements after storing in a container kept at a temperature below 90 K. For the PXRD measurements, the hydrate sample was finely powdered in a nitrogen atmosphere at a temperature below 100 K. The fine-powdered hydrate samples were toploaded on a copper (Cu) specimen holder. PXRD measurements were done using Cu Kα radiation by a θ/2θ step scan mode with a step width of 0.02° (40 kV, 40 mA; Rigaku model Ultima III). Here, the parallel beam optics was employed for the PXRD measurement to avoid shift in the 2θ angles of diffraction. To eliminate further systematical errors, LaB6 as an external standard reference (NIST Standard Reference Material 660a) was also measured. Analysis of the lattice constants was done by a full-pattern fitting method using the Rietveld program RIETAN-FP.13 For the computations, the initial coordinates of the water oxygen atoms in the sH clathrate unit cell for the simulations were taken from clathrate X-ray crystallography.14 The initial positions of the water hydrogen atoms about the oxygen atoms were assigned to be consistent with the ice rules while simultaneously minimizing the total energy and dipole moment of the unit cell.15 In the simulations, a 3 × 3 × 3 replica of the sH hexagonal unit cell with 918 water molecules is used. The center of mass of each large guest (from the optimized gasphase structure) and CH4 helper gas molecules were initially placed in the center of the respective cages. Guest structures and positions in the cages equilibrate during the simulation. The initial, gas phase structures of the large guest molecules were determined by quantum chemical structural optimization using the Gaussian 03 suite16 and the B3LYP/6-311++G(d,p) level of theory. The intermolecular van der Waals potentials between atoms i and j on different molecules are modeled as sum of Lennard-Jones (LJ) and electrostatic point charge interactions

In the present study, we quantitatively analyze the dependence of the lattice constants and thermodynamics stability on the molecular sizes of LMGSs in greater detail allowing for guest molecular flexibility. Lattice constants for six sH clathrate hydrates are measured and calculated by PXRD measurements and Parrinello−Rahman molecular dynamics (MD) simulations. We determine conformation changes of the LMGSs as a result of encapsulation and the orientation of the LMGSs with respect to the long axis of the sH large cages from the simulations. These factors are related to the lattice constant and hydrate stability. Furthermore, a general discussion is given on the effect of the structural variability on the experimental characterization of the sH clathrate hydrates using X-ray diffraction and NMR methods.



EXPERIMENTAL AND COMPUTATIONAL METHODS Six hydrate crystal samples were synthesized in a high-pressure vessel from mixtures of distilled water, pure CH4 (0.9999 mole fraction certified purity; Takachiho Chemical Industrial, Co., Ltd.), and each of the following LMGSs: 2-methylbutane (C5H12, also known as isopentane); 2,2-dimethylbutane (C6H14, also known as neohexane); 2,3-dimethylbutane; 3,3dimethyl-2-butanone (C6H12O, also known as pinacolone); methylcyclopentane (C6H12); and methylcyclohexane (C7H14). All LMGSs were purchased from Sigma-Aldrich with minimum purity of 97 mol %. The molecular structures of these LMGSs are shown in Figure 1. Hydrate crystal samples for PXRD measurements were prepared with liquid water, a LMGS liquid, and CH4 gas using a high-pressure cylinder. The main component of the apparatus was a stainless steel cylindrical vessel with inner volume of 200 cm3 which is equipped with a magnetic stirrer to agitate the fluids and hydrate crystals inside the vessel. The vessel was immersed in a bath filled with aqueous ethylene glycol solution. Each experimental run is commenced by placing 35 cm3 of liquid water and 20 cm3 of each LMGS liquid in the vessel. This ratio of water to LMGS is determined through trial and error. Hydrate sample made from this ratio contained less ice and solid LMGS. The vessel was then immersed in the temperaturecontrolled bath equipped with a PID-controlled heater and a cooler set. The temperature inside the vessel was kept at 275 K. The air was purged from the vessel by repeating the procedure of pressurization with CH4 to 0.4 MPa and depressurization to atmospheric pressure three times. CH4 gas was supplied from a high-pressure cylinder through a pressure-regulating valve until the pressure in vessel increase to 3.0 MPa to avoid the formation of simple CH4 sI hydrate. During hydrate formation in the vessel, a line connecting the vessel and CH4 supply cylinder was intermittently closed. Upon hydrate formation, when the pressure decreased to the equilibrium pressure of sH hydrate formed with CH4 and each LMGS, the vessel was recharged with CH4 gas to 3.0 MPa and repeated until no further pressure reduction was observed. Inside of the vessel was continuously agitated at 400 rpm after hydrate nucleation. Nearly complete conversion of water to hydrate was obtained when no further pressure reduction was observed. After the pressure stabilized, the vessel was subsequently removed from the temperature-controlled bath and immediately immersed in liquid nitrogen bath. After the temperature in the vessel decreased below 200 K, the vessel was removed from the liquid nitrogen pool and quickly disassembled. Then, the lower part of the vessel containing the hydrate sample was again placed in a

⎡⎛ ⎞12 ⎛ ⎞6 ⎤ qiqj σij σij V (rij) = ∑ 4εij⎢⎢⎜⎜ ⎟⎟ − ⎜⎜ ⎟⎟ ⎥⎥ + r ⎝ rij ⎠ ⎦ 4πε0rij i,j ⎣⎝ ij ⎠

(1)

where σij and εij are the distance and energy parameters of the ij pair separated by a distance of rij and qi and qj are the electrostatic point charges on the atoms. Water molecules in the clathrate were modeled using the TIP4P four-charge model,17 while the large cage guest molecules are allowed flexibility and their intra- and intermolecular potential interactions were modeled with the general AMBER force field.18 The Tse−Klein−McDonald potential was chosen for methane.19 The values for the parameters σii and εii for selected atom types are given in Table S1with potentials between unlike atoms calculated using the standard combination rules, εij = (εiiεjj)1/2 and σij = (σii + σjj)/2. Partial electrostatic charges on the atoms of the structurally optimized guest molecules were calculated using the “charges from electrostatic potential grid” (CHELPG) method21 implemented in the Gaussian 03. The point charges qi on the water and optimized guest molecules are given in Table S2. Nonisotropic constant pressure/temperature NpT (Parrinello−Rahman)22 molecular dynamics simulations on periodic simulation cells were performed using the DL_POLY software program version 2.2023 with pressure and temperature regulated using the modified Nosé−Hoover thermostat/ barostat algorithm22 with thermostat and barostat relaxation C

dx.doi.org/10.1021/jp4005899 | J. Phys. Chem. C XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry C

Article

Figure 2. Powder X-ray diffraction patterns of the sH hydrates at 100 K: (a) CH4 + 2-methylbutane hydrate; (b) CH4 + 2,2-dimethylbutane hydrate; (c) CH4 + 2,3-dimethylbutane hydrate; (d) CH4 + pinacolone hydrate; (e) CH4 + methylcyclopentane hydrate; (f) CH4 + methylcyclohexane hydrate. The curves in each pattern represent the observed intensities. The arrows show the patterns subjected to solid LMGSs. The asterisks show the calculated peak positions for hexagonal ice.

Table 1. Lattice Constants (Å) for Binary sH Clathrate Hydrates with Methane Measured by PXRD Measurements at 100 K and Calculated from Nonisotropic Constant Pressure−Temperature Molecular Dynamics Simulations at 100 K with 24.5 bar Pressurea lattice constants measured by PXRD/Å

lattice constants calculated by MD/Å

hydrate formers

a/Å

c/Å

a/c

a/Å

c/Å

a/c

2-methylbutane 2,2-dimethylbutane 2,3-dimethylbutane pinacolone methylcyclopentane methylcyclohexane

12.135 12.184 12.190 12.170 12.134 12.174

10.018 9.971 9.997 9.992 10.031 10.019

1.2113 1.2219 1.2193 1.2180 1.2097 1.2151

12.01 12.02 12.05 12.01 12.02 11.98

9.92 9.90 9.88 9.90 9.93 9.92

1.2109 1.2195 1.2201 1.2144 1.2060 1.2122

a The uncertainty in the experimentally measured lattice constants is 0.001 Å. The uncertainty in the lattice constants calculated by the simulation is between 0.1 and 0.2 Å.

times of 0.2 and 1.0 ps, respectively. The Verlet leapfrog algorithm was used with a time step of 1 fs. Long-range electrostatic interactions were calculated using the Ewald summation method,22 and all intermolecular interactions in the simulation box were calculated within a cutoff distance of Rcutoff = 13.0 Å. All systems were initialized by running an isotropic NpT simulation for 500 ps at the temperature and pressure of interest. Afterward, a nonisotropic NpT simulation performed for an additional time of 500 ps with the first 200 ps used for temperature scaled equilibration. Trajectory averaging was performed on additional simulations run for 200 ps. Simulations were run at 100 and 274 K temperatures and 0.245 kbar pressure. The 274 K simulations where performed with methane in the small cages. Unit cell lattice constants and the average largest (heavy atom)−(heavy atom) distance of the

guests in the large clathrate hydrate cages were determined after the simulations. We calculated the average tilt angle of the largest heavy atom distance in the molecule with respect to the c-axis of the large cages, which corresponds to the largest cage dimension. For guest molecules with a butane backbone, the CCCC dihedral angles are additionally determined for the guest molecules. These dihedral angles can be seen in the Newman projections of the molecules shown in Figure 1 and Table S5.



RESULTS AND DISCUSSION

The results of PXRD measurements performed at 100 K are shown in Figure 2. The crystal structure of all the hydrates tested in the present study is identified to be sH. The values of D

dx.doi.org/10.1021/jp4005899 | J. Phys. Chem. C XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry C

Article

the lattice constants a and c of the sH hydrates were measured and are given in Table 1 along with corresponding values from MD simulations. To quantitatively study the dependence of the lattice constants on molecular sizes and orientations of LMGSs in the sH large cages, at a time t in the simulation, we calculate all heavy atom (carbon···carbon and carbon···oxygen) distances in the encapsulated LMGSs from MD simulations and use the maximum length, dmax(t), of each molecule i of the 27 LMGS molecules in the large cages of the simulation cell to determine the average largest heavy atom distance, ⟨dclathrate(t)⟩ ⟨dclathrate(t )⟩ =

1 27

high activation energies of 20−30 kJ/mol to transform to each other (see Figure S1). The different conformers of the cycloalkane guests have considerably different effective lengths, but the energy difference between these conformers and the substantial activation energy required to convert the conformers to each other means that only the two conformers shown in Figure 1 will be observed in the hydrates. See below and ref 11 for further discussions. The tilt angle, θi(t), of the vector dmax,i(t) with respect to the c-axis of the simulation cell (the long axis of the sH large cages) from a snapshot at time t as defined in Figure 3 is also

27

∑ dmax ,i(t )

(2)

i=1

The calculated values of ⟨dclathrate(t)⟩ are given in Table 2. Table 2. Average Value of ⟨dclathrate(t)⟩ (Å), Tilt Angle ⟨θ(t)⟩ (in deg), and the Maximum Heavy Atom Distance of the Molecules in the Gas Phase, dgas (Å), in Binary sH Clathrate Hydrate of LMGS + Methane for a Simulation Snapshot at Time t Calculated from Parrinello−Raman Constant Pressure−Temperature Molecular Dynamics Simulations at 274 K and 24.5 bar guest

⟨dclathrate(t)⟩/Å

dgas/Å

⟨θ(t)⟩/deg

2-methylbutane 2,2-dimethylbutane 2,3-dimethylbutane methylcyclohexane methylcyclopentane pinacolone

3.826 3.923 3.955 4.396 3.758 3.860

3.923 3.928 3.938 4.410 3.844 3.927

23.7 25.5 28.2 17.8 20.5 26.2

Figure 3. Maximum heavy atom length, dmax(t), and the tilt angle, θ(t), of a pinacolone guest configuration with respect to the c-axis (the long axis of the sH large cage) at a time t in the simulation trajectory.

The 2-methylbutane molecule has two conformers shown in Figure 1 and characterized in Table S5. The more stable conformer (“anti”) in the gas phase has a dmax = 3.92 Å, while the second conformer (“gauche”) which is 3.7 kJ/mol higher in energy in the gas phase is smaller with dmax = 3.72 Å. There are two anti conformers and one gauche conformer for this molecule. The 2,2-dimethylbutane molecule has a single conformer. Rotating the dihedral angle in this molecule gives rise to the similarly structured conformers in these guest molecules. The 2,3-dimethylbutane molecule has two conformers (see Figure 1 and Table S5) which are close in energy. The one “anti, anti” conformer in the gas phase has a dmax = 3.94 Å, while the two “anti, gauche” conformers have dmax = 3.91 Å. Methylcyclopentane has two stable conformers; the conformer shown in Figure 1 with dmax = 3.84 Å is more stable by 4.8 kJ/mol than the conformer formed by flipping the C1 carbon to the other side of the ring. This second conformer, however, has a smaller dmax = 3.19 Å. Finally, the chair conformation of methylcyclohexane (shown in Figure 1) is more stable than the boat by 8.8 kJ/mol. The chair conformation has dmax = 4.41 Å while that of the boat is dmax = 3.82 Å. From the geometric characterizations, we observe that the greatest change in the effective molecule length of the guests upon conformational change for the substituted butanes occurs in 2-methylbutane. Conformational changes in 2,2-dimethyllbutane give identical structure. Conformational changes in 2,3dimethylbutane give structures with the approximately the same maximum length. The different conformers may be energetically accessible at normal temperatures but require relatively

calculated and the average values, ⟨θ(t)⟩, for the LMGSs in each hydrate are given in Table 2. We discuss trajectory averages ⟨dmax⟩ and ⟨θi⟩ below. The values for maximum length of the LMGS molecules as determined from gas phase structural optimization of these molecules are also given in Table 2 for comparison. Most encapsulated molecules become shorter than the gas phase structures. Figure 4 shows the correlation between the experimental and calculated lattice constants with ⟨dclathrate(t)⟩ for the sH hydrates. While the absolute values of the predicted lattice constants are smaller than experimental values, the trends of the a and c lattice constants measured by PXRD measurements are well reproduced with MD simulations. These results confirm usefulness of the ⟨dclathrate(t)⟩ in predicting the changes in lattice constant. With further refining of the force field used in the simulations, quantitative agreement between experimental and simulated lattice constants should be achievable. The values of ⟨dclathrate(t)⟩ of methylcyclopentane and methylcyclohexane LMGSs with cyclic structure are 3.76 and 4.40 Å, respectively. Within a family of LMGSs with four carbon main chain structures, 2-methylbutane, 2,2-dimethylbutane, 2,3-dimethylbutane, and pinacolone, the values of ⟨dclathrate(t)⟩ fall within the narrow range of 3.82−3.96 Å. As ⟨dclathrate(t)⟩ increases, within each family of LMGSs, the lattice constants are greater along the a-axis while smaller along the caxis. The lattice constants a and c expand or contract within E

dx.doi.org/10.1021/jp4005899 | J. Phys. Chem. C XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry C

Article

cyclohexane with larger ⟨dclathrate(t)⟩ has a lower equilibrium pressure. Within the alkyl family of LMGSs, the equilibrium CH4 pressure of the sH hydrate of 2-methylbutane is much higher than those of sH hydrate formed with other LMGSs, even though its ⟨dclathrate(t)⟩ is only slightly shorter than the other guests. There is no clear trend for the dependence of equilibrium pressure on ⟨dclathrate(t)⟩ and it is not possible to solely use ⟨dclathrate(t)⟩ to correlate with the thermodynamic stability of sH hydrates. To better understand the relationship between the molecular structure of LMGSs and the thermodynamic stability, it is necessary to consider the molecular structure and motion of LMGSs in the large cages, which are discussed below. Sample snapshots of 2-methylbutane and methylcyclohexane LMGSs in sH large cages at 274 K are shown in Figures 6 and

Figure 4. Lattice constants for sH hydrates with CH4 and various LMGSs at 100 K: (a) a-axis; (b) c-axis. The filled symbols are lattice constants measured by PXRD, and empty symbols are calculated from MD simulations: (▲) 2-methylbutane; (▼) 2,2-dimethylbutane; (◀) 2,3-dimethylbutane; (▶) pinacolone; (■) methylcyclopentane; (●) methylcyclohexane.

Figure 6. Details from a snapshot the simulation of the 2methylbutane sH clathrate hydrate at 274 K. The guest molecules are shown by the side of each cage for clarity. Both anti- and gaucheconformers of this 2-methylbutane are observed in the cages.

(0.5 ± 0.1)% in each family of hydrates. The crystal lattice expands along the a-axis but contracts along the c-axis, and hence the ratios of a/c vary to within (1.0 ± 0.2)% of each other, which is larger than the variation range of the separate a and c lattice vectors. The correlation between experimentally measured equilibrium methane pressures24−27 with ⟨dclathrate(t)⟩ is shown in Figure 5. Within the cycloalkane family of LMGSs, methyl-

7, and snapshots of 2,3-dimethylbutane configuration in the cages are shown in Figure S2. It is clear from these figures that guest molecule conformations and orientations inside the different large cages can differ significantly. The combination of the intramolecular and confinement forces of the guest in the cage can lead to guest molecules adapting conformations

Figure 5. Equilibrium pressure of sH hydrates formed with CH4 + various LMGSs in terms of ⟨dclathrate⟩ calculated from MD simulations: (▲) 2-methylbutane;23 (▼) 2,2-dimethylbutane;24 (◀) 2,3-dimethylbutane;25 (▶) pinacolone;24 (■) methylcyclopentane;25 (●) methylcyclohexane.26

Figure 7. Details from a snapshot of the simulation of the methylcyclohexane sH clathrate hydrate. The guest molecules are shown by the side of each cage for clarity. The guests do not show great change in shape but can be oriented in different directions and tilted by different amounts with respect to the long axis of the cages. F

dx.doi.org/10.1021/jp4005899 | J. Phys. Chem. C XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry C

Article

Figure 8. Variation of the dmax,i(t) and the tilt angle θ(t) (in deg) from a snapshot of the 27 2-methylbutane (●) and methylcyclohexane (□) guest molecules at time t in the simulation cell. The variation in the maximum lengths of the 2-methylbutane is greater than for the rigid methylcyclohexane. The variation in size of 2-methylbutane guests is related to the presence of both anti- and gauche-conformers of this molecule in the simulation.

different from those of the most stable structure in the gas phase. The tilt angles of the guest molecules with respect to the long axis of the large sH cages are also determined by a combination of the particular guest shape and the thermal energies of the guests and cage water molecules. The values of ⟨dclathrate(t)⟩ and dgas, which is the maximum heavy atom distance of LMGSs in gas phase, are given in Table 2. It is clear that the average largest heavy atom distance of the LMGSs in the sH hydrates are different from those in gas phase. The difference between ⟨dclathrate(t)⟩ and dgas of 2methylbutane is particularly large, which suggests that 2methylbutane can undergo the largest structural changes in the cages. This can be understood with reference to the dmax values for the different conformers of the LMGS molecules given in Table S4. The different conformers of 2,2-dimethylbutane, 2,3dimethylbutane, and pinacolone have identical or similar dmax values, while the activation barrier required to access the less stable conformers of methylcyclopentane and methylcyclohexane are so high that only the low energy conformers will be seen in the cages. In the discussions that follow, we will mostly focus on the variable-shape 2-methylbutane and the rigid fixedshape methylcyclohexane for further analysis. Variation in dmax,i(t) for 2-methylbutane and methylcyclohexane molecules from a snapshot of a simulation at 274 K are seen in Figures 6 and 7. The values of the dmax,i(t) for all 27 guest molecules for these two LMGSs from a simulation snapshot at time t are shown in Figure 8. We discuss the behavior of the trajectory-averaged quantities in the context of Figure 9. Anti and gauche conformations, which correspond to large and small values of dmax,i(t), respectively, are seen in 2methylbutane guest molecules in Figure 6. The distribution of dmax,i(t) in Figure 8 is seen to be bimodal around the stable anti (dmax = 3.9 Å) and gauche (dmax = 3.2 Å) conformations. In contrast, for methylcyclohexane, the dmax,i(t) values shown in Figure 8 are distributed in a fairly narrow range around 4.4 Å. To minimize repulsions from the cage, the encapsulated 2methylbutane molecules change conformation, which leads to smaller maximum heavy atom distances compared to the isolate guests. The simulation prediction regarding the stability of gauche conformers has been verified for guest species in sII clathrate

Figure 9. (a) Distribution of dihedral angles ϕi(t) among the 2methylbutane guest molecules at a time t in the 274 K simulation. The dihedral angles at ∼60° and ∼180° represent gauche- and anticonformers of the molecules in the clathrate hydrate. (b) Distribution of average dihedral angles ⟨ϕi⟩ of the guest molecule from a 200 ps simulation. The greater spread of ⟨ϕi⟩ shows that during the simulation the dihedral angles (guest configurations) change and the structure of the guests is dynamic.

hydrates. For example, the gauche conformation is observed for the n-butane guest31 and the boat conformation of tetrahyG

dx.doi.org/10.1021/jp4005899 | J. Phys. Chem. C XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry C

Article

dropyran32 guest as determined by single crystal X-ray diffraction. Neither of these is the most stable conformation of the sII guest in the gas phase. In the sH guests with the fourcarbon backbone, the energy difference between the anti and gauche forms of the guests in the gas phase is ∼4 kJ/mol or less (see Table S5).11 This value gives an indication of the magnitude of the guest−cage interaction energies which may favor the formation of the gauche forms in the cages. The differences observed between ⟨dclathrate(t)⟩ and the gas phase values for 2-methylbutane in Table 2 are due to changes in the dihedral angles for a fraction of these guest molecules in the clathrate hydrate sample. The differences between ⟨dclathrate(t)⟩ and the gas phase dmax are smaller in the other guests in conformational changes either do not change the dmax by a great amount, or the energy barrier to producing the other conformer is too high. As mentioned previously, the changes in dmax values in the guest molecules are related to changes in the dihedral angle ϕ in the molecules. The distribution of dihedral angles among the 27 2-methylbutane guests in a snapshot at time t, ⟨ϕi(t)⟩, and averaged over 200 ps, ⟨ϕi⟩, of the simulation trajectory are shown in Figure 9a. Similar to the dmax,i(t) values, the distributions of the dihedral angles ⟨ϕi(t)⟩ in the 2methylbutane molecule in Figure 9a are bimodal near 60° and 180° angles which show the presence of gauche- and anticonformers in the butyl chain of this molecule. The distribution of the trajectory averaged ⟨ϕi⟩ for each guest is shown in Figure 9b and is more disperse. The bimodal distribution of the dihedral angles in ⟨ϕi⟩ is not clear with 2-methylbutane molecules showing a spread of average dihedral angles between 60° and 180°. The reason for the spread in dihedral angles is that some 2-methylbutane molecules undergo conformation changes during the simulation and switch between the two antiand gauche-forms, and the average dihedral angle can thus show values in between 60° and 180°. Our previous calculations11 show that the energy of the staggered form of the 2-methylbutane molecule is high, and the 120° dihedral angle calculated for ⟨ϕi⟩ of some guests is the result of averaging between anti- and gauche-forms and does not indicate the presence of unstable staggered conformation among the guests. The transformation of the guests between stable anti and gauche conformers during the simulation trajectory is why we chose to do the analyses on geometrical quantities at a time t and not over the entire trajectory. However, we do not expect the other large molecules of the present study, i.e. 2,2dimethylbutane, 2,3-dimethylbutane, methylcyclohexane, methylcyclopentane, and pinacolone, to have more than one dmax value, and so the averages of this properties at a time t and over the entire trajectory should yield similar results for these guests. The distributions of tilt angles, θi(t), with respect to the long axis of the cages for the 27 2-methylbutane and methylcyclohexane guests are shown in Figure 8. The tilt angles for the rigid methylcyclohexane are distributed over a somewhat smaller range than 2-methylbutane but still show considerable variability. As given in Table 2, the average tilt angles of all LMGS molecules at any time ⟨θ(t)⟩ show some correlation to their corresponding ⟨dclathrate(t)⟩ values; however, it is not generally true that the molecules with smaller ⟨dclathrate(t)⟩ have the smaller tilt angles ⟨θ(t)⟩. The presence of bulky side chains on the LMGS in directions perpendicular to dmax of Figure 3 also affects the tilt angle. The bulkier (“wider”) methylcyclohexane molecules, for example, have a smaller average tilt angle,

which may be related to the requirement of incorporating the wide molecule in the sH large cages. The large sH cages are partially tapered near both poles. Values of dmax,i(t) and θi(t) for the 27 2-methylbutane and methylcyclohexane guests in the simulation are given in Table S4 and are plotted in Figure 10. There is considerable scatter in

Figure 10. Tilt angle θi(t) (in deg) against dmax,i(t) for the 27 2methylbutane (●) and methylcyclohexane (□) guest molecules in the simulation cell.

the relation between molecule sizes and tilt angles in the sample. This partly shows why average properties like ⟨dclathrate⟩ and ⟨θ⟩ alone do not directly correlate to the hydrate phase stability, as shown in Figure 5. The observations on the structural variation of the guests in the cages in the trajectory are significant from the point of view of modeling sH hydrate phases with statistical mechanical models such as the van der Waals−Platteeuw theory. A single assumed guest structure (usually based on gas-phase optimized structure) may not give a good description of state of the guests in the hydrate phase. If a single guest structure or size variable is used in such a model, the description may be untransferable or inconsistent with other simulations and force fields. Any model for predicting sH stabilities should be able to incorporate the guest flexibility in its framework. The observations of variability in guest conformations and positioning in the cages can also be important in interpreting and solving for guest positions in the sH clathrate hydrates from PXRD28 or single-crystal X-ray diffraction.29 As the simulations show, in addition to the large range of spatial disorder of the guests in the cages, there may be conformation variance in the structure of the guests, as well, which makes fitting the atomic positions of the guests rather difficult. The range of guest structures may also be reflected in the chemical shifts11 or relaxation times30 in the NMR spectra of large guests. A further intriguing aspect of the MD simulations is the dynamic conformational changes of some guests in the large sH cages. The anti−gauche transformation of guests like 2methylbutane should lead to other spectroscopic signatures in the sH clathrate hydrates.



CONCLUSIONS We performed PXRD measurements and MD simulations for different LMGS to analyze in detail the effect of guest molecular size and conformational changes on the crystal H

dx.doi.org/10.1021/jp4005899 | J. Phys. Chem. C XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry C

Article

(3) Ripmeester, J. A.; Ratcliffe, C. I. Xenon-129 NMR Studies of Clathrate Hydrates: New Guests for Structure II and Structure H. J. Phys. Chem. 1990, 94, 8773−8776. (4) Mori, Y. H. Recent Advances in Hydrate-Based Technologies for Natural Gas Storage − a Review. J. Chem. Ind. Eng. (China) 2003, 54 (Suppl.), 1−17. (5) Mao, W. L.; Mao, H. K. Hydrogen Storage in Molecular Compounds. Proc. Natl. Acad. Sci. U. S. A. 2004, 101, 708−710. (6) Duarte, A. R. C.; Shariati, A.; Rovetto, L. J.; Peters, C. J. Water Cavities of sH Clathrate Hydrate Stabilized by Molecular Hydrogen: Phase Equilibrium Measurements. J. Phys. Chem. B 2008, 112, 1888− 1889. (7) Brewer, P. G.; Friederich, G.; Peltzer, E. T.; Orr, F. M., Jr. Direct Experiments on the Ocean Disposal of Fossil Fuel CO2. Science 1999, 284, 943−945. (8) Ogawa, T.; Itoh, T.; Watanabe, K.; Tahara, K.; Hiraoka, R.; Ochiai, R.; Ohmura, R.; Mori, Y. H. Development of a Novel HydrateBased Refrigeration System: A Preliminary Overview. Appl. Therm. Eng. 2006, 26, 2157−2167. (9) Tezuka, K.; Taguchi, T.; Alavi, S.; Sum, A. K.; Ohmura, R. Thermodynamic Stability of Structure H Hydrates Based on the Molecular Properties of Large Guest Molecules. Energies 2012, 5, 459−465. (10) Udachin, K. A.; Ratcliffe, C. I.; Enright, G. D.; Ripmeester, J. A. Transformation of the Hexagonal-Structure Clathrate Hydrate of Cyclooctane to a Low-Symmetry from Below 167 K. Angew. Chem., Int. Ed. 2008, 47, 9704−9707. (11) Lee, J.-W.; Lu, H.; Moudrakovski, I. L.; Ratcliffe, C. I.; Ohmura, R.; Alavi, S.; Ripmeester, J. A. C-13 NMR Studies of Hydrocarbon Guests in Synthetic Structure H Gas Hydrates: Experiment and Computation. J. Phys. Chem. A 2011, 115, 1650−1657. (12) Takeya, S.; Hori, A.; Uchida, T.; Ohmura, R. Synthesis and Characterization of a Structure H Hydrate Formed with Carbon Dioxide and 3,3-Dimethyl-2-butanone. J. Phys. Chem. B 2006, 110, 12943−12947. (13) Izumi, F.; Momma, K. Rietvelt-Analysis Program RIETAN-98 and Its Applications to Zeolites. Solid State Phenom. 2007, 130, 15−20. (14) Udachin, K. A.; Ratcliffe, C. I.; Enright, G. D.; Ripmeester, J. A. Structure H Hydrate: A Single Crystal Diffraction Study of 2,2Dimethylpentane·5(Xe,H2S)·34H(2)O. Supramol. Chem. 1997, 8, 2732. (15) Okano, Y.; Yasuoka, K. Free-Energy Calculation of Structure-H Hydrates. J. Chem. Phys. 2006, 124, 024510. (16) Gaussian 03, Revision C.02: Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Montgomery, Jr., J. A.; Vreven, T.; Kudin, K. N.; Burant, J. C.; Millam, J. M.; Iyengar, S. S.; Tomasi, J.; Barone, V.; Mennucci, B.; Cossi, M.; Scalmani, G.; Rega, N.; Petersson, G. A.; Nakatsuji, H.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Klene, M.; Li, X.; Knox, J. E.; Hratchian, H. P.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Ayala, P. Y.; Morokuma, K.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Zakrzewski, V. G.; Dapprich, S.; Daniels, A. D.; Strain, M. C.; Farkas, O.; Malick, D. K.; Rabuck, A. D.; Raghavachari, K.; Foresman, J. B.; Ortiz, J. V.; Cui, Q.; Baboul, A. G.; Clifford, S.; Cioslowski, J.; Stefanov, B. B.; Liu, G.; Liashenko, A.; Piskorz, P.; Komaromi, I.; Martin, R. L.; Fox, D. J.; Keith, T.; AlLaham, M. A.; Peng, C. Y.; Nanayakkara, A.; Challacombe, M.; Gill, P. M. W.; Johnson, B.; Chen, W.; Wong, M. W.; Gonzalez, C.; Pople, J. A. Gaussian, Inc., Wallingford, CT, 2004. (17) Abascal, J. L. F.; Sanz, E.; Garcia Fernandez, R.; Vega, C. A Potential Model for the Study of Ices and Amorphous Water: TIP4P/ Ice. J. Chem. Phys. 2005, 122, 234511. (18) Cornell, W. D.; Cieplak, P.; Bayly, C. L.; Gould, I. R.; Merz, K. M., Jr.; 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 Acid, and Organic Molecules. J. Am. Chem. Soc. 1995, 117, 5179−5197.

lattice parameters and thermodynamic stability of sH hydrates. The sH hydrates in this work have 2-methylbutane, 2,2dimethylbutane, 2,3-dimethylbutane, 3,3-dimethyl-2-butanone, methylcyclopentane, and methylcyclohexane as LMGSs with methane as help gas. The results of PXRD measurements and MD simulations showed that the lattice constants vary with the maximum length of the guests, ⟨dclathrate(t)⟩ at a time t. The trends in the expansion/contraction of the a and c lattice constants as measured by PXRD measurements are well reproduced by MD simulations; however, the TIP4P potential for water and the GAFF force fields for the guests give a systematic underestimation of the these lattice constants. An increase in ⟨dclathrate(t)⟩ correlates well with a larger a lattice constant and smaller c lattice constant. In addition to effects of molecular size, we studied the average tilt angle of LMGSs in the large cages and determined that generally the larger molecules tend to have a smaller tilt angle in the cages. The tilt angles of the LMGSs also depend on the width of the molecule and not only the maximum length of the molecule. We observed that the 2methylbutane guest molecules show a bimodal distribution of dmax,i(t) and tilt angles θi(t) in the large cages which is related to the different possible dihedral conformers of the guests in the large sH cages. A calculation of trajectory averages of these quantities shows that they undergo changes of conformation between anti and gauche forms during the simulation. Our results show that a clear dependence of thermodynamic stability (based on methane pressure needed to form the hydrate) cannot be obtained exclusively with the single factor geometric factor ⟨dclathrate⟩. Other factors, including the tilt angle distribution and the dmax,i distribution (or equivalently the dihedral angle distribution), may affect the accommodation of the guests in the large sH cages and the stability of the sH phase. These factors must also be considered when understanding the effect of the guest on the lattice constants of the sH phase. The structure and stability of the sH clathrate hydrates depend on many geometrical and chemical details of the guest molecules in the hydrate phase. As we advance in our understanding of the hydrates, more of these factors need to be included in the microscopic and statistical mechanical analyses of these phases.



ASSOCIATED CONTENT

S Supporting Information *

Tables S1−S5; Figures S1 and S2. This material is available free of charge via the Internet at http://pubs.acs.org.



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected] (R.O.); [email protected] (S.A.). Notes

The authors declare no competing financial interest.



REFERENCES

(1) Sloan, E. D.; Koh, C. A. Clathrate Hydrates of Natural Gases, 3rd ed.; CRC Press: Boca Raton, FL, 2008. (2) Alavi, S.; Udachin, K. A.; Ratcliffe, C. I.; Ripmeester, J. A. Clathrate Hydrates. In Supramolecular Chemistry: From Molecules to Nanomaterials; John Wiley & Sons: New York, 2012. I

dx.doi.org/10.1021/jp4005899 | J. Phys. Chem. C XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry C

Article

(19) 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. (20) Hirschfelder, J. O.; Curtiss, C. F.; Bird, R. B. Molecular Theory of Gases and Liquids; Wiley: New York, 1967. (21) Breneman, C. M.; Wiberg, K. B. Determining Atom-Centered Monopoles from Molecular Electrostatic Potentials. The Need for High Sampling Density in Formamide Conformational Analysis. J. Comput. Chem. 1990, 11, 361−373. (22) (a) Frenkel, D.; Smit, B. Understanding Molecular Simulation; Academic Press: San Diego, 2000. (b) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Oxford Science Publications: Oxford, 1987. (23) DL_POLY 2.20; Forester, T. R., Smith, W., Todorov, I. T., Eds.; CCLRC, Daresbury Laboratory: Daresbury, UK, 2008. (24) Mehta, A. P.; Sloan, E. D. Structure H Hydrate Phase Equilibria of Methane + Liquid Hydrocarbon Mixtures. J. Chem. Eng. Data 1993, 38, 580−582. (25) Ohmura, R.; Uchida, T.; Takeya, S.; Nagao, J.; Minagawa, H.; Ebinuma, T.; Narita, H. J. Phase Equilibrium for Structure-H Hydrates Formed with Methane and either Pinacolone (3,3-Dimethyl-2butanone) or Pinacolyl Alcohol (3,3-Dimethyl-2-butanol). J. Chem. Eng. Data 2004, 48, 1337−1340. (26) Mehta, A. P.; Sloan, E. D. Structure H Hydrate Phase Equilibria of Paraffins Naphthenes, and Olefins with Methane. J. Chem. Eng. Data 1994, 39, 887−890. (27) Ohmura, R.; Uchida, T.; Takeya, S.; Nagao, J.; MinagaRwa, H.; Ebinuma, T.; Narita, H. J. Clathrate Hydrate Formation in (Methane + Water + Methylcyclohexanone) Systems: the First Phase Equilibrium Data. J. Chem. Thermdyn. 2003, 35, 2045−2054. (28) Takeya, S.; Udachin, K. A.; Moudrakovski, I. L.; Susilo, R.; Ripmeester, J. A. Direct Space Methods for Powder X-ray Diffraction for Guest-Host Materials: Applications to Cage Occupancies and Guest Distributions in Clathrate Hydrates. J. Am. Chem. Soc. 2010, 132, 524−531. (29) Alavi, S.; Udachin, K. A.; Ripmeester, J. A. Effect of Guest-Host Hydrogen Bonding on the Structures and Properties of Clathrate Hydrates. Chem.Eur. J. 2010, 16, 1017−1025. (30) 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. (31) Luzi, M.; Schicks, J. M.; Naumann, R.; Erzinger, J.; Udachin, K.; Moudrakovski, I.; Ripmeester, J. A.; Ludwig, R. Proc. 6th Inter. Conf. Gas Hydrates (ICGH 2008), Vancouver, British Columbia, Canada, July 6−10, 2008. (32) Udachin, K. A.; Ratcliffe, C. I.; Ripmeester, J. A. Single Crystal Diffraction Studies of Structure I, II and H Hydrates: Structure, Cage Occupancy and Composition. J. Supramol. Chem. 2002, 2, 405−408.

J

dx.doi.org/10.1021/jp4005899 | J. Phys. Chem. C XXXX, XXX, XXX−XXX