Transformation of hydrogen bond network during the CO2 clathrate hydrate dissociation
Kirill Gets, Vladimir Belosludov, Ravil Zhdanov, Yulia Bozhko, Rodion Belosludov, Oleg Subbotin, Nikita Marasanov, Yoshiyuki Kawazoe
TTransformation of hydrogen bond network during the CO clathrate hydrate dissociation Kirill Gets *, Vladimir Belosludov , Ravil Zhdanov , Yulia Bozhko , Rodion Belosludov , Oleg Subbotin , Nikita Marasanov , Yoshiyuki Kawazoe Novosibirsk State University, 1, Pirogova str., Novosibirsk, 630090, Russia Nikolaev Institute of Inorganic Chemistry SB RAS, 3, Acad. Lavrentiev Ave., Novosibirsk, 630090, Russia Institute for Materials Research, Tohoku University, 2-1-1 Katahira, Aoba-ku, Sendai, 980–8577, Japan New Industry Creation Hatchery Center, Tohoku University, 4-4-6 Aramaki aza Aoba, Aoba-ku,Sendai, 980-8579, Japan SRM institute of Science and Technology, SRM Nagar, Kattankulathur, Kancheepuram District,Tamil Nadu, 603 203, India
Abstract
The object of this study is the kinetic process of solid-liquid first-order phase transition –melting of carbon dioxide CS-I hydrate with various cavity occupation ratios. The work was donewithin a framework of study on the local structure of water molecules. These include the timedepending change of the short-range order at temperatures close to the melting point and comparisonwith hexagonal ice structure. Using molecular dynamics method, dependencies of the internal energyof the studied systems on the time of heating were found. Jumps in the internal energy of solids in therange at 275-300 K indicate a phase transition. The study of oxygen-oxygen radial distribution andhydrogen-oxygen-oxygen mutual orientation angles between molecules detached at no more than 3.2Å allowed to find the H-bond coordination number of all molecules and full number of H-bonds andshowed the instant (less than 1 ns) reorganization of short-range order of all molecules. The structureanalysis of every neighbor water molecules pairs showed the ~10-15% decrease of H-bond numberafter the melting whereas all molecules form single long-range hydrogen bond network. The analysisof hydrogen bond network showed the minor changes in the H-bond interaction energy at solid-liquidphase transition.
Keywords
Gas hydrate, phase transition, molecular dynamics, hydrogen bond network, short-range order.
1. Introduction
Gas hydrates are crystalline water-based compounds. Their structure consists of water as thehost molecules and trapped inside gas as guest molecules. The most common structures are: CS-I,CS-II and HS-III. The cavities in host lattice consist of hydrogen bonded tetragons, pentagons andhexagons. For example, the CS-I hydrate consists of large cavities (twelve pentagons and twohexagons, 5 ) and small cavities (twelve pentagons, 5 ). Each cavity can hold a number of variousguest molecules depends on their size [1]. Thermodynamic stability of clathrate hydrates depends on p-T conditions as well as type of guest molecules inside of hydrate cavities. [2, 3].Gas hydrates are considered as attractive media to store and sequester carbon dioxide (CO ) inthe solid form [4-6]. This problem is considered as part of a task of replacing methane in methanehydrate deposits, which are associated with a high threat of methane in the greenhouse effect in thecase of methane hydrate melting and gas release [7-9].Development of technologies for carbon dioxide utilization and storage in hydrate formrequires a clear understanding of the processes of formation and dissociation of CO hydrates.Recently, the molecular dynamic (MD) studies of gas hydrate described the dissociation processes byradual destruction of hydrate structure at hydrate-liquid water interface (dissolution) [10-12]. Theother possibility of the melting caused by hydrate heating only gives an idea of the structural changesin the host lattice close to the solid-liquid phase transition. This is an attempt to find the answer to thefundamental question: whether the melting of solid is a first-order phase transition or a continuousprocess [13]. Melting of a hydrate can be considered as a heterogeneous/homogeneous process in apresence/absence of imperfections [14], such as, for example, empty cavities. The results of MDsimulations of ideal crystals shows the 20-30% increase in value of melting temperature incomparison to experimental data [15, 16]. However it can diminish the high value of meltingtemperature by the increasing the system size [14, 17].Molecular dynamics is a powerful and precise theoretical method to study both macroscopicproperties and microscopic structure and to describe a relationship between them. However, the studyof the kinetic processes that taken place at gas hydrate formation and dissociation are limited [10,11]. The most common water models for describing hydrates melting are: TIP4P [18], SPC/E [19],TIP4P-Ew [20], TIP4P/Ice [21], TIP4P/2005 [22]. Any of these potentials are able to give aqualitative description of the processes of formation and decomposition of gas hydrates. The TIP4Pand SPC/E potentials are preferred for gas hydrate dissociation kinetic studying [12, 23, 24], and therest of potentials for the studying of hydrate formation [17, 25]. A more detailed description can befound in the work [26].Melting of gas hydrates was studied in a series of works [27-32] which are devoted to thestudy of thermodynamic properties, such as free energy or chemical potential of hydrate systems [33,34], effects of “help” gases or inhibitors on hydrate stability [35-37], but at the moment there is nodescription of the melting of hydrates within the framework of the hydrogen bond networkreorganization and changes of energy and structure properties of every close H O pair (short-rangeorder of every molecule). Our work is devoted to study the kinetics processes of melting CO hydrate in the absence ofliquid water and with constant uniform heating of the hydrate and hence to study the phase transitionon molecular level. The local reorganization can be generalized and described in the terms ofhydrogen bond network. The hydrogen bond network unites all water molecules. Its local structure[38] and collective dynamics [39] are responsible for water properties. The process of clathratehydrate melting rises a problem of study the transformation of hydrogen bond network on the liquidwater – clathrate hydrate phase transition [40, 41] as well as on the surface of the resultant gasbubbles.In this paper, we propose a description of melting process via the structural analysis ofneighboring water molecule ensemble (R O−O < 3.2 Å) and their interaction energy. The study isconducted on a series of snapshots (instantaneous configurations) obtained by MD during themodeling of crystalline hexagonal ice (I h ), CS-I hydrate with 0%, 25% and 50% ratio of cavityoccupation by the CO molecules close to the melting point and. The aim is to study the changing ofstructural (short-range order) and energy properties of CO hydrate and ice during the meltingprocess under the heating.
2. Calculation details
To model the CO hydrate melting and compare its structure characteristics with melting Ihice, the MD method within the approximation of a rigid molecule was used (LAMMPS softwarepackage, [42]). The interaction between water molecules was described by the modified simple pointcharge extended (SPC/E) potential [19]. It was chosen because it qualitatively describes the kineticsof dissociation [12, 23, 24]. The modified version was used according to good description ofthermodynamics for these phases in our previous studies [43, 44]. The short-range interaction wasdescribed by Lennard-Jones potential with parameter σ = 3.1556 Å and the energy parameter ε =0.65063 kJ mol − . The charges on hydrogen (q H = +0.4238|e|) and on oxygen (q O = −0.8476|e|)atoms were obtained to describe the long-range Coulomb interaction. Other parameters were setriginal [19]. The CO molecules were considered as single particle and their parameters were σ = 4.0Å and the energy parameter ε = 1.5798 kJ mol−1 [45-47] . The supercell cell of 10 x 10 x 10 unit cells of CS-I hydrate with 0%, 25% and 50% rate ofsingle carbon dioxide occupation of cavity was chosen for the calculations. The procedure forobtaining the structure of crystalline ice I h has been described previously in detail [48].The structure changing during the melting process the hydrate and the ice were modeledusing an NPT ensemble. The pressure was set to 1 bar, the initial temperature was set to 275 K andthe heating rate differs from ~0.05 K/ns for CS-I structures, and ~0.1 K/ns for I h ice. During thesimulation, the coordinates of all atoms within the model were recorded at each 5 fs for empty CS-Ihydrate, 100 fs for 25% and 50% cavity filled CS-I hydrate, and 200 fs for I h ice, determining theinstantaneous configurations of the systems (snapshots). The temperature and pressure of the modelsystems are controlled by Nose-Hoover style non-Hamiltonian equations of molecular motion [49].To study hydrogen bond network structure characteristics the analysis of a series of snapshotsbefore the melting point and after it has been performed. This includes:a) determining the dependence of the distribution of the number of closest O···O pairs (R O−O < 3.2 Å)on the distance between them, on the mutual orientation angles ∠ H−O···O. Defining the H-bondformation between the O···O pair using the geometric criteria (R
O−O < 3.2 Å, ∠ H−O···O < 30 o ) [50].This criteria were chosen due to slight differences between energy, geometric and hybrid criteria [51,52];b) determining the dependence of the distribution of the number of CO ···CO on the distancebetween them;c) determining the dependence of the number of H-bonds per molecule as a function of time andtemperature and the average number of molecules linked by H-bonds in the nearest environment(short-range order) of all molecules;d) calculating the interaction energies between neighbor water molecules as a sum of Coulomb andvan der Waals (Lennard-Jones) potentials.
3. Results and discussion
The time dependence of caloric curves for ice and CO hydrate are presented in Fig. 1. The caloriccurve represents the internal energy of the system i.e. the sum of potential and kinetic energies. Theempty CS-I structure collapsed almost at the beginning of simulation at 275 K. The hydrate with 25%occupation melted after 15 ns at 275 K and at 50% of occupation, the hydrate melted at 300 K afterthe 487 ns of heating. As the solids melt the internal energy rises in absolute value by ~10%, exceptthe melting of CS-I hydrate with 50% cavity occupation ratio, where the energy increasing value is~14% that could be influenced by the presence of high number of carbon dioxide molecules as wellas higher temperature. These results agree well with the prediction calculated from the Clapeyronequation, the experimental data on CO clathrate hydrate dissociation [53] and the simulation resultsof melting of water in solid phase [54, 55]. The observed jump of internal energy is a sign of first-order phase transition: the meltingprocess take less than one nanosecond in the comparison to the surface melting, i.e in the presence ofliquid water on the surface. In that case, the dissociation could take 20–50 ns [10] for the cell of thesame size and more for larger cell model, because first-order phase transition taking place throughthe full volume and should characterizes by the intermittent processes that are continuous at thesurface melting. A short time of phase transition (structure transformation time) has collective natureand it takes place due to the fact that the network of hydrogen bonds even in the liquid phase isdynamic [39] and its structure changes have collective nature as was found [56-59].Transition instantaneity will also be clearly seen at the following figures proving the first-order phase transition.ig.1. Caloric curves for CS-I hydrates with different occupancy rates at 275 K and 300 K and I h iceat 279 K. Dependence of kinetic and potential energy sum on time. Caloric curve of ice wasnormalized to the values as they consist of 46 molecules.Fig. 2. Dependence of H-bond number N H-bonds normalized to number of molecules N mol as a functionof time for CS-I hydrates with different occupancy rates at 275 K and 300 K and I h ice at 279 K.The study of molecule pairs, which mutual geometry satisfies the hydrogen bond geometricalcriteria, made it possible to calculate the number of hydrogen bonds in hydrate and ice before andafter melting point. These pairs form the H-bond network and the dependence of their normalizednumber on time is presented in Fig. 2. Before the simulation started (without thermal motion) everywater molecule of I h ice and hydrates has exactly 4 molecules in the short-range order. It correspondsto N H-bonds / N mol exactly equals 2. However, due to the thermal motion some H-bonds in the solidstructure are broken and N
H-bonds / N mol is lower than 2. It can be seen that the number of H-bonds varies with time, but is near the average value. Thechange in the number of hydrogen bonds with time is connected with the constant rebuilding of theocal order of each water molecule — the change in the angles of mutual orientation and distancesbetween molecules due to thermal motion. As it was shown earlier, the transition from the solid to theliquid phase takes less than one nanosecond: ubiquitous rebuilding of short-range order structureclearly indicates a collective process. For melted ice and hydrate average value of N
H-bonds / N mol depends on temperature and lowers with temperature increase, however the effect of CO moleculepresence improves the stability of hydrate structure with increase of cavity occupancy rate. Anotherpoint is the melted state when CO molecules are wedged in H O molecule short-range order,preventing the formation of H-bonds that can be seen for melted 50% occupied hydrate. Empty and25% occupied CS-I and I h ice structures are in good agreement with results obtained by Pauling [60]which states that ice melting leads to the breaking of ~15% of H-bonds. We found that in every snapshot (at any time) of water system, all water molecules form adynamic H-bond network. Thus, the H-bond network is preserved when solid water melts.Fig. 3. Dependence of H-bond number N H-bonds normalized to number of O···O pairs (R
O−O < 3.2 Å)N
O-O as a function of time for CS-I hydrates with different occupancy rates at 275 K and 300 K and I h ice at 279 K.Geometry changes that H-bond network undergoes upon melting of CO CS-I hydrate and I h ice can be estimated by changes in the H-bond lengths (R O−O ) and the mutual orientation angles( ∠ H−O···O) between neighbor H O molecules (R
O−O < 3.2 Å). Figure 3 shows the time dependenceof H-bond number normalized to the number of neighbor molecule pair, which may not satisfy theangular criterion of the hydrogen bond, in CS-I hydrate with different occupancy rate, I h ice. Eachwater molecule in the crystal structure of the CS-I hydrate and I h ice has a tetrahedral orderedenvironment (short-range order). At high temperatures, the spread in the lengths of hydrogen bondsand the angles of mutual orientation between the neighboring molecules increases due to the thermalmotion. However, the number of neighboring molecule pairs, which internal mutual orientation angleis above 30 degrees, was close to 1%. This means that almost all pairs satisfy the H-bond criterion inangle and N H-bonds /N O-O relation. The comparison of data presented in Fig. 2 and Fig. 3, allows to conclude that N
O-O normalized to N mol before and after melting fluctuates around value of 2. Thus, the main change in theH-bond geometry upon melting is that 10–15% of pairs do not satisfy the angular criterion of thehydrogen bond and the violation of the tetrahedral short-range order takes place. Nevertheless themost of the O···O pairs satisfy the angle H-bond criterion. ig. 4. Time dependence of average interaction energy between O···O pair, H-bonded and non H-bonded pairs in CS-I hydrates with 0%, 25% and 50% occupancy rates at 275 K and 300 K and I h iceat 279 K.The sharp change in the geometry of the H-bonds does not lead to a significant degradation inthe stability of the H-bond network. Figure 4 shows the time dependence of average interactionenergy between every neighbor molecules, molecules that form H-bond and neighbor moleculeswhich satisfy H-bond criterion on distance and do not satisfy the mutual orientation angle H-bondcriterion (non H-bond) in CS-I hydrate and I h ice before and after melting. The obtained resultsindicates the completeness of the melting process.The solid phase is characterized by the short-range ordering of each molecule; therefore,before reaching the melting point, the average value of the interaction energy for the H-bonded pairscoincides with the average value for all pairs of neighbor molecules. After melting, the number ofnon H-bonded pairs significantly increases that leads to a decrease in the average values fluctuationscatter. For the same reason, the average energy value for all pairs increases by ~15%. However,during the melting the binding energy between the H-bonded molecules changes only by 5-7%. Thismeans that the H-bond geometry does not undergo significant changes in the melting process andindicates the preservation of tetrahedral short-range order in the H-bond network after melting.The insignificant changing in the interaction energy between the H-bonded molecules couldbe the origin of H-bond network stability due to the fact that the energy of thermal motion in thetemperature region of 275 – 300 K is lower than 1 kcal/mol. On average, the interaction energy ofany nearby molecules exceeds the value of thermal energy.The strong fluctuation of the average energy value for pairs of non H-bonded pairs can beexplained by thermal motion and ~1% number. Before melting this energy fluctuates close to -2 kcal/mol values, but after melting it fluctuates in the region [-1; -0.8] kcal/mol. In the solid phase, thisenergy is significantly lower than thermal energy and is responsible for the structure stability in thesolid phase. In the liquid phase, this energy is comparable to thermal energy, which explains theability of its structure to reorganize. The increase of average interaction energy value after melting isdue to the further disordering of the short-range order of molecules in the liquid phase. ig. 5. Dependence of molecule fraction having 1, 2, 3 and 4 H-bonds as a function of time of CS-Ihydrates with different occupancy rates at 275 K and 300 K and I h ice at 279 K.To determine the short-range order of each molecule in the terms of H-bonding and to studythe preservation of tetrahedral type of hydrogen bond network during the solid melting, the timedependencies of the fraction of molecules participating in 1, 2, 3 and 4 H-bonds were obtained andpresented in Fig. 5. It is clear that in the solid state the overwhelming majority of molecules hastetrahedral surrounding. Less than 10% of molecules have only 3 water molecules in their short-rangeorder. It can be expected than that thermal motion leads to breaking of the number of H-bonds insolid state, but the presence of CO molecules makes hydrate structure stable at higher temperature.In the case of water, only about 50% of the water molecules form 4 hydrogen bonds with theirneighbors. Upon melting the fraction of molecules taking part in the formation of 4 H-bonds instantly decreasing by 40-50%, which leads to significant increasing the fraction of molecules that forming 3H-bonds. This makes tetrahedral H-bonding weaker (with rising of absolute value of full energy) andleads to the local structure changing. The fraction of molecules that taking part in the formation of 3and 4 H-bonds make up over 85% in total that allows the H-bond network to unite all molecules. Thefraction of molecules participating in 2 H-bonds at given temperatures could reach 10%. The numberof molecules forming 1 H-bond is insignificantly low. Thus before and after melting H-bond networkmostly consists of molecules that form 4 hydrogen bonds – H O molecules pairs, which geometrysatisfies to H-bond criteria and has a small scatter of lengths and angles. Fraction of moleculesforming 4 H-bonds in water decreasing with temperature (and thermal motion) increase.Figure 2 and Figure 5 show that the presence of small number of CO molecules almost doesnot affect the tetrahedral ordering of H-bond network of solid hydrate. However the presence ofcavities leads to instability of ~3-5% of hydrogen bonds at comparing I h ice and CS-I structures. a) (b)(c) 485 ns (d) 487 ns (e) 487.25 ns (f) 489 nsFig. 6. Distribution P d of CO ···CO pair lengths d O-O in the CS-I hydrates with 25% (a) and 50% (b)cavity occupancy rate with a 0.01 Å step at the initial (ideal) state, before melting and after melting.Spatial distribution of CO molecules before melting (c), at melting point (d) and right after meltingpoint (e) and after 1.75 ns (f) of 50% occupied CS-I hydrate with corresponding time mark.Figures 6a,b show the distribution of CO ···CO pair lengths of initial and final configurationsand before and after the melting point, respectively. The sharp peaks of the initial state correspond tothe regular arrangement of CO molecules in hydrate structure. These peaks blur during the heatinghowever CO molecules remain inside the deformed cavities. After the melting point the CO molecules are distributed evenly throughout the whole volume (Fig. 6e). The small peak in theinterval of 4-6 Å as shown in Fig. 6b represent the formation of small CO gas bubbles .
4. Conclusions
The investigation of the kinetic process of carbon dioxide CS-I hydrate and I h ice melting atthe molecular level showed the instantaneous change of hydrogen bond network as well as short-range order of every H O molecule during the solid-liquid transition. This indicates that the nature ofdissociation process is a collective phenomenon for ice and carbon dioxide hydrate. The calculatedcaloric curves clearly showed the jump in the internal energy of the ice, hydrate of CO gas with 0%,25% and 50% cavity occupation ratio in the region 275 K - 300 K. The step by step conductedanalysis of sequential structures revealed the instant (less than 1 nanosecond) first-order solid-liquidphase transition at 275, 279 and 300 K.The analysis showed that every H O molecule taking part in the formation of dynamichydrogen bond network. It was shown that integrity of network was not affected by the phasetransition and it persists after the melting point was reached when ~15-17% of H-bonds are broken.The interaction energies of the most of H-bonded pairs are greatly exceed the thermal energy eithern solid and liquid phase. The corresponding energies of non H-bonded pairs slightly exceeds andcomparable with thermal energies in solid and liquid phase, respectively. Overwhelming majority ofO···O pairs in solid state are bound by H-bonds defining the solid properties. However, the formationof a hydrogen bond between molecules strongly depends on the distance between the oxygen atomsand the angle of mutual orientation, which can vary greatly with time even in the liquid state.The presence of cavities in crystalline structure decreases its stability however the increase ofCO molecules occupied the cavities stabilizes the hydrate structure. Acknowledgements
This work of K.G., V.B., R.Z., Y.B., O.S., N.M. and R.G. is supported by the Russian ScienceFoundation under grant No 18-19-00124 (caloric curves, hydrogen bond structure and energyanalysis) that held at the Novosibirsk State University. R.B. and Y.K. are thankful to the Ministry ofEducation, Culture, Sports, Science, and Technology of Japan (Grant No. 17H03122) for financialsupport (large-scale molecular dynamics simulation) and are grateful for the continuous support ofthe crew at the Center for Computer Materials Science at the Institute for Materials Research,Tohoku University, Sendai.
References
1. N.J. English, J.M.D. MacElroy, Perspectives on molecular simulation of clathrate hydrates: Progress, prospects and challenges, Chem. Eng. Sci. 121 (2015) 133–156. https://doi.org/10.1016/j.ces.2014.07.047.[2] E.D. Sloan, C. Koh, Clathrate Hydrates of Natural Gases, third Edition, CRC Press, Boca Raton, 2007. [3] M.R. Walsh, C.A. Koh, D.E. Sloan, A.K. Sum, D.T. Wu, Microsecond simulations of spontaneous methane hydrate nucleation and growth, Science 326 (2009) 1095–1098. https://doi.org/10.1126/science.1174010.[4] C.A. Rochelle, A.P. Camps, D Long, A. Milodowski, K. Bateman, D. Gunn, P. Jackson, M.A. Lovell, J. Rees, Can CO hydrate assist in the underground storage of carbon dioxide? Geol. Soc. Spec. Publ. 319 (2009) 171–183. https://doi.org/10.1144/SP319.14.[5] B. Tohidi, J. Yang, M. Salehabadi, R. Anderson, A. Chapoy, CO Hydrates Could Provide Secondary Safety Factor in Subsurface Sequestration of CO , Environ. Sci. Technol. 44 (2010) 1509–1514. https://doi.org/10.1021/es902450j.[6] N. Goel, In situ methane hydrate dissociation with carbon dioxide sequestration: Current knowledge and issues, J. Pet. Sci. Eng. 51 (2006) 169. https://doi.org/10.1016/j.petrol.2006.01.005.[7] W.F. Waite, L.A. Stern, S.H. Kirby, W.J. Winters, D.H. Mason, Simultaneous determination of thermal conductivity, thermal diffusivity and specific heat in sI methane hydrate, Geophys. J. Int. 169(2007) 767-774. https://doi.org/10.1111/j.1365-246X.2007.03382.x.[8] R. Boswell, Resource potential of methane hydrate coming into focus, J. Pet. Sci. Eng. 56 (2007) 9-13. https://doi.org/10.1016/j.petrol.2006.09.002.[9] R. Boswell, T.S. Collett, Current perspectives on gas hydrate resources, Energy Environ. Sci. 4 (2011) 1206-1215. https://doi.org/ 10.1039/C0EE00203H .[10] S. Sarupria, P.G. Debenedetti, Molecular dynamics study of carbon dioxide hydrate dissociation,
J. Phys. Chem. A
115 (23) (2011)
Piñeiro, M.M., Vega, C., 2015. Molecular dynamics simulation of CO2 hydrates: Prediction of three phase coexistence line.
J. Chem. Phys. (12) 124505. https://doi.org/10.1063/1.4916119.[12] Alavi, S., Ripmeester, J. A., 2010. Nonequilibrium adiabatic molecular dynamics simulations of methane clathrate hydrate decomposition.
J. Chem. Phys. (14), 144703. https://doi.org/10.1063/1.3382341.[13] Dash, J. G., 1999. History of the search for continuous melting.
Rev. Mod. Phys. (5), 1737. https://doi.org/10.1103/RevModPhys.71.1737.14] S. Liang, L. Yi, D. Liang, Molecular insights into the homogeneous melting of methane hydrates, J. Phys. Chem. C (49) (2014) 28542-28547. https://doi.org/10.1021/jp511362s.[15] Jin, Z. H., Gumbsch, P., Lu, K., Ma, E., 2001. Melting mechanisms at the limit of superheating.
Phys. Rev. Lett. (5), 055703. https://doi.org/10.1103/PhysRevLett.87.055703.[16] P.M. Agrawal, B.M. Rice, D.L. Thompson, Molecular dynamics study of the effects of voids and pressure in defect-nucleated melting simulations, J. Chem. Phys. (21) (2003) 9680-9688. https://doi.org/10.1063/1.1570815.[17] Conde, M. M., Vega, C., 2010. Determining the three-phase coexistence line in methane hydrates using computer simulations.
J. Chem. Phys. (6), 064507. https://doi.org/10.1063/1.3466751.[18] W.L. Jorgensen, J. Chandrasekhar, J.D. Madura, R.W. Impey, M.L. Klein, Comparison of simple potential functions for simulating liquid water,
J. Chem. Phys. (2) (1983) 926-935. https://doi.org/10.1063/1.445869.[19] H.J.C. Berendsen, J.R. Grigera, T.P. Straatsma, The missing term in effective pair potentials, J. Phys. Chem. 91 (24) (1987) 6269-6271. https://doi.org/10.1021/j100308a038.[20] H.W. Horn, W.C. Swope, J.W. Pitera, J.D. Madura, T.J. Dick, G.L. Hura, T. Head-Gordon, Development of an improved four-site water model for biomolecular simulations: TIP4P-Ew,
J. Chem. Phys. (20) (2004) 9665-9678. https://doi.org/10.1063/1.1683075.[21] Abascal, J. L. F., Sanz, E., García Fernández, R., Vega, C., 2005. A potential model for the study of ices and amorphous water: TIP4P/Ice.
J. Chem. Phys. (23), 234511. https://doi.org/10.1063/1.1931662.[22] Abascal, J. L., Vega, C., 2005. A general purpose model for the condensed phases of water: TIP4P/2005.
J. Chem. Phys. (23), 234505. https://doi.org/10.1063/1.2121687.[23] N. Choudhary, S. Das, S. Roy, R. Kumar, Effect of polyvinylpyrrolidone at methane hydrate-liquid water interface. Application in flow assurance and natural gas hydrate exploitation,
Fuel ,
186 (2016)
J. Phys. Chem. C (23) (2013) 12172-12182.https://doi.org/10.1021/jp4023772.[25] Y.T. Tung, L.J. Chen, Y.P. Chen, S.T. Lin, The growth of structure I methane hydrate from molecular dynamics simulations,
J. Phys. Chem. B (33) (2010) 10804-10813. https://doi.org/10.1021/jp102874s.[26] N. Choudhary, S. Chakrabarty, S. Roy, R. Kumar, A comparison of different water models for melting point calculation of methane hydrate using molecular dynamics simulations,
Chem. Phys.
516 (2019)
J. Chem. Phys. (5), 056101. https://doi.org/10.1063/1.4790647.[28] Michalis, V. K., Costandy, J., Tsimpanogiannis, I. N., Stubos, A. K., Economou, I. G., 2015. Prediction of the phase equilibria of methane hydrates using the direct phase coexistence methodology.
J. Chem. Phys. (4), 044501. https://doi.org/10.1063/1.4905572.[29] D. Jin, B. Coasne, Molecular Simulation of the Phase Diagram of Methane Hydrate: Free Energy Calculations, Direct Coexistence Method, and Hyperparallel Tempering,
Langmuir , (42) (2017) 11217-11230. https://doi.org/10.1021/acs.langmuir.7b02238.[30] M.H. Waage, T.J. Vlugt, S. Kjelstrup, Phase Diagram of Methane and Carbon Dioxide Hydrates Computed by Monte Carlo Simulations, J. Phys. Chem. B (30) (2017) 7336-7350. https://doi.org/10.1021/acs.jpcb.7b03071.[31] L. Jensen, K. Thomsen, N. von Solms, S. Wierzchowski, M.R. Walsh, C.A. Koh, E.D. Sloan, D.T. Wu, A.K. Sum, Calculation of liquid water− hydrate− methane vapor phase equilibria from molecular simulations,
J. Phys. Chem. B (17) (2010) 5775-5782. https://doi.org/10.1021/jp911032q.32] Smirnov, G. S., Stegailov, V. V., 2012. Melting and superheating of sI methane hydrate: Molecular dynamics study.
J. Chem. Phys. (4), 044523. https://doi.org/10.1063/1.3679860.[33] C.Y. Geng, H Wen, H. Zhou, Molecular Simulation of the Potential of Methane Reoccupation during the Replacement of Methane Hydrate by CO , J. Phys. Chem. A 113 (2009) 5463–5469. https://doi.org/10.1021/jp811474m.[34] E.M. Yezdimer, P.T. Cummings, A.A. Chialvo, Determination of the Gibbs Free Energy of Gas Replacement in SI Clathrate Hydrates by Molecular Simulation, J. Phys. Chem. A 106 (2002) 7982–7987. https://doi.org/10.1021/jp020795r.[35] N.I. Papadimitriou, I.N. Tsimpanogiannis, A.K. Stubos, A. Martin, L.J. Rovetto, C.J. Peters, Unexpected Behavior of Helium as Guest Gas in sII Binary Hydrates, J. Phys. Chem. Lett. 1 (2010) 1014–1017. https://doi.org/10.1021/jz9004625.[36] Susilo, R., Alavi, S., Ripmeester, J. A., Englezos, P., 2008. Molecular dynamics study of structure H clathrate hydrates of methane and large guest molecules. J. Chem. Phys. 128, 194505. https://doi.org/10.1063/1.2908074[37] B. Wathen, P. Kwan, Z. Jia, V. Walker, Modeling the Interactions between Poly(N-Vinylpyrrolidone) and Gas Hydrates: Factors Involved in Suppressing and Accelerating Hydrate Growth, In: D.J..K. Mewhort, N.M. Cann, G.W. Slater, T.J. Naughton (Eds.), High Performance Computing Systems and Applications, Springer, Berlin/Heidelberg, 2010, pp. 117-133.[38] Nilsson, A., Pettersson, L.G.M., 2015. The structural origin of anomalous properties of liquid water. Nat. Commun. 6, 8998. https://doi.org/10.1038/ncomms9998.[39] D. Eisenberg, W. Kauzmann, The Structure and Properties of Water, Oxford University Press, New York, 1969.[40] S. Alavi, K. Shin, J.A. Ripmeester, Molecular dynamics simulations of hydrogen bonding in clathrate hydrates with ammonia and methanol guest molecules J. Chem. Eng. Data 60 (2014) 389-397. https://doi.org/10.1021/je5006517[41] English, N.J., Clarke, E.T., 2013. Molecular dynamics study of CO hydrate dissociation: Fluctuation-dissipation and non-equilibrium analysis. J. Chem. Phys. 139, 094701. https://doi.org/10.1063/1.4819269.[42] S. Plimpton, Fast parallel algorithms for short-range molecular dynamics, J. Comput. Phys. 117 (1995) 1-19. https://doi.org/10.1006/jcph.1995.1039.[43] Subbotin, O.S., Adamova, T.P., Belosludov, R.V., Mizuseki, H., Kawazoe, Y., Rodger, P.M., Belosludov, V.R., 2009, Theoretical study of phase transitions in Kr and Ar clathrate hydrates from structure II to structure I under pressure. J. Chem. Phys. 131, 114507. https://doi.org/10.1063/1.3212965.[44] Belosludov, V.R., Bozhko, Y.Y., Gets, K.V., Subbotin, O.S., Kawazoe, Y., 2018. Clathratehydrates for energy storage and transportation. JPCS
Phys. Chem. Chem. Phys. (18) (2018) 12637-12641, https://doi.org/10.1039/C8CP01595C[47] Belosludov, V.R., Bozhko, Y.Y., Subbotin, O.S., Belosludov, R.V., Zhdanov, R.K., Gets, K.V.,Kawazoe, Y., 2018. Influence of N on Formation Conditions and Guest Distribution of Mixed CO + ₂ ₂ CH Gas Hydrates. ₄ Molecules (12), 3336. https://doi.org/10.3390/molecules23123336.[48] R.V. Belosludov, K.V. Gets, O.S. Subbotin, R.K. Zhdanov, Y.Y. Bozhko, V.R. Belosludov, J.I. Kudoh, Modeling the polymorphic transformations in amorphous solid ice, J. Alloy. Compd. 707 (2017) 108-113. https://doi.org/10.1016/j.jallcom.2016.12.197.[49] V.K. Malinovsky, R.K. Zhdanov, K.V. Gets, V.R. Belosludov, Y.Y. Bozhko, V.A. Zykova,N.V. Surovtsev, Origin of the anomaly in the behavior of the viscosity of water near 0° C, JETPletters (11) (2015) 732-736. https://doi.org/10.1134/S0021364015230095.50] M.H. Mir, J.J. Vittal, Phase transition accompanied by transformation of an elusive discrete cyclic water heptamer to a bicyclic (H O) cluster, Angew. Chem. Int. Ed. 46 (2007) 5925–5928. https://doi.org/10.1002/ange.200701779.[51] Matsumoto, M., 2007. Relevance of hydrogen bond definitions in liquid water. J. Chem. Phys. 126, 054503-6. https://doi.org/10.1063/1.2431168.[52] Kumar, R., Schmidt, J.R., Skinner, J.L., 2007. Hydrogen bonding definitions and dynamics in liquid water. J. Chem. Phys. 126, 204107. https://doi.org/10.1063/1.2742385.[53] G.K. Anderson, Enthalpy of dissociation and hydration number of carbon dioxide hydrate from the Clapeyron equation, The Journal of Chemical Thermodynamics , (7) (2003) 1171-1183. https://doi.org/10.1016/S0021-9614(03)00093-4[54] García Fernández, R., Abascal, J. L., Vega, C., 2006. The melting point of ice Ih for common water models calculated from direct coexistence of the solid-liquid interface. J. Chem. Phys. (14),144506. https://doi.org/10.1063/1.2183308. [55] E.N. Brodskaya, Molecular dynamics simulation of ice nanocluster in supercooled water,
Molecular Physics , (17-18) (2007) 2211-2216. https://doi.org/10.1080/00268970701516404.[56] A. Tokmakoff, Shining light on the rapidly evolving structure of water, Science 317 (2007) 54-55. https://doi.org/ 10.1126/science.1144515 . [57] Iwashita, T., Wu, B., Chen, W.R., Tsutsui, S., Baron, A.Q., Egami, T., 2017. Seeing real-space dynamics of liquid water through inelastic x-ray scattering. Sci. Adv. 3, e1603079. https://doi.org/ 10.1126/sciadv.1603079 .[58] S. Perticaroli, B. Mostofian, G. Ehlers, J.C. Neuefeind, S.O. Diallo, C.B. Stanley, L. Daemen, T.Egami, J. Katsaras, X. Cheng, J.D. Nickels, Structural relaxation, viscosity, and network connectivityin a hydrogen bonding liquid, Phys. Chem. Chem. Phys. 19 (2017) 25859-25869.https://doi.org/ 10.1039/c7cp04013j .[59] Elton, D.C., Fernández-Serra, M., 2016. The hydrogen-bond network of water supports propagating optical phonon-like modes. Nat. Commun. 7, 10193. h ttps://doi.org/10.1038/ncomms10193 .[60] L. Pauling, The structure of water. In: D. Hadzi and H. W. Thompson (Eds.), Hydrogen bonding,Pergamon Press, New York, 1959, pp. 1-6. Author contributions
K.G., Y.B. and V.B. designed the research idea. N. M., O.S. and R. Z. developed software program.R. Z. and K.G. analyzed calculation results. R.B. and Y.K. performed MD simulations on SC inJapan. All authors discussed the results and prepared the manuscript.
Corresponding authors
Dr. Kirill V. GetsE-mail address: [email protected] postal address: Nikolaev Institute of Inorganic Chemistry Siberian Branch of Russian Academyof Sciences, 3, Acad. Lavrentiev Ave., Novosibirsk, 630090, Russia