Molecular dynamics simulation of structure H clathrate-hydrates of binary guest molecules

Journal of Natural Gas Chemistry 20(2011)577–584

Molecular dynamics simulation of structure H clathrate-hydrates of binary guest molecules Hamid Erfan-Niya,

Hamid Modarress∗

Department of Chemical Engineering, Amirkabir University of Technology, Hafez Avenue 15914, Tehran, Iran [ Manuscript received March 14, 2011; revised July 27, 2011 ]

Abstract Molecular dynamics (MD) simulations are performed to study the stability of structure H clathrate-hydrates of methane+large-molecule guest substance (LMGS) at temperatures of 270, 273, 278 and 280 K under canonical (NVT-) ensemble condition in a 3×3×3 structure H unit cell replica with 918 TIP4P water molecules. The studied LMGS are 2-methylbutane (2-MB), 2,3-dimethylbutane (2,3-DMB), neohexane (NH), methylcyclohexane (MCH), adamantane and tert-butyl methyl ether (TBME). In the process of MD simulation, achieving equilibrium of the studied system is recognized by stability in calculated pressure for NVT- ensemble. So, for the accuracy of MD simulations, the obtained pressures are compared with the experimental phase diagrams. Therefore, the obtained equilibrium pressures by MD simulations are presented for studying the structure H clathrate-hydrates. The results show that the calculated temperature and pressure conditions by MD simulations are consistent with the experimental phase diagrams. Also, the radial distribution functions (RDFs) of host-host, host-guest and guest-guest molecules are used to analysis the characteristic configurations of the structure H clathrate-hydrate. Key words structure H clathrate-hydrate; methane+large-molecule guest substance; molecular dynamics; stability; radial distribution function

1. Introduction Clathrate-hydrates are nonstoichiometric inclusion compounds consisting of a three-dimensional host lattice of water molecules, in which guest molecules are encaged in polyhedral cells formed by hydrogen-bonded water molecules. Hydrates formed by natural gases and water, either in nature or in industrial processes, are commonly known as gas hydrates [1−3]. Gas hydrates naturally occur under conditions of high pressure and low temperature [4]. The clathrate-hydrate stability is due to the interactions between the water molecules forming a host lattice and the gas molecules occupying the cavities of the lattice as a guest molecule [1,5]. These interactions play an important role for clathrate-hydrates in being stabilized at suitable pressures and temperatures [6,7]. It is essential to study gas hydrates not only because of the great potential they have as a future energy source but also because of the problems they create to the petroleum industry during the production, transportation and processing of natural gas and oil [7]. During oil and gas production and transport in a cold climate, clathrate-hydrates can be easily formed in pipelines and production equipment. The traditional remedy for this problem is to inject large volumes (up to 40%) of ∗

alcohols and/or glycols into production streams [8,9]. Depending on the size and nature of guest molecules and the thermodynamic conditions, clathrate-hydrates are formed in three main types of structure I, structure II and structure H [1,10]. Structure I is a cubic structure with a lattice parameter ˚ [11] containing 46 H2 O molecules. The structure of 12.05 A I unit cell consists of two pentagonal dodecahedron (512) and six tetrakaidecahedron (512 62 ) [1]. Structure II is also a cu˚ [12] containbic structure with a lattice parameter of 17.3 A ing 136 water molecules and its unit cell contains eight large cages (512 64 ) and 16 small cages (512). Structure H is a hexag˚ and onal structure [13] with lattice parameters of a = 12.2 A ˚ There are 34 water molecules per hexagonal cell c = 10.1 A. with 3 small cages (512 ), 2 medium cages (43 56 63 ) and one large cage (51268 ) [14]. Structure H hydrate, as the newest structure of the clathrate-hydrates, was first discovered by Ripmeester et al. [13]. Stability of structure H hydrate requires the presence of a small help gas such as methane or xenon together with a large-molecule guest substance (typically having diameters ˚ like adamantane and methylcyclohexane between 7 and 9 A) (MCH). Whereas hydrates are a widely known problem in the natural gas industry due to their tendency to plug pipelines,

Corresponding author. Tel: +98-21-64543176; Fax: +98-21-6405847; E-mail: [email protected] (H.Modarress)

Copyright©2011, Dalian Institute of Chemical Physics, Chinese Academy of Sciences. All rights reserved. doi:10.1016/S1003-9953(10)60242-3

578

Hamid Erfan-Niya et al./ Journal of Natural Gas Chemistry Vol. 20 No. 6 2011

the existence of structure H hydrate causes a potential problem for the petroleum industry since several structure H formers, such as isopentane, neohexane, and methylcyclopentane (MCP) compose a small fraction of crude oil [15]. Phase equilibrium data of the structure H hydrates of methane and large-molecules guest substances (LMGSs) were investigated and reviewed by Sloan et al. [1]. With few exceptions, all structure H equilibria data were obtained with methane as the small component. Also, because there were always four phases present (LW -H(sH)-V-LHC ) with three components (including water), the Gibbs Phase Rule provides for only one composition of each phase which satisfies the equilibrium conditions. Consequently, all data were taken without measurement of any phase composition, and only the equilibrium temperature and pressure of the four phases were reported [1]. Initial phase equilibrium data of structure H hydrate was reported for methane+adamantane by Lederhos et al. [16]. Consequently, several authors published phase equilibrium data of structure H hydrate for a wide range of largemolecule guest substance (LMGS) in the presence of methane. Becke et al. [17] measured the high-pressure structure H equilibria of methane+methylcyclohexane (MCH). Some other authors [18−27] measured phase equilibrium data for a series of methane+LMGS systems with methane as a help guest. The temperature and pressure conditions of structure H hydrate stability were found to be consistent with those of hydrocarbon production, processing and transportation facilities [15]. The statistical thermodynamics model for hydrate equilibria has been explained formally by van der Waals and Platteeuw (vdWP) [28]. This theory is based on there being equilibrium between guest molecules in gaseous state and those confined in hydrate cage. It has been used successfully to predict phase behavior of many hydrates, both simple and complex (more than one type of guest molecule). However, some empirical parameters are always required in the calculations. It would be attractive to argue the thermodynamic stability of the gas hydrates only in terms of intermolecular interactions between water and guest molecules, thereby avoiding having to introduce the empirical parameters [1]. A direct application of vdWP theory sometimes leads to an incorrect phase diagram, in particular for large guest species [6]. Computer simulations, particularly molecular dynamics (MD) simulation, are a powerful tool for studying thermodynamic properties of water and gas hydrates [29]. More recently, computational methods have been used to study gas hydrates and in many of these cases, methane hydrates receive the most attention to validate the published experimental data on thermodynamics of hydrate formation [30−39]. Most molecular simulation studies have focused on the simple energy force fields, in which treating water molecules as rigid bodies and methane molecules as spherical particles [30−39]. Several molecular simulation studies on clathratehydrates have been reported regarding the stability of clathrate-hydrates. Rodger [40] applied the molecular simulation to investigate the van der Waals and Platteeuw statistical mechanical theory. Their studies discuss the importance of repulsive forces from the guest molecules on stabi-

lizing the hydrate water lattice [1]. Tanaka et al. [41−44] investigated the thermodynamic stability of the structure I and structure II hydrates using molecular simulations. They proposed a generalized version of vdWP theory, which was developed for each situation, and discussed the thermodynamic stability of clathrate-hydrates only in terms of intermolecular interactions. Their approach, however, requires an equation based on the statistical thermodynamics. Okano and Yasuoka [45] investigated the thermodynamic stability of the structure H hydrates of methane and various largemolecule guest substances under constant pressure and temperature (NPT-) ensemble using spherical one-site LennardJones (LJ) potential model for the guest substances. They calculated the free-energy difference changing on σ and ε only of the LMGS molecules. The use of one-site model for large-molecule guest substances may be useful as one can readily test various chemical species of large-molecule guest substances by changing the parameter set in the potential model. But, it is reported that the shape of the large-molecule guest substances is an important factor in determining the thermodynamics properties of the hydrates [1]. Miyoshi et al. [46] performed MD simulations on structure H hydrates formed with methane+2,2-dimethylpentane and methane+3,3dimethylpentane under NPT- ensemble conditions. They obtained the trajectories and the distribution profile of the sites of the dimethylpentane molecules. Their results demonstrated the distinct molecular motions of the two dimethylpentanes in the cages of hydrates. Also, Susilo et al. [47] employed MD simulations on structure H clathrate-hydrates of methane and large guest molecules to provide further insight regarding the dependence of methane occupancy on the type of LMGS and pressure. They found that thermodynamically, methane occupancy depends on the pressure but not the nature of LMGS. Their simulations confirm the instability of hydrate when small and medium cages are empty. It is concluded from literature review that not enough attention has been paid on the comparison of calculated equilibrium data by MD simulations with phase diagrams. In this research, the available phase equilibrium data for structure H clathrate-hydrate of methane+LMGSs are utilized to investigate the effects of molecular dynamics simulations for the successful prediction of structure H equilibria; where the large-molecule guest substances (LMGSs) are 2methylbutane (2-MB), 2,3-dimethylbutane (2,3-DMB), neohexane (NH), methylcyclohexane (MCH), adamantane and tert-butyl methyl ether (TBME). 2. Computational details 2.1. Interatomic potential model The force field for water interactions is transferable intermolecular potential 4 point (TIP4P) model [48], which is a four site water model with rigid atomic positions. The intermolecular potential for the interaction of a pair of molecules is of the following form

Journal of Natural Gas Chemistry Vol. 20 No. 6 2011

V (rO−O ) =

A r12 O−O

−

C r6O−O

qi qj +∑ i,j rij

(1)

where, i and j label the charge sites, rO−O is the distance between two oxygen atoms, A and C are the potential parameters. Interaction parameters of TIP4P model are shown in Table 1. Table 1. Parameters of TIP4P potential model [48] TIP4P 0.957 104.52 600.0 610.0 0 0.52 −1.04 0.150

˚ r(OH), A ∠HOH, deg A × 10−3 , kcal·A˚ 12 ·mol−1 C, kcal·A˚ 6 ·mol−1 q(O) q(H) q(M) ˚ r(OM), A

The van der Waals interactions between guest-guest and host-guest molecules are based on the Lennard-Jones (12−6) ˚ Coulombic interacpotential with a cutoff distance of 15 A. tions between point charges qi and qj located on the atomic nuclei i and j are used to model the electrostatic intermolecular interactions. The standard Lorentz-Berthelot combination rules [49], εij = (εii εjj )1/2 and σij = (σii +σjj )/2 are used to derive the Lennard-Jones potential parameters between unlike atom-type force centers i and j from the values of the parameters between similar atom types. Therefore, the force fields for guest-guest and host-guest interactions are expressed by the following intermolecular potential U (r) =

N −1 N

∑∑

i=1 j>i



 4εij

σij rij

12



σij − rij

6 

qi qj + 4πε0 rij



(2) The interaction potential parameters of studied LMGSs are available in Refs. [47,50−52]. 2.2. Simulation method The model system of methane+LMGS and water are built for the study of structure H clathrate-hydrate stability. The methane molecules are enclosed in both small and medium cages, while LMGS molecules are only enclathrated in large cage. Molecular dynamics simulations of structure H clathrate-hydrates are performed using LAMMPS code [53]. LAMMPS integrates Newton’s motion equations for collections of atoms, molecules, or macroscopic particles that interact via short- or long-range forces with a variety of initial and/or boundary conditions. In the first step of simulation, the initial positions of all atoms and molecules are determined. The initial configurations of the guest molecules and oxygen atoms of water molecules in the unit cell of structure H hydrate are obtained from the experimental data [45,54,55]. The hexagonal simulation cells for structure H clathrates consists of

579

3×3×3 replicas of the unit cell, with initial lattice parame˚ and c = 10.1 A. ˚ Guest molecules (methane ters of a = 12.2 A and LMGS) are initially placed in the center of the appropriate cages and their positions are allowed to equilibrate during the simulation. The motion equations are integrated using Verlet algorithm [56], and rigid water molecule constraints are implemented with SHAKE algorithm [57]. Also, the temperature is controlled by a Nos´e-Hoover thermostat [58,59]. The trajectories are generated in the canonical (NVT-) ensemble, with the integration time step of 1 fs, and the periodic boundary conditions are applied in all the directions. The temperatures are set at T = 270, 273, 278 and 280 K to simulate the thermodynamically stable hydrates. The total simulation time for each case is set on 300 ps. This simulation time is sufficient for obtaining converged simulation energies and pressures. Ewald summation method [60] is used for long-range electrostatic interactions with a precision of 1×10−6 and all intermolecular interactions in the simulation box are calculated ˚ One unit cell of structure within a cutoff distance of 15.0 A. H clathrate-hydrate contains 34 water molecules. In the 918 oxygen atoms of 3×3×3 MD cells, there are a total of 81 small cages, 54 medium cages and 27 large cages. We fill different composition of guest molecules (methane+LMGS) in each appropriate cage in all the simulations and after 24 simulations in order to obtain the thermodynamic properties, energies of systems and radial distribution functions. 3. Results and discussion In this work, hydrate systems of methane+LMGS with various compositions are simulated under canonical (NVT) ensemble conditions. The unit cells are placed inside of hexagonal box of simulation with size equal to the lattice parameters of structure H. The simulation box is held fixed while the structure evolved with time and equilibrated at the desired temperature. The main goal in NVT- calculation is to calculate equilibrium pressures of hydrates with different guests at different temperatures. The obtained data should be compared with experimental phase diagrams to investigate the existence of the obtained P-T data in the hydrate phase region. In Table 2, the equilibrium pressures and energies of systems are calculated by MD simulation at temperatures of 270, 273, 278 and 280 K, respectively. The obtained energies are kinetic energy, potential energy, total energy (sum of kinetic and potential energies), van der Waals energy, electrostatic energy and configuration energy (sum of van der Waals and electrostatic energies). All MD simulations reached an equilibrium state when the obtained energies of the system fluctuate over constant amount. These MD calculations are performed for various gas compositions of methane and LMGSs; since, the equilibrium hydrate structure is sensitive to the large-molecule guest substances and methane as well as their sizes. At equilibrium state, the structure of clathrate-hydrate is visualized to see whether it is in the stable form. Figure 1 is the network

580

Hamid Erfan-Niya et al./ Journal of Natural Gas Chemistry Vol. 20 No. 6 2011

of H2 O for clathrate-hydrate of neohexane (NH)+methane at the initial (t = 0 ps) and end of simulation (t = 300 ps). This figure is shown for whole calculation cell, where a unit cell has methane and large-molecule guest substances. As

seen in Figure 1, the hydrate structure is in stable condition. For other cases of methane+LMGS clathrate-hydrates, the results are similar.

Table 2. Thermodynamic properties of structure H hydrate of methane and large-molecule guest substance under canonical (NVT-) ensemble conditions Guest molecules Adamantane+methane

MCH+methane

TBME+methane

NH+methane

2-Methylbutane+methane

2,3-Dimethylbutane+methane

Temperature (K) 270 273 278 280 270 273 278 280 270 273 278 280 270 273 278 280 270 273 278 280 270 273 278 280

EvdW (kJ/mol) 2.1002 2.1031 2.0982 2.1037 1.8295 1.8128 1.8161 1.8227 1.7626 1.7612 1.6002 1.5708 1.8831 1.8622 1.8784 1.8718 2.1465 2.1416 2.1250 2.1402 2.0246 2.0398 2.0336 2.0331

Eelect (kJ/mol) −15.0595 −15.0397 −14.9926 −14.9893 −13.0243 −12.9725 −12.9165 −12.9099 −12.9893 −12.9922 −12.6381 −12.4937 −13.0929 −13.0479 −13.0330 −13.0119 −15.4808 −15.4391 −15.3765 −15.3683 −15.2941 −15.3008 −15.2130 −15.2014

Econfig (kJ/mol) −12.9593 −12.9366 −12.8944 −12.8856 −11.1948 −11.1597 −11.1004 −11.0872 −11.2267 −11.2310 −11.0379 −10.9229 −11.2098 −11.1857 −11.1546 −11.1401 −13.3343 −13.2975 −13.2515 −13.2281 −13.2695 −13.2610 −13.1794 −13.1683

Ekin (kJ/mol) 2.1720 2.1904 2.2295 2.2348 1.9696 1.9917 2.0353 2.0525 2.0324 2.0345 2.0658 2.1058 1.9986 2.0097 2.0414 2.0535 2.2623 2.2914 2.3207 2.3418 2.2442 2.2611 2.2854 2.2997

Epot (kJ/mol) −75.6888 −75.6655 −75.6229 −75.6142 −83.6354 −83.5991 −83.5384 −83.5249 −87.5349 −87.5380 −87.3391 −87.2199 −85.1974 −85.1727 −85.1415 −85.1264 −78.8630 −78.8259 −78.7789 −78.7557 −78.2187 −78.2100 −78.1278 −78.1167

Etot (kJ/mol) −73.5168 −73.4751 −73.3934 −73.3794 −81.6659 −81.6074 −81.5031 −81.4724 −85.5026 −85.5034 −85.2733 −85.1141 −83.1989 −83.1630 −83.1000 −83.0729 −76.6007 −76.5345 −76.4582 −76.4139 −75.9744 −75.9489 −75.8424 −75.8170

Pressure (kPa) 1450 1627 2605 3071 1237 1538 2733 3819 1405 1802 2545 3779 1074 1336 2509 3532 2007 2353 3880 4837 1185 1817 2930 3964

Figure 1. Time variation of network of water at the time of (a) 0 ps (initial condition) and (b) 300 ps (end of simulation) for NH+CH4 hydrate at T = 273 K, where green, blue and red colors are O atoms of water, NH molecules and C atoms of methane, respectively

According to the experimental data [1,15−27] for structure H clathrate-hydrates, the phase equilibria data for binaryguest mixtures containing methane and LMGSs are described. To indicate the consistency of the obtained data by MD calculations with phase equilibria data for each binary guest molecules, phase diagrams of pressure versus temperature are plotted in Figure 2. However, in most references [1,15−27], phase diagrams are represented in semi-logarithmic plot of

pressure versus temperature, where there is more than one data set. This kind of plot is useful when one of the plotted variables covers a large range of values and the other has only a restricted range. Since, we simulate only the hydrate phase; therefore, the obtained equilibrium pressures should be higher than the phase boundary and accurately be existed in (sH) hydrate phase. In order to better indicate the obtained equilibrium pressures (by MD simulations) in the hydrate phase, it is

Journal of Natural Gas Chemistry Vol. 20 No. 6 2011

reasonable to use simple linear plot. Figure 2 demonstrates that the obtained equilibrium P-T data exist in the (structure H) hydrate phase which is compatible with experimental phase diagrams. To evaluate the characteristic configurations of the crystalline structure H clathrate-hydrates, the radial distribution functions (RDFs) are used. The pair radial distribution functions (gij (r)) give the probability of finding a pair of particles i and j with distance r apart, relative to the probability expected for a completely random distribution of particles at the same density [61,62]. Therefore, the microstructure of

581

the structure H hydrate is described by the host-host, hostguest and guest-guest radial distribution functions gOO (r), gOG (r), gOM (r) and gGM (r) (Figure 3) for clathrate-hydrates including methane (M) and LMGS as a guest molecules in NVT- ensemble, where G denotes the large-molecule guest substances such as 2-methylbutane (2-MB), neohexane (NH), 2,3-dimethylbutane, adamantane, methylcyclohexane (MCH) and tert-butyl methyl ether (TBME). The obtained RDFs are at temperature of 280 K. These RDFs have been calculated to ˚ for the studied structure H clathratethe length of r = 15 A hydrates.

Figure 2. Pressure versus temperature plots showing a comparison between MD simulation results and experimental phase diagrams for structure H clathratehydrates. (a) Methane+2-methylbutane, (b) Methane+NH, (c) Methane+2,3-dimethylbutane, (d) Methane+MCH, (e) Methane+adamantane, (f) Methane+TBME

582

Hamid Erfan-Niya et al./ Journal of Natural Gas Chemistry Vol. 20 No. 6 2011

It may be seen that the host lattice of clathrate-hydrate structures has similar structural characteristics and is reasonably similar to previous computer simulation results [47,54]. The gOO (r) function shown in Figure 3(a) is defined with respect to the oxygen-oxygen distance between two water molecules. The gOO (r) is the most informative, since it invites comparison with the known pair occurrences which have regularly spaced oxygens radial distribution functions in the water and clathrate hydrates. At short distances (less than atomic diameter) g(r) is zero. The first peak location of oxygen-oxygen

˚ for all studied atoms of water molecules was at 2.775 A methane+LMGS hydrates; and this means that all the neighboring oxygen atoms are at the obtained distances. Therefore, the polyhedrons composed of oxygen atoms have the equal edges in hydrate structure [61,62]. The second peak located ˚ for all six studied cases. The third peak location at 4.575 A ˚ for 2-MB+methane and 2,3-DMB+methane was at 6.525 A ˚ for other studied methane+LMGS hyhydrates and at 6.375 A drates. The results for other host-host radial distribution functions (gOH (r) and gHH (r)) were similar.

Figure 3. Comparison of radial distribution functions obtained at 280 K. (a) gOO (r), (b) gOG (r), (c) gOM (r), (d) gGM (r); where G is 2-methylbutane (2-MB), 2,3-dimethylbutane (2,3-DMB), adamantane, methyl-cyclohexane (MCH), neohexane (NH) and tert-butyl methyl ether (TBME), respectively

The host-host radial distribution functions give us useful information about host lattice structure. As we know, in a crystal which referred to as a long-range order system, RDF has an infinite number of sharp peaks whose heights and separations are characteristic of the lattice structure [63]. This state cannot be seen in a gas or liquid which has disordered or shortrange order system. In a liquid, molecules are in an approximately close packed arrangement in which first two or three shells around given molecule are identified but well-defined positional relationships disappear gradually over longer range [64]. In real gases, molecules are essentially in random relative positioning and there is a single peak in g(r) at the point which potential energy has a minimum [65]. These points demonstrate that the host lattice structure preserves its crystalline stability during the simulation time.

The RDFs of oxygen atoms (water)-guest molecules are shown in Figure 3(b) and 3(c). Since the RDFs of oxygen atoms (water)-guest molecules are similar to those of hydrogen atoms (water)-guest molecules; we focus only the RDFs of oxygen atoms-guest molecules. The first peak position of oxygen-LMGS (OG) molecules oc˚ for 2MB+methane 2,3-DMB+methane and curred at 5.025 A ˚ for adamantane+methane NH+methane hydrates, at 4.575 A ˚ hydrate, at 3.825 A for MCH+methane hydrate, and at ˚ for TBME+methane hydrate. This peak is due to the 5.325 A water molecules that are composed of 512 68 cages (large cage) containing large-molecule guest substances. Also, the first peak in oxygen-methane radial distribution function gOM (r) ˚ for 2MB+methane, adamantine+methane located at 3.825 A ˚ for 2,3and MCH+methane hydrates, and at 3.675 A

Journal of Natural Gas Chemistry Vol. 20 No. 6 2011

DMB+methane, NH+methane and TBME+methane hydrates, which is represented in Figure 3(c) due to the distance between oxygen and methane in the structure. This peak is due to the water molecules that are composed of small (512 )and medium (43 56 63 ) cages containing small methane molecules. The RDF of guest-guest molecule gGM (r) is shown in Figure 3(d). The first peak position of LMGS-methane (GM) ˚ for 2-MB+methane, 2,3-DMB+methane and was at 7.575 A ˚ for adamantane+methane NH+methane hydrates, at 6.825 A ˚ hydrate, at 6.975 A for MCH+methane hydrate, and at ˚ for TBME+methane hydrate, which is due to nearest 7.425 A distance between LMGS molecule and methane in the small, medium and large cages of structure H clathrate-hydrate. Similar results can be obtained for other MD simulated guestguest radial distribution functions of gGG (r) and gMM (r). It is concluded from radial distribution functions (RDFs) that the selected methane+LMGS hydrate structures in this study are indeed stable over the length of our MD simulations under suitable thermophysical conditions. Moreover, the presence of the methane and large-molecule guest substances in the structure H clathrate-hydrate helps to more stability of overall structure. 4. Conclusions Molecular dynamics simulations have been employed to study the structure H methane+LMGS clathrate-hydrates. The effects of methane+LMGS binary guest molecules on the thermodynamic stability for the structure H hydrates have been clarified by evaluating the thermodynamic properties of the hydrates with various binary guest molecules of methane+LMGS. In this work, MD simulations of structure H clathrate-hydrates under NVT- conditions are performed. In these studies, a set of simulations with different binary guest compositions are carried out to obtain equilibrium pressures of structure H clathrate-hydrates. Then, the pressures at equilibrium of MD simulations are compared with experimental phase diagrams which indicates that the MD calculated P-T equilibrium data exist in the hydrate phase region. To analyze the structure of clathrate-hydrate system, we used the radial distribution functions and our results indicate that the distribution functions are shown as a clear evidence for stability of structure H clathrate-hydrate. As a result, the presence of methane+LMGS binary guest molecules in the clathratehydrate helps to increase the stability of structure H clathratehydrate. It is concluded that hydrate structures are stable in the simulated cases under appropriate thermo-physical conditions and the obtained equilibrium data by MD simulations are reliable according to previous experimental and theoretical results. References [1] Sloan E D, Koh C A. Clathrate Hydrates of Natural Gases. 3rd Ed. Boca Raton: CRC Press, 2008 [2] Wang L, Dong S L. J Nat Gas Chem, 2010, 19(1): 43

583

[3] Wang L, Dong S L. J Nat Gas Chem, 2007, 16(4): 423 [4] Lang X M, Fan S S, Wang Y H. J Nat Gas Chem, 2010, 19(3): 203 [5] Hao W F, Wang J Q, Fan S S, Hao W B. Energy Conv Manag, 2008, 49(10): 2546 [6] Tanaka H. Fluid Phase Equilib, 1998, 144(1-2): 361 [7] Gaudette J, Servio P. J Chem Eng Data, 2007, 52(4): 1449 [8] Kvamme B, Kuznetsova T, Aasoldsen K. Mol Simul, 2005, 31(14-15): 1083 [9] Kvamme B, Kuznestova T, Aasoldsen K. J Mol Graph, 2005, 23(6): 524 [10] Chihaia V, Adams S, Kuhs W F. Chem Phys, 2005, 317(2-3): 208 [11] McMullan R K, Jeffrey G A. J Chem Phys, 1965, 42(8): 2725 [12] Mak T C W, McMullen R K. J Chem Phys, 1965, 42(8): 2732 [13] Ripmeester J A, Tse J S, Ratcliffe C I, Powell B M. Nature, 1987, 325(6100): 135 [14] Klauda J B, Sandler S I. Chem Eng Sci, 2003, 58(1): 27 [15] Mehta A P, Sloan E D. AIChE J, 1996, 42(7): 2036 [16] Lederhos J P, Mehta A P, Nyberg G B, Warn K J, Sloan E D. AIChE J, 1992, 38(7): 1045 [17] Becke P, Kessel D, Rahmanian I. In: Proceedings of the European Petroleum Conference, 1992. 159 [18] Mehta A P, Sloan E D. J Chem Eng Data, 1993, 38(4): 580 [19] Mehta A P, Sloan E D. J Chem Eng Data, 1994, 39(4): 887 [20] Mehta A P, Sloan E D. AIChE J, 1994, 40(2): 312 [21] Thomas M, Behar E. In: Proceedings of the 73rd GPA Convention, 1994. 100 [22] Tohidi B, Danesh A, Burgass R, Todd A. In: Proceedings of the Second International Conference on Natural Gas Hydrates, Toulouse, France, 1996. 109 [23] Mooijer-van den Heuvel M M, Peters C J, de Swaan Arons J. Fluid Phase Equilib, 2000, 172(1): 73 [24] Nakamura T, Makino T, Sugahara T, Ohgaki K. Chem Eng Sci, 2003, 58(2): 269 [25] Ohmura R, Matsuda S, Uchida T, Ebinuma T, Narita H. J Chem Eng Data, 2005, 50(3): 993 [26] H¨utz U, Englezos P. Fluid Phase Equilib, 1996, 117(1-2): 178 [27] Makogon T Y, Mehta A P, Sloan E D. J Chem Eng Data, 1996, 41(2): 315 [28] van der Waals J H, Platteeuw J C. Adv Chem Phys, 1959, 2: 1 [29] Cheng W, Wu H C, Ye X Q, Zhou H Y. Prog Nat Sci, 2004, 14(11): 1015 [30] Alavi S, Ripmeester J A, Klug D D. J Chem Phys, 2005, 123(2): 024507 [31] English N J, MacElroy J M D. J Comput Chem, 2003, 24(13): 1569 [32] Moon C, Taylor P C, Rodger P M. J Am Chem Soc, 2003, 125(16): 4706 [33] Storr M T, Taylor P C, Monfort J P, Rodger P M. J Am Chem Soc, 2004, 126(5): 1569 [34] Chialvo A A, Houssa M, Cummings P T. J Phys Chem B, 2002, 106(2): 442 [35] Cao Z T, Tester J W, Sparks K A, Trout B L. J Phys Chem B, 2001, 105(44): 10950 [36] Zele S R, Lee S Y, Holder G D. J Phys Chem B, 1999, 103(46): 10250 [37] Forrisdahl O K, Kvamme B, Haymet A D J. Mol Phys, 1996, 89(3): 819 [38] Klauda J B, Sandler S I. Chem Eng Sci, 2003, 58(1): 27 [39] Klauda J B, Sandler S I. J Phys Chem B, 2002, 106(22): 5722 [40] Rodger P M. J Phys Chem, 1990, 94(15): 6080

584 [41] [42] [43] [44] [45] [46] [47] [48] [49] [50] [51] [52]

Hamid Erfan-Niya et al./ Journal of Natural Gas Chemistry Vol. 20 No. 6 2011

Tanaka H, Kiyohara K. J Chem Phys, 1993, 98(5): 4098 Tanaka H, Kiyohara K. J Chem Phys, 1993, 98(10): 8110 Tanaka H, J Chem Phys, 1994, 101(12): 10833 Tanaka H, Nakatsuka T, Koga K. J Chem Phys, 2004, 121(11): 5488 Okano Y, Yasuoka K. J Chem Phys, 2006, 124(2): 024510 Miyoshi T, Ohmura R, Yasuoka K. Mol Simul, 2007, 33(1-2): 65 Susilo R, Alavi S, Ripmeester J A, Englezos P. J Chem Phys, 2008, 128(19): 194505 Jorgensen W L, Chandrasekhar J, Madura J D, Impey R W, Klein M L. J Chem Phys, 1983, 79(2): 926 Poling B E, Prausnitz J M, O’Connell J P. The Properties of Gases and Liquids. 5th Ed. New York: McGraw-Hill, 2001 Papadimitriou N I, Tsimpanogiannis I N, Peters C J, Papaioannou A T, Stubos A K. J Phys Chem B, 2008, 112(45): 14206 Alavi S, Ripmeester J A, Klug D D. J Chem Phys, 2006, 124(20): 204707 Chang J, Sandler S I. J Chem Phys, 2004, 121(15): 7474

[53] Plimpton S J. J Comput Phys, 1995, 117(1): 1 [54] Pratt R M, Mei D H, Guo T M, Sloan E D. J Chem Phys, 1997, 106(10): 4187 [55] Susilo R, Alavi S, Ripmeester J, Englezos P. Fluid Phase Equilib, 2008, 263(1): 6 [56] Verlet L. Phys Rev, 1967, 159(1): 98 [57] Ryckaert J P, Ciccotti G, Berendsen H J C. J Comput Phys, 1977, 23(3): 327 [58] Nose S. J Chem Phys, 1984, 81(1): 511 [59] Hoover W G. Phys Rev A, 1985, 31(3): 1695 [60] Allen M P, Tildesley D J. Computer Simulation of Liquids. Oxford: Clarendon Press, 1987 [61] Ferdows M, Ota M. JKAU: Eng Sci, 2005, 16: 131 [62] Ferdows M, Ota M. Chem Eng Technol, 2005, 28(2): 168 [63] Konnert J H, Karle J, D’Antonio P. In: ASM Handbook, Vol. 10, Materials Characterization, Ohio, USA, 1986. 393 [64] Soper A K. Chem Phys, 2000, 258(2-3): 121 [65] Ben-Naim A. Molecular Theory of Solutions. New York: Oxford University Press, 2006