Formation of Methane Hydrate in the Presence of Natural and Synthetic Nanoparticles

Natural gas hydrates occur widely on the ocean-bed and in permafrost regions, and have potential as an untapped energy resource. Their formation and growth, however, poses major problems for the energy sector due to their tendency to block oil and gas pipelines, whereas their melting is viewed as a potential contributor to climate change. Although recent advances have been made in understanding bulk methane hydrate formation, the effect of impurity particles, which are always present under conditions relevant to industry and the environment, remains an open question. Here we present results from neutron scattering experiments and molecular dynamics simulations that show that the formation of methane hydrate is insensitive to the addition of a wide range of impurity particles. Our analysis shows that this is due to the different chemical natures of methane and water, with methane generally excluded from the volume surrounding the nanoparticles. This has important consequences for our understanding of the mechanism of hydrate nucleation and the design of new inhibitor molecules.


Neutron scattering experiments
Our neutron scattering experiments were conducted on the NIMROD 1 and SANDALS 2 time-of-flight neutron diffractometers at the ISIS pulsed neutron source, STFC Rutherford Appleton Laboratory (Didcot, UK). These instruments are optimized for studies of liquids and amorphous materials containing a high proportion of hydrogen ( 1 H), and they provide continuous access to a momentum transfer range 0.02 < Q < 50Å −1 where; and θ is the scattering angle, λ the wavelength, and d the spacing between the relevant crystal planes.
The hydrate samples were prepared in situ on the beam-line in a cylindrical geometry null scattering titanium/zirconium alloy pressure cell, of the type originally used by Buchanan et al. 3 This cell has height 40 mm, inner diameter 15 mm and wall thickness 3 mm (see Fig. S1). The bottom of this cell has a dead-volume that contains a 10 mm steel ball bearing, and the entire cell system can be inverted with a frequency of ca. 0.5 Hz in the neutron beam to allow mixing of the sample and the pressurizing methane gas. Temperature was controlled to within ±0.05 o C via a circulating water-glycol heat bath. For each experiment the sealed/evacuated sample cell was first loaded with 8.3 cm 3 of liquid (D 2 O or D 2 O + clay/silica) via a bleed-in pipe at the cell base. Pressurized methane (CH 4 ) was then introduced over the liquid using a pressurized gas hand pump. Standard working conditions were 180 bar methane and 278 K. At this working pressure the sI methane hydrate is stable below 293 K. 3 Methane pressure was maintained during the experiment by top-up from the hand-pump.
To ensure thorough and reproducible mixing of the methane gas and solution, we em-S3 ployed two agitation regimes, referred to as 'standard' and 'short'. The first stage in both regimes was 'pre-production', which consisted of 15 min data collection, followed by 15 min shaking, then another 15 min data collection, before cooling over a 30 min interval from 298 K to 278 K. In the standard regime, we then performed the following: 'stage 0', 15  The 'short agitation time', t a,sh , is defined as the time immediately after stage 1. Both of these agitation regimes are shown schematically in Fig. S1.
For data correction and calibration, neutron scattering patterns were measured from each sample and also collected from: (i) the empty instrument, (ii) the empty instrument with the sample insert (but no sample container), (iii) the empty sample container, and; (iv) a 10 mm diameter cylindrical rod of vanadium that acts as a known standard scatterer. 4 Background, multiple scattering, absorption, and normalization correction procedures were implemented by the in-house "Gudrun" suite of programs. 4 After corrections and normalization, the quantity that is obtained from our neutron diffraction experiment is the differential cross-section (DCS), which is defined as follows: 4 where c α and σ α = 4πb 2 α are, respectively, the atomic fraction and scattering cross-section of species α, and N is the number of atoms per unit volume. F (Q) is the interference function S4 (INT), which contains only information about the structure within the sample: where b α is the coherent scattering length of species α, S αβ (Q) is the partial structure factor, and the sum is over all different pairs of particle types. Table S1 contains values of σ α and b α for relevant nuclides. Fig. S2 illustrates the DCS and INT functions for the control CH 4 -D 2 O before and after agitation. The DCS function oscillates about the mean scattering level, α c α b 2 α , which is proportional to the average neutron scattering cross-section per atom. By reference to Table S1 we see therefore that, due to the very strong neutron scattering of hydrogen (  to broaden further, but this effect was not observed in our data. As mentioned in the main article, all of the systems are low-viscosity liquids, with the exception of 2 wt% Laponite R B which is a strong gel former, and 2 wt% Laponite R RD which forms a thixotropic gel over ca. 6 hours. In the case of Laponite R RD and Laponite R B, increasing the concentration of nanoparticles to 2 wt % leads to the formation of a gel phase, in which the clay particles are fully dispersed due to interparticle electrostatic double-layer repulsion. 5 For Laponite R RD this process takes place over several hours, and so we were able to compare hydrate formation in the liquid and gel phases, whereas for Laponite R B the gel forms over minutes. Fig. S4 shows that the presence of a gel phase strongly inhibits the formation of hydrate. Indeed, in the case of 2 wt % Laponite R B, only trace amounts of hydrate were formed even after S6 very extended agitation (4 hours of total agitation). This observation is consistent with an excluded volume for hydrate around each clay particle, which is maximized in the gel-phase due to the more uniform distribution of particles than in the liquid. We note from Table 1 that the average particle-particle separation is around 25 nm at this concentration, which would imply that methane hydrate formation is inhibited up to a distance of at least 10 nm from each particle.

Neutron scattering samples
The aim of this study is to understand the effects on methane hydrate formation of the presence of clay and silica particles of varying hydrophilicity/hydrophobicity and at varying concentrations. The samples studied are summarized in Table 1.
The clays we used are in the 2:1 family, for which the end members are talc and pyrophyllite (uncharged, hydrophobic) and mica (highly charged, hydrophilic). To ensure dispersion of the clays we typically prepared them in sodium substituted form, with the exception of the high charge vermiculite which was prepared with propylammonium (C 3 H 9 NH + 3 ).
Two synthetic clays were studied, Laponite R RD hectorite with dry formula unit: 6 and Laponite R B fluorohectorite with dry formula unit: These synthetic clays have disc-like particles of thickness 0.92 nm and average diameter 25 nm. The surface charge density of these materials is −0.12 C m −2 giving ζ-potentials of around 40 mV.
We also studied two natural clays. Wyoming montmorillonite (SWy-2) was prepared in S7 its sodium exchanged form to give dry formula unit: The layer charge density is 0.10 C m −2 and ζ-potential around −37 mV. This natural clay has greater particle-to-particle compositional variations than the synthetic Laponites, including the possibility of more hydrophobic regions. Secondly of the natural clays we prepared high charge Eucatex vermiculite with interlayer propylammonium ions of composition: This material occurs naturally as macroscopic flakes, and it was exfoliated by exchanging with propylammonium counterions, dried and ground in a rotary mill to give typical particle All measurements were made in heavy water D 2 O as this provides a strong coherent signal and avoids the high background resulting from incoherent scattering from H 2 O. In the case 2 wt % Laponite R B we conducted agitation over around 4 hours in total, with data collection of around an hour after each step. The lack of hydrate Bragg peaks after this extended regime pointed strongly towards hydrate inhibition around the solid particles. In the case of 2 wt % Laponite R RD we conducted two complete cycles of agitation. In cycle #1 we observed Bragg peak formation after the usual 15 minutes final agitation. The hydrate crystals were then melted at 298 K. After this aging, allowing a gel to form, the process was repeated in cycle #2. As shown in Fig. S4, this second cycle showed much reduced Bragg S8 peak intensities.
In addition, to examine the effect of possible impurities and cell surfaces, particularly after initial mixing, we also conducted experiments on the following samples, using the 'short' agitation regime: control (no additives), 0.5 wt % Laponite R RD, excess 360 mesh Fe powder (Sigma Aldrich), 0.5 wt % graphene oxide (prepared by standard Hummer's method), and 5 wt % C 12 E 6 hexaethylene glycol monododecyl ether surfactant (Sigma). Results are shown in Fig. S5 and show no positive effect on the growth or rate of hydrate formation.

Molecular dynamics simulations
Coarse grained MD simulations contained 6846 water and 1154 methane molecules. Water was modeled using the mW model, 7 which represents a water molecule as a single site that interacts with its neighbors through both two-and three-body potentials. Despite its simplicity, mW gives a good description of water's structural and thermodynamic quantities, 7 and has been successfully applied in studies of heterogeneous ice nucleation (see e.g. Refs. [8][9][10][11][12]. It has also been used in studies of homogeneous hydrate nucleation. 13,14 Methane-methane and methane-water interactions were described by the potential outlined by Jacobson and surface hydrophilicity, a 100 ns isothermal-isobaric simulation was performed at 250 K and 900 atm, which resulted in a phase separated mixture of methane and water. To quantify the amount of dissolved methane, we considered the largest cluster of connected methane molecules (using a neighbor cutoff of 0.5 nm) to be the gas phase. 18 Methane molecules not part of this cluster were considered to be in solution. Typical initial water-to-methane ratios for the fcc surface varied from 25:1 to 41:1, significantly higher than 6:1 that one obtains S10 by only considering the total number of methane molecules in the simulation. Five initial configurations were selected from the second half of this trajectory, and the velocities were randomized according to the Maxwell-Boltzmann distribution for a temperature of 250 K, and the target temperature of the thermostat was decreased at a rate of 0.1 K/ns. The equations of motion were integrated until nucleation was observed, which took on the order of 100 ns. We note that although these conditions result in a high driving force for nucleation, similar protocols have been used previously to successfully investigate heterogeneous ice nucleation. [8][9][10][11][12] This suggests that if the surfaces were to act as catalysts for methane hydrate formation, then this would be observed with the simulation techniques used in this study.
We have also conducted simulations with an all-atom model of methane, water and the clay mineral kaolinite. The force field of Tse, Klein and McDonald was used to describe the methane, 19 TIP4P/ICE was used for the water 20 and the CLAYFF potential was used for the kaolinite. 21 The total number of water and methane molecules was 2944 and 512, respectively, and the surface area of the kaolinite was approximately 4.6×5.4 nm 2 . The water and methane were constrained using the SHAKE algorithm. 22 Particle-mesh Ewald 23,24 with a grid-spacing of 0.1 nm and an interpolation order of 4 was used to treat the long-range electrostatics.
A cut-off of 0.9 nm was used for the real-space interactions. Dynamics were propagated with a leap-frog integrator as implemented in the GROMACS 4.5.5 simulation package 25 with a time step of 2 fs. The system was prepared by melting a hydrate crystal at 425 K and 400 bar (394.8 atm) for 20 ns, resulting in a phase separated system with the methane at the silicate terminated face, and a planar interface separating the methane and water. Initial configurations were then drawn from this trajectory, with the velocities randomized with a target temperature of 245 K. The temperature was maintained with a Nosè-Hoover chain (length 10) and pressure of 500 bar was maintained using a Parrinello-Rahman barostat, allowing the direction normal to the kaolinite surface to fluctuate. S11

cleation occurs
The gray shaded regions in Fig. 5 provide estimates for the range of heights above the surface at which nucleation was observed. To provide this estimate, we first identified which water molecules were 'hydrate-like' according to the CHILL+ algorithm. 26  In Table S2 we show the measured non-equilibrium freezing temperatures T f (temperatures at which nucleation of gas hydrate was observed). Especially at the fcc surfaces, it is clear that there is some variability in T f , and it appears that as the surface hydrophilicity S12 increases, so too does T f . While one might interpret this as evidence of the surface promoting nucleation, inspection of Figs. 5, S7 and S9 shows that this increase in T f appears to be correlated with an increase in the concentration of dissolved methane. The number densities of methane in the bulk regions are also given in Table S2. In fact, as the hydrophilicity of the surface increases, there is a driving force for more water to bind to the surface; this causes the curvature of the methane-water interface to increase, which in turn increases the amount of dissolved methane in solution (due to Laplace pressure effects). For an in-depth discussion of the effects of methane-water interfacial curvature on hydrate nucleation see Ref. 18. Due to the finite size of our simulations, it is difficult to evaluate how the surfaces investigated would affect the dissolved methane concentration in the thermodynamic limit.
We can conclude from this analysis, however, that the presence of particles can affect the rate of hydrate nucleation indirectly, if they have a significant effect on uptake of methane into solution.   Figure S1: Schematic of the 'standard' and 'short' agitation regimes used in the neutron scattering experiments. The 'pre-production' stage consisted of: 15 min data collection, then 15 min shaking, followed by another 15 min data collection. This was done at a temperature of 298 K. The sample was then cooled to 278 K over 30 min. We then followed either the 'standard' or 'short' agitation regime, as outlined in the schematic. We define t = 0 immediately after stage 0 in both regimes. The "standard agitation time", t a,st , is the time elapsed after stage 4 in the standard regime. Similarly, the "short agitation time", t a,sh is the time elapsed after stage 1 in the short regime. After the agitation time, data was collected in 15 min intervals. The pressure was 180 bar CH 4 throughout, for both agitation regimes. The inset (bottom right) shows a picture of the shaker cell used, along with a schematic of the rocking motion. In the case of Laponite R RD, hydrate forms after the first cycle of agitation in the liquid phase, but is strongly inhibited after the second cycle of agitation and gel formation. In the case of the Laponite R B fluoro-clay sample, gel formation is rapid and hydrate formation is strongly inhibited even after more than 5 hours in the first cycle (the sample was agitated for 4 hours in total). These data are strong evidence for an excluded volume to hydrate formation around each of the clay particles.