Formation of Clathrate Hydrates of Water-Soluble Guest Molecules

Sep 6, 2016 - We examine the formation mechanism of EO and THF hydrates using ... rates of THF hydrate largely deviate from the Wilson–Frenkel model...
0 downloads 0 Views 1MB Size
Subscriber access provided by UNIV OF CALIFORNIA SAN DIEGO LIBRARIES

Article

Formation of Clathrate Hydrates of Water-Soluble Guest Molecules Takuma Yagasaki, Masakazu Matsumoto, and Hideki Tanaka J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.6b06498 • Publication Date (Web): 06 Sep 2016 Downloaded from http://pubs.acs.org on September 11, 2016

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

The Journal of Physical Chemistry C is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Formation of Clathrate Hydrates of Water-Soluble Guest Molecules

Takuma Yagasaki, Masakazu Matsumoto, and Hideki Tanaka* Research Institute for Interdisciplinary Science, Okayama University, Okayama, 700-8530, Japan

1

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abstract Clathrate hydrates of water-soluble guest molecules, such as ethylene oxide (EO) and tetrahydrofuran (THF), have been often investigated in experimental studies instead of gas hydrates because their dissociation temperatures are higher than the ice point under ambient pressure. We examine the formation mechanism of EO and THF hydrates using molecular dynamics simulations. The crystal growth rates are determined by the simulations of the hydrate/solution two-phase coexistence. It is found that the growth rate of EO hydrate is an order of magnitude higher than that of THF hydrate. The growth rates of THF hydrate largely deviate from the Wilson-Frenkel model while the model well approximates the growth rates of EO hydrate, indicating that trapping of guest molecules on the hydrate surface, which causes the slowing down of crystal growth of THF hydrate, is insignificant for EO hydrate. We also perform long-time simulations of aqueous EO and THF solutions to examine nucleation of clathrate hydrate. Spontaneous nucleation occurs only in the EO solution within the simulation time. Similar to previous studies on methane hydrate, the obtained solid structure exhibits no long-ranged order. It is found that the 512 hydrate cage, which is the most dominant cage type in the early stage of the nucleation of methane hydrate, is not a major cage type in the nucleation process of EO hydrate.

2

ACS Paragon Plus Environment

Page 2 of 33

Page 3 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Introduction Clathrate hydrates are ice-like solids in which water cages entrap guest molecules.1 There are several types of hydrate structures.1-4 In many cases, small guest molecules such as H2 and large molecules like propane form cubic structure II (CS-II) hydrates while medium-size molecules such as CO2 form cubic structure I (CS-I) hydrates. There are vast quantities of natural gas hydrates in ocean sediments, and they are considered as a potential energy resource.5-8 Other than a resource of natural gas, clathrate hydrates have versatile utilities such as gas separation,9,10 gas storage,5,11,12 and cold storage.13,14 Thus, understanding the mechanism of formation of clathrate hydrates is crucial for their utilization. Although clathrate hydrate formers are typically small gas molecules, some water-soluble molecules also form clathrate hydrates. Hydrates of small gas molecules form only at high pressures when the temperature is higher than the ice point. In contrast, clathrate hydrates of water-soluble molecules form under ambient pressure at temperatures somewhat higher than the ice point, and it is easy to obtain a large single crystal.15 Moreover, it is possible to estimate the intrinsic formation rate of clathrate hydrates, which does not include effects of slow transfer of guest molecules from the gas phase to the growing hydrate surface.16 Therefore, this type of clathrate hydrates has been examined instead of gas hydrates in experimental studies on growth kinetics, morphology, and effects of additives.1,15-19 Tetrahydrofuran (THF) is the most popular water-soluble clathrate hydrate former. THF hydrate is a CS-II clathrate hydrate and it is produced at 4.4 ºC under ambient pressure.16 Natural gases are mixtures of methane and a small amount of larger molecules such as propane, and they form CS-II hydrates.1 Therefore, THF hydrate has been used as an analogue of natural gas hydrates. For

3

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

example, low dosage hydrate inhibitors (LDHI), which are used to prohibit formation of hydrate plugs in gas pipelines, are often tested with THF hydrate before field trials.20-22 Molecular dynamics (MD) simulations have provided a wealth of information on formation of clathrate hydrates.23-54 However, most of the studies focused only on gas hydrates. The number of MD studies on water-soluble hydrate formers is small despite their experimental importance.44,52-59 An important question arises as to whether the formation mechanism of clathrate hydrates of water-soluble guests is similar to that of gas hydrates. If there are substantial differences, water-soluble hydrate formers would be inadequate as analogues of gas hydrate formers. Recently, we investigated the crystal growth of THF hydrate in its melt using MD simulations.58 It was found that the growth of THF hydrate is slow compared with ice Ih and gas hydrates under highly supersaturated conditions. In THF hydrate, the guest molecules are encapsulated only in 51264 cages and all 512 cages are left empty because of the large molecular size of THF. There are open 512 cages at the hydrate surface and an open 512 cage can partially include a THF molecule. We calculated the potential of mean force for a THF molecule. It was shown that there are free energy barriers for desorption from the open 512 cage at the surface, and the barriers are higher than thermal energy. These free energy barriers cause the slow growth of THF hydrate since the THF molecules in the open 512 cages must be eliminated in the crystal growth process to form the correct crystal structure. Such a mechanism does not work for natural gas hydrates because the main component of natural gas, methane, is small enough to be encaged in the 512 cages. We also found that the mechanism of the adsorption of THF on the hydrate surface is essentially the same as that of a type of LDHI.58,60 It was suggested that use of THF hydrate instead of natural gas hydrates might be problematic because of 4

ACS Paragon Plus Environment

Page 4 of 33

Page 5 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

possible competitive effects between the adsorption of THF and that of LDHI. In this paper, we examine formation of clathrate hydrate of a different water-soluble hydrate former, ethylene oxide (EO). EO is smaller than THF in molecular size. EO hydrate was used as an analogue of pure methane hydrate in experimental studies because both of them are CS-I clathrate hydrates.17,19,61,62 To the best of our knowledge, there is no MD study on formation of EO hydrate. We modify a force field of EO so as to reproduce the dissociation temperature of EO hydrate. MD simulations of the hydrate/solution two-phase coexistence are performed to evaluate the growth rate. It is shown that the crystal growth of EO hydrate is much faster than that of THF hydrate. The growth rates are analyzed by a kinetic model for crystal growth known as the Wilson-Frenkel model.63,64 We also perform long simulations of aqueous EO and THF solutions. Spontaneous hydrate nucleation is observed only in the EO solution. The obtained solid structures are analyzed using the algorithm developed by Jacobson et al.65 We finally discuss the difference in the nucleation rate between EO and THF hydrates.

Computational Details MD simulations are performed using the GROMACS 4.6 package.66,67 The pressure is kept by the Parrinello-Rahman method,68,69 and the temperature is maintained by the Nosé-Hoover method.70,71 The particle mesh Ewald (PME) method is employed to calculate the electrostatic interactions.72,73 The crystal growth rate is determined by the direct coexistence method.48,74,75 First, a 6 × 6 × 12 unit cell replica of CS-I hydrate is prepared. The initial crystalline structure obeys the Bernal–Fowler’s ice rule and it has no net dipole moment.76 Each 51262 cage

5

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

is occupied by an EO molecule but all 512 cages are left empty in the initial configuration. The numbers of water molecules and EO molecules are 19872 and 2592, respectively. This composition corresponds to the EO concentration of 11.5 mole %. An NVT simulation is performed for 100 ps at T = 900 K for a half of the system to melt to an aqueous solution while the remaining molecules are retained in the initial hydrate configuration applying appropriate constraints. Then, an NPT simulation is performed for 40 ps at T = 200 K and P = 1 bar to relax the hydrate structure. The resultant structure, which is shown in Figure 1a, is used as the initial configuration of the crystal growth/dissociation simulations. The product NPT run is performed at P = 1 bar for temperatures ranging from 220 K to 286 K.

Figure 1 Snapshots of the hydrate/solution two-phase coexistence system at (a) t = 0 ns and (b) t = 250 ns. The simulation is performed at T = 278 K. EO molecules in the 51262 cages, those in the 512 cages, and those dissolved in the aqueous solution are colored gray, red, and blue, respectively. The hydrogen bond network of water in the hydrate is represented by lines while the water molecules in the solution are hidden for clarity. The green spheres in panel (b) show vacant 51262 cages in the grown hydrate.

6

ACS Paragon Plus Environment

Page 6 of 33

Page 7 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Water molecules are described by the TIP4P/Ice model.77 This model well reproduces the melting temperature of ice Ih and clathrate hydrates.48,58,77,78 The geometric parameters of EO are taken from Ref. 79. Lennard-Jones (LJ) parameters of EO are the same as the united-atom OPLS force field.80 The Lorentz-Berthelot rule is applied to any unlike LJ interactions. In the OPLS force field, the partial charge on the oxygen, qO, is -0.50 e and the charge on CH2, qC, is 0.25 e. We perform an MD simulation at 280 K starting from the structure shown in Figure 1a with the original OPLS force field. As shown in Figure 2, the potential energy of the system increases with time when qC = 0.25 e, indicating that the hydrate slab dissociates at this temperature. The experimental two-phase coexistence temperature for the 11.5 mole % aqueous EO solution is ~11 ºC17,61 and the reported ice point of the TIP4P/Ice model is 270 K with an uncertainty of 3 K.75 Thus, the hydrate slab should remain almost unchanged at T = 280 K. We therefore modify the force field of EO by changing the partial charge. Figure 2 shows that the slope of the potential energy curve decreases with decreasing qC, i.e., the hydrate/solution two-phase coexistence temperature becomes higher as the hydrophilicity of guest decreases. We perform single-phase simulations (the system consists of 1242 water and 162 EO molecules), and find that the potential energy of the solution system with qC = 0.2500 e is lower than that with qC = 0.2375 e by 0.37 kJ per one mole of water while the corresponding result for the hydrate phase is only 0.03 kJ mol-1. Guest molecules form hydrogen bonds with surrounding water molecules in the aqueous solution while not in the hydrate crystal. Thus, the increase in the hydrophilicity of guest enhances the stability of the aqueous phase relative to the hydrate phase, and it leads to lower hydrate/solution two-phase coexistence temperature. As shown in Figure 2, the hydrate slab remains unchanged for 7

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 33

long time as expected when qC = 0.2425 e. Therefore, all the following simulations are performed using this value. The force field parameters are summarized in Table 1.

Figure 2 Time evolution of the potential energy of the hydrate/solution two-phase coexistence system at 280 K.

The red, orange, purple, black, cyan, and green curves

are the results of simulations of qC = 0.2500 e, 0.2475 e, 0.2450 e, 0.2425 e, 0.2400 e, and 0.2375 e, respectively. The potential energy is given in kJ per one mole of water.

Table 1 Force field parameters (a) Geometric parameters bond length / Å water

EO

THF

angle / deg

O-H

0.9572

H-O-H

104.52

O-M

0.1577

O-C

1.4360

C-C

1.4720

O-Cα

1.4100

Cα-O-Cα

111.00

Cα-Cβ

1.5300

O-Cα-Cβ

109.40

σ/Å

ε / kJ mol-1

q/e

(b) Pair potential parameters

8

ACS Paragon Plus Environment

Page 9 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

water

EO

THF

O

3.1668

0.88217

0.0

H

0.0

0.0

0.5897

M

0.0

0.0

-1.1794

O

3.0

0.71128

-0.4850

C

3.8

0.49371

0.2425

O

3.0

0.71128

-0.568



3.8

0.49371

0.284



3.8

0.49371

0.0

We compare the partial charge for EO, qC = 0.2425 e determined based on the dissociation temperature, with the electrostatic potential (ESP) derived charge at the B3LYP/aug-cc-pVTZ level of theory. The ESP charge is calculated from the Merz-Singh-Kollman scheme using the Gaussian 09 package and the solvent effect is taken into account with the polarizable continuum model.81-83 The obtained ESP charge for the CH2 group of EO (the sum of the partial charges on a C atom and two H atoms) is 0.1635 e, and this value is much smaller than qC = 0.2425 e. This difference mainly arises from the fact that the EO molecules are described by a united-atom model in our MD simulations. When hydrogen atoms of the CH2 group are explicitly considered, the coulombic interaction energy between the CH2 group and a surrounding water molecule is dominated by the interaction between the positive charge on an H atom of the CH2 group and the negative charge site on the water molecule, M. The interaction between the negative charge on the C atom of the CH2 group and the M site on the water molecule is insignificant because the C-M distance is longer than the H-M distance. Therefore, the charge on the CH2 group is effectively more positive than the simple sum 9

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

of the partial charges of the C and H atoms, 0.1635 e. The united-atom model includes this effect and thus its charge, qC = 0.2425 e, is larger than 0.1635 e. We also perform MD simulations of aqueous EO and THF solutions to observe hydrate nucleation. We use a force field model for THF revised based on the OPLS model in our previous paper.58 A THF molecule is represented by a rigid plane pentagon consisting of five interaction sites (-Cβ-Cα-O-Cα-Cβ-). The force filed parameters are given in Table 1. The dissociation temperature of THF hydrate described by the combination of this model with the TIP4P/Ice model is 275 K, which is very close to the experimental value of 277.6 K. The EO and THF solutions are prepared by melting a 4 × 4 × 4 unit cell replica of EO hydrate crystal (2944 water and 512 EO molecules) and a 3 × 3 × 3 unit cell replica of THF hydrate crystal (3672 water and 216 THF molecules), respectively. An NVT simulation is performed for 100 ps at T = 900 K to melt the hydrate, which is followed by a 100 ps NPT run at T = 250 K and P = 1 bar to relax the liquid configuration and the volume of the system. The product NPT run is then performed at several temperatures in the deeply supercooled region at P = 1 bar. Guest molecules are classified according to the tetrahedrality of surrounding water molecules. The F3 parameter is calculated for water molecules in the first hydration shell (r < 5.8 Å) of a guest molecule.87 The guest molecule is categorized as a molecules dissolved in the aqueous solution if the number of surrounding water molecules with a highly tetrahedral structure (F3 < 0.3) is less than 18. Otherwise, the guest molecule is considered to come under a solid-like molecule.

Results and Discussion Figure 3 presents the time evolution of the potential energy of the hydrate/solution

10

ACS Paragon Plus Environment

Page 10 of 33

Page 11 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

two-phase coexistence system at several temperatures. The potential energy increases with time at T = 286 K owing to dissociation of the hydrate slab. The sharp rise at t = 80 ns is caused by the collapse of the thin hydrate slab. The slope of the potential energy decreases with decreasing temperature, and it becomes zero at T = 280 K, which is identified with the dissociation temperature of EO hydrate of this model. The slope becomes negative at T = 278 and 276 K, indicating that the hydrate slab grows as time evolves at these temperatures.

Figure 3 Time evolution of the potential energy of the hydrate/solution two-phase coexistence system at several temperatures.

Figure 1 presents snapshots of the simulation performed at T = 278 K.

In this

figure, EO molecules in the 51262 and 512 cages are colored gray and red, respectively. Figure 1b shows that most of the 51262 cages are occupied by EO molecules in the grown hydrate crystal. As shown by green spheres, there are only two empty 51262 cages. A certain amount of 512 cages are also occupied by guest molecules. In order to evaluate the composition of the grown hydrate, we define that a water molecule is in the hydrate slab when the water molecule is located within a distance of 5.8 Å from any one of solid-like EO molecules (this threshold length is determined by the first minimum of the EO-water radial distribution function in the hydrate phase). The number of water 11

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 33

molecules in the hydrate increases from 10454 at t = 0 ns to 15808 at t = 250 ns and the number of EO molecules in the hydrate increases from 1254 to 2047 in the same period. Hence, the number of EO molecules per unit cell in the grown crystal is  



 6.81, where  is the number of water molecules in a unit cell

of CS-I hydrate, 46. This value is in agreement with the experimental values at similar conditions of 6.72 to 6.82.17 Figure 4a presents the distribution functions of the water molecules in the hydrate slab at t = 0 and 250 ns obtained from the simulation performed at T = 278 K. We see that the left interface is located at z ~ −38 Å at t = 0 ns and z ~ −60 Å at t = 250 ns. To determine the position of the hydrate/solution interface, z0, we fit the distribution function with    tanh

$$% &

,

(1)

where ρh is the density of water molecules in the hydrate and w is the width of the interface. The red curve in Figure 4a is the fitting result at t = 250 ns. Figure 4b plots the time evolution of the surface position, z0(t), at T = 278 K. The growth rate of the interface is determined from the slope of z0(t).

Figure 4 (a) Distribution functions of the water molecules in the hydrate slab at t = 0 ns 12

ACS Paragon Plus Environment

Page 13 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(black) and t = 250 ns (purple). The red curve is obtained from the fitting to Eq. (1). (b) Time evolution of the location of the left hydrate/solution interface at T = 278 K.

Figure 5a presents the temperature dependence of the crystal growth rate of EO hydrate. Two points are plotted at each temperature because of the presence of two interfaces. The thermodynamic driving force for crystal growth, i.e., the difference in the chemical potential between the solution and hydrate phases, increases with decreasing temperature. However, thermal molecular motions are highly suppressed in the deeply supercooled region. Thus, the growth rate first increases and then decreases with decreasing temperature. A similar behavior has been observed in MD simulations of crystal growth of ice Ih.58,88,89

Figure 5 Growth rates of (a) EO hydrate and (b) THF hydrate.

The solid curve is

obtained from the fitting to the Wilson-Frenkel equation, Eq. (2).

We could not find experimental data of the growth rate of EO hydrate in literature. Alternatively, we compare the present result with growth rates of similar crystals calculated from MD simulations. Figure 5a shows that the growth rate is maximized at ~260 K, which is 20 K lower than the melting point of 280 K. The maximum growth 13

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 33

rate is roughly 0.4 Å/ns. This value is close to the maximum growth rate of ice Ih for the TIP4P/Ice model of ~0.45 Å/ns at T ~ 260 K.58 The maximum growth rate of EO hydrate is also close to the growth rate of H2S hydrate in the hydrate/solution/gas three-phase coexistence simulations at T = 275 K and P = 10 MPa reported by Liang and Kusalik, 0.39 Å/ns (the temperature of 275 K corresponds to ~15 K lower than the dissociation temperature of H2S hydrate at the same pressure).29 Figure 5b presents the growth rates of THF hydrate calculated in our previous study.58 The maximum growth rate of THF hydrate is ~0.03 Å/ns at T ~ 270 K. This value is an order of magnitude lower than the maximum growth rate of EO hydrate at T ~ 260 K. As mentioned in the Introduction, the slow growth of THF hydrate was attributed to the trapping of guest molecules in the open 512 cages on the growing surface.58 It is possible to evaluate the significance of the surface trapping effect using the Wilson-Frenkel equation given by

'(  )*( +1 − exp 0

123 45 4 67 445

89,

(2)

where kB, D(T), ∆Hf, Tm, and C are the Boltzmann constant, diffusion coefficient, enthalpy of fusion, dissociation temperature, and a constant, respectively.63,64 Only the thermal driving force and the diffusion coefficient are taken into account in this model. Thus, deviations from Eq. (2) become larger when the surface trapping effect is remarkable. We fit the growth rates of EO and THF hydrates with Eq. (2) using D(T), ∆Hf, and Tm as inputs. The enthalpy of fusion and the dissociation temperature are 5.2 kJ mol-1 and 280 K for EO hydrate and 4.6 kJ mol-1 and 275 K for THF hydrate. The diffusion coefficients of the aqueous EO solution (2944 water and 512 EO molecules) and the THF solution (3672 water and 216 THF molecules) are shown in Figure 6. The diffusion coefficient of the guest species is used for the fitting (the use of the diffusion 14

ACS Paragon Plus Environment

Page 15 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

coefficient of water yields almost the same result because there is no substantial difference in the temperature dependence of the diffusion coefficient between water and guest both for the EO and THF solutions). Figure 5b demonstrates that the growth rates of THF hydrate largely deviate from the Wilson-Frenkel model. In contrast, as shown in Figure5a, the Wilson-Frenkel equation well approximates the temperature dependence of the growth rate of EO hydrate, indicating that the surface trapping effect is insignificant for the crystal growth of EO hydrate. This is because the molecular size of EO is small and an EO molecule can be included in a 512 cage in the bulk crystal. This result and the similarity in the growth rate of EO and H2S hydrate suggest that the growth process of EO hydrate is similar to that of gas hydrates under high pressures at which gas molecules well dissolve in the liquid phase.

Figure 6 Self diffusion coefficients of EO (red) and water (black) in the aqueous EO solution and those of THF (blue) and water (green) in the aqueous THF solution. The diffusion coefficients of the EO solution cannot be evaluated at temperatures lower than 220 K because of the nucleation of hydrate.

In the hydrate/solution two-phase coexistence simulations, the {001} surfaces are exposed to the solution. This condition is chosen because early simulation studies 15

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

mainly examined this surface. However, it was reported that the growth rates of EO and THF hydrates depend on the crystallographic plane.15,16,19,52,90 Experimental studies demonstrated that slowly grown EO hydrate becomes a rhombic dodecahedral crystal, indicating that the growth of the {011} surface is slower than other surfaces.15,19 The shape of a single crystal of THF hydrate is octahedral, indicating that the {111} surface is the slowest growing surface for THF hydrate.15 We perform two-phase coexistence simulations for the {011} surface of EO hydrate and the {111} surface of THF hydrate at T = 260 K. The numbers of guest and water molecules are 2304 and 17664 for EO hydrate and 1536 and 26112 for THF hydrate. Snapshots of these simulations are provided in Supporting Information. We find that the growth rate of the {011} surface of EO hydrate, 0.32 Å ns-1, is lower than that of the {001} surface, 0.38 Å ns-1. The growth rate of the {111} surface of THF hydrate, 0.01 ns-1, is also lower than the corresponding result of the {001} surface, 0.02 Å ns-1. These results are in agreement with the experimentally observed crystal morphologies of EO and THF hydrates. The origin of the surface dependence is unclear because there are several possible mechanisms, such as the differences in the hydrogen bonding topology, the surface trapping effect, and the energy required to create a step at the surface.15,52,58,91 Further research is needed on this issue. Next, we discuss nucleation of clathrate hydrates. We perform MD simulations of the aqueous EO solution in the deeply supercooled region. Figure 7a shows the time evolution of the potential energy. The potential energy remains constant at higher two temperatures, 228 K and 224 K. In contrast, a decrease in the potential energy is observed at T = 220 K and 216 K, indicating that nucleation of clathrate hydrate occurs at these temperatures. The temperature, 220 K, is 60 K lower than the dissociation 16

ACS Paragon Plus Environment

Page 16 of 33

Page 17 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

temperature of EO hydrate. We plot the number of solid-like EO molecules against time in Figure 7b. This figure also shows that nucleation indeed occurs at T = 220 K and 216 K.

Figure 7 (a) Time evolution of the potential energy of the EO/water mixture. The simulations are started from a liquid configuration. (b) The number of solid-like guest molecules. The inset magnifies the initial part of the result at T = 220 K. The orange curve in the inset shows the number of EO molecules in the largest solid cluster at T = 220 K.

The inset in Figure 7b magnifies the initial part of the time evolution of the number of solid-like molecules at T = 220 K. The inset also shows the number of EO molecules in the largest cluster of solid-like EO molecules. We here assume that two solid-like EO molecules are in the same cluster when the distance between the centers of mass of these molecules is less than 9 Å (this threshold is determined from the first minimum of the radial distribution function of guest molecules in CS-I hydrate). In the trajectory at T 17

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

= 220 K, a cluster forms at t ~ 70 ns and grows to a size of Ns = 18 at t = 90 ns. However, this cluster completely disappears at t ~ 200 ns. Formation and dissociation of a relatively large cluster (Ns = 10) is also found at t ~ 250 K. A cluster forms at t ~ 290 ns and this cluster reaches the critical size. Figure 8 shows the growth behavior of this cluster. The cluster grows beyond Ns = 30 at t ~ 450 ns but then shrinks to Ns = 13 at t ~ 560 ns, i.e., the cluster at t = 460 ns shown in Figure 8 is still smaller than the critical size. The cluster continuously grows with time for t > ~560 ns. The final hydrate structure seems to be amorphous since it has no long-ranged order. The final solid structure at T = 216 K also does not have long-ranged order. Amorphous clathrate is commonly observed in nucleation simulations of clathrate hydrates of gas molecules, and it has been assumed as a precursor of crystalline hydrate.24,28,30,34-38,40,43-45 The present study demonstrates that amorphous clathrate forms irrespective of the solubility of guest molecules. A similar result was obtained in a simulation study using a coarse-grained force field model.43

18

ACS Paragon Plus Environment

Page 18 of 33

Page 19 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 8 Snapshots of the nucleation simulation of EO hydrate at T = 220 K. Solid-like and liquid-like EO molecules are shown as red and blue triangles, respectively. Hydrogen bonds of water molecules around solid-like EO molecules are represented by lines. Water molecules in the liquid phase are not shown.

We count the numbers of four regular hydrate cages, 512, 51262, 51263, and 51264 using the algorithm developed by Jacobson et al.65 Early MD studies showed that the dominant cage species is 512 in the amorphous methane hydrate.28,36,38,39 MD simulations of Zhang and coworkers demonstrated that there exist multiple formation pathways for methane hydrate, and direct formation of crystalline CS-I structure is one of them.37 Even in this case, 512 is the most abundant cage type in the early stage of the nucleation (300 < t < 500 ns in Figure 5b of ref. 37). In contrast, as shown in Figure 9, 19

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

the number of 51262 cages is larger than that of 512 cages in the amorphous EO hydrate. This difference is attributed to the size of guest molecules. Most hydrate formers can be categorized as one of three types: small molecules which favor 512 cages, medium-size molecules which fit with 51262 cages, and large guest molecules which can be included only in 51264 cages.1 CS-I hydrate consists of 512 and 51262 cages and CS-II hydrate is composed of 512 and 51264 cages. Therefore, the medium and large-size guest molecules form CS-I and CS-II hydrates, respectively. The small guest molecules form CS-II hydrate because the number density of 512 cages in CS-II is higher than that in CS-I. Both EO and methane form CS-I hydrate and thus they are categorized as medium-size hydrate formers. However, the size of methane is close to the boundary between the medium and small-size guest species (O2 forms CS-II although its size, 4.2 Å, is very close to the size of methane, 4.36 Å).1 Thus, the ratio of 512 cages to 51262 cages is larger for methane hydrate than for EO hydrate.

Figure 9 Time evolution of the numbers of 512 (red), 51262 (green), 51263 (black), and 51264 cages (blue) at (a) T = 220 K and (b) T = 216 K. The cyan curve shows the number 20

ACS Paragon Plus Environment

Page 20 of 33

Page 21 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

of empty 512 cages (the red curve is the sum of the numbers of filled and empty 512 cages).

It is worthwhile to mention the difference between all-atom and coarse-grained MD simulations. Jacobson et al. performed MD simulations of nucleation of water-soluble CS-I former using the coarse-grained mW model.43 The guest molecule was designed to mimic trimethylene oxide, which is somewhat larger than EO, and the concentration of the guest in the initial aqueous solution was half of the present study. The final solid structure of their study is amorphous. They reported that the most abundant cage type in the final amorphous structure is 512, and 85 % of them are empty cages. In the present all-atom MD study using the TIP4P/Ice model, in contrast, the number of 512 cages is smaller than that of 51262 cages, and roughly half of 512 cages are empty (Figure 9). The difference between the two studies suggests that the stability of empty 512 cages is sensitive to the type of force field, although there might be some effects arising from the differences in the guest size and the concentration. We examine the stability of amorphous EO clathrate. A spherical cluster of amorphous EO clathrate is obtained from a configuration along the trajectory at T = 220 K by fixing the positions of the water and guest molecules within 15 Å from a guest molecule and melting the rest of the molecules. MD simulations are carried out at several temperatures from this configuration to determine the dissociation temperature of the cluster. The simulations are also performed for a spherical cluster of crystalline EO hydrate with the same radius. Figure 10a presents the time evolution of the potential energy of the system including a spherical crystalline cluster. The potential energy decreases at T = 240 K while it increases at T = 245 K, indicating that the dissociation 21

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

temperature of the crystalline cluster is ~242.5 K. This temperature is lower than the dissociation temperature of the crystalline hydrate slab, 280 K, because of the Gibbs-Thomson effect (the dissociation temperature decreases with increasing the curvature of the interface).92-94 Figure 10b shows the corresponding result for the spherical amorphous cluster. The dissociation temperature of the amorphous cluster, ~232.5 K, is lower than that of the crystalline cluster of ~242.5 K, showing that the amorphous structure is less stable than the crystalline structure. A similar result was reported for methane hydrate described by the mW model.45

Figure 10 Time evolution of the potential energy of the system including (a) a crystalline EO hydrate cluster and (b) an amorphous EO clathrate cluster.

We also carry out MD simulations of the aqueous THF solution for 2000 ns at T = 228, 224, 220, 216, 208, and 200 K. The lowest temperature, T = 200 K, is lower than the dissociation temperature of THF hydrate by 75 K. In addition, the simulation time for the THF solution is twice longer than that for the EO solution of 1000 ns. However, no nucleation is observed in the THF solution. The nucleation rate is proportional to the diffusion coefficient in classical nucleation theory.95 As shown in Figure 6, the diffusion dynamics of the aqueous THF 22

ACS Paragon Plus Environment

Page 22 of 33

Page 23 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

solution is slower than that of the EO solution and the difference increases with decreasing temperature. The ratio of diffusion coefficient of EO to that of THF is 16 at T = 224 K. This is one reason for the slow nucleation of THF hydrate compared with EO hydrate in the deeply supercooled region. There might be other mechanisms for the slow nucleation rate of THF hydrate. For example, formation of a hydrate cluster of critical size is expected to be inhibited by the surface trapping effect. However, we cannot examine this mechanism because no nucleation occurs in this study. The nucleation rate depends not only on the diffusion coefficient but also on thermodynamic properties such as the surface free energy,95 however, the estimation of them is not a simple task.47 Further investigations are required for understanding the difference in nucleation mechanism between THF and EO hydrates.

Conclusions We have investigated formation mechanism of EO and THF hydrates using MD simulations. A force field model for EO is revised based on the united-atom OPLS model to reproduce the dissociation temperature of EO hydrate. The growth rate of EO hydrate is close to that of H2S hydrate while it is much higher than that of THF hydrate. The temperature dependence of the growth rate of EO hydrate is well approximated by the Wilson-Frenkel equation, showing that the surface trapping effect, which results in the slow growth of THF hydrate, does not affect the growth process of EO hydrate. The present study suggests that the growth mechanism of EO hydrate is similar to that of gas hydrates under high pressures. Spontaneous hydrate nucleation is observed in the aqueous EO solution at supercooling of ~60 K. The obtained solid structure has no long-ranged order.

23

ACS Paragon Plus Environment

The 512

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

hydrate cage, which is the most dominant cage type in the nucleation of methane hydrate, is not a major cage type in the nucleation process of EO hydrate. We have also performed simulations of the THF solution. However, hydrate nucleation is not observed within the simulation time of 2000 ns due to the slow diffusion dynamics of the aqueous THF solution in the deeply supercooled region. Clathrate hydrates of water-soluble guest molecules have been frequently used in experimental studies as analogues of gas hydrates. We believe that water-soluble hydrate formers are also useful in MD simulations. There are several difficulties in simulations of gas hydrate formation: i) Transfer of gas molecules from the gas phase to the growing surface is a slow process even at high pressures. Relatively small simulation boxes have been employed to reduce the simulation time. ii) The shape of the gas phase is often spherical or cylindrical (this can be avoided by performing simulations starting from the two-phase coexistence with planar liquid/gas interfaces and the gas phase region is large enough). In this case, the curvature of the liquid/gas interface, and thus the Young-Laplace pressure of the gas bubble, changes as the formation of gas hydrate proceeds. As a result, the thermodynamic driving force for hydrate formation becomes time-dependent. This would make a kinetic analysis complex. iii) It is difficult to examine the effects of LDHIs on the hydrate formation because LDHIs are amphiphilic and they are expected to favor the liquid/gas interface rather than the growing hydrate surface (this is not problematic in real systems because the distance between the liquid/gas interface and the hydrate surface is much longer).20-22 iv) The shape of the nucleus is an important factor in terms of nucleation kinetics. It is difficult to assess the intrinsic shape when there are bubbles near the nucleus. v) Simulations have been often performed for highly supersaturated solutions 24

ACS Paragon Plus Environment

Page 24 of 33

Page 25 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

without the initial gas phase to investigate homogeneous nucleation. In this case, the system size must be small because the free energy barrier for bubble nucleation decreases with increasing system size and thus bubble quickly form in large simulation boxes.96 The use of water-soluble hydrate formers instead of gas molecules allows us to avoid all of these problems (a simulation with a constant driving force, which is mentioned in ii), can be performed when the occupancy of the hydrate is known and the concentration of the guest in the solution is set so that it is identical to that in the hydrate. This technique has been used in experimental studies15,16). EO is more appropriate than THF because of the absence of the surface trapping effect and the higher nucleation and growth rates. MD simulations of EO hydrate would yield a wealth of information on the formation mechanism of clathrate hydrates.

Supporting Information Snapshots of the crystal growth simulations for the {011} surface of EO hydrate and the {111} surface of THF hydrate.

Author information Corresponding Author Phone: +81-(0)86-251-7769 e-mail: [email protected]

Acknowledgments The present work was supported by FLAGSHIP2020, MEXT within the priority 25

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

study5 (Development of new fundamental technologies for high-efficiency energy creation, conversion/storage and use) using computational resources of the K computer provided by the RIKEN Advanced Institute for Computational Science through the HPCI System Research project (Project ID: hp150267) and JSPS KAKENHI Grant Number JP16K17857. A part of calculations were performed on the computers at Research Center for Computational Science, Okazaki, Japan.

26

ACS Paragon Plus Environment

Page 26 of 33

Page 27 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

References (1) Sloan, E. D.; Koh, C. A. Clathrate Hydrates of Natural Gases; CPC Press: Boca Raton, 2008. (2) Ripmeester, J. A.; Tse, J. S.; Ratcliffe, C. I.; Powell, B. M. A New Clathrate Hydrate Structure. Nature 1987, 325, 135-136. (3) Matsumoto, M.; Tanaka, H. On the Structure Selectivity of Clathrate Hydrates. J. Phys. Chem. B 2011, 115, 8257-8265. (4) Huang, Y.; Zhu, C.; Wang, L.; Cao, X.; Su, Y.; Jiang, X.; Meng, S.; Zhao, J.; Zeng, X. C. A New Phase Diagram of Water under Negative Pressure: The Rise of the Lowest-Density Clathrate S-III. Sci. Adv. 2016, 2, e1501010. (5) Sloan, E. D. Fundamental Principles and Applications of Natural Gas Hydrates. Nature 2003, 426, 353-359. (6) Kvenvolden, K. A. Gas Hydrate and Humans. Ann. N.Y. Acad. Sci. 2000, 912, 17-22. (7) Boswell, R. Resource Potential of Methane Hydrate Coming into Focus. J. Pet. Sci. Eng. 2007, 56, 9-13. (8) Kida, M.; Suzuki, K.; Kawamura, T.; Oyama, H.; Nagao, J.; Ebinuma, T.; Narita, H.; Suzuki, H.; Sakagami, H.; Takahashi, N. Characteristics of Natural Gas Hydrates Occurring in Pore-Spaces of Marine Sediments Collected from the Eastern Nankai Trough, Off Japan. Energy Fuels 2009, 23, 5580-5586. (9) Kang, S.-P.; Lee, H. Recovery of Co2 from Flue Gas Using Gas Hydrate:  Thermodynamic Verification through Phase Equilibrium Measurements. Environ. Sci. Technol. 2000, 34, 4397-4400. (10) Linga, P.; Kumar, R.; Englezos, P. The Clathrate Hydrate Process for Post and Pre-Combustion Capture of Carbon Dioxide. J. Hazard. Mater. 2007, 149, 625-629. (11) Mao, W. L.; Mao, H. K. Hydrogen Storage in Molecular Compounds. Proc. Natl. Acad. Sci. U.S.A. 2004, 101, 708-710. (12) Mao, W. L.; Mao, H.-k.; Goncharov, A. F.; Struzhkin, V. V.; Guo, Q.; Hu, J.; Shu, J.; Hemley, R. J.; Somayazulu, M.; Zhao, Y. Hydrogen Clusters in Clathrate Hydrate. Science 2002, 297, 2247-2249. (13) Imai, S.; Okutani, K.; Ohmura, R.; Mori, Y. H. Phase Equilibrium for Clathrate Hydrates Formed with Difluoromethane + Either Cyclopentane or Tetra-N-Butylammonium Bromide. J. Chem. Eng. Data 2005, 50, 1783-1786. (14) Delahaye, A.; Fournaison, L.; Marinhas, S.; Chatti, I.; Petitet, J.-P.; Dalmazzone, D.; Fürst, W. Effect of THF on Equilibrium Pressure and Dissociation Enthalpy of CO2 Hydrates Applied to Secondary Refrigeration. Ind. Eng. Chem. Res. 2006, 45, 391-397. (15) Larsen, R.; Knight, C. A.; Sloan Jr, E. D. Clathrate Hydrate Growth and Inhibition. Fluid Phase Equilib. 1998, 150–151, 353-360. (16) Bollavaram, P.; Devarakonda, S.; Selim, M. S.; Sloan, E. D. Growth Kinetics of Single Crystal SII Hydrates: Elimination of Mass and Heat Transfer Effects. Ann. N.Y. Acad. Sci. 2000, 912, 533-543. (17) Huo, Z.; Jager, M. D.; Miller, K. T.; Sloan Jr, E. D. Ethylene Oxide Hydrate Non-Stoichiometry:

27

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Measurements and Implications. Chem. Eng. Sci. 2002, 57, 705-713. (18) McLaurin, G. E.; Whalley, E. Negative Octahedral Snowflakes or Tyndall Figures in Tetrahydrofuran Clathrate Hydrate. Nature 1988, 332, 711-712. (19) Larsen, R.; Knight, C. A.; Rider, K. T.; Sloan, E. D. Growth and Inhibition of Ethylene Oxide Clathrate Hydrate. Ann. N.Y. Acad. Sci. 2000, 912, 441-451. (20) Kelland, M. A. History of the Development of Low Dosage Hydrate Inhibitors. Energy Fuels 2006, 20, 825-847. (21) Kelland, M. A. Production Chemicals for the Oil and Gas Industry; CPC Press: Boca Raton, 2014. (22) Perrin, A.; Musa, O. M.; Steed, J. W. The Chemistry of Low Dosage Clathrate Hydrate Inhibitors. Chem. Soc. Rev. 2013, 42, 1996-2015. (23) Moon, C.; Taylor, P. C.; Rodger, P. M. Molecular Dynamics Study of Gas Hydrate Formation. J. Am. Chem. Soc. 2003, 125, 4706-4707. (24) Hawtin, R. W.; Quigley, D.; Rodger, P. M. Gas Hydrate Nucleation and Cage Formation at a Water/Methane Interface. Phys. Chem. Chem. Phys. 2008, 10, 4853-4864. (25) Vatamanu, J.; Kusalik, P. G. Unusual Crystalline and Polycrystalline Structures in Methane Hydrates. J. Am. Chem. Soc. 2006, 128, 15588-15589. (26) Vatamanu, J.; Kusalik, P. G. Molecular Insights into the Heterogeneous Crystal Growth of SI Methane Hydrate. J. Phys. Chem. B 2006, 110, 15896-15904. (27) Vatamanu, J.; Kusalik, P. G. Heterogeneous Crystal Growth of Methane Hydrate on Its SII [001] Crystallographic Face. J. Phys. Chem. B 2008, 112, 2399-2404. (28) Vatamanu, J.; Kusalik, P. G. Observation of Two-Step Nucleation in Methane Hydrates. Phys. Chem. Chem. Phys. 2010, 12, 15065-15072. (29) Liang, S.; Kusalik, P. G. Crystal Growth Simulations of H2S Hydrate. J. Phys. Chem. B 2010, 114, 9563-9571. (30) Liang, S.; Kusalik, P. G. Exploring Nucleation of H2S Hydrates. Chem. Sci. 2011, 2, 1286-1292. (31) Liang, S.; Kusalik, P. G. Nucleation of Gas Hydrates within Constant Energy Systems. J. Phys. Chem. B 2013, 117, 1403-1410. (32) Pirzadeh, P.; Kusalik, P. G. Molecular Insights into Clathrate Hydrate Nucleation at an Ice-Solution Interface. J. Am. Chem. Soc. 2013, 135, 7278-7287. (33) Lauricella, M.; Meloni, S.; Liang, S.; English, N. J.; Kusalik, P. G.; Ciccotti, G. Clathrate Structure-Type Recognition: Application to Hydrate Nucleation and Crystallisation. J. Chem. Phys. 2015, 142, 244503. (34) 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. (35) Walsh, M. R.; Beckham, G. T.; Koh, C. A.; Sloan, E. D.; Wu, D. T.; Sum, A. K. Methane Hydrate Nucleation Rates from Molecular Dynamics Simulations: Effects of Aqueous Methane Concentration,

28

ACS Paragon Plus Environment

Page 28 of 33

Page 29 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Interfacial Curvature, and System Size. J. Phys. Chem. C 2011, 115, 21241-21248. (36) Walsh, M. R.; Rainey, J. D.; Lafond, P. G.; Park, D.-H.; Beckham, G. T.; Jones, M. D.; Lee, K.-H.; Koh, C. A.; Sloan, E. D.; Wu, D. T.et al. The Cages, Dynamics, and Structuring of Incipient Methane Clathrate Hydrates. Phys. Chem. Chem. Phys. 2011, 13, 19951-19959. (37) Zhang, Z.; Walsh, M. R.; Guo, G. J. Microcanonical Molecular Simulations of Methane Hydrate Nucleation and Growth: Evidence That Direct Nucleation to SI Hydrate Is among the Multiple Nucleation Pathways. Phys. Chem. Chem. Phys. 2015, 17, 8870-8876. (38) Sarupria, S.; Debenedetti, P. G. Homogeneous Nucleation of Methane Hydrate in Microsecond Molecular Dynamics Simulations. J. Phys. Chem. Lett. 2012, 3, 2942-2947. (39) Jiménez-Ángeles, F.; Firoozabadi, A. Nucleation of Methane Hydrates at Moderate Subcooling by Molecular Dynamics Simulations. J. Phys. Chem. C 2014, 118, 11310-11318. (40) Jiménez-Ángeles, F.; Firoozabadi, A. Enhanced Hydrate Nucleation near the Limit of Stability. J. Phys. Chem. C 2015, 119, 8798-8804. (41) English, N. J.; Lauricella, M.; Meloni, S. Massively Parallel Molecular Dynamics Simulation of Formation of Clathrate-Hydrate Precursors at Planar Water-Methane Interfaces: Insights into Heterogeneous Nucleation. J. Chem. Phys. 2014, 140, 204714. (42) English, N. J.; MacElroy, J. M. D. Theoretical Studies of the Kinetics of Methane Hydrate Crystallization in External Electromagnetic Fields. J. Chem. Phys. 2004, 120, 10247-10256. (43) Jacobson, L. C.; Hujo, W.; Molinero, V. Nucleation Pathways of Clathrate Hydrates: Effect of Guest Size and Solubility. J. Phys. Chem. B 2010, 114, 13796-13807. (44) Jacobson, L. C.; Hujo, W.; Molinero, V. Amorphous Precursors in the Nucleation of Clathrate Hydrates. J. Am. Chem. Soc. 2010, 132, 11806-11811. (45) Jacobson, L. C.; Molinero, V. Can Amorphous Nuclei Grow Crystalline Clathrates? The Size and Crystallinity of Critical Clathrate Nuclei. J. Am. Chem. Soc. 2011, 133, 6458-6463. (46) Nguyen, A. H.; Jacobson, L. C.; Molinero, V. Structure of the Clathrate/Solution Interface and Mechanism of Cross-Nucleation of Clathrate Hydrates. J. Phys. Chem. C 2012, 116, 19828-19838. (47) Knott, B. C.; Molinero, V.; Doherty, M. F.; Peters, B. Homogeneous Nucleation of Methane Hydrates: Unrealistic under Realistic Conditions. J. Am. Chem. Soc. 2012, 134, 19544-19547. (48) Conde, M. M.; Vega, C. Determining the Three-Phase Coexistence Line in Methane Hydrates Using Computer Simulations. J. Chem. Phys. 2010, 133, 064507. (49) Michalis, V. K.; Costandy, J.; Tsimpanogiannis, I. N.; Stubos, A. K.; Economou, I. G. Prediction of the Phase Equilibria of Methane Hydrates Using the Direct Phase Coexistence Methodology. J. Chem. Phys. 2015, 142, 044501. (50) Nada, H. Growth Mechanism of a Gas Clathrate Hydrate from a Dilute Aqueous Gas Solution: A Molecular Dynamics Simulation of a Three-Phase System. J. Phys. Chem. B 2006, 110, 16526-16534. (51) Tung, Y.-T.; Chen, L.-J.; Chen, Y.-P.; Lin, S.-T. Molecular Dynamics Study on the Growth of

29

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Structure I Methane Hydrate in Aqueous Solution of Sodium Chloride. J. Phys. Chem. B 2012, 116, 14115-14125. (52) Nada, H. Anisotropy in Growth Kinetics of Tetrahydrofuran Clathrate Hydrate: A Molecular Dynamics Study. J. Phys. Chem. B 2009, 113, 4790-4798. (53) Wu, J.-Y.; Chen, L.-J.; Chen, Y.-P.; Lin, S.-T. Molecular Dynamics Study on the Equilibrium and Kinetic Properties of Tetrahydrofuran Clathrate Hydrates. J. Phys. Chem. C 2015, 119, 1400-1409. (54) Wu, J.-Y.; Chen, L.-J.; Chen, Y.-P.; Lin, S.-T. Molecular Dynamics Study on the Growth Mechanism of Methane Plus Tetrahydrofuran Mixed Hydrates. J. Phys. Chem. C 2015, 119, 19883-19890. (55) 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. (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) Alavi, S.; Ripmeester, J. A.; Klug, D. D. Molecular-Dynamics Simulations of Binary Structure II Hydrogen and Tetrahydrofurane Clathrates. J. Chem. Phys. 2006, 124, 14704. (58) Yagasaki, T.; Matsumoto, M.; Tanaka, H. Mechanism of Slow Crystal Growth of Tetrahydrofuran Clathrate Hydrate. J. Phys. Chem. C 2016, 120, 3305-3313. (59) Song, B.; Nguyen, A. H.; Molinero, V. Can Guest Occupancy in Binary Clathrate Hydrates Be Tuned through Control of the Growth Temperature? J. Phys. Chem. C 2014, 118, 23022-23031. (60) Yagasaki, T.; Matsumoto, M.; Tanaka, H. Adsorption Mechanism of Inhibitor and Guest Molecules on the Surface of Gas Hydrates. J. Am. Chem. Soc. 2015, 137, 12079-12085. (61) Glew, D. N.; Rath, N. S. Aqueous Nonelectrolyte Solutions. Part XIII. Ice and Hydrate Freezing Points of Aqueous Ethylene Oxide Solutions and the Formula of Congruent Ethylene Oxide Hydrate. Can. J. Chem. 1995, 73, 788-796. (62) Hollander, F.; Jeffrey, G. A. Neutron Diffraction Study of the Crystal Structure of Ethylene Oxide Deuterohydrate at 80°K. J. Chem. Phys. 1977, 66, 4699-4705. (63) Kirkpatrick, R. Crystal Growth from the Melt: A Review. Am. Mineral. 1975, 60, 798-814. (64) Jackson, K. The Interface Kinetics of Crystal Growth Processes. Interface Sci. 2002, 10, 159-169. (65) Jacobson, L. C.; Hujo, W.; Molinero, V. Thermodynamic Stability and Growth of Guest-Free Clathrate Hydrates: A Low-Density Crystal Phase of Water. J. Phys. Chem. B 2009, 113, 10298-10307. (66) Hess, B.; Kutzner, C.; van der Spoel, D.; Lindahl, E. GROMACS 4: Algorithms for Highly Efficient, Load-Balanced, and Scalable Molecular Simulation. J. Chem. Theory Comput. 2008, 4, 435-447. (67) Van der Spoel, D.; Lindahl, E.; Hess, B.; Groenhof, G.; Mark, A. E.; Berendsen, H. J. C. GROMACS: Fast, Flexible, and Free. J. Comput. Chem. 2005, 26, 1701-1718.

30

ACS Paragon Plus Environment

Page 30 of 33

Page 31 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(68) Nosé, S.; Klein, M. L. Constant Pressure Molecular Dynamics for Molecular Systems. Mol. Phys. 1983, 50, 1055-1076. (69) Parrinello, M.; Rahman, A. Polymorphic Transitions in Single Crystals: A New Molecular Dynamics Method. J. Appl. Phys. 1981, 52, 7182-7190. (70) Nosé, S. A Molecular Dynamics Method for Simulations in the Canonical Ensemble. Mol. Phys. 1984, 52, 255-268. (71) Hoover, W. G. Canonical Dynamics: Equilibrium Phase-Space Distributions. Phys. Rev. A 1985, 31, 1695-1697. (72) Darden, T.; York, D.; Pedersen, L. Particle Mesh Ewald: An N⋅Log(N) Method for Ewald Sums in Large Systems. J. Chem. Phys. 1993, 98, 10089-10092. (73) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G. A Smooth Particle Mesh Ewald Method. J. Chem. Phys. 1995, 103, 8577-8593. (74) Ladd, A. J. C.; Woodcock, L. V. Triple-Point Coexistence Properties of the Lennard-Jones System. Chem. Phys. Lett. 1977, 51, 155-159. (75) García Fernández, R.; Abascal, J. L. F.; Vega, C. The Melting Point of Ice Ih for Common Water Models Calculated from Direct Coexistence of the Solid-Liquid Interface. J. Chem. Phys. 2006, 124, 144506. (76) Petrenko, V. F.; Whitworth, R. W. Physics of Ice; Oxford University Press: Oxford, 1999. (77) Abascal, J. L. F.; Sanz, E.; García Fernández, R.; Vega, C. A Potential Model for the Study of Ices and Amorphous Water: Tip4p/Ice. J. Chem. Phys. 2005, 122, 234511. (78) Míguez, J. M.; Conde, M. M.; Torré, J.-P.; Blas, F. J.; Piñeiro, M. M.; Vega, C. Molecular Dynamics Simulation of CO2 Hydrates: Prediction of Three Phase Coexistence Line. J. Chem. Phys. 2015, 142, 124505. (79) Cunningham, G. L.; Boyd, A. W.; Myers, R. J.; Gwinn, W. D.; Le Van, W. I. The Microwave Spectra, Structure, and Dipole Moments of Ethylene Oxide and Ethylene Sulfide. J. Chem. Phys. 1951, 19, 676-685. (80) Briggs, J. M.; Matsui, T.; Jorgensen, W. L. Monte Carlo Simulations of Liquid Alkyl Ethers with the Opls Potential Functions. J. Comput. Chem. 1990, 11, 958-971. (81) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.et al: Gaussian 09. Gaussian, Inc.: Wallingford, CT, USA, 2009. (82) Singh, U. C.; Kollman, P. A. An Approach to Computing Electrostatic Charges for Molecules. J. Comput. Chem. 1984, 5, 129-145. (83) Besler, B. H.; Merz, K. M.; Kollman, P. A. Atomic Charges Derived from Semiempirical Methods. J. Comput. Chem. 1990, 11, 431-439. (84) Abascal, J. L. F.; Vega, C. A General Purpose Model for the Condensed Phases of Water:

31

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

TIP4P/2005. J. Chem. Phys. 2005, 123, 234505. (85) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. Comparison of Simple Potential Functions for Simulating Liquid Water. J. Chem. Phys. 1983, 79, 926-935. (86) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. The Missing Term in Effective Pair Potentials. J. Phys. Chem. 1987, 91, 6269-6271. (87) Baez, L. A.; Clancy, P. Computer-Simulation of the Crystal-Growth and Dissolution of Natural-Gas Hydrates. Ann. N.Y. Acad. Sci. 1994, 715, 177-186. (88) Rozmanov, D.; Kusalik, P. G. Temperature Dependence of Crystal Growth of Hexagonal Ice (Ih). Phys. Chem. Chem. Phys. 2011, 13, 15501-15511. (89) Weiss, V. C.; Rullich, M.; Köhler, C.; Frauenheim, T. Kinetic Aspects of the Thermostatted Growth of Ice from Supercooled Water in Simulations. J. Chem. Phys. 2011, 135, 034701. (90) Carver, T. J.; Drew, M. G. B.; Rodger, P. M. Characterisation of the {111} Growth Planes of a Type II Gas Hydrate and Study of the Mechanism of Kinetic Inhibition by Poly(Vinylpyrrolidone). J. Chem. Soc., Faraday Trans. 1996, 92, 5029-5033. (91) Liu, X. Y.; Boek, E. S.; Briels, W. J.; Bennema, P. Prediction of Crystal Growth Morphology Based on Structural Analysis of the Solid-Fluid Interface. Nature 1995, 374, 342-345. (92) Johnson, C. A. Generalization of the Gibbs-Thomson Equation. Surf Sci. 1965, 3, 429-444. (93) Uchida, T.; Ebinuma, T.; Ishizaki, T. Dissociation Condition Measurements of Methane Hydrate in Confined Small Pores of Porous Glass. J. Phys. Chem. B 1999, 103, 3659-3662. (94) Yagasaki, T.; Matsumoto, M.; Andoh, Y.; Okazaki, S.; Tanaka, H. Effect of Bubble Formation on the Dissociation of Methane Hydrate in Water: A Molecular Dynamics Study. J. Phys. Chem. B 2014, 118, 1900-1906. (95) Auer, S.; Frenkel, D. Numerical Prediction of Absolute Crystallization Rates in Hard-Sphere Colloids. J. Chem. Phys. 2004, 120, 3015-3029. (96) Wedekind, J.; Reguera, D.; Strey, R. Finite-Size Effects in Simulations of Nucleation. J. Chem. Phys. 2006, 125, 214505.

32

ACS Paragon Plus Environment

Page 32 of 33

Page 33 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

TOC graphic

33

ACS Paragon Plus Environment