Mechanical Properties of Methane Hydrate - American Chemical Society

Oct 10, 2018 - Faculty of Engineering, China University of Geosciences, Wuhan, Hubei ... NTNU Nanomechanical Lab, Norwegian University of Science and ...
0 downloads 0 Views 4MB Size
Subscriber access provided by University of Sunderland

C: Energy Conversion and Storage; Energy and Charge Transport

Mechanical Properties of Methane Hydrate: Intrinsic Differences from Ice Pinqiang Cao, Jianyang Wu, Zhisen Zhang, Bin Fang, and Fulong Ning J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.8b06002 • Publication Date (Web): 10 Oct 2018 Downloaded from http://pubs.acs.org on October 13, 2018

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 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 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.

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 36 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

Mechanical Properties of Methane Hydrate: Intrinsic Differences from Ice Pinqiang Cao,1, 2 Jianyang Wu,2, 3,* Zhisen Zhang,2 Bin Fang,1 and Fulong Ning,1, 4,* 1

2

Faculty of Engineering, China University of Geosciences, Wuhan, Hubei 430074, PR China

Department of Physics, Research Institute for Biomimetics and Soft Matter, Jiujiang Research

Institute and Fujian Provincial Key Laboratory for Soft Functional Materials Research, Xiamen University, Xiamen 361005, PR China 3

4

NTNU Nanomechanical Lab, Norwegian University of Science and Technology (NTNU), Trondheim 7491, Norway

Laboratory for Marine Mineral Resources, Qingdao National Laboratory for Marine Science and Technology, Qingdao 266237, China

Abstract: Understanding fundamental mechanical behaviors of ice-like crystals is of importance in many engineering aspects. Herein, mechanical characteristics of monocrystalline methane hydrate (MMH) and hexagonal ice (Ih) under mechanical loads are contrasted by atomistic simulations. Effects of engineering strain rate, temperature, crystal-orientation and occupancy of guest molecule on the mechanical properties of MMH are investigated. Results show that engineering strain rate, temperature and occupancy of guest molecules in 51262 cages greatly affect the mechanical strength and failure strain of MMH, whereas the effect of crystal-orientation on the tensile response of MMH such as along the

[1 0 0] and [1 1 0] directions is negligible. Particularly, occupancy of guest molecule in 51262 cage primarily governs the mechanical strength and elastic limits of MMH. For Ih, it is tensile stiffer than that of MMH at 263.15 K and 10 MPa, and shows unique mechanical characteristics such as tension-induced stiffening and compression-induced remarkable softening under the [0 0 0 1] directional load. Both

*

Email: [email protected]; [email protected] 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

Page 2 of 36

crystals demonstrate brittle fracture behavior but different plasticity with dislocation-free in MMH yet dislocation activities in Ih. The intrinsic differences in mechanical properties of MMH and monocrystalline Ih mainly result from the host-guest molecule interactions and relative angles which tetrahedral hydrogen bonds make to the loading direction. These mechanical characteristics present microscopic insights to understand the mechanical responses of natural occurring and artificial synthetic gas hydrates.

2

ACS Paragon Plus Environment

Page 3 of 36 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

1 Introduction Natural gas hydrates are ice-like crystalline compounds, in which natural gas as guest molecules, e.g., methane, ethane, carbon dioxide, and so on, are physically encapsulated within polyhedral water cages of hydrogen-bonding framework

1-5

. The guest molecules play a significant role in the structural and

thermodynamic stability of latticed gas hydrates. Natural gas hydrates occur extensively in various environmental settings, such as the earth’s deep ocean floor continental margins,6-7 permafrost regions2, 8-9

, and petrochemical production lines10. Methane hydrates are the most common and abundant in

nature. Because of huge global amount of carbon embedded within hydrates, they are regarded as a huge untapped future energy resource

7, 11-14

. On the other hand, natural gas hydrates play a key role in

geo-mechanical stability of hydrate-bearing sediments, which may have implications for geological hazards15-17, climate change16-19. Therefore, a large body of both experimental and theoretical works has been performed to investigate the physico-chemical properties of hydrate-bearing sediments in the past few decades

15, 20-23

. However, due to the complexity in hydrate occurrence habits, compositions and

sedimentary heterogeneity, mechanical properties of hydrate-bearing sediments are still poorly understood. Above of all, the mechanical behaviors of pure gas hydrates are important data for the assessment of the mechanical stability of hydrate-bearing sediments. To date, a number of investigations on elastic properties of pure gas hydrates have been performed experimentally and theoretically. Using in-situ high-pressure Brillouin spectroscopy in a diamond-anvil cell at 296 K, Shimizu et al.24 reported that elastic modulus of monocrystalline methane hydrate (MMH) ranges 3.5 to 17 GPa, and their results indicated that MMH is slightly more compressible than hexagonal 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

Page 4 of 36

ice (Ih). Using the first-principles calculations, Jendi et al.25 showed that Young’s modulus of MMH is around 16.0- 17.8 GPa, and its uniaxial tensile and compressive strength ranges from 1.41-1.70 GPa and 3.01-4.00 GPa for different crystalline directions, respectively. Jia et al.26 via molecular dynamics (MD) simulations reported that the maximum strength of MMH and Ih are approximately 1.20 and 0.8 GPa at 250 K and 40 MPa, respectively. Limited by preparation of high-quality gas hydrate samples in the laboratory and recovered from natural settings, mechanical properties of pure gas hydrates, especially the response mechanism under mechanical loads, are still largely insufficient. In this work, mechanical properties of MMH are comprehensively investigated by means of MD simulations, including the effects of strain rate, temperature and occupancy of guest molecules on the mechanical properties of MMH. Specifically, because appearance of MMH is analogous to ice, mechanical properties of Ih are also studied for comparison.

2 Methodology MMH is cubic structure I (SI) hydrate in this study. For molecular model of MMH structure, the prime positions of the oxygen atoms of water molecules were taken from the X-ray diffraction analysis data27. Methane guest molecules were placed at the center of the water cages. SI MMH is composed of basic repeating two types of water cages, namely small CH4@512 and large CH4@51262, as displayed in Figures 1a-1c. MMH with an 8×8×8 surpercell having dimensions of 96.24 × 96.24 × 96.24 Å3 was generated for mechanical tests (Figure 1g). To reveal the effect of crystal-orientation on the mechanical

4

ACS Paragon Plus Environment

Page 5 of 36 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

behaviors, cubic MMH was intersected to achieve other orthogonal structures. The structures with

[1 1 0] and [1 1 1] orthogonal directions possess the dimensions of 104.58 × 98.60 × 98.60 Å3 and 98.60 × 100.0 × 115.20 Å3, respectively. To uncover the effect of occupancy of guest molecules on the mechanical properties of MMH, the guest molecules in cages of 512, 51262, and 512&51262 were removed, respectively. With regards to ice, initial atomic structure of monocrystalline Ih was obtained by the X-ray diffraction analysis data28. One unit cell of Ih is composed of 16 water molecules. Figures 1d-1f present the three geometrical configurations of Ih with different crystal planes which are yellow-highlighted. The molecular model of Ih with 13 ×14 × 11 unit cells containing 32,032 water molecules was constructed for mechanical tests (Figure 1h). The dimensions of molecular model of Ih are 101.66 × 103.04 × 99.33 Å3, which is close to that of molecular model of MMH. Similarly, effect of crystalline orientation on the mechanical properties of Ih was considered. A monatomic model was adopted to describe both water and methane molecules in MMH and water molecules in Ih. Each molecule is represented as a single sphere in the monatomic model. The coarsening level of methane molecule is identical to the OPLS united-atom (OPLSUA) methane model29. The tetrahedral short-ranged interactions between monatomic water and water-methane and methane-methane30-31 are mimicked by the available Stillinger-Weber (SW) potential as follows:32.

E  2 (rij )   3 (rij , rik ,ijk ) i

j i

i

 2 (rij )  A [ B( )4  1]exp( rij

(1)

j i k  j

 rij  a

5

ACS Paragon Plus Environment

)

(2)

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

3 (rij , rik ,ijk )   [cos ijk  cos 0 ]2 exp(

 rij  a

Page 6 of 36

) exp(

 rij  a

)

(3)

where rij and rik are distances between particles i and j, and between particles i and k, respectively. θijk is the angle subtended by the position vectors of rij and rik . The constant A, B, γ, a and θ0 parameters are taken as 7.049556277, 0.6022245584, 1.2, 1.8, and 109.5º, respectively. The value of the parameter λ determines the strength of the tetrahedral interactions. Other parameter details can be found in Table 1. Such coarse-grained model30-31 of water shows over 2-orders of magnitude more efficient than fully atomistic models in reproducing a range of properties of liquid- and solid-stated water30. Additionally, it also has been successfully utilized to predict the mechanical properties of methane hydrates 33.

Prior to the mechanical loads, molecular structures were quasi-statically optimized to a local minimum energy configuration using conjugate gradient method with an energy tolerance of 1.0× 0-4eV and a force tolerance of 1.0×10-4 eV/Å34. Relaxations were conducted with 1.0×106 timesteps under NPT (constant number of particles, constant pressure, and constant temperature) ensemble, and then with another 1.0×106 timesteps under NVT (constant number of particles, constant volume, and constant temperature) ensemble. The Nosé-Hoover thermostat and barostat techniques with damping times of 0.1 and 1 ps were utilized to control the temperature and pressure, respectively. The velocity-verlet integration algorithm was adopted to integrate the equations with timestep of 1 fs. Periodic boundary conditions (PBC) were applied in three orthogonal directions. In this study, the definition of tensile (compressive) strain is similar to the engineering strain that is defined as: 6

ACS Paragon Plus Environment

Page 7 of 36 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



l  l0 l0



l l0

(4)

where ε, l0 and l are the engineering strain, the original length and the final length of simulation box along the deformational direction, respectively. Young’s modulus is determined from the slope in the linear elastic regime of the loading curve as:

E

 ( ) 

(5)

where E, σ and ε are the Young's modulus, the loading stress and the engineering strain, respectively. Here, the Young’s modulus of MMH is calculated by linearly fitting the mechanical stress-strain curves in the elastic regime of 0-0.05. The fracture strength of MMH is defined by the maximum stress in the loading curve, and a specimen fails by dissociation of cages at the maximum stress. With regard to mechanical tests, uniaxial loads were achieved by using deformation control method. The procedures were carried out on the relaxed specimens with a constant strain rate of 1.0×108/s by uniformly rescaling particle-coordinates in every 1000 timesteps. MD loading simulations were conducted under a modified NPT ensemble, specifically, NVT in the loading direction yet NPT in the lateral directions. The deformational MD steps were carried out in NLzpxpyT ensemble, where Lz represents the length of simulation box along the loading direction. Nosé-Hoover anisotropic barostat and thermostat methods were adopted to maintain pressure only in lateral directions independently and temperature, respectively. Such simulation setup allows the Poisson effect during the mechanical loading. Stress of particles was collected on the basis of definition of virial stress during the MD runs. Both stress 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 36

and potential energy of particles in the systems were averaged every 5000 timesteps to remove the thermal fluctuations. All the MD calculations were carried out using the Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) software package35. 3 Results and discussion 3.1 Effect of engineering strain rate The mechanical stress-strain curves of MMH at 263.15 K and confining pressure of 10 MPa under three different strain rates along the [1 0 0] directional deformation are plotted in Figures 2a and 2b. Apparently, the mechanical strength and failure strain of MMH are sensitive to the applied engineering strain rates. Upon tension, MMH show fracture strengths of approximately 0.94, 0.90 and 0.85 GPa as the applied strain rate is 1.0×109/s, 1.0×108/s and 1.0×107/s, respectively. Similarly, the fracture strains also increases with increasing strain rate, with the values of around 0.145, 0.134, and 0.125 for the strain rates of 1.0×109/s, 1.0×108/s and 1.0×107/s, respectively. The relationship between fracture strength and strain rate can be expressed by an Arrhenius relation36. 1 m

  A exp(

Q ) RT

(6)

where  and  are the engineering strain rate and the fracture strength, respectively. Q , R, T and m are the activation energy, the universal gas constant, the deformation temperature and sensitivity of the strain-rate, respectively. m and A are constant parameters determined from the method of least squares. Both tensile and compressive strengths are well-described by the above equation. As an example, the tension simulation data yields the following relation:

8

ACS Paragon Plus Environment

Page 9 of 36 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

ln( )  ln(exp(35.316)) 

1 78300 ln( )   0.0207  ln( )  0.493 48.309 8.314  263

(7)

where the units of  ,  and T are s-1, GPa and K, respectively. Q and R are chosen as 78300 J·mol-1 37

and 8.314 J·mol-1·K-1, respectively. Prior to fracture, pronounced strain-induced softening behavior is

identified for compressive loads as indicated in Figure 2b, in contrast to the case of tension loads. Moreover, MMH along the [1 0 0] direction exhibits higher tensile strengths than the compressive ones. For example, for a given strain rate of 1.0×109/s, it yields tensile strength of 0.94 GPa under tension, whereas under compression the fracture strength is approximately 0.75 GPa. As is known, one unit-cell of MMH is composed of 2·512+6·51262 cages formed by water molecules, and the cages are occupied by guest molecules. Both 512 and 51262 cages are non-spherical but faceted morphology. Moreover, the anharmonic responses of local hydrogen bonds and bond angles may occur under mechanical loads. The asymmetry of tension-compression mechanical responses for the faceted molecular cages of 512 and 51262 is mainly responsible for the different mechanical properties such as the fracture strength, failure strain for tension and compression loads. However, Young’s modulus is insensitive to the strain rate and loading mode. Young’s moduli by linearly fitting the loading curves in the range of ε ≤ 0.05 are around 7.78 GPa and 7.42 GPa for tension and compression, respectively, in agreement with the previous studies

24, 26, 33

, yet lower than those from first-principles25 at zero K. According to our results obtained

from the MD simulations, no significant differences in the elastic regions are observed from Figures 2a-2b with the three strain rates. Owing to length and time scale limitations in MD simulations, the strain rate of 1.0×108/s was selected to investigate the mechanical behaviors of many nanomaterials38-42. Therefore, the considerably engineering strain rate of 1.0×108/s was employed to obtain the reasonable 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

Page 10 of 36

mechanical behaviors of MMH and Ih and saved the computational resource in our investigations. To further reveal the behaviors observed from Figures 2a and 2b, a series of snapshots of molecular structures subjected to both tension and compression at a specific strain rate of 1.0×108/s were captured and shown in Figures 2c and 2d. Prior to fracture, both tensile and compressive elasticities are structurally characterized by misshaping the 512 and 51262 cages formed by water molecules, as indicated by snapshot #1 and #2 in Figures 2c and d. Similar morphology of misshaped water cages at the same elastic strain is identified under both tension and compression. This explains the similar values of tensile and compressive Young’s moduli of MMH. However, the snapshots in Figures 2c-#2 and 2d-#2 reveal that both faceted molecular cages of 512 and 51262 under compression is able to be more deformed than that under tension, illustrating the higher compressive failure strain. The high compressive strain and the anharmonic deformational response of faceted cages mainly explain the pronounced strain-induced softening behavior of MMH subjected to compression. Snapshots #3 and #4 in Figures 2c and 2d illustrate the failure patterns of MMH subjected to tension and compression loads, respectively. Upon tension, MMH fails by direct cleavage along the (1 0 0) crystallographic plane that is perpendicular to the tensile direction, resulting in release of methane molecules. This type of fracture pattern results from that part of hydrogen bonds in hexagonal ring of 51262 cage that form zero angle to the tensile direction are uniformly dissociated. Under compression, however, structure fails by fracturing along the 1 2 0  crystallographic plane that makes finite angle to the compressive direction, as indicated by the purple-highlighted structure, immediately followed by crystal-rotation to relieve strain energy. This results in formation of fresh grain boundary (GB). Differing from the case of tension, 10

ACS Paragon Plus Environment

Page 11 of 36 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

compression-induced fracture mainly results from the larger hydrogen-bonded angular deformation. Posterior to the failures, misshapen cage structures are fully recovered.

3.2 Effect of temperature A series of MD simulations were carried out to investigate the effect of temperature on the mechanical properties of MMH along the [1 0 0] direction at 10 MPa under both tension and compression. Figures 3a

and

3b

show

the

resulting

stress-strain

curves.

Obviously,

MMH

exhibits

greatly

temperature-dependent mechanical responses. MMH yields an apparent strain-induced softening phenomenon for all measured temperatures under compression, in contrast to those subjected to tension. Remarkably, at lower temperature conditions, MMH shows both higher yield stresses and strains. Young’s moduli determined from the stress-strain curves at strains ranging from 0 to 0.05 are plotted in Figures 3c-3d. Apparently, Young’s modulus of MMH increases as the temperature decreases under both tension and compression. As an example, tensile Young’s modulus of MMH varies from 7.23 to 8.28 GPa with temperature decreasing from 283.15 to 203.15 K. This temperature-induced softening of MMH can be explained by strong thermal vibrations of molecules at high temperature. Such strong thermal vibrations of molecules lead to mechanical weakening in the single-crystal. Also, due to the thermal agitation, the hydrogen bonds more easily reach the critical molecular distance and dissociate. For a given temperature, tensile Young’s modulus is slightly higher than compressive Young’s modulus. Our MD findings are consistent with the previous experimental and theoretical results

24, 26, 33

slightly lower than those reported by Jendi et al. using density functional theory (DFT)

11

ACS Paragon Plus Environment

25

, but

. The

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 36

development of molecular structures reveals that temperature dependence of failure pattern under both tension and compression is negligible.

3.3 Effect of crystal-orientation Uniaxial deformations were imposed along the [1 0 0] , [1 1 0] , and [1 1 1] directions to study effect of crystal-orientation on the mechanical properties of MMH at 263.15 K and 10 MPa by means of MD simulations. Figure 4 shows the resulting mechanical stress-strain curves. Upon tension, it is observed from Figure 4a that MMH shows almost identical mechanical responses under [1 0 0] and [1 1 0] directional loads. However, robust mechanical responses of MMH subjected to [1 1 1] direction are identified. With regard to the cases under compression as shown in Figure 4b, the mechanical behaviors of MMH present pronounced loading direction-dependence. In sharp contrast to the cases under tension, there are significant differences in the mechanical responses between under [1 0 0] and [1 1 0] directional loads. In terms of yield strength, they can be sorted as [1 1 1] > [1 1 0] > [1 0 0] . Intriguingly, single-crystal, under [1 0 0] compression, shows obvious strain-induced softening behavior prior to fracturing, in contrast to those under other two directional loads. Young’s moduli of MMH obtained by fitting the elastic regime in the range of ε ≤ 0.05 from the mechanical curves are also calculated. Subjected to elongation along the three directions, MMH shows almost identical Young’s modulus of around 7.84 GPa. On the contrary, upon compression, Young’s moduli of MMH obtained from [1 0 0] , [1 1 0] and [1 1 1] directional mechanical curves are around 7.42, 7.80 and 7.98 GPa, respectively. The crystal-orientation dependence of mechanical properties of MMH originates from

12

ACS Paragon Plus Environment

Page 13 of 36 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

directional dependence of deformation behavior of unique geometrical 512, 51262 water cages within MMH. Figures 5a-5c show the side-viewed configurations of different oriented cages as the loads are imposed along the [1 0 0] , [1 1 0] and [1 1 1] directions, respectively. Apparently, differences in combined orientations of cages within MMH between the three directional loads are responsible for the direction-dependent mechanical properties. Notably, a flatwise compressive load along the [1 0 0] direction that is perpendicular to the hexagonal rings of polyhedral 51262 water cages results in the remarkable strain-induced softening (see Figure 5a-3). However, slanted polyhedral 512, 51262 water cages subjected to the [1 1 1] directional compression loads exhibit higher stiffness and higher strength than those under other directional loadings. That can be attributed to the configurations of hexagonal rings from the loading directional views. We also found that the effect of crystal orientation on the tensile response of MMH such as along the [1 0 0] and [1 1 0] directions is negligible, but not for the compressive response. This is primarily attributed to the asymmetry of tension-compression mechanical responses in the locally molecular structures of the crystal. As is known, MMH is structurally characterized by the host water cage-based framework connected by the hydrogen bonds between water molecules. Upon different directional loads, there exist different variations in the local hydrogen bond length and hydrogen bond angle at different strains. For example, the localized hydrogen bonds making small angles to the tensile-loading direction are significantly stretched with increasing tensile strain, whereas subjected to compression they are contracted with increasing strain. And the particle interactions are described by the anharmonic many-body Stillinger-Weber potential30-31. Therefore, different mechanical responses in locally deformed molecular cages are expected. This mainly explains 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 36

different effects of crystal orientation between the tensile and compressive loads. The factor of the bonds making finite angles to the loading direction plays a vital role in the mechanical properties of crystalline materials43-45. For example, diamond exhibits the weakest strength under tensile load along the [1 1 1] direction, originating from that the C-C bonds are parallel to the constraint43. The strongly directional hydrogen bonds in hydrate mainly produce highly anisotropic stress response.

3.4 Effect of occupancy of guest molecules in polyhedral cages SI MMH is composed of 512 and 51262 cages, and those cages are stabilized by guest molecules. Commonly, polyhedral cages in naturally occurring methane hydrates are not fully occupied. In this section, we thus try to figure out the effect of occupancy of the guest molecules on the mechanical responses of MMH. Actually, it also reflects the vital effect of host-guest interaction on the mechanical response of MMH. Four extreme occupancies of guest molecules in the polyhedral cages of MMH are taken into account. (a) Both 512 and 51262 water cages in MMH are fully occupied by methane molecules, named as pristine structure. (b) 512 cage is empty, while the 51262 cage is fully occupied by methane molecules, named as CH4-free@512 MMH. (c) 512 cage is fully occupied by methane molecules, while the 51262 cage is empty, named as CH4-free@51262 MMH. (d) Both 512 and 51262 water cages are empty, named as CH4-free@512&CH4-free@51262 MMH. Figure 6 displays the resulting stress-strain curves of MMH subjected to tension and compression along the [1 0 0] direction. Obviously, methane-free in polyhedral cages has an influence on the mechanical properties of MMH. For both tensile and compressive loads (Figures 6a-6b), CH4-free@51262 greatly weakens the hydrate structure. Mechanical

14

ACS Paragon Plus Environment

Page 15 of 36 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

strengths of MMH decrease by approximately 18% and 34%, corresponding to degradations from around 900 MPa to 736 MPa and from 743 MPa to 493 MPa under tension and compression, respectively. Moreover, failure strains of MMH also exhibit reductions of approximately 11% and 33%, corresponding to reductions from around 0.135 to 0.12 and from 0.15 to 0.10, under tension and compression, respectively. By comparison, CH4-free@51262 MMH exhibits a lower reduction in both failure strength and failure strain under tension than those under compression. The Young’s modulus of MMH under tension decreases monotonically from 7.78 to 7.05 GPa in the order of pristine structure, CH4-free@512 MMH, CH4-free@51262 MMH, and CH4-free@512&CH4-free@51262 MMH, respectively. Also, the compressive Young’s modulus exhibits the similar tendency, but with lower values deceasing from 7.42 to 6.47 GPa. It can be concluded that the guest molecule-occupancy in the 51262 cages determines the mechanical strength of SI MMH. This mainly results from the fact that the ratio of small (512 cages) to large cages (51262 cages) is 2 to 6 in SI MMH. Our findings also support the fact that large 51262 cages are always occupied to about 100% in reality, while small 512 cages can have varying occupancies depending on the pressure and temperature conditions46-50. Snapshots of MMH under uniaxial straining at specific strains are captured and shown in Figure 7 to reveal the deformation mechanisms of MMH from molecular insights. Subjected to tensile loadings, CH4-free@512 MMH fails in brittle manner (Figure 7a), similar to pristine structure. With regard to the CH4-free@51262 and CH4-free@512&CH4-free@51262 MMH, similar structural evolutions are identified prior to fracture (Figures 7b-7c). Strikingly, for CH4-free@51262 MMH, guest molecules in small 512 water cages are able to diffuse into large 51262 water cages at the fracture region to form 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

Page 16 of 36

metastable crystalline structures (Figure 7b). This process involves the crystallization and reformation of hydrate structures. Such phenomenon results in the steady recovery of mechanical stress with the increasing strain in mechanical stress-strain curves (see Figure 6a). However, the crack of CH4-free@512&CH4-free@51262 MMH rapidly propagates to achieve complete rupture, leading to that water cages are collapsed to form amorphous water structure (Figure 7c). As for the case under compression (Figures 7d-7f), either CH4-free@512 or CH4-free@51262 MMH fails by forming a bicrystal with GBs (Figures 7d-7e). Subsequently, GB sliding, crystallization and reformation collectively dominate the plastic deformation at high strain levels. CH4-free@512&CH4-free@51262 MMH is directly collapsed to achieve amorphous water structure (Figure 7f), similar to the case under tension.

3.5 Mechanical properties of Ih for comparison Mechanical properties of Ih subjected to uniaxial tension and compression are also studied in comparison with MMH to reveal the nature origin of the similarities and differences in mechanical properties of both crystals. Figures 8a and 8b present the simulated stress-strain curves of Ih along the three different crystallographic directions together with aforementioned mechanical responses of fully-occupied SI MMH. By comparison, it is found that there are distinct differences in the mechanical responses between MMH and Ih. Depending on the crystalline directions and loading modes, Ih yields Young’s modulus varying from 8.9 - 9.9 GPa, which falls in the range of available experimental measurements (8.6 - 12 GPa)

51-53

, yet is slightly higher than that (around 7.2 GPa) from previous MD simulations with fully

atomistic TIP4P/ice water model at 250 K and 40 MPa

26

, but lower than that from the first-principles

16

ACS Paragon Plus Environment

Page 17 of 36 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

calculations25. Apparently, Ih is stiffer than MMH under tension at 263.15K and 10MPa. This mainly comes from the different mechanical rigidity in locally molecular structures; More localized hydrogen bond and bond angle in ice are slightly tensile stiffer than those in MMH. Furthermore, previous studies via first-principles calculations showed a higher Young’s modulus of ice than that of MMH25, and all of work started with the relaxed atomic structures with lattice vectors achieving a zero stress state25. But a lower Young’s modulus of ice than that of MMH at 250 K and 40 MPa was found in the MD results by Jia et al.26 Upon tension, Ih shows a maximum tensile strength for the [0 0 0 1] directional load, which is comparable to those of MMH. Under the [1 1 2 0] and [0 1 1 0] directional tensions, however, strengths of Ih are approximately 14% ~ 16% lower than those of MMH. Also, fracture strains of MMH are around 12.5% - 45% higher than those of Ih along the [1 1 2 0] , [0 0 0 1] , and [0 1 1 0] directions. Subjected to uniaxial compression, remarkably strain-induced softening is identified for load along the [0 0 0 1] direction, similar to that of MMH along [1 0 0] directional compressive load. Moreover, the

smoothly nonlinear reduction in compressive strength of Ih is recognized with strain varying from around 0.12-0.158, resulting in negative compressive stiffness (see Figure 8b). This strain-induced softening behavior is in sharp contrast to the case of tension load where strain-induced stiffening behavior with strain ranging from around 0.075-0.10 is observed. Subsequently, sudden rise in compressive strength occurs as compressive strain is increased from around 0.158 - 0.183. This indicates occurrence of large-scale structural transformation during the compressive process. In terms of strength, Ih subjected to the [0 0 0 1] direction load shows the strongest compressive strength of around 1.06 GPa. For the case of [0 1 1 0] directional load as shown in Figure 8b, Ih presents the weakest mechanical 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

Page 18 of 36

strength of around 0.68 GPa before the fracture occurs at 263.15 K and 10 MPa. With regard to the elastic strain limits, MMH also demonstrates higher fracture strains than Ih for [1 1 2 0] and [0 1 1 0] uniaxial loads. With regard to the case of [0 0 0 1] direction load (see Figure 8b), however, large elastic strain limits that is comparable to MMH results from the long-range strain-induced softening phenomenon. As shown in Figures 8c-8d, Ih shows temperature-dependence, in accordance with that of MMM in our simulations. Yet it exhibits pressure-independence in the confining pressure range of 0.1MPa-10MPa. Figure 9 shows a succession of localized molecular snapshots that capture the microstructural developments of Ih at 223.15 K and ambient pressure along various crystallographic directions to understand the deformational mechanisms. Within elastic deformation regime, polyhedrons made up of water molecules in Ih structure are arranged for all three directional. For [1 1 2 0] and [0 1 1 0] uniaxial loads, all the tetrahedral hydrogen bonds in Ih make finite angles to the straining directions, whereas under [0 0 0 1] load one part of tetrahedral hydrogen bonds is parallel to the loading direction. This mainly results in the different mechanical responses. Similarly, using first-principle calculations, Li et al.44 found that the strongly directional C-C bonds in cubic diamond determine highly anisotropic mechanical stress response. For example, the significant strain-induced stiffening in cubic diamond subjected to [1 1 1] directional compression results from that one part of tetrahedral C-C bonds makes zero angle to the loading direction. Upon tension, as the applied strain reaches critical values, Ih fails by crystallographic slips as visually demonstrated by the unidentified water particles (Figures 9a-9c), which is corresponding to the deep drops of stresses in the curves of Figure 8. Such failure pattern is 18

ACS Paragon Plus Environment

Page 19 of 36 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

insensitive to the crystalline loading directions and has been also identified by experimental observations54. One type of crystallographic slips under tension is identified and listed in Table 2. As a result of the crystalline slips, bicrystals with GBs form. Mismatched lattice in newly formed GBs as displayed by the broken blue-colored water particles in Figures 9a-3, 9b-3 and 9c-3 clearly suggests debonding and sliding behavior associated with the brittle fracturing in Ih. The development of dislocation activities during the plastic deformation is also observed. In contrast, dislocation-free is observed in MMH during the whole deformation. The dislocations occurring in Ih are identified by the means of the Dislocation Extraction Algorithm (DXA) method55. Posterior to the failure, dislocation motion, GB sliding, amorphization and recrystallization dominate the following plastic deformations, resulting in the oscillation of tensile loads in the stress-strain curves.  1 2 1 0 dislocation is identified in Ih for the three directional loading tests as shown in Figures 9a-9c. This type of dislocation has been also previously detected in Ih and hexagonal GaN56-57. For the cases under compression, however, two distinct deformation behaviors are observed after the first peak strength. For [1 1 2 0] and

[0 1 1 0] directional compression (Figure 9d and Figure 9f), the mechanisms of failure and plastic deformation are similar to those under tensile loadings. For [0 0 0 1] directional compression, large-scale structural transition initiates as the strain reaches critical value of around 0.13 (Figures 9e-2 and -3), explaining the remarkably strain-induced softening behavior as presented in Figure 8. Surprisingly, at a critical compression strain of around 0.16, a new solid phase of water ice forms from Ih, strikingly differing from that under tension. This new crystal lattice has not been identified by available experiments and theories. Viewing from [0 0 0 1] direction, it still shows hexagonal lattice, and here it 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

Page 20 of 36

is defined as iso-hexagonal ice (Iiso-h). Significantly uprising compressive strength followed by the strain-induced softening is attributed to compression of Iiso-h. Lastly, Iiso-h loses structural stability by crystallographic slips as the compressive strain exceeds a critical value of around 0.2. This destabilization greatly releases the strain energy, resulting in large-scale structural recovery from compressed Iiso-h. Crystalline orientation of recovered Ih slightly differs from that of original counterpart, indicating rotation of crystal during plastic deformation (see Figure 9e). In contrast to other two cases along the [1 1 2 0] and [0 1 1 0] directions, two types of crystallographic slips are identified (see Table 2). Likewise, recrystallization, dislocations and new GB sliding mainly govern the later plastic deformation. Apparently, the plastic deformation of the two crystals in mechanical behaviors is accommodated by the competing dynamic processes of crystal decomposition and reformation, agreement with the previous findings26, 33, 58-60. GBs in ice crystals are composed of five-, six- and seven-membered rings, and disordered water structures, whereas those of hydrate crystals are mainly represented by metastable pentakaidecahedral and hexadecahedral cages, and disordered water structures.

4. Conclusions In summary, classical MD simulations with coarse-grained model were carried out to determine mechanical properties of MMH and Ih, including the mechanical strength, failure strain, and Young’s modulus, under various loading conditions. MD results exhibit that Young’s modulus of MMH is insensitive to the engineering strain rate, and CH4-free@512, whereas it shows temperature, and

20

ACS Paragon Plus Environment

Page 21 of 36 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

CH4-free@51262 dependence. For Ih, tension and compression simulations show that ice crystal, under

[0 0 0 1] crystalline directional straining, exhibits superior mechanical properties than those under [1 0 1 0] and [1 1 2 0] ones. By comparison, MMH is less tensile stiffer than ice structure at 263.15 K and 10 MPa. Remarkably, under [0 0 0 1] directional compression, one new solid phase of water ice forms from hexagonal ice, resulting in surprising transition from strain-softening to strain-hardening. The plastic deformation of both crystals under compression is facilitated by GB formation, GB sliding, amorphization and recrystallization. Within plastic strain regime, ice presents dislocation structures, while dislocation-free is detected in MMH. The host-guest molecule interactions and finite angles which tetrahedral hydrogen bonds make to the loading direction can mainly be responsible for the intrinsic differences in mechanical properties of MMH and monocrystalline Ih. The present study provides molecule insight into the intrinsic mechanical properties of MMH and Ih.

Acknowledgments This work was financially supported by the National Natural Science Foundation of China (Grant No. 41672367, 11772278, 11502221 and 51274177), Qingdao National Laboratory for Marine Science and Technology Open Fund (QNLM2016ORP0203), China Geological Survey Project (DD20160216, DD20189320), the Fundamental Research Funds for the Central Universities (Xiamen University: Grant No. 20720180014, 20720180018 and 20720160088), Fujian Provincial Department of Science & Technology (2017J05028), the project sponsored by the Scientific Research Foundation for the Returned Overseas Chinese Scholars from State Education Ministry, Doctoral Fund of the Ministry of Education

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

Page 22 of 36

(20130121110018), ―111‖ Project (B16029) and the 1000 Talents Program from Xiamen University. The computational resources were provided by Information & Network Center for Computational Science at Xiamen University and High Performance Computing Center of CUG.

References 1.

Sloan, E. D. Fundamental Principles and Applications of Natural Gas Hydrates. Nature 2003, 426, 353-363.

2.

Kvenvolden, K. A. Gas Hydrates-Geological Perspective and Global Change. Rev. Geophys. 1993, 31,

173-187. 3.

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. 4.

Jacobson, L. C.; Hujo, W.; Molinero, V. Amorphous Precursors in the Nucleation of Clathrate Hydrates. J.

Am. Chem. Soc. 2010, 132, 11806-11811. 5.

Henley, H.; Thomas, E.; Lucia, A. Density and phase equilibrium for ice and structure I hydrates using the

Gibbs–Helmholtz constrained equation of state. Chem. Eng. Res. Des. 2014, 92, 2977-2991. 6.

Max, M. D.; Lowrie, A. Oceanic Methane Hydrates: A ―Frontier‖ Gas Resource. J. Pet. Geol. 2007, 19,

41-56. 7.

Makogon, Y. F.; Holditch, S. A.; Makogon, T. Y. Natural gas-hydrates — A potential energy source for the

21st Century. J. Pet. Sci. Eng. 2007, 56, 14-31. 8.

Buffett, B.; Archer, D. Global inventory of methane clathrate: sensitivity to changes in the deep ocean. Earth

Planet. Sci. Lett. 2004, 227, 185-199. 9.

Max, M. D.; Clifford, S. M. The state, potential distribution, and biological implications of methane in the

Martian crust. J. Geophys. Res.-Planets. 2000, 105, 4165-4171. 10. Hammerschmidt, E. G. Formation of Gas Hydrates in Natural Gas Transmission Lines. Ind. Eng. Chem. 1934, 26, 851-855. 11. Kvenvolden, K. A. Methane hydrate — A major reservoir of carbon in the shallow geosphere? Chem. Geol. 1988, 71, 41-51. 22

ACS Paragon Plus Environment

Page 23 of 36 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

12. MacDonald, G. J. The Future of Methane as an Energy Resource. Annu. Rev. Energy. 1990, 15, 53-83. 13. Boswell, R. Is gas hydrate energy within reach? Science 2009, 325, 957-958. 14. Klauda, J. B.; Sandler, S. I. Global Distribution of Methane Hydrate in Ocean Sediment. Energy Fuels 2005, 19, 459-470. 15. Ning, F.; Yu, Y.; Kjelstrup, S.; Vlugt, T. J. H.; Glavatskiy, K. Mechanical properties of clathrate hydrates: status and perspectives. Energy Environ. Sci. 2012, 5, 6779-6795. 16. Kvenvolden, K. A. Potential effects of gas hydrate on human welfare. Proc. Natl. Acad. Sci. 1999, 96, 3420-3426. 17. Maslin, M.; Owen, M.; Betts, R.; Day, S.; Dunkley Jones, T.; Ridgwell, A. Gas hydrates: past and future geohazard? Philos. Trans. R. Soc. A-Math. Phys. Eng. Sci. 2010, 368, 2369-2393. 18. Lambeck, K.; Esat, T. M.; Potter, E. K. Links between climate and sea levels for the past three million years. Nature 2002, 419, 199-206. 19. Xu, W.; Lowell, R. P.; Peltzer, E. T. Effect of seafloor temperature and pressure variations on methane flux from a gas hydrate layer: Comparison between current and late Paleocene climate conditions. J. Geophys. Res.-Solid Earth 2001, 106, 26413-26423. 20. Liu, W.; Zhao, J.; Luo, Y.; Song, Y.; Li, Y.; Yang, M.; Zhang, Y.; Liu, Y.; Wang, D. Experimental measurements of mechanical properties of carbon dioxide hydrate-bearing sediments. Mar. Pet. Geol. 2013, 46, 201-209. 21. Casco, M. E.; Silvestrealbero, J.; Rey, F.; Bansode, A.; Urakawa, A.; Peral, I.; Kaneko, K. Methane hydrate formation in confined nanospace can surpass nature. Nat. Commun. 2015, 6, 6432. 22. Uchida, T.; Takeya, S.; Chuvilin, E. M.; Ohmura, R.; Nagao, J.; Yakushev, V. S.; Istomin, V. A.; Minagawa, H.; Ebinuma, T.; Narita, H. Decomposition of methane hydrates in sand, sandstone, clays, and glass beads. J. Geophys. Res.-Solid Earth 2004, 109, B05206. 23. Waite, W. F.; Santamarina, J. C.; Cortes, D. D.; Dugan, B.; Espinoza, D. N.; Germaine, J.; Jang, J.; Jung, J. W.; Kneafsey, T. J.; Shin, H.; Soga, K.; Winters, W. J.; Yun, T. S. Physical properties of hydrate-bearing sediments. Rev. Geophys. 2009, 47, RG4003. 23

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 24 of 36

24. Shimizu, H.; Kumazaki, T.; Kume, T.; Sasaki, S. Elasticity of single-crystal methane hydrate at high pressure. Phys. Rev. B 2002, 65, 212102. 25. Jendi, Z. M.; Servio, P.; Rey, A. D. Ideal Strength of Methane Hydrate and Ice Ih from First-principles. Cryst. Growth Des. 2015, 15, 5301-5309. 26. Jia, J.; Liang, Y.; Tsuji, T.; Murata, S.; Matsuoka, T. Microscopic Origin of Strain Hardening in Methane Hydrate. Sci. Rep. 2016, 6, 23548. 27. Mcmullan, R. K.; Jeffrey, G. A. Polyhedral Clathrate Hydrates. IX. Structure of Ethylene Oxide Hydrate. J. Chem. Phys. 1965, 42, 2725-2732. 28. Bernal, J.; Fowler, R. A theory of water and ionic solution, with particular reference to hydrogen and hydroxyl ions. J. Chem. Phys. 1933, 1, 515-548. 29. Jorgensen, W. L.; Madura, J. D.; Swenson, C. J. Optimized Intermolecular Potential Functions for Liquid Hydrocarbons. J. Am. Chem. Soc. 1984, 106, 6638-6646. 30. Molinero, V.; Moore, E. B. Water Modeled as an Intermediate Element between Carbon and Silicon. J. Phys. Chem. B 2009, 113, 4008-4016. 31. Jacobson, L. C.; Molinero, V. A Methane-Water Model for Coarse-Grained Simulations of Solutions and Clathrate Hydrates. J. Phys. Chem. B 2010, 114, 7302-7311. 32. Stillinger, F. H.; Weber, T. A. Computer Simulation of Local Order in Condensed Phases of Silicon. Phys. Rev. B 1985, 31, 5262-5271. 33. Wu, J.; Ning, F.; Trinh, T. T.; Kjelstrup, S.; Vlugt, T. J. H.; He, J.; Skallerud, B. H.; Zhang, Z. Mechanical instability of monocrystalline and polycrystalline methane hydrates. Nat. Commun. 2015, 6, 8743. 34. Polak, E.; Ribière, G. Note sur la convergence de methodes de directions conjugées. Rev. Francaise Inf. Rech. Opér. 1969, 16, 35-43. 35. Plimpton, S. Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comput. Phys. 1995, 117, 1-19. 36. Sellars, C. M.; McTegart, W. J. On the mechanism of hot deformation. Acta Metall. 1966, 14, 1136-1138.

24

ACS Paragon Plus Environment

Page 25 of 36 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

37. Kim, H. C.; Bishnoi, P. R.; Heidemann, R. A.; Rizvi, S. S. H. Kinetics of methane hydrate decomposition. Chem. Eng. Sci. 1987, 42, 1645-1653. 38. Zhang, Y. Y.; Pei, Q. X.; Wang, C. M. Mechanical properties of graphynes under tension: A molecular dynamics study. Appl. Phys. Lett. 2012, 101, 081909. 39. Qing-Xiang, P.; Yong-Wei, Z.; Vivek, B. S. Mechanical properties of methyl functionalized graphene: a molecular dynamics study. Nanotechnology 2010, 21, 115709. 40. Pei, Q. X.; Zhang, Y. W.; Shenoy, V. B. A molecular dynamics study of the mechanical properties of hydrogen functionalized graphene. Carbon 2010, 48, 898-904. 41. He, L.; Guo, S.; Lei, J.; Sha, Z.; Liu, Z. The effect of Stone–Thrower–Wales defects on mechanical properties of graphene sheets – A molecular dynamics study. Carbon 2014, 75, 124-132. 42. Bucholz, E. W.; Sinnott, S. B. Mechanical behavior of MoS2 nanotubes under compression, tension, and torsion from molecular dynamics simulations. J. Appl. Phys. 2012, 112, 123510. 43. Blase, X.; Gillet, P.; San Miguel, A.; Mélinon, P. Exceptional Ideal Strength of Carbon Clathrates. Phys. Rev. Lett. 2004, 92, 215505. 44. Li, B.; Sun, H.; Chen, C. Extreme Mechanics of Probing the Ultimate Strength of Nanotwinned Diamond. Phys. Rev. Lett. 2016, 117, 116103. 45. Zhang, M.; Liu, H.; Li, Q.; Gao, B.; Wang, Y.; Li, H.; Chen, C.; Ma, Y. Superhard BC3 in Cubic Diamond Structure. Phys. Rev. Lett. 2015, 114, 015502. 46. Takeya, S.; Udachin, K. A.; Moudrakovski, I. L.; Susilo, R.; Ripmeester, J. A. Direct Space Methods for Powder X-ray Diffraction for Guest−Host Materials: Applications to Cage Occupancies and Guest Distributions in Clathrate Hydrates. J. Am. Chem. Soc. 2010, 132, 524-531. 47. Davidson, D. W.; Ratcliffe, C. I.; Ripmeester, J. A. 2H and13C NMR Study of Guest Molecule Orientation in Clathrate Hydrates. J. Inclusion Phenom. 1984, 2, 239-247. 48. Udachin, K. A.; Ratcliffe, C. I.; Ripmeester, J. A. Single Crystal Diffraction Studies of Structure I, II and H Hydrates: Structure, Cage Occupancy and Composition. J. Supramol. Chem. 2002, 2, 405-408.

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

Page 26 of 36

49. Udachin, K. A.; Ratcliffe, C. I.; Ripmeester, J. A. Structure, Composition, and Thermal Expansion of CO 2 Hydrate from Single Crystal X-ray Diffraction Measurements. J. Phys. Chem. B 2001, 105, 4200-4204. 50. Alavi, S.; Dornan, P.; Woo, T. K. Determination of NMR Lineshape Anisotropy of Guest Molecules within Inclusion Complexes from Molecular Dynamics Simulations. ChemPhysChem 2008, 9, 911-919. 51. Petrovic, J. Review mechanical properties of ice and snow. J. Mater. Sci. 2003, 38, 1-6. 52. Fletcher, N. H.; Hasted, J. B. The Chemical Physics of Ice. Phys. Today 1971, 24, 49-49. 53. Godio, A.; Rege, R. B. The mechanical properties of snow and ice of an alpine glacier inferred by integrating seismic and GPR methods. J. Appl. Geophys. 2015, 115, 92-99. 54. Trickett, Y.; Baker, I.; Pradhan, P. The orientation dependence of the strength of ice single crystals. J. Glaciol. 2000, 46, 41-44. 55. Alexander, S.; Vasily, V. B.; Athanasios, A. Automated identification and indexing of dislocations in crystal interfaces. Modell. Simul. Mater. Sci. Eng. 2012, 20, 085007. 56. Hayes, C.; Webb, W. Dislocations in ice. Science 1965, 147, 44-45. 57. Blumenau, A. T.; Elsner, J.; Jones, R.; Heggie, M. I.; Öberg, S.; Frauenheim, T.; Briddon, P. R. Dislocations in hexagonal and cubic GaN. J. Phys.: Condens. Matter 2000, 12, 10223. 58. Stern, L. A.; Durham, W. B.; Kirby, S. H. Grain-size-induced weakening of H2O ices I and II and associated anisotropic recrystallization. J. Geophys. Res.-Solid Earth 1997, 102, 5313-5325. 59. Stern, L. A.; Kirby, S. H.; Durham, W. B. Peculiarities of methane clathrate hydrate formation and solid-state deformation, including possible superheating of water Ice. Science 1996, 273, 1843-1848. 60. Stern, L. A.; Kirby, S. H.; Durham, W. B. Polycrystalline Methane Hydrate:  Synthesis from Superheated Ice, and Low-Temperature Mechanical Properties. Energy Fuels 1998, 12, 201-211.

26

ACS Paragon Plus Environment

Page 27 of 36 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

Tables Table 1. Interaction parameters for the coarse-grained model

ε (kcal/mol)

σ (Å)

λ

Water

6.189

2.3925

23.15

Methane

0.340

4.08

0

Water-Methane

0.180

4.00

0

Table 2. Crystallographic slips of monocrystalline Ih subjected to uinaxial mechanical loads along the [1 1 2 0] ,

[0 0 0 1] and [0 1 1 0] directions.

Load

Load directions

Slip plane

Slip type

Modes

Tension

[1 1 2 0]

112   2

[0 0 0 1]

0112  

[0 1 1 0]

0112  

[1 1 2 0]

112   0

Pyramidal slip

Prismatic slip

112   2 Compression

[0 0 0 1]

Pyramidal slip

0112  

0110  

[0 1 1 0]

27

ACS Paragon Plus Environment

Prismatic slip

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 28 of 36

Figures and Figures legends

Figure 1 Molecular models of SI MMH and Ih. (a) Samll (512) water cages composed of 12 pentagonal faces. (b) Large (51262) water cages composed of 12 pentagonal and 2 hexagonal faces. (c) SI MMH unit cell. (d)-(f) Molecular structures of Ih with specific crystal planar surfaces colored by yellow. The oxygen atoms and hydrogen atoms in the interior are rendered by red and light pink, respectively. (g)-(h) Molecular atomic structures of SI MMH and Ih.

28

ACS Paragon Plus Environment

Page 29 of 36 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 2 Mechanical properties of MMH for various engineering strain rates. (a) Tensile stress versus strain curves along the [1 0 0] direction. (b) Compressive stress versus strain curves along the [1 0 0] direction. Snapshots of MMH at the engineering strain rate of 1.0×108/s. (c) Tensile deformation test conditions along the [1 0 0] direction. (d) Compressive deformation test conditions along the [1 0 0] direction. Molecule particles are colored according to the potential energy. Specifically, uniformly selected molecule particles are colored by purple for clearly monitoring structural changes during deformation.

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

Page 30 of 36

Figure 3 Mechanical properties of MMH for various temperatures under mechanical loads. (a) and (b) Mechanical stress versus strain curves along the [1 0 0] direction under tension and compression, respectively. (c) Young’s modulus versus temperature under tension along the [1 0 0] direction. (d) Young’s modulus versus temperature under compression along the [1 0 0] direction.

30

ACS Paragon Plus Environment

Page 31 of 36 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 4 Mechanical properties of MMH along various directions under mechanical loads. (a) Mechanical stress versus strain curves under tension. (b) Mechanical stress versus strain curves under compression.

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

Page 32 of 36

Figure 5 Water cage configurations of MMH along various directions under mechanical loads. (a) Water cage configurations along the [1 0 0] direction. (b) Water cage configurations along the [1 1 0] direction. (c) Water cage configurations along the [1 1 1] direction. Black arrows indicate the mechanical loading directions. The water particles are colored by red. The carbon atoms in water cage of CH4@512 and CH4@51262 are colored by brown and blue, respectively. The snapshots of number 1 are water cage of CH4@512, whereas the snapshots of numbers 2-3 are water cage of CH4@51262.

Figure 6 Mechanical properties of MMH under mechanical loads. (a) Tensile stress-strain curves of MH along the

[1 0 0] direction. (b) Compressive stress-strain curves of MH along the [1 0 0] direction.

32

ACS Paragon Plus Environment

Page 33 of 36 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 7 Snapshots of MMH under mechanical loads. (a)-(c) Tensile deformation test conditions along the [1 0 0] direction for guest molecule CH4-free@512, CH4-free@51262, and CH4-free@512&CH4-free@51262, respectively. (d)-(f) Compressive deformation test conditions along the [1 0 0] direction for guest molecule CH4-free@512, CH4-free@51262, and CH4-free@512&CH4-free@51262, respectively. Molecule particles are colored according to the

33

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 34 of 36

potential energy. Specifically, uniformly selected molecule particles are colored by yellow for clearly monitoring structural changes during deformation.

Figure 8 Mechanical properties of MMH and Ih under mechanical loads. (a) Tensile stress-strain curves of MMH and Ih under various directional loads at 263.15K and 10MPa. (b) Compressive stress-strain curves of MMH and Ih under various directional loads at 263.15K and 10MPa. Solid line represents MMH, and dash dot line represents hexagonal ice in Figures 8a-8b. (c) Tensile stress-strain curves of Ih under various directional loads. (d) Compressive stress-strain 34

ACS Paragon Plus Environment

Page 35 of 36 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

curves of Ih under various directional. Upon mechanical loads along [0 0 0 1] , Ih exhibits unique elastic responses that are in strikingly contrast to that under [1 1 2 0] and [0 1 1 0] .

Figure 9 Localized snapshots of Ih for various crystallographic directions at 223.15 K and ambient pressure along various crystallographic directions. (a)-(c) Molecular snapshots of monocrystals subjected to [1 1 2 0] , [0 0 0 1] and

35

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 36 of 36

[0 1 1 0] directional tensions at different strains, respectively. (d)-(f) Localized snapshots of monocrystals for [1 1 2 0] , [0 0 0 1] and [0 1 1 0] directional compression loadings at different strains, respectively. Black arrows indicate the mechanical loading directions. Water particles are colored according to the identified type of structural phase; Hexagonal, cubic and unidentified water structures are rendered by orange, light blue and purple for clarification, respectively. Specifically, uniformly selected water particles are colored by blue for clearly monitoring structural changes during deformation. Strain-induced dislocations in Ih are illustrated by solid segments. Different colored segments represent different type of dislocation; 1/3 1 2 1 0 ,  0 0 0 1 , 1 1 0 0 , 1/3 1 1 0 0 and unidentified dislocations are rendered by green, yellow, pink, sky-blue and red for clarification, respectively.

TOC Graphic

36

ACS Paragon Plus Environment