Identification of Clathrate Hydrates, Hexagonal Ice, Cubic Ice, and

Article pubs.acs.org/JPCB

Identification of Clathrate Hydrates, Hexagonal Ice, Cubic Ice, and Liquid Water in Simulations: the CHILL+ Algorithm Andrew H. Nguyen and Valeria Molinero* Department of Chemistry, The University of Utah, 315 South 1400 East, Salt Lake City, Utah 84112-0850, United States ABSTRACT: Clathrate hydrates and ice I are the most abundant crystals of water. The study of their nucleation, growth, and decomposition using molecular simulations requires an accurate and efficient algorithm that distinguishes water molecules that belong to each of these crystals and the liquid phase. Existing algorithms identify ice or clathrates, but not both. This poses a challenge for cases in which ice and hydrate coexist, such as in the synthesis of clathrates from ice and the formation of ice from clathrates during self-preservation of methane hydrates. Here we present an efficient algorithm for the identification of clathrate hydrates, hexagonal ice, cubic ice, and liquid water in molecular simulations. CHILL+ uses the number of staggered and eclipsed water−water bonds to identify water molecules in cubic ice, hexagonal ice, and clathrate hydrate. CHILL+ is an extension of CHILL (Moore et al. Phys. Chem. Chem. Phys. 2010, 12, 4124−4134), which identifies hexagonal and cubic ice but not clathrates. In addition to the identification of hydrates, CHILL+ significantly improves the detection of hexagonal ice up to its melting point. We validate the use of CHILL+ for the identification of stacking faults in ice and the nucleation and growth of clathrate hydrates. To our knowledge, this is the first algorithm that allows for the simultaneous identification of ice and clathrate hydrates, and it does so in a way that is competitive with respect to existing methods used to identify any of these crystals. cules,41−44 a tetrahedral order parameter for water combined with the number of guest neighbors around a central one,30,45 Voronoi tessellation of the space around the guest molecules,46 and a “mutually coordinated guest” order parameter that considers the order of guest and water molecules.47 These algorithms provide order parameters that are useful for distinguishing clathrates from liquid phases andin the case of the Vertex order parameter37to identify clathrate hydrate polymorphs. Existing order parameters are able to identify either clathrates or ice but not both, although a combination of order parameters can be used to achieve that goal. Rodger and coworkers used a combination of three order parameters to distinguish molecules in hexagonal ice and clathrate from those in solution.43,48 The three order parameters are F3 (deviation from the tetrahedral angle in the hydrogen-bonded network; similar to the tetrahedral order parameter in ref 49), F4t (average triple product of vectors in the dihedral of the hydrogen bonding system), and F4φ (torsional/dihedral in the hydrogen bonding network). Ice and clathrate were only identified when the three order parameters were used simultaneously. Similarly, Radhakrishnan and Trout proposed the use of a combination of bond orientational order parameters and a tetrahedral order parameter for water, plus guest order parameters to bias the formation of clathrate45 or

1. INTRODUCTION Water forms a wealth of crystalline phases. At ambient and moderate pressures, these include hexagonal ice, cubic ice (which has never been observed pure in experiments), stacking disordered ice I with cubic and hexagonal ice layers,1−9 and clathrate hydrates of structure sI, sII, TS, and HS-I.10−12 When methane clathrate hydrate is taken outside its region of stability, at temperatures just below the melting point of ice, it selfpreserves from further decomposition by the formation of an ice layer on its surface.13−20 This phenomenon is important for storage and transportation of gas as solid hydrate. On the other hand, gas clathrates are customarily synthesized in the laboratory from ice and natural gases such as methane and propane.15,21−25 To study through molecular simulations these processes in which ice and clathrate coexist, usually in the presence of a liquid phase, it is necessary to count with a computationally efficient and easy to implement algorithm that is able to identify water molecules in the different crystals and the liquid phase. There are several order parameters that have been successful in distinguishing crystals from liquids.4,26−31 Order parameters based on the local bond order ql based on spherical harmonics of order l = 3, 4, or 6 can distinguish ice from liquid.4,27,32−34 A neural network algorithm has been recently proposed and trained to identify several ice polymorphs.35 Algorithms to identify clathrate hydrates are limited to a few based on identifying the characteristic clathrate polyhedral cages (512, 51262, 51263, and 51264)36−39 and half polyhedral cages (56, 5661, and 5662),40 a vertex order parameter that can identify clathrate polymorphs based on the types of polyhedral cages to which each water molecule in the hydrate belongs,37 order parameters that evaluate dihedral angles between neighboring mole© XXXX American Chemical Society

Special Issue: Branka M. Ladanyi Festschrift Received: October 12, 2014 Revised: November 10, 2014

A

dx.doi.org/10.1021/jp510289t | J. Phys. Chem. B XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry B

Article

hexagonal ice50,51 from liquid water in umbrella sampling simulations. To monitor or bias the nucleation and growth of clathrate hydrates from ice, or vice versa, it is convenient to identify the different crystal phases and the liquid with a single efficient order parameter. The CHILL+ algorithm introduced in this work achieves this goal. CHILL+ is an extension of the CHILL algorithm,4 which identified cubic ice, hexagonal ice, interfacial ice, and liquid water. CHILL+ keeps the previous functionality, improves the detection of hexagonal ice, and also identifies clathrate hydrates. This paper is organized as follows: Section 2 presents the simulations used to parametrize and test CHILL+, section 3 the CHILL+ algorithm and parameters, and section 4 the validation of CHILL+ for the analysis of liquid water, stacking disorder in ice, and the nucleation and growth of clathrate hydrates. Section 5 presents the main conclusions of this work.

64. The simulation starts with a two-phase system (aqueous solution and a slab of small (S) guest fluid) composed of 5258 mW water molecules and 2742 S guest molecules. The simulation was evolved at 210 K and 500 atm for 50 ns. More details about the simulations of clathrate growth and nucleation can be found in refs 65 and 64, respectively.

3. CHILL+ ALGORITHM AND PARAMETERS CHILL+ identifies water molecules belonging to either liquid or crystals using the correlation of orientational order of a water molecule with its four closest neighbors. CHILL+ is an extension of the CHILL algorithm introduced in ref 4, which identifies ice but not clathrates. CHILL+ uses the number of staggered and eclipsed water−water bonds to distinguish cubic ice, hexagonal ice, and clathrate hydrates from liquid. Each of the four water−water bonds made by each molecule in cubic ice (ice Ic) is a staggered bond, while each molecule in hexagonal ice (ice Ih) has three water−water staggered bonds (the bonds in the basal plane) and one eclipsed bond (the bond parallel to the c-axis of the crystal).4,79 Each water molecule in clathrate hydrates, on the other hand, forms four eclipsed bonds. Figure 1 illustrates eclipsed and staggered bonds in

2. SIMULATIONS Molecular dynamics simulations were carried out using LAMMPS.52 Water was modeled using the mW water model, a monatomic model with short-range three-body tetrahedral “hydrogen-bonding” interactions.53 The mW model has been extensively used to study liquid water and solutions and the nucleation and growth of ice and clathrates.4,5,33,36,37,40,53−78 The equations of motion were integrated with the velocity Verlet algorithm using 10 fs time steps.53 All simulations were carried out using the isothermal−isobaric (NpT) ensemble with periodic boundary conditions, using a Nosé−Hoover thermostat and barostat with damping constants of 1 and 5 ps, respectively. Five sets of simulations were used to evaluate CHILL+. (1) Equilibrium simulations of crystal phases were used to parametrize the CHILL+ algorithm. Simulation cells of the crystal phases (cubic ice, hexagonal ice, sI, TS, and sII clathrate) were built from their respective unit cells to contain approximately 4000 water molecules. Equilibration simulations of the crystal phases were performed at 1 atm and 250 K for cubic and hexagonal ices and for sII hydrates, at 240 K for sI hydrate, at 235 K for TS and HS-I hydrates, and also at 270 K for hexagonal ice. These temperatures were selected to ensure the empty hydrates were just below their melting points.69 Each of these simulations was evolved for 10 ns. (2) Equilibrium simulations of liquid water consisting of 4096 mW water molecules evolved at p = 1 atm and T = 250, 270, and 298 K, for 20 ns at each temperature. (3) Growth simulations of hexagonal ice were used to assess the ability of CHILL+ to identify hexagonal ice and the formation of stacking faults. Two simulations consisting of 13824 mW water molecules were prepared by exposing (a) the prismatic plane of hexagonal ice in a 137 Å × 53 Å × 57 Å cell and (b) the basal plane of hexagonal ice in a 46 Å × 53 Å × 173 Å cell. Two-phase systems were prepared exposing each ice plane: liquid water in coexistence with a hexagonal ice slab that was one-third the length of the long axis. Simulations of these systems for 5 ns were performed at 250 K and 1 atm. (4) Simulation of growth of empty sI clathrate was utilized to compare the performance of CHILL+ and the cage recognition algorithm of ref 36. The simulation trajectory was obtained from ref 65 and consists of a simulation cell with 19872 mW water molecules evolved for 120 ns at 240 K and 1 atm from an initial configuration that contains a slab of empty sI clathrate surrounded by liquid water. (5) Simulation of nucleation of gas hydrate obtained from ref

Figure 1. Staggered (S) and eclipsed (E) bonds in sI clathrate and hexagonal ice. From left to right, the panels show an eclipsed bond (dark blue line) in sI clathrate, an eclipsed bond in ice Ih (dark blue line, parallel to the c-axis), and a staggered bond in ice Ih (green line). The red balls represent the four closest neighbors to the pair of molecules in the bond.

hexagonal ice and sI clathrate. Table 1 indicates the correspondence between the number of eclipsed and staggered bonds and the assignment made by the CHILL+ algorithm. The staggered and eclipsed bonds in CHILL and CHILL+ are identified with the correlation of orientational bond order parameter developed by ten Wolde et al.28 on the basis of the bond-order parameters of Steinhardt et al.27 The local order of rank l around each water molecule i is defined by a local orientational bond order parameter vector ql(i) with 2l + 1 complex components around the four nearest neighbors qlm(i) =

1 4

4

∑ Ylm( rij)̂ j=1

(1)

qlm(i) projects the orientational structure of the four closest neighbors of molecules based on spherical harmonics. The correlation, c(i, j), of the local order is measured by the dot product of ql(i) between a pair of neighboring molecules B



Article

Table 1. Identification of Ice, Clathrate, and Liquid in CHILL+ Using the Number of Staggered and Eclipsed Bonds structure cubic ice hexagonal ice interfacial ice clathrate interfacial clathrate liquid a

E

S

neigh.a

comment

0 1 any 0 4 3 N/A

4 3 2 3 0 any N/A

4 4 4 4 4 4 any

no change was made to staggered bond from ref 4 wider eclipsed range than in ref 4 identifies 99% hexagonal ice up to 270 K must have at least one first neighbor water with more than two staggered bonds4 must have at least one first neighbor water with more than one staggered bond4 bulk clathrate and part of the interface of clathrates can be identified with four eclipsed bonds partial clathrate cages and threads of clathrate-like order not cubic, hexagonal, interfacial ice, clathrate, or interfacial clathrate4

Number of water neighbors within 3.5 Å (first coordination shell).

c(i , j) =

=

Table 2. Definition of Staggered (S) and Eclipsed (E) Bonds in CHILL and CHILL+

ql(i)·ql(j) |ql(i)||ql(j)| l * (j) ∑m =−1 qlm(i)qlm l

l

* (i))1/2 (∑ * (j))1/2 (∑m =−1 qlm(i)qlm q (j)qlm m =−1 lm

a

(2)

It was found that l = 3 or 4 are best to distinguish between ice and liquid.4 The CHILL algorithm uses a correlation of q3 to distinguish crystals (cubic and hexagonal ice) from liquid and interfacial ice. CHILL identifies a staggered bond when c(i, j) ≤ −0.8 and an eclipsed bond when −0.05 ≥ c(i, j) ≥ −0.2.4 Figure 2 shows the distribution of the correlation of local order c(i, j) with l = 3 for ice Ih, ice Ic, liquid water, and sII

algorithm

eclipsed bonds (E)

staggered bonds (S)

CHILL+ CHILLa

0.25 ≥ c(i, j) ≥ −0.35 −0.05 ≥ c(i, j) ≥ −0.2

c(i, j) ≤ −0.8 c(i, j) ≤ −0.8

From ref 4.

the upper boundary of E from 0.20 to 0.25 only increased the number of water molecules identified as hydrate by less than 0.1%. The wider range for the eclipsed bond is required to identify water molecules that belong to the six-member ring boat conformation present in ice Ih as well as in the flat pentagonal and hexagonal rings of the clathrate hydrates. CHILL was developed for the identification of ice nucleated from deeply supercooled water. For a simulation of hexagonal ice at 250 K and 1 atm, CHILL identifies approximately 67% of ice as hexagonal, whereas CHILL+ recognizes 99% of ice as hexagonal.

4. VALIDATION OF CHILL+ 4.1. Liquid Water. Increasing the range of c(i, j) for eclipsed bonds produces a higher overlap of correlation of individual bond order in the liquid and the hexagonal ice and clathrate crystals; hence, we scrutinize the structure of liquid water with CHILL+. We did not observe polyhedral cages such as dodecahedra (512, 20 water molecules) and other clathrate cages within the 20 ns simulations of 4096 water molecules at T = 250, 270, and 298 K. Less than 1% of liquid water is identified by CHILL+ as either cubic or hexagonal ice. CHILL+ identified approximately 2−3% of water molecules in the liquid as clathrate when defining as part of the clathrate a water molecule with four eclipsed bonds (4E). These molecules with four eclipsed bonds form very small clusters in the liquid: the average sizes and standard deviations of water clusters with 4E in the liquid are 4 ± 1 at 250, 270, and 298 K. Throughout this work, we use a clustering cutoff distance of 3.5 Å that ensures all first neighbors are accounted for. The largest cluster size with 4E contains 8, 6, and 8 water molecules in the simulations at 250, 270, and 298 K. If we extend the definition of clathrate as molecules having 3E or 4E, CHILL+ identified approximately 3−5% of water molecules as clathrate. These molecules identified as clathrate were also scattered, forming small clusters: the largest 3E+4E clusters contained at most 13 molecules and averages 6 ± 2 for liquid water at 250, 270, and 298 K. Less than 1% of these 3E+4E water molecules belonged to clusters identified as clathrate with more than 10 molecules, which is half the size of the smallest clathrate cage. Although 3E +4E clusters represent a small fraction of the liquid, to ensure that the molecules in the liquid are not mislabeled as crystal, it

Figure 2. Distribution of correlation of orientational order parameters from spherical harmonics with l = 3 for hexagonal ice, cubic ice, liquid water, and sII clathrate at 250 K and 1 atm. The range that defines a staggered bond is identical in CHILL and CHILL+ and is shown with a violet bar. The cyan bar shows the CHILL+ range for the definition of an eclipsed bond (−0.35 to 0.25). The brown bar shows the narrower range for an eclipsed bond in CHILL. The wider range for an eclipsed bond in CHILL+ identifies more hexagonal ice at higher temperature than CHILL (99% vs 67% at 250 K). Although we only show the distribution of correlation of bond orientational order for the sII clathrate, the results are the same for the other clathrate polymorphs (sI, TS, and HS-I).

clathrate, all of them at 250 K. The violet and brown bars indicate the range for staggered and eclipsed bonds, respectively, in CHILL. The range for staggered bonds remains unchanged in CHILL+, but the range for eclipsed bonds is extended to 0.25 ≥ c(i, j) ≥ −0.35 (cyan bar) to allow the identification of clathrate and hexagonal ice at temperatures up to their melting points. Table 2 summarizes the ranges for eclipsed and staggered bonds in CHILL and CHILL+. The boundary of the E bond can be shifted slightly without significantly impacting the identification. For example, shifting C



Article

is important to filter the results by selecting the largest clusters of ice and clathrate. 4.2. Stacking Disorder in Ice. CHILL+ identifies at least 99% of hexagonal ice for equilibrated crystals of hexagonal ice in the simulations at 1 atm and T = 250 and 270 K. To demonstrate the ability of CHILL+ to identify the ice I polymorphs, liquid and interfacial ice, we performed simulations of growth of hexagonal ice at 250 K and 1 atm for two systems: one exposing the prismatic plane and the other exposing the basal plane to the supercooled liquid. The basal plane of hexagonal ice grows stacking disordered ice I (Figure 3), with layers of hexagonal and cubic ice, as shown in previous Figure 4. Growth of empty sI clathrate at 240 K and 1 atm shown by the time evolution of the number of water molecules in clathrate hydrate. Red represents clathrate cages (although it is possible to identify half cages with the cage code algorithm,36 we only considered complete cages), blue represents the water molecules with four eclipsed bonds (4E), and green represents the water molecules with at least three eclipsed bonds (3E+4E). 4E undercounts the water molecules in the cages by 2.5%; however, 3E+4E performs as well as the cage code in identifying all the polyhedral cages as well as some partial order at the interface, known as the clathrate halo (see Figure 5).

is reproduced very well by CHILL+ counting of molecules with at least three eclipsed bonds (3E+4E), i.e., interfacial clathrate plus clathrate of Table 1. Figure 5 displays a snapshot of the

Figure 3. Growth of hexagonal ice from liquid at 250 K and 1 atm results in a stacking disordered ice I crystal. The upper snapshot shows the initial block of hexagonal ice (cyan sticks), with the interfacial ice (purple balls; on the two sides of the periodic cell) in contact with liquid water (not shown). Ice Ih grows into a stacking disordered ice I crystal, with layers of cubic ice (yellow sticks) and hexagonal ice.

simulations with mW water,4,5,73−76,80 atomistic models,2,7,81−84 and experiments.1,3,6,8 The wider range of the eclipsed bond in CHILL+ ensures that less than 1% of the molecules within the crystal are identified as intermediate or interfacial ice. The growth of ice Ih exposing the prismatic plane does not produce stacking faults (as expected, as these would run perpendicular to the growing plane, making their nucleation more difficult). CHILL+ accurately monitors the growth of hexagonal ice, of cubic ice, or of a crystallite containing both ice I polymorphs. 4.3. Growth and Nucleation of Clathrate Hydrates. CHILL+ identifies that 100% of the water molecules in bulk clathrate crystals sI, sII, HS-I, and TS at 0 K and that 97.5−99% of the water molecules at a temperature just below their melting points have four eclipsed bonds. All the hydrate molecules are accounted for if clathrates are identified as those with either 4E or 3E. The definition of clathrate as a molecule with at least three eclipsed bonds recognizes not only full polyhedral cages but also partial and half cages. This allows for the observation of clathrate ordering surrounding the clathrate nucleus, i.e., the clathrate “halo”,43 as well as for the identification of small clathrate nuclei. In what follows, we demonstrate the use of CHILL+ for the analysis of simulations of growth and nucleation of clathrate hydrates. The cage identification code36 identifies the polyhedral cages 12 12 2 12 3 (5 , 5 6 , 5 6 , and 51264) that form Frank−Kasper clathrates, and it has been extensively used to track the progress of hydrate nucleation and growth in several studies.36−38,59,65−67,77,78,85,86 Figure 4 shows the number of water molecules in the hydrate phase as empty sI hydrate grows from supercooled water at 240 K and 1 atm.36 The growth of the hydrate crystal, as accounted for by the complete clathrate cages identified by the cage code,

Figure 5. CHILL+ recognizes clathrate cages, plus the interfacial clathrate halo. The figure shows a snapshot of the growth simulation of empty sI clathrate at 240 K and 1 atm. Water molecules in the liquid phase (between the two growing interfaces) are not shown here to ease the visualization. The polyhedral clathrate cages (identified using the cage algorithm of ref 36) are shown with red sticks. The clathrate crystal identified by CHILL+ as molecules in the largest cluster with at least three eclipsed bonds is shown with blue (4E) and green (3E) beads. Gray sticks connect the clathrate molecules identified by CHILL+ but not recognized by the cage code. Except for the few interfacial hydrate molecules, these are perfectly overlaid with the identification by the cage code. CHILL+ identifies half and partial cages at the clathrate surface.

growing crystal: there is an excellent correspondence between molecules identified as hydrates by CHILL+ and the cage code. CHILL+ with 3E+4E identifies all clathrate polyhedral clathrate cages, as well as half cages and partial hydrate order, the clathrate halo,43 present at the clathrate−solution interface. CHILL+ with the stricter 4E criterion identifies 97.5% of the water molecules found by the cage code. As evinced by the difference between complete cages and 4E when the crystallization is complete, the 2.5% difference mostly arises from thermal fluctuations within the bulk clathrate crystal that result in correlation values, c(i, j), outside the defined region. D



Article

5. CONCLUSIONS We introduced a new algorithm, CHILL+, which identifies water molecules that belong to hexagonal ice, cubic ice, clathrate hydrates, and liquid water based on the number of eclipsed and staggered water−water bonds. CHILL+ uses the positions of tetrahedrally coordinated water molecules to identify crystalline order. We illustrated its use for the monatomic water model mW, but the algorithm can also be used with atomistic models of water, and models of other tetrahedral substances, such as carbon, silicon, germanium, tin, silica, and germania, that also form clathrate and/or diamond crystals. CHILL+ is an extension of CHILL,4 which identifies hexagonal and cubic ice but not clathrates and uses a narrower range to identify eclipsed bonds. By increasing the range of the eclipsed bonds, CHILL+ significantly improves the detection of hexagonal ice up to its melting point, from 67% with CHILL to 99% with CHILL+. In this sense, CHILL+ supersedes CHILL for the analysis of ice nucleation and growth. Identification of hydrates during their nucleation and growth with CHILL+ and the cage identification code of ref 36 is similarly accurate. The cage code identifies not only the number but also the type of polyhedral cages in amorphous and crystalline clathrate (which allows for the identification of hydrate polymorphs with the vertex order parameter37); however, the recognition of hydrates with the cage code becomes expensive for large systems. On the other hand, CHILL+ recognizes the hydrate cages and the partial hydrate order at the clathrate surface (the clathrate halo) and is less computationally expensive, but it does not distinguish among the different types of clathrate cages and cannot identify hydrate polymorphs. CHILL+ is at least 5 times faster than the cage code when analyzing a simulation cell with ∼20000 water molecules. The efficiency of both codes can be improved by coupling them with a molecular dynamics or Monte Carlo code from which they can borrow the neighbor list, because the calculation of the neighbor list is one of the most expensive tasks in CHILL+. Even with the improvement of the neighbor list, the cage code would still be limited by the cost of the

We conclude that, with the parameters of CHILL+ of Table 2, the largest cluster of 3E+4E is best for identification of clathrate hydrates. The identification of the interfacial clathrate halo by CHILL+ allows for the characterization of the initial steps in the nucleation of hydrates. The halo reveals local order at a subcage resolution, providing insights on how clathrate cages are nucleated. Figure 6 shows the number of water molecules in the

Figure 6. Number of water molecules in the large hydrate cluster during nucleation of S hydrate. Representative simulation of nucleation at 210 K and 500 atm.64 The color code is the same as in Figure 4. Time t0 marks the formation of the critical clathrate nucleus, which grows to form an amorphous clathrate. The inset displays the number of water molecules identified as hydrate during the induction period.

hydrate phase along a simulation of the nucleation of clathrate hydrate with S guests at 210 K and 500 atm.64 Note that 3E+4E tracks almost exactly the molecules identified with the cage code when the nucleus has over ∼700 water molecules. Smaller hydrate nuclei are quite disordered,37,64,66,67,87 as can be seen in Figure 7, which illustrates the incipient evolution of clathrate nuclei toward a clathrate crystallite. Same as discussed above for the growth simulations, the order parameter 3E+4E reveals some hydrate-like ordering beyond that of the full polyhedral cages, which is useful to monitor or bias a simulation toward the formation of clathrates without requiring full formation of cages.

Figure 7. Incipient clathrate hydrate nuclei. Snapshots show the largest clathrate nuclei identified with CHILL+ and cage code. The largest clathrate nuclei were determined using a clustering cutoff distance of 3.5 Å. The color code in Figure 5 is also used here. CHILL+ shows the emerging nucleus at t0 (see Figure 6) which the cage code was not able to identify. Nuclei at subsequent times show the growth of the largest clathrate nucleus with CHILL+ showing the ordering surrounding the cages and the cages themselves. The S guest molecules are not shown here, as they are not considered by CHILL+ in the identification of the clathrate nuclei. E



Article

algorithm that searches for the connectivity in all possible fivemember rings and six-member rings before combining those rings to form polyhedral cages. The CHILL+ algorithm allows for an efficient identification when ice and clathrates are both present, for example, in the synthesis of clathrates from ice, the nucleation of ice from clathrates in self-preservation, and the competing growth of ice and clathrates from supercooled solutions. The figure that accompanies the Abstract illustrates such an application of CHILL+ to a simulation cell that comprises ice in contact with both a methane fluid and a saturated aqueous solution containing a clathrate nucleus. We conclude that CHILL+ is promising to monitor the evolution of clathrate hydrates or ice, even in the presence of the other crystal.

■

(14) Kuhs, W. F.; Genov, G.; Staykova, D. K.; Hansen, T. Ice Perfection and Onset of Anomalous Preservation of Gas Hydrates. Phys. Chem. Chem. Phys. 2004, 6 (21), 4917. (15) Kuhs, W. F.; Staykova, D. K.; Salamatin, A. N. Formation of Methane Hydrate from Polydisperse Ice Powders. J. Phys. Chem. B 2006, 110 (26), 13283−13295. (16) Falenty, A.; Kuhs, W. F. “Self-Preservation” of Co2 Gas Hydrates–Surface Microstructure and Ice Perfection. J. Phys. Chem. B 2009, 113 (49), 15975−15988. (17) Sun, D.; Shimono, Y.; Takeya, S.; Ohmura, R. Preservation of Carbon Dioxide Clathrate Hydrate at Temperatures Below the Water Freezing Point under Atmospheric Pressure. Ind. Eng. Chem. Res. 2011, 50 (24), 13854−13858. (18) Takeya, S.; Ripmeester, J. A. Dissociation Behavior of Clathrate Hydrates to Ice and Dependence on Guest Molecules. Angew. Chem., Int. Ed. Engl. 2008, 47 (7), 1276−1279. (19) Takeya, S.; Shimada, W.; Kamata, Y.; Ebinuma, T.; Uchida, T.; Nagao, J.; Narita, H. In Situ X-Ray Diffraction Measurements of the Self-Preservation Effect of Ch4hydrate. J. Phys. Chem. A 2001, 105 (42), 9756−9759. (20) Takeya, S.; Yoneyama, A.; Ueda, K.; Mimachi, H.; Takahashi, M.; Sano, K.; Hyodo, K.; Takeda, T.; Gotoh, Y. Anomalously Preserved Clathrate Hydrate of Natural Gas in Pellet Form at 253 K. J. Phys. Chem. C 2012, 116 (26), 13842−13848. (21) Rivera, J. J.; Janda, K. C. Ice Particle Size and Temperature Dependence of the Kinetics of Propane Clathrate Hydrate Formation. J. Phys. Chem. C 2012, 116 (36), 19062−19072. (22) Kuhs, W. F.; Klapproth, A.; Chazallon, B. In Chemical Physics of Air Clathrate Hydrates, Physics of Ice Core Records, Hokkaido University Press: 2000, pp 373−392. (23) Staykova, D. K.; Kuhs, W. F.; Salamatin, A. N.; Hansen, T. Formation of Porous Gas Hydrates from Ice Powders: Diffraction Experiments and Multistage Model. J. Phys. Chem. B 2003, 107 (37), 10299−10311. (24) Zylyftari, G.; Ahuja, A.; Morris, J. F. Nucleation of Cyclopentane Hydrate by Ice Studied by Morphology and Rheology. Chem. Eng. Sci. 2014, 116, 497−507. (25) Moudrakovski, I. L.; Sanchez, A. A.; Ratcliffe, C. I.; Ripmeester, J. A. Nucleation and Growth of Hydrates on Ice Surfaces: New Insights From129xe Nmr Experiments with Hyperpolarized Xenon. J. Phys. Chem. B 2001, 105 (49), 12338−12347. (26) Nelson, D. R. Order, Frustration, and Defects in Liquids and Glasses. Phys. Rev. B: Condens. Matter Mater. Phys. 1983, 28 (10), 5515−5535. (27) Steinhardt, P. J.; Nelson, D. R.; Ronchetti, M. BondOrientational Order in Liquids and Glasses. Phys. Rev. B 1983, 28 (2), 784−805. (28) ten Wolde, P. R.; Ruiz-Montero, M. J.; Frenkel, D. Numerical Calculation of the Rate of Crystal Nucleation in a Lennard-Jones System at Moderate Undercooling. J. Chem. Phys. 1996, 104 (24), 9932. (29) Mason, P. E.; Brady, J. W. “Tetrahedrality” and the Relationship between Collective Structure and Radial Distribution Functions in Liquid Water. J. Phys. Chem. B 2007, 111 (20), 5669−79. (30) Radhakrishnan, R.; Trout, B. L. Order Parameter Approach to Understanding and Quantifying the Physico-Chemical Behavior of Complex Systems. In Handbook of Materials Modeling; Yip, S., Ed.; Springer: Dordrecht, The Netherlands, 2005; pp 1613−1626. (31) Santiso, E. E.; Trout, B. L. A General Set of Order Parameters for Molecular Crystals. J. Chem. Phys. 2011, 134 (6), 064109. (32) Tanaka, H. Bond Orientational Order in Liquids: Towards a Unified Description of Water-Like Anomalies, Liquid-Liquid Transition, Glass Transition, and Crystallization: Bond Orientational Order in Liquids. Eur. Phys. J. E: Soft Matter Biol. Phys. 2012, 35 (10), 113. (33) Li, T.; Donadio, D.; Russo, G.; Galli, G. Homogeneous Ice Nucleation from Supercooled Water. Phys. Chem. Chem. Phys. 2011, 13 (44), 19807−19813.

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS We gratefully acknowledge the support by the National Science Foundation through awards CHE-1012651 and CHE-1309601. We thank the Center of High Performance Computing of the University of Utah for allocation of computing time and technical support.

■

REFERENCES

(1) Kuhs, W.; Bliss, D.; Finney, J. High-Resolution Neutron Powder Diffraction Study of Ice Ic. J. Phys. Colloques 1987, 48 (C1), C1-631− C1-636. (2) Carignano, M. A. Formation of Stacking Faults During Ice Growth on Hexagonal and Cubic Substrates. J. Phys. Chem. C 2007, 111 (2), 501−504. (3) Hansen, T. C.; Koza, M. M.; Kuhs, W. F. Formation and Annealing of Cubic Ice: I. Modelling of Stacking Faults. J. Phys.: Condens. Matter 2008, 20, 285104. (4) Moore, E. B.; de la Llave, E.; Welke, K.; Scherlis, D. A.; Molinero, V. Freezing, Melting and Structure of Ice in a Hydrophilic Nanopore. Phys. Chem. Chem. Phys. 2010, 12 (16), 4124−4134. (5) Moore, E. B.; Molinero, V. Is It Cubic? Ice Crystallization from Deeply Supercooled Water. Phys. Chem. Chem. Phys. 2011, 13 (44), 20008−20016. (6) Malkin, T. L.; Murray, B. J.; Brukhno, A. V.; Anwar, J.; Salzmann, C. G. Structure of Ice Crystallized from Supercooled Water. Proc. Natl. Acad. Sci. U.S.A. 2012, 109 (4), 1041−1045. (7) Pirzadeh, P.; Kusalik, P. G. On Understanding Stacking Fault Formation in Ice. J. Am. Chem. Soc. 2011, 133 (4), 704−707. (8) Kuhs, W. F.; Sippel, C.; Falenty, A.; Hansen, T. C. Extent and Relevance of Stacking Disorder in “Ice I(C)”. Proc. Natl. Acad. Sci. U.S.A. 2012, 109 (52), 21259−21264. (9) Morishige, K.; Uematsu, H. The Proper Structure of Cubic Ice Confined in Mesopores. J. Chem. Phys. 2005, 122 (4), 44711. (10) Sloan, E. D.; Koh, C. A. Clathrate Hydrates of Natural Gases, 3rd ed.; CRC Press/Taylor-Francis: Boca Raton, FL, 2007. (11) Udachin, K. A.; Enright, G. D.; Ratcliffe, C. I.; Ripmeester, J. A. Clathrate Hydrates: Some New Structural Information. In ACS Division of Fuel Chemistry, Preprints, 1997, 42 (2), 467−469. (12) Allen, K. W.; Jeffrey, G. A. On the Structure of Bromine Hydrate. J. Chem. Phys. 1963, 38 (9), 2304−2305. (13) Takeya, S.; Ripmeester, J. A. Anomalous Preservation of Ch4 Hydrate and Its Dependence on the Morphology of Hexagonal Ice. ChemPhysChem 2010, 11 (1), 70−73. F



Article

(56) Moore, E. B.; Molinero, V. Ice Crystallization in Water’s “NoMan’s Land”. J. Chem. Phys. 2010, 132 (24), 244504. (57) Moore, E. B.; Molinero, V. Growing Correlation Length in Supercooled Water. J. Chem. Phys. 2009, 130 (24), 244505. (58) Moore, E. B.; Allen, J. T.; Molinero, V. Liquid-Ice Coexistence Below the Melting Temperature for Water Confined in Hydrophilic and Hydrophobic Nanopores. J. Phys. Chem. C 2012, 116 (13), 7507− 7514. (59) Jacobson, L. C.; Molinero, V. A Methane-Water Model for Coarse-Grained Simulations of Solutions and Clathrate Hydrates. J. Phys. Chem. B 2010, 114 (21), 7302−7311. (60) Limmer, D. T.; Chandler, D. Phase Diagram of Supercooled Water Confined to Hydrophilic Nanopores. J. Chem. Phys. 2012, 137 (4), 044509. (61) Limmer, D. T.; Chandler, D. Premelting, Fluctuations and Coarse-Graining of Water-Ice Interfaces. arXiv preprint arXiv:1407.3514, 2014. (62) Holten, V.; Limmer, D. T.; Molinero, V.; Anisimov, M. A. Nature of the Anomalies in the Supercooled Liquid State of the Mw Model of Water. J. Chem. Phys. 2013, 138 (17), 174501. (63) Li, T.; Donadio, D.; Galli, G. Ice Nucleation at the Nanoscale Probes No Man’s Land of Water. Nat. Commun. 2013, 4, 1887. (64) Jacobson, L. C.; Hujo, W.; Molinero, V. Nucleation Pathways of Clathrate Hydrates: Effect of Guest Size and Solubility. J. Phys. Chem. B 2010, 114 (43), 13796−13807. (65) 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 (37), 19828−19838. (66) 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 (16), 6458−6463. (67) Jacobson, L. C.; Hujo, W.; Molinero, V. Amorphous Precursors in the Nucleation of Clathrate Hydrates. J. Am. Chem. Soc. 2010, 132 (33), 11806−11811. (68) Nguyen, A. H.; Molinero, V. Cross-Nucleation between Clathrate Hydrate Polymorphs: Assessing the Role of Stability, Growth Rate, and Structure Matching. J. Chem. Phys. 2014, 140 (8), 084506. (69) Nguyen, A. H.; Molinero, V. Stability and Metastability of Bromine Clathrate Polymorphs. J. Phys. Chem. B 2013, 117 (20), 6330−6338. (70) Shepherd, T. D.; Koc, M. A.; Molinero, V. The Quasi-Liquid Layer of Ice under Conditions of Methane Clathrate Formation. J. Phys. Chem. C 2012, 116 (22), 12172−12180. (71) Lupi, L.; Kastelowitz, N.; Molinero, V. Vapor Deposition of Water on Graphitic Surfaces: Formation of Amorphous Ice, Bilayer Ice, Ice I, and Liquid Water. J. Chem. Phys. 2014, 141, 18C508. (72) Lupi, L.; Molinero, V. Does Hydrophilicity of Carbon Particles Improve Their Ice Nucleation Ability? J. Phys. Chem. A 2014, 118 (35), 7330−7337. (73) Johnston, J. C.; Molinero, V. Crystallization, Melting, and Structure of Water Nanoparticles at Atmospherically Relevant Temperatures. J. Am. Chem. Soc. 2012, 134 (15), 6650−6659. (74) Hudait, A.; Molinero, V. Ice Crystallization in Ultrafine WaterSalt Aerosols: Nucleation, Ice-Solution Equilibrium, and Internal Structure. J. Am. Chem. Soc. 2014, 136 (22), 8081−8093. (75) Lupi, L.; Hudait, A.; Molinero, V. Heterogeneous Nucleation of Ice on Carbon Surfaces. J. Am. Chem. Soc. 2014, 136 (8), 3156−3164. (76) Bullock, G.; Molinero, V. Low-Density Liquid Water Is the Mother of Ice: On the Relation between Mesostructure, Thermodynamics and Ice Crystallization in Solutions. Faraday Discuss. 2013, 167, 371−388. (77) 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 (40), 23022−23031. (78) 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 (48), 19544−19547.

(34) Brukhno, A. V.; Anwar, J.; Davidchack, R.; Handel, R. Challenges in Molecular Simulation of Homogeneous Ice Nucleation. J. Phys.: Condens. Matter 2008, 20 (49), 494243. (35) Geiger, P.; Dellago, C. Neural Networks for Local Structure Detection in Polymorphic Systems. J. Chem. Phys. 2013, 139 (16), 164105. (36) Jacobson, L. C.; Hujo, W.; Molinero, V. Thermodynamic Stability and Growth of Guest-Free Clathrate Hydrates: A LowDensity Crystal Phase of Water. J. Phys. Chem. B 2009, 113 (30), 10298−10307. (37) Jacobson, L. C.; Matsumoto, M.; Molinero, V. Order Parameters for the Multistep Crystallization of Clathrate Hydrates. J. Chem. Phys. 2011, 135 (7), 074501. (38) 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.; Sum, A. K. The Cages, Dynamics, and Structuring of Incipient Methane Clathrate Hydrates. Phys. Chem. Chem. Phys. 2011, 13 (44), 19951. (39) Guo, G.-J.; Zhang, Y.-G.; Liu, C.-J.; Li, K.-H. Using the FaceSaturated Incomplete Cage Analysis to Quantify the Cage Compositions and Cage Linking Structures of Amorphous Phase Hydrates. Phys. Chem. Chem. Phys. 2011, 13 (25), 12048. (40) Bi, Y.; Li, T. Probing Methane Hydrate Nucleation through the Forward Flux Sampling Method. J. Phys. Chem. B 2014, DOI: 10.1021/jp503000u. (41) Baez, L. A.; Clancy, P. Computer Simulation of the Crystal Growth and Dissolution of Natural Gas Hydratesa. Ann. N. Y. Acad. Sci. 1994, 715 (1), 177−186. (42) Moon, C.; Taylor, P. C.; Rodger, P. M. Molecular Dynamics Study of Gas Hydrate Formation. J. Am. Chem. Soc. 2003, 125 (16), 4706−4707. (43) Moon, C.; Hawtin, R. W.; Rodger, P. M. Nucleation and Control of Clathrate Hydrates: Insights from Simulation. Faraday Discuss. 2007, 136, 367−382. (44) Rodger, P. M.; Forester, T. R.; Smith, W. Simulations of the Methane Hydrate/Methane Gas Interface near Hydrate Forming Conditions Conditions. Fluid Phase Equilib. 1996, 116 (1−2), 326− 332. (45) Radhakrishnan, R.; Trout, B. L. A New Approach for Studying Nucleation Phenomena Using Molecular Simulations: Application to Co[Sub 2] Hydrate Clathrates. J. Chem. Phys. 2002, 117 (4), 1786. (46) Chakraborty, S. N.; Grzelak, E. M.; Barnes, B. C.; Wu, D. T.; Sum, A. K. Voronoi Tessellation Analysis of Clathrate Hydrates. J. Phys. Chem. C 2012, 116 (37), 20040−20046. (47) Barnes, B. C.; Beckham, G. T.; Wu, D. T.; Sum, A. K. TwoComponent Order Parameter for Quantifying Clathrate Hydrate Nucleation and Growth. J. Chem. Phys. 2014, 140 (16), 164506. (48) 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 (32), 4853−4864. (49) Errington, J. R.; Debenedetti, P. G. Relationship between Structural Order and the Anomalies of Liquid Water. Nature 2001, 409 (6818), 318−321. (50) Radhakrishnan, R.; Trout, B. L. Nucleation of Hexagonal Ice (Ih) in Liquid Water. J. Am. Chem. Soc. 2003, 125 (25), 7743−7747. (51) Radhakrishnan, R.; Trout, B. L. Nucleation of Crystalline Phases of Water in Homogeneous and Inhomogeneous Environments. Phys. Rev. Lett. 2003, 90 (15), 158301/1−158301/4. (52) Plimpton, S. J. Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comput. Phys. 1995, 117, 1−19. (53) Molinero, V.; Moore, E. B. Water Modeled as an Intermediate Element between Carbon and Silicon. J. Phys. Chem. B 2009, 113 (13), 4008−4016. (54) DeMille, R. C.; Molinero, V. Coarse-Grained Ions without Charges: Reproducing the Solvation Structure of Nacl in Water Using Short-Ranged Potentials. J. Chem. Phys. 2009, 131 (3), 034107. (55) Moore, E. B.; Molinero, V. Structural Transformation in Supercooled Water Controls the Crystallization Rate of Ice. Nature 2011, 479 (7374), 506−508. G



Article

(79) Svishchev, I.; Kusalik, P. Crystallization of Liquid Water in a Molecular Dynamics Simulation. Phys. Rev. Lett. 1994, 73 (7), 975− 978. (80) Solveyra, E. G.; de la Llave, E.; Scherlis, D. A.; Molinero, V. Melting and Crystallization of Ice in Partially Filled Nanopores. J. Phys. Chem. B 2011, 115 (48), 14196−14204. (81) Choi, S.; Jang, E.; Kim, J. S. In-Layer Stacking Competition During Ice Growth. J. Chem. Phys. 2014, 140 (1), 014701. (82) Benet, J.; Macdowell, L. G.; Sanz, E. A Study of the Ice−Water Interface Using the Tip4p/2005 Water Model. Phys. Chem. Chem. Phys. 2014, 16, 22159−22166. (83) Seo, M.; Jang, E.; Kim, K.; Choi, S.; Kim, J. S. Understanding Anisotropic Growth Behavior of Hexagonal Ice on a Molecular Scale: A Molecular Dynamics Simulation Study. J. Chem. Phys. 2012, 137 (15), 154503. (84) Rozmanov, D.; Kusalik, P. G. Anisotropy in the Crystal Growth of Hexagonal Ice, I(H). J. Chem. Phys. 2012, 137 (9), 094702. (85) Sarupria, S.; Debenedetti, P. G. Homogeneous Nucleation of Methane Hydrate in Microsecond Molecular Dynamics Simulations. J. Phys. Chem. Lett. 2012, 3 (20), 2942−2947. (86) Walsh, M.; Koh, C.; Sloan, E.; Sum, A.; Wu, D. Microsecond Simulations of Spontaneous Methane Hydrate Nucleation and Growth. Science 2009, 326 (5956), 1095. (87) Vatamanu, J.; Kusalik, P. G. Observation of Two-Step Nucleation in Methane Hydrates. Phys. Chem. Chem. Phys. 2010, 12 (45), 15065−15072.

H


Identification of Clathrate Hydrates, Hexagonal Ice, Cubic Ice, and

Recommend Documents