Calculations of NMR properties for sI and sII clathrate hydrates of carbon dioxide

Chemical Physics 433 (2014) 31–41

Contents lists available at ScienceDirect

Chemical Physics journal homepage: www.elsevier.com/locate/chemphys

Calculations of NMR properties for sI and sII clathrate hydrates of carbon dioxide Paweł Siuda, Joanna Sadlej ⇑ Faculty of Chemistry, University of Warsaw, Pasteura 1, 02-093 Warsaw, Poland

a r t i c l e

i n f o

Article history: Received 8 July 2013 In final form 15 January 2014 Available online 31 January 2014 Keywords: NMR DFT Clathrate hydrates

a b s t r a c t Nuclear shielding and spin–spin coupling constants (intra- and intermolecular) have been calculated for cages forming sI and sII clathrate hydrates of carbon dioxide (for all atoms of host and guest molecules). Structures of 512, 51262 and 51264 cages have been constructed using neutronographic data and DFT/ B3LYP calculations conducted with HuzIII-su3 basis set for NMR parameters determination. Based on those results it is possible to discriminate between CO2 molecules residing in each type of the cage. The analysis of NMR parameters calculated for water molecules is focused on their dependence on geometry of the molecular environment. It is possible to connect changes in NMR parameters with types of Hbond patterns present in cages of hydrates and the strength of H-bonds formed. Moreover, our results show that topologically differentiable water molecules forming cages are characterized by distinct NMR parameters, for example 17O shielding constants for water molecules of different topologies differ by 1.6 and 2.1 ppm for cages 51262 and 51264, respectively. This observation could be confirmed experimentally. Ó 2014 Elsevier B.V. All rights reserved.

1. Introduction Nuclear Magnetic Resonance spectroscopy is a widespread and well established method of structural analysis. Due to the sensitivity of measured parameters to the electronic structure, conformational changes of a molecule and its environment NMR spectra are a very valuable source of information on the molecular structure. In particular, NMR spectroscopy plays an important role in detection and characterization of hydrogen bonds, therefore this technique is being widely used now to study structures, conformations and properties of molecules [1]. Usually, magnetic shielding constants and indirect spin–spin coupling constants (SSCC) have been exploited for the structural NMR analysis. The size and changes of NMR parameters are very sensitive measure of the geometries as these parameters strongly depend on intermolecular bond distances, angles and environment. NMR chemical shifts are widely used to characterize the structure, to understand the nature of the guest–cage interactions, and dynamics of guest molecules in clathrate hydrates [2–4]. Clathrate hydrates (CH) are a class of crystalline non-stoichiometric inclusion compounds with cages formed by water molecules. There are three main clathrate hydrate’s structures sI, sII and sH [5,6].

⇑ Corresponding author. E-mail address: [email protected] (J. Sadlej). http://dx.doi.org/10.1016/j.chemphys.2014.01.007 0301-0104/Ó 2014 Elsevier B.V. All rights reserved.

The type of hydrate formed depends on the size and nature of the guest molecules. Unit cell of type sI CH is formed by 46 water molecules, made up of 8 polyhedral cages: two pentagonal dodecahedral, noted as 512 and six large 14-sided cages noted as 51262 [7]. The base of the notation designates the type of face, while the exponent – the number of faces of the same type. The ratio of large to small cages in sI structure is 3:1. Small guest molecules form structure sI, generally. Larger guest molecules form structure sII, with unit cell composed of 136 water molecules. They are made up of 24 polyhedral cages: 8 large 16-sided cages 51264 and 16 pentagonal small 12-sided cages 512, with the ratio of large to small cages in sII structure 1:2. The sH type contains six polyhedral cages: one large 51268, three medium 4351263, and two small 512 [8]. To support this hydrate structure, a mixture of molecules that differ in size is needed (small reside in smaller cavities, while bigger occupy big cavity). Carbon dioxide may form clathrate hydrates of structures denoted sI consisting of cages 512 and 51262, and, when mixed with other gases (N2, CH4, H2 [9,10]), sII – consisting of cages 512 and 51264. Faces of all cages are therefore pentameric or hexameric rings of water molecules. In three dimensional clathrate hydrate structure, every water molecule is located in a junction of four cages and is a member six rings of water molecules. Due to this fact, it is possible to distinguish water molecules forming only pentameric rings (denoted 56), five pentameric and one hexameric

32

P. Siuda, J. Sadlej / Chemical Physics 433 (2014) 31–41

(5561) and four pentameric and two hexameric (5462). In sI structure water molecules of all three topologies are present with 8:12:3 ratio, while in the structure sII waters of only 56 and 5561 topologies are found with the ratio of 5:12. Going from three dimensional crystal to single cage all 56 water molecules are becoming 53, all 5462 are becoming 5261 while half of 5561 molecules is attaining 53 topology and second half 5261. Therefore in single cages water molecules may have only two types of topology, that is 53 and 5261 (for more details please see paragraph 1 in Supplementary Materials). Fragments of three-dimensional crystalline CHs of sI and sII type showing topologically distinguishable water molecules are depicted on Figs. 1 and 2, respectively (see also movies presenting those structures in Supplementary Materials). In principle, atoms in molecules of different topologies, being crystalographically distinguishable [11], are surrounded by different electronic densities and therefore should have different values of NMR

parameters; however, the differences may be small and difficult to observe experimentally. It has been demonstrated both experimentally (employing macroscopic and microscopic techniques) and theoretically (ab initio, molecular dynamics and Monte Carlo simulations) that there is a great number of guest molecules stabilizing CH structures. Among them, noble gases (Ne, Ar, Kr, Xe) [3,12,13], small homo- and heteronuclear molecules (H2, N2, CO, CO2) [13–16] and hydrocarbons (methane, ethane, propane) [17,18], that lack strong interaction with water. On the other hand, there are many dipolar, organic (formaldehyde, tetrahydrofuran, ethylene oxide) [19,20] and inorganic molecules (HCN) [21], that are capable to form hydrogen bonds with water molecules. Natural ‘‘gas hydrates’’ are considered to be an unexploited fuel source for the future, while synthetic hydrates are recognized as novel materials (hydrogen storage, cool energy storage) [22,23] and may play a role in the

Fig. 1. Topologically distinguishable water molecules in sI clathrate hydrate. Green – 56, red – 5561, blue – 5462. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Fig. 2. Topologically distinguishable water molecules in sII clathrate hydrate. Green – 56, red – 5561. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

P. Siuda, J. Sadlej / Chemical Physics 433 (2014) 31–41

effort to reduce greenhouse effect (carbon dioxide sequestration). Methane hydrate represents the most important species, as a potential future energy resource and geo-hazard for climate stability caused by gas hydrate dissociation [7,24]. The second most important guest gas forming CH is carbon dioxide. It is a major greenhouse gas and much effort is devoted currently to lower its anthropogenic emission to the atmosphere. One of the options is storage of CO2 in the form of clathrate hydrate. According to molecular dynamics calculations performed by Alavi and Woo [25] the sH CO2 are even more stable than the methane clathrate under the same conditions. Another study shows that methane could be replaced by CO2 in naturally occurring hydrates [26], what would enable both methane extraction and CO2 sequestration at the same time. Due to those reasons CO2 hydrates are extensively studied experimentally [14] and theoretically [27]. Spectroscopy (IR or Raman vibrational spectra and NMR spectra) is widely used to trace the growth of hydrates in environments close to those found in natural occurrences [28], as well as to study weak interactions of the guest with the cage walls. Previously, we have demonstrated the use of calculated NMR parameters for analysis of the molecular interactions of the methane and host water molecules [2]. The calculations of the chemical shifts and, for the first time, the indirect intra- and intermolecular spin–spin coupling constants (SSCCs) show that the environment of the encapsulating cage noticeably affects the parameters of the methane molecule in the frame of the static model. In the case of CH with CO2 the small 512 cage is common to all sI, sII and sH structures. Using 13C NMR chemical shifts it was found that the non-spherical CO2 molecule could be a guest molecule in the small cages in the three structures [3]. Some studies have found that the interactions between CO2 and water molecules affect the dynamics of the host molecules and the stability of the clathrate considerably [29,30]. In addition, it is known that experimental NMR spectroscopy display strong guest–water interaction arising from the rotational motions of guest CO2 molecules [31]. Further experimental investigations has proven, that the anisotropy of NMR signals of CO2, which is a consequence of non-uniform motion of CO2 inside the cavities of the hydrate, could be successfully modeled by means of molecular dynamics technique [4]. This work on guest molecule in CHs is a continuation of our previous study of NMR parameters in water clusters [32] and methane clathrates [2]. The calculations of the isotropic and anisotropic chemical shifts of the CO2 and H2O molecules and SSCCs for cages found in sI and sII structures (small 512 found in both sI and sII structures and larger cages 51262 of sI and 51264 of sII structure) are presented. Our main interest is in the basic aspects of CHs, most notably the understanding of the relation between NMR parameters of the guest molecule and the structure of the cage. There are two possibilities to analyze our results, namely the cages may be perceived simply as clusters of CO2 interacting with water molecules and as models of real carbon dioxide CHs. When they are seen as clusters the analysis of our results will be performed in terms of DAA (one donor–two acceptor type) and DDA type (double donor–one acceptor of proton) water molecules. As for the second view few remarks should be made. In real crystals almost all water molecules (except for defects) are four coordinated of DDAA type (double donor–double acceptor type water molecules), while in our models all water molecules are three coordinated (of DAA and DDA types) [32]. As was already stated, crystals of sI and sII hydrates are built of three types of topologically different water molecules denoted 56, 5561 and 5462. For 512 cage all water molecules are of 53 topology, while for 51262 and 51264 cages water molecules of both 53 and 5261 topologies are present. Those differences in H-bond frameworks and topological characteristics of water molecules between single cages and three dimensional structures should be considered to properly

33

compare our results with experimental ones (see Supplementary Materials). Moreover, external cages present in the three dimensional crystalline structure of real clathrate hydrates may also affect the NMR parameters of guest molecule. However, it was found recently for methane clathrate hydrate, that the presence of external cages influences NMR parameters negligibly [33]. As the structures of CO2 and methane clathrate hydrates are analogous, we assume that the above observation should be true for CO2 clathrate hydrate too. As the influence of external cages do not have important influence on NMR properties of guest molecules, we have also decided not to include analysis of the influence of different proton configurations on those properties, as it may be assumed, that it will also be limited. Specific configurations could influence NMR properties for water molecules involved, but in the bulk phase most of those specific features would be lost due to Boltzmann averaging of all possible proton configurations present. Moreover, our systems are static, in the sense, that we are using single geometries for all cages. In real CHs crystals, due to rotational freedom of CO2 molecules inside cages [34], many orientations of CO2 are possible. Experimental NMR parameters are therefore averages of all orientations present in the real system weighted by the Boltzmann’s factors of energy. In principle, to reconstruct those experimental values, many possible structures should be generated (with the use of, for example, molecular dynamics), NMR parameters for them calculated and summed up with proper weights. It should also be mentioned, that many hydrogen configurations in the cages are possible and they should be averaged as well in order to represent real hydrate in detail. The last comment on the limitations of our model addresses rovibrational effects. We do not calculate them in this work, but wherever experimental or computational data exist, they are included in the discussion of the results. It has been demonstrated by us previously for methane clathrates [2], that despite of all limitations of this model, NMR parameters obtained are very close to experimental ones for crystalline CH structures. Our aim here is to expand this knowledge to CO2 behavior in clathrate hydrate. The structure of this paper is as follows: first, the methods employed for geometry optimization and the calculations of NMR parameters are described. Section 3 presents the calculations and the discussion of the results. A brief summary is presented in the last section.

2. Computational details 2.1. Geometry optimization As the model systems we have chosen structures of cage 512 found in sI and sII crystalline hydrates, 51262 cage of sI hydrate and 51264 cage characteristic for sII hydrate. Locations of water’s oxygen atoms residing in vertices of the cages were based on the neutronographic data [7]. Then, the hydrogen atoms were added according to Bernal and Fowler [35] ice rules. In the last step CO2 molecules were inserted into each cage. During the geometry optimization positions of water oxygen atoms were frozen to preserve the overall structure characteristic for clathrate hydrates. DensityFunctional Theory (DFT), using the hybrid three-parameter Becke– Lee–Yang–Parr (B3LYP) functional [36,37] with the basis set aug-cc-pVDZ [38] was employed. This method was shown to give reliable results for similar systems [2]. During the geometry optimization no counterpoise correction were made for the basis set superposition error (BSSE). At each stationary point, vibrational frequencies were calculated in order to confirm the nature of the stationary point on the potential energy surface. This study showed, that 512 cage is the minimum, while 51262 and 51264

34

P. Siuda, J. Sadlej / Chemical Physics 433 (2014) 31–41

structures had 4 and 6 negative frequencies respectively, being saddle points. As our cages were designed to represent structures present in bulk crystal, we were not interested in finding the global minima. We wanted to obtain structures similar to most abundant ones in real clathrate crystals, so we have limited ourselves to obtain stationary geometries. Gaussian 03 package [39] was used to conduct geometry optimization calculations. 2.2. NMR parameters calculations To accurately predict NMR parameters computationally, proper description of electronic density close to the nuclei is required [1]. On the other side, presence of H-bonds enforces usage of diffuse functions to accurately describe the influence of strong interactions on the electronic density. According to the results of our former studies on similar systems [2,40] which have proven to be reliable, DFT/B3LYP approach and HuzIII-su3 [41] basis set dedicated to NMR parameters were chosen. It is well known that conventional DFT functionals as B3LYP, produce values that are too deshielded relative to experiment and to the best ab initio calculated values [42], therefore we decided to discuss the relative changes under complexation only. On the other hand, this methodology is inexpensive. Calculations were conducted using the Dalton [43] package. 3. Results and discussion 3.1. Structures Stationary structures obtained for all three cages are illustrated in Fig. 3. The smallest one, denoted 512, consists of 20 water molecules forming 12 pentagonal faces. The cages, denoted 51262 and 51264, are built of 24 and 28 water molecules, respectively. In addition to 12 pentagonal faces they consist of two and four hexagonal faces. Every water molecule building the cages forms 3 hydrogen bonds. Half of them act as double donor and single acceptor of proton (we denote them DDA), while the other half donate one and accept two protons (we denote them DAA) [32]. This arrangement of hydrogen bond network is typical to polyhedral water clusters (PWC’s) [44]. The DAA-type water molecules have one OH bond not involved in hydrogen bond (dangling OH bond). The DDA-type water molecules, on the other hand, have a lone electron pair not accepting a hydrogen bond. In real crystals of clathrate hydrates dangling OH bonds and oxygen lone electron pairs are not abundant – they exist only in defects of the ordered periodic structure. The number of such defects (and, therefore, the number of dangling OH bonds and lone electron pairs) in real crystals is negligible compared to the total number of host water molecules. In our structures it is however n/2 (where n is the total number of host water molecules). This feature needs to be taken into account while treating our structures as models of crystalline clathrate hydrates. Moreover, cage 51262 has additionally two substantially weaken H-bond between host water molecules, which could result in formation of Bjerrum L-defects [45], in which water molecules have one OH bond turned toward the oxygen of carbon dioxide molecule. This defect is caused by a strong interaction between guest carbon dioxide molecule and host water molecule. This attraction is strong enough to weaken H-bond between neighboring water molecules and eventually break it. In our case two such situations are clearly observed. One of them is depicted on Fig. 4. This kind of interaction may in principle led also to chemical reaction of carbonic acid or carbonate formation, as showed by Stirling et al. [46]. Comparison of average H-bond parameters between water molecules and corresponding parameters of ‘defected’ H-bonds is pre-

Fig. 3. Structures of cages 512, 51262 and 51264 with enclathrated carbon dioxide molecule: (a) 512, found in sI and sII hydrates (b) 51262, found in sI hydrates (c) 51264, found in sII hydrates.

sented in Table 1. Significant elongation of the hydrogen bond length and the change of the OAHAO angle are clear signs of substantial weakening of H-bonds caused by the interaction with CO2

35

P. Siuda, J. Sadlej / Chemical Physics 433 (2014) 31–41

found in three dimensional crystalline structure than it is for Hbond characteristics. Unfortunately, there is still no experimental data available for topologically different water molecules forming CHs. 3.2. Shielding constants of CO2 Shielding constants for the 13C and 17O for three cages and CO2 monomer are presented in Table 2 together with the available experimental results.

Fig. 4. Fragment of the cage 51262. Interaction of CO2 oxygen atom with hydrogen atom of water molecule bends and as a consequence weakens H-bond between water molecules containing oxygen atoms marked as O1 and O2. Water molecules containing atoms O1 and O3 form typical, strong H-bond.

Table 1 Comparison of average and defected H-bonds for 51262 cage.

Average H-bond Defect 1 Defect 12

rOH (Å)

OAHAO angle (°)

1.781 1.850 1.953

166 144 138

Table 2 Comparison of calculated and experimental (see details in text) 13C and 17O shielding constants (in ppm) for carbon dioxide molecules in 512, 51262 and 51264 cages and gaseous carbon dioxide (geometry optimized at B3LYP/aug-cc-pVDZ, NMR at B3LYP/ HuzIII-su3). Atom type

In:

dexp TMS

C

Cage 512 62 Cage 512 62 Cage 512 62 Gas phase

123.1 127.8 124.6 125.5

O

Cage 512 62 Cage 512 62 Cage 512 62 Gas phase

rcalc [10] [10] [10] [47–49]

225.8 [50]

dexp CO2

45.2 48.3 46.8 46.3

2.4 2.3 0.9 0.0

186.6 199.3 203.1 205.5

dcalc CO2 1.1 2.0 0.5 0.0 18.9 6.2 2.4 0.0

molecule. Detailed structural parameters of the cages are given in Tables 3 and 4. Even more important for comparison of calculated constants with the experimental data is the analysis based on topological characteristics of water molecules. Topology of water molecules forming single cages may be more directly connected with those

3.2.1. Carbon atoms The absolute value of the shielding constant for CO2 in the gas phase (monomeric) that we have obtained (46.3 ppm) is close to the best value for DFT/B3LYP (48.0 ppm) published by Auer et al. [51] with larger basis set and slightly different interatomic distances. Benchmark computations at CCSD(T) level with very large basis set (see [51] for details) gave 60.2 ppm. According to our calculations, the shifts in the 13C shielding constants upon complexation (with respect to the CO2 monomer) are equal 0.9, 2.0 and 0.5 ppm for 512, 51262, and 51264 cages, respectively. The greatest change is found for the intermediate cage 51262. The 13C shielding constants for the 512 cage is even smaller that value for monomeric carbon dioxide. Thus, the interaction of carbon dioxide with different cages can have opposite influence on 13C shielding constants. Transfer from the gas phase to the cage of a clathrate may weaken or enhance shielding effect, depending on the cage type. Similar experimental observation was made for propane sII hydrate [18]. Propane carbons in clathrate hydrate show inverted carbon shielding constants when compared to the gas phase. Hence, it is possible to observe the inversion of the shielding constant changes for both hydrophobic and hydrophilic guest molecules. As opposed to methane and propane, CO2 molecule is not hydrophobic and does not tend to minimize the interaction with host lattice. In all three cages, CO2 molecule is interacting favorably with surrounding water molecules. This interaction is even promoting the formation of L-defects in hydrogen bond network in the cage 51262. This may lead to a conclusion, that hydrogen bonds network is stronger in the 512 than in 51262 and 51264 cages. It may be attributed to the presence of six-membered water rings in the bigger cages. The H-bonds in flat six-membered water rings are weaker than in five membered rings because of higher tension. The angle between H-bonds equals to ca. 120° and 108°, respectively. The second value is much closer to tetrahedral, what makes 512 structure more stable than 51262 cage. Let us compare now the obtained results with the available experimental data. It needs to be stressed, that no vibrational correction is included in our calculated values. According to Auer et al. [51], zero-point vibrational correction (ZPE) for 13C shift at CCSD(T)/quadruple-zeta level of theory amounts to 1.34 ppm. This correction is of systematic nature, so it will not affect the rel-

Table 3 Average absolute 17O and 1H shielding constants (in ppm) and covalent OAH bond lengths (in Å) for water molecules forming cages 512, 51262 and 51264, divided according to type of H-bond pattern and topological criteria. Star denotes dangling hydrogens (geometry optimized at B3LYP/aug-cc-pVDZ, NMR at B3LYP/HuzIII-su3). Standard deviation in % of average value in parentheses. Molecule type

512

51262

51264

rO

rH

rOH

rO

rH

rOH

rO

rH

rOH

DAA

271.3 (2)

22.3 (5) 29.5⁄ (1)

1.015 0.979⁄

279.4 (1)

23.7 (5) 30.1⁄ (1)

0.996 0.970⁄

281.7 (1)

23.9 (5) 30.4⁄ (1)

0.990 0.965⁄

DDA 53 52 61

277.0 (2) 274.1 (2)

24.7 (4) 23.9 (6)

0.995 1.002

281.7 (2) 279.8 (1) 281.4 (2)

25.8 (5) 24.8 (7) 25.3 (6)

0.982 0.987 0.986

287.6 (1) 282.9 (2) 285.0 (1)

25.9 (5) 24.8 (5) 25.3 (6)

0.977 0.979 0.983

Values for monomer water, all in ppm: calculated (this work): rO = 325.00; rH = 31.34 exp.: [52] rO = 322.81; rH = 30.05.

36

P. Siuda, J. Sadlej / Chemical Physics 433 (2014) 31–41

Table 4 Calculated and empirical JXY (in Hz) and O  O and CAO distances (in Å) for CO2 molecules enclathrated by 512, 51262 and 51264 cages and for gas phase CO2 molecule (geometry optimized at B3LYP/aug-cc-pVDZ, NMR at B3LYP/HuzIII-su3). 512 1

JCO

51262 calc

24.63 32.87 2

JOOcalc 1.35

rCO

1

1.155 1.193

20.61 24.06

rOO 2.348

JCO

2

calc

JOOcalc 1.55

51264 1

rCO

1

1.152 1.170

25.95 25.70

1.164 1.171

27.86 emp 3.47 [53] calc 22.27 [53]

rOO 2.322

ative changes of shielding constants. What may affect the comparison, however, is the difference in conditions used for different experiments. For a more detailed discussion of this issue, see Supplementary Materials. Experimental solid state 13C NMR spectra for CO2 clathrate hydrates of sI and sII types were measured by Seo and Lee [10] at 243 K, leading to 123.1, 127.8 and 124.6 ppm, for 512, 51262 and 51264 cages, respectively. The empirical 13C chemical shift of 125.5 ppm for gas phase molecule presented in Table 2 was estimated from the experimental value [47] (referenced to TMS) measured at 300 K by including temperature corrections for both CO2 and TMS molecules. According to Morin et al. [49] TMS signal would shift going from 300 to 240 K by 0.6 ppm. No analogous data exists for CO2, but there is a paper describing OCS [48], for which temperature dependence of 13C shielding constant may be expected to be of the same order of magnitude. This correction amounts to around 0.1 ppm, what in total gives empirical 13C chemical shift of 125.5 ppm (with respect to TMS) for gas phase CO2 at 240 K. Thus, the experimental shift induced by encapsulation for the cage 51262 is very close to the calculated one (2.3 as compared to 2.0 ppm). Two other values for the cages 512 and 51264 differ slightly more, ca. 1.3 ppm (2.4 vs. 1.1 ppm and 0.9 vs. 0.5 ppm, respectively), but the overall trend is well preserved. Experimental anisotropies of NMR signals for CO2 clathrate hydrates were reported by Ripmeester and Ratcliffe [3] and Seo and Lee [10]. Comparison of those data with our results could serve as another test for validity of approximations of our model. According to IUPAC recommendation [56], shielding anisotropy equals f = rZZ  riso. As cages 512 of sI and 51264 of sII provide cavities of roughly spherical symmetry no orientation of guest molecule is favored, therefore isotropic NMR signals are expected. Cage 512 found in sII structure has slightly distorted symmetry that results in anisotropy of the NMR signal. However, as explained in the Introduction, we use sI 512 cage to model both sI and sII structures as we are concentrated on shielding constants and spin–spin coupling constants in this work. Therefore only for sI 51262 cage our result may be confronted with experimental result. As the CO2 molecule is rotating in equatorial plane (perpendicular to sixfold symmetry axis) of 51262 cage, static value obtained in our calculation should be adjusted according to relation

fexp ¼ 0:5ð3cos2 h  1Þfcalc before comparison to experimental result. Theta is the angle between the principal component of the shielding tensor and rotation axis, and in real system is averaged over the motion of the guest. In our case the theta equals 65.1° corresponding to fexp = 54.6, what is very close to 63.1 [3] and 55.3 [10] reported in the literature for spectra obtained around 270 K. It may be stated then, that cage 51262 represents the situation found in hydrate sI around this temperature (for more data on anisotropy see Table 1 of Supplementary Materials).

JCOcalc

Gas phase CO2

rCO

2

JOOcalc 1.80

rOO 2.335

JCO

2

JOO 1.97 calc 2.36 [55]

rCO 1.169 1.160 [54] 1.160 [54]

exp exp

rOO 2.338 exp 2.324 [54]

3.2.2. Oxygen atoms Oxygen shielding constants presented in Table 2 do not show the pattern observed for the carbon atom. The values of chemical shift (with respect to monomer value – last column of the table) form a decreasing trend starting from 512 encaged molecule through 51262 to 51264 encaged molecule. The induced changes (18.9, 6.2 and 2.4 ppm, respectively) can be explained in terms of size of the cages: in 51264, CO2 molecule has more space and is able to rotate almost freely. It is not possible in the smallest cage, where the radius of available cavity equals to 3.95 Å [7] and is comparable to van der Waals size of CO2, what precludes free rotation. Therefore, it should be expected that the interaction with the host lattice 51262 and 51264 is weaker and parameters such as nuclear shielding constants are less affected with respect to gas phase. Experimental 17O chemical shift are available only for gaseous CO2, therefore it is not possible to confront calculated effect of enclathration of CO2 molecules with experimental value. However, it may be stated that DFT/B3LYP is giving reliable results for 17O shielding constants, underestimating exact values by 20 ppm. As a comparison, best to date results at CCSD(T)/pz3d2f level overestimate the actual shielding constant by around 4.7 ppm [57].

3.3. Shielding constants of H2O It was already mentioned, that water molecules forming cages may be divided into two groups (DAA and DDA), depending on their H-bonds characteristics. Table 3 presents values of shielding constants averaged in each group of water molecules for oxygen and hydrogen atoms. Hydrogen atoms of DAA waters were additionally divided into two groups: one, which contain hydrogen atoms involved in H-bonding, and the other, so called dangling OAH bonds (dangling hydrogens). Standard deviations are also provided as a percent of average value of the shielding constant.

3.3.1. Oxygen atoms The experimental absolute shielding constant for oxygen in gaseous water equals to 322.81 Hz [52] DFT/B3LYP/HuzIIIsu3 result of 325.0 Hz (this work) is very close to above cited value. Let us now discuss the environment-induced changes in the 17O shielding constant of water molecules in CHs. Firstly, the calculated oxygen shielding constants, presented in Table 3, differ significantly both between cages and different H-bonding patterns. However, with growing size of the cage absolute shielding constants are monotonically growing toward the experimental value (322.81 ppm) for both types of water molecules. Secondly, the shielding constants decrease at two types of water molecules, DAA and DDA, interacting with CO2 in three investigated cages: 53.7 and 48.0 ppm for 512, 45.6 and 43.3 for 51262, 43.3 and 37.4 ppm for 51264 cages (the reference state is the calculated r(17O) = 325.0 ppm), respectively. Thus, the interaction-induced

P. Siuda, J. Sadlej / Chemical Physics 433 (2014) 31–41

shifts of the oxygen nuclei are dependent on environment and are smaller for the DDA than for DAA-type water molecules. This trend may be interpreted as a result of decreasing strength of interaction between water molecules as the cage gets larger. There are two observations supporting this statement: (i) larger cages contain more hexagonal faces which are geometrically less comfortable for the formation of H-bonds between water molecules. In pentamers angle between neighboring bonds is around 108°, while in hexamers goes up to around 120°. The tetrahedral angle of water is close to 105°, so pentamers are much less strained that hexamers. That steric inconvenience is leading to decreasing strength of H-bonds in hexamers. (ii) Another evidence of this tendency is clearly seen in the changes of average length of covalent OAH bonds involved in the formation of H-bonds between water molecules. Table 3 summarizes average OH bond lengths for both types of water molecules found in all three cages. For DAA molecules OH bond length is decreasing from 1.015 Å in cage 512 to 0.990 Å in cage 51264. Similarly, for DDA molecules it is lowering from 0.995 Å in cage 512 to 0.977 Å in cage 51262. Decreasing of the intramolecular OH distance is concerted with the weakening of the H-bond strength as the cage size grows. As mentioned in the Introduction, topology of the cages enables another distinction between water molecules forming cages. As each molecule resides in a vertex of the polyhedron built of pentamers (in cages 512, 51262, 51264) and hexamers (in cages 51262, 51264), some of them are neighboured by 3 pentamers (denoted 53, Fig. 5), while others by 2 pentamers and one hexamer (denoted 5261, Fig. 6). Analysis of OH bond lengths for molecules of 53 and 5261 types, presented in Table 3, leads to the conclusion consistent with those drawn in former paragraphs. As more hexamers are added to the structure and the size of the cage grows, OH bond are getting shorter and H-bonds formed between water molecules are getting weaker. Absolute shielding constants for oxygen atoms grouped according to topological criterion show monotonic growth toward value characteristic for H2O monomer, as seen from Table 3. Average values for 53 molecules are close to those observed for DAA waters, while 5261 are close to DDA. However, differences between 5261 and 53 are smaller than differences between values for DDA and DAA molecules.

Fig. 5. Fragment of the cage 512 depicting water molecule of 53 topology. Central water molecule is located in the vertex of three pentamers.

37

Fig. 6. Fragment of the cage 51262 depicting water molecule of 5261 topology. Central water molecule is located in the vertex of two pentamers and one hexamer.

3.3.2. Hydrogen atoms Absolute shielding constants of protons are presented in Table 3. Three groups of atoms are distinguished for each cage: protons involved in H-bonding, protons belonging to DDA or DAA water molecules and the so called dangling protons. Dangling proton’s values of SC are the closest to the monomer value (differences are 1.8, 1.2 and 0.9 ppm for cages 512, 51262 and 51264, respectively). The H-bonded protons in DAA water molecules exhibit greater changes in proton’s SC (9.0, 7.7 and 7.4 ppm) than those in DDA molecules (6.6, 5.5 and 5.4 ppm). Within each group of protons, differences between cages are small – much smaller than analogous differences between shielding constants for oxygen nuclei. It would not be possible to distinguish cages on the basis of those values, especially for 51262 and 51264 cages, as corresponding shielding constants are almost equal. However, differences between average shielding constants of groups of protons are more significant. Protons involved in H-bond formation and belonging to DDA water molecules are characterized by greater shielding constants than analogous protons of DAA water molecules by at least 2 ppm for all three cages. Dangling protons have shielding constants greater than those of DDA hydrogen atoms by not less than 4 ppm. Differences between dangling and H-bonded protons are clear also from comparison of the relevant average of shielding constants for those two groups. Non H-bonded protons show very small deviation in all cages (around 1%), what can be treated as a result of homogeneous structural and electronic environment within this group. Standard deviation for both groups of H-bonded protons is of the order of 5%. This is an evidence of much greater diversity of structural and electronic environments in those groups than in non H-bonded protons. Protons could also be divided according to topological criteria into groups belonging to 53 and 5261 water molecules. As no dangling protons are present in the crystal structure (apart from defects) only shielding constants of H-bonded protons were averaged to obtain the data presented in Table 3. Although differences between values for protons of 53 and 5261 water molecules are small (0.5 ppm for both cages), they are within the accuracy of spectrometers. It should be possible to observe distinct signals

38

P. Siuda, J. Sadlej / Chemical Physics 433 (2014) 31–41

from protons of topologically different water molecules forming crystal hydrates (more on that in Supplementary Materials).

A second intramolecular coupling constant, which could be interesting for characterization of the guest molecule in CH is the 2 JOO. The induced changes in the 2JOO upon enclathration (with respect to the CO2 monomer) are much smaller than for 1JCO, ranging from ca. 0.6 to ca. 0.2 Hz. There is no published experimental data for the 2JOO in CO2. San Fabian et al. [55] report a value of 2.36 Hz, calculated at the MCSCF/BS2 level, which agrees well with 1.97 Hz, obtained in this work. The predominant contribution to the 1JCO coupling is the FC term, which is responsible for around 90% of the total value (see Table 2 of Supplementary Materials). Among other three components, only PSO and SD are of importance, but as they are of opposite size, their influence on the total coupling is limited. For 2JOO coupling, FC and SD components are of similar size but opposite signs. Therefore, final value is very close to PSO term (much smaller than FC and SD), as DSO is negligible. Differences in CAO bond lengths and SSCC in the enclathrated CO2 molecule can be correlated with the changes of molecular environment. Analysis of the influence of the closest O  O and O  H contacts on SSCC in CO2 molecule shows that interactions with oxygens of water are much more important for SSCC (and CAO bond lengths), even though O  H distances are in all cases shorter. Generally, interaction with few closest water oxygens is more important than interaction with a single, closest water oxygen. Those two correlations are not surprising, and data they are based on can be found in Table 3 of Supplementary Materials. Intramolecular OH and HH coupling constants of water molecules are listed in Table 5. First, we discuss the effects caused by the formation of hydrogen bonding on 1JOH. There are three different types of OH bonds present in both cages: first category are OH bonds of DDA water molecules, both involved in DDA H-bonding, with SSCC denoted as 1J DDA OH and next two are present in DAA water molecules – the OH bonds involved in H-bonding and not involved in DAA H-bonding (having dangling hydrogens), with SSCC’s de1 DAA noted as 1J DDA OH and J OH , respectively. Total values of coupling constants are for all groups within a range of 2.3 Hz only. The formation of the hydrogen bond in 512 cage influences strongest the 1J DAA OH by 8.19 Hz (from 76.04 in 1 DDA monomer to 84.23 Hz). The changes of 1J DAA are OH and of J OH slightly smaller (7.44 and 8.03 Hz, respectively). Similar changes one can find for 51262 cage (7.19, 6.40 and 6.20 Hz 1 DAA 1 DDA 6 for 1J DAA OH , J OH and J OH , respectively). In the cage 5 the trend is reversed (7.04, 7.83 and 8.50 Hz), but the changes are of the

3.4. Intramolecular spin–spin coupling constants of CO2 and H2O molecules The calculated 1JCO and 2JOO SSCCs for CO2 molecule are collected in Table 4. Empirical data and results of other authors calculations at DFT/B3LYP level are also presented together with the values of CAO bond lengths and O  O distances. Intramolecular 1 JOH, 2JHH coupling constants of water molecules are listed in Table 5. Although oxygen atoms of CO2 molecule are chemically indistinguishable, their interactions with water molecules building cages are non-symmetrical. Because of that, CAO bond lengths and the 1JCO are not equal for any of three cages. However, deviations in these values from monomeric ones are getting smaller with increasing size of the cage. All 1JCO values are positive. Upon enclathration reductions are observed in all cases except one. The changes are ca. 3, 7 and 2 Hz for the shorter CAO bond length, and ca. 5, 4, and 2 Hz for the longer CAO bond length, for 512, 51262 and 51264 cages, respectively. Variability of those changes clearly shows, that specific interactions present in different cages and depending on CO2 position and orientation are crucial for the values obtained. Therefore a more detailed study investigating wide spectrum of possible orientations of CO2 molecule within cages would be needed to quantitatively describe the changes of 1JCO upon enclathration. Based on our data, it can be stated qualitatively, that enclathration reduces 1JCO value. To the best of our knowledge, there is no published experimental SSCC for CO2 molecule enclathrated in CH. Experimental result for gas phase CO2 gives 1JCO of 16.1 Hz [58]. Inclusion of the ZPE correction (2.63 Hz [53]) gives the value 13.47 Hz, which can be compared to our static result more directly. In fact, it is twice smaller than 27.86 Hz, obtained by us. Complete basis set limit calculations at the DFT/B3LYP level and SOPPA approach (best results to date) both tend to overestimate the 1JCO for CO2 by at least 50% [53], so our result is not surprising. Moreover, CAO bond length obtained in our calculations is longer by 0.009 Å than experimental value of 1.160 Å used in reference calculations, what additionally enhances the difference in obtained values.

Table 5 Average values of intramolecular SSCCs of XJYZ type and Y–Z distances for 512, 51262 and 51264 cages (geometry optimized at B3LYP/aug-cc-pVDZ, NMR at B3LYP/HuzIII-su3). Standard deviation in % of average value in parentheses. Coupled atoms

J [Hz]

Distance [Å]

Cage

J [Hz]

Distance [Å]

H  H DAA

3.40 (38) 7.72 (18) 7.91 (13)

1.637 (1) 1.569 (1) 1.559 (0)

512 51262 51264

7.68 (16) 7.69 (11)

1.567 (1) 1.559 (1)

H  H DDA

3.81 (54) 7.95 (16) 7.66 (11)

1.625 (1) 1.563 (1) 1.560 (1)

512 51262 51264

3.81 (54) 8.32 (15) 7.38a

1.625 (1) 1.558 (1) 1.564a

OAH⁄ DAA

84.23 (3) 83.23 (4) 83.08 (3)

0.979 (1) 0.970 (0) 0.965 (0)

512 51262 51264

OAH DAA

83.48 (4) 82.64 (5) 83.87 (3)

1.015 (1) 0.996 (1) 0.990 (1)

512 51262 51264

83.35 (5) 84.28 (3)

0.985 (1) 0.981 (1)

84.07 (4) 82.24 (4) 84.54 (3)

0.995 (1) 0.982 (1) 0.977 (1)

512 51262 51264

83.88 (4) 81.29 (4) 84.58 (2)

1.002 (1) 0.988 (1) 0.982 (1)

OAH DDA

Coupled atoms H  H 5261

H  H 53

OAH 5261

OAH 53

Values for monomer (in Hz): calculated (this work): 1JOH = 76.04, 2JHH = 7.72 calculated (CCSD/pz3d2f) [59]: 1JOH = 78.85, 2JHH = 7.84 exp.: 1JOH = 80.6 [60], 2JHH = 6.89 [61] a Only one water molecule in this category.

39

P. Siuda, J. Sadlej / Chemical Physics 433 (2014) 31–41

same range as in two former cases. Thus, independent on the type of water molecule, the hydrogen bonding formation appears to decrease the absolute value of 1JOH. Dividing water molecules according to topological criteria gives further insight into influence of the environment on intramolecular SSCCs. Within our approach, all water molecules forming cage 512 are of 53 topology, while for larger 51262 and 51264 cages both water topologies are present. Again H-bonded OAH groups are used as in real crystals dangling OAH bonds are found only in defects. For 51262 cage there is a clear difference of 2 Hz between 1JOH values for 53 and 5261 (81.29 and 83.35 Hz, respectively). In the 51264 cage, both values are very close to each other (84.58 and 84.28 Hz, respectively). This result is however based on only 5 representatives of 53 topology present in 51264 cage, as compared to 37 representatives of 5261 topology. Based on data for 51262 cage it can be stated that different topologies give distinguishable values of 1JOH. Independent of the scheme of partition, all 1JOH presented in Table 5 are dominated by FC term giving at least 85% of total value (see Table 4 of Supplementary Materials). Among other three components only PSO is of importance, adding not more than 15% of the final value, while DSO and SD influence on the final value of coupling together on average is less than 1%. Values of intramolecular coupling constants 2JHH for water are presented in the same Table 5. The changes in 2JHH SSCC in the DAA molecules are 4.32, 0.00 and 0.19 Hz (with respect to 7.72 Hz obtained for a single molecule) for 512, 51262 and 51264 cages, respectively. The DDA-type water molecules show very similar trend giving 3.91, 0.23 and 0.06 Hz. Changes observed for both types of water molecules are strongly correlated with the changes of H  H distance. This intramolecular distance is much longer in 512 cage (by 0.07 Å for both type of water molecules, which is 5% of the H  H distance), than in two bigger cages. Subsequently, 2JHH is twice smaller in 512 cage than in bigger cages, where H  H distance is close to monomeric value. It can be stated that final values of 2JHH are determined by the type of water cage, as this determines the length of the H  H distance. The type of H-bonding network water molecules are involved in is of minor importance for 2JHH values. To compare 2JHH values with data for crystalline hydrates it is more natural to divide results according to topological criteria. For 512 cage only 53 water molecules are present and for bigger 51262 and 51264 cages both water topologies are found. In 51262 cage the difference between two topologies equals 0.7 Hz (8.32 and 7.68 Hz for 53 and 5261, respectively). In the case of 51264 cage there is only one representative of 53 topology, so reliable comparison with average for 5261 topology is not possible. However, it may be expected that also for 51264 cage visible difference in 2JHH could be observed, given more representatives (a larger system). H-bonds in 51264 cage are generally weaker than in 51262 cage so their influence on water geometries is smaller. On the other hand, 2JHH depend very strongly on water geometry, so even those smaller changes in geometry should be visible in 2JHH values. To conclude, it may be stated that different topology of water molecule is reflected in distinguishable values of 2JHH coupling constant. FC contributions (see Table 4 of Supplementary Materials) to 2 JHH are predominant, while other contributions are comparable for all cages. It is clear that the crucial part of the 2JHH coupling can be described by FC mechanism. 3.5. Intermolecular

2h

JOO and

Table 6 Average values of intermolecular SSCCs of 2hJOO type and rOO for 512, 51262 and 51264 cages (optimized at B3LYP/aug-cc-pVDZ, NMR at B3LYP/HuzIII-su3). Standard deviation in % of average value in parentheses. Coupled atoms

Cage

J [Hz]

rOO

DAA–DAA

512 51262 51264 512 51262 51264 512 51262 51264

4.81 4.72 4.79 2.93 3.00 2.86 4.40 3.36 4.16

2.759 2.748 2.742 2.752 2.743 2.745 2.749 2.752 2.757

DAA–DDA

DDA–DDA

(0) (1) (0) (0) (1) (0) (1) (0) (0)

sponding four components to them. Three groups of the coupling constants can be distinguished for three 512, 51262 and 51264 cages, namely DAA–DAA, DAA–DDA and DDA–DDA. All of couplings are transferred through hydrogen bond, so one of the oxygens is always donating hydrogen atom, while the other acts as an acceptor. Distances between coupled oxygens are very similar and independent both of types of atoms involved, and of the cage size. However, that is due to geometry optimization procedure in which oxygens were kept frozen to preserve the character of crystalline structure. Within three distinguished groups of constants both total values and four components are very similar and independent of cage size. Interestingly, values for DAA–DDA pairs do not lie between DAA–DAA and DDA–DDA, as could be expected. While coupling constants for DAA–DAA pairs are the greatest, constants for DAA–DDA are the smallest. This counterintuitive pattern was also observed before for methane clathrate hydrate [2]. It is due to FC component, which dominates the total value of the coupling. Values of three other components SD, PSO and DSO for DAA–DDA coupling lie between those for DAA–DAA and DDA–DDA pairs. However, their absolute values are much smaller than FC component and they differ in sign, partially canceling each other, so the sum of their values amount to less than 15% of FC. Now, we would like to discuss the other type of hydrogen bond transmitted coupling constants 1hJOH. Table 7 contains values of intermolecular O  H coupling constants transmitted through Hbond. Four distinct groups of such couplings could be distinguished for all three cages, as both O and H atoms can belong to DDA or DAA water molecules. Fig. 7 depicts different types of possible HTable 7 Average values of intermolecular SSCCs of 1hJOH type, average values of inter- and intramolecular distances OAH, O  H and O  O for 512, 51262 and 51264 cages used to describe relative strength of H-bonds formed by different pairs of water molecules (geometry optimized at B3LYP/aug-cc-pVDZ, NMR at B3LYP/HuzIII-su3). Coupled atoms (O)DAA2  (H)DAA

Cage

J [Hz]

rOO

rOH

rOH (bond)

12

5 51262 12 4 5 6

7.11 7.44 7.28 7.27

2.759 2.748 2.742 2.750

1.755 1.767 1.760 1.761

1.004 0.981 0.982 0.989

512 51262 51264

7.95 6.94 7.02 7.30

2.752 2.743 2.745 2.747

1.767 1.804 1.790 1.787

0.985 0.938 0.955 0.959

512 51262 51264

6.59 6.83 6.26 6.56

2.752 2.743 2.745 2.747

1.727 1.736 1.747 1.737

1.025 1.007 0.998 1.010

512 51262 51264

7.68 6.93 7.33 7.31

2.749 2.752 2.747 2.749

1.760 1.770 1.765 1.765

0.989 0.982 0.982 0.984

Average (O)DAA  (H)DDA

Average (O)DDA  (H)DAA

1h

JOH coupling constants Average

Let us discuss first the intermolecular H-bond transmitted O  O SSCCs denoted 2hJOO. Values of these constants between oxygen nuclei of neighboring water molecules are presented in Table 6, while Table 5 in Supplementary Materials contain values of corre-

(17) (15) (18) (40) (61) (50) (11) (39) (18)

(O)DDA  (H)DDA

Average

40

P. Siuda, J. Sadlej / Chemical Physics 433 (2014) 31–41

Fig. 7. Depiction of possible types of H-bonds observed in all three cages. Differentiation based on type of water molecules taking part in H-bond formation (DDA and DAA). O(DDA)  H(DAA) form the strongest, while O(DAA)  H(DDA) weakest H-bond.

bonds as found in cage 51262. All four combinations of coupled atoms share characteristic scheme of contributions to final coupling. The average 1hJOH are substantial, as shown in Table 7. They are within the range of 6–8 Hz. These couplings are dominated by the FC term for all intermolecular coupling types (see Table 6 of Supplementary Materials). Moreover, they all lie within 2 Hz margin. Despite of that it is possible to observe correlations of the 1hJOH value with the strength of H-bonds involving different types of O and H atoms and with changes of distances between those atoms. Detailed analysis of the strength of H-bonds formed is presented in Supplementary Materials. Based on this analysis can be stated that the average strength of hydrogen bonds formed between water molecules is decreasing from (O)DDA  (H)DAA through (O)DAA  (H)DAA and (O)DDA  (H)DDA to (O)DAA  (H)DDA. Moreover, while (O)DDA  (H)DAA pair of water molecules form the strongest H-bonds, the relative strengths of three other types of H-bonds are very similar.

4. Conclusions Calculations of NMR parameters for CO2 molecule enclathrated by 512, 51262 and 51264 cages of sI and sII hydrates have been performed at DFT/B3LYP level with HuzIIIsu3 basis set. Detailed analysis of the shielding and the indirect spin–spin intra- and intermolecular coupling constants values have been presented and discussed. An attempt was made to construct the connection between the NMR parameters and the topology of cages and types of hydrogen bonding in water molecules. The results are summed up in the following points:

1. Based on the observation of the complexation-induced shifts in 13C shielding constants of the carbon atom in CO2 it is possible to discriminate both between cages and structures sI and sII of the clathrate hydrates. Transfer from gas phase to clathrate cage may weaken or enhance shielding effect, depending on the cage type. 2. The 17O shielding constants of CO2 molecule form a decreasing trend starting from 512 encaged molecule through 51262 encaged and 51264 encaged molecule to monomeric one. The complexation-induced changes are: 18.9, 6.2 and 2.4 ppm, respectively. These changes can be explained in terms of size of the cages. 3. The interaction-induced shifts of the shielding constants of oxygen nuclei 17O of water molecules are dependent on environment and are smaller for the DDA then for DAAtype water molecules. This trend may be interpreted as a result of decreasing strength of interaction between water molecules as the cage gets larger. 4. Intramolecular indirect 1JCO spin–spin coupling constants in CO2 molecule can also serve to distinguish between gaseous and enclathrated molecules. Similarly, to the carbon nucleus shielding constant, depending on the cage, the coupling constant value may both increase or decrease with respect to isolated CO2 molecule. 5. The complexation-induced changes in the intramolecular proton–oxygen 1JOH and proton–proton coupling constants 2 JHH do not vary significantly in DAA and DDA type of water in cages. Thus generally, independent of the type of water molecule, the hydrogen bonding formation appears to increase the absolute values of 1JOH and 2JHH (with one exception, for 2JHH in cage 512).

P. Siuda, J. Sadlej / Chemical Physics 433 (2014) 31–41

6. On the other hand, the differences between intermolecular 1h JOH transmitted through hydrogen bonds H  O are more substantial for different combinations of water molecules and concerted with the strength of H-bonds formed. The increase of their values is connected with the shortening of the intramolecular, covalent OAH bond and the elongation of the intermolecular O  H distance. The analysis of data suggests: hydrogen bonds between DDA type water molecules acting as a proton acceptor from DAA type water molecules are stronger than hydrogen bonds formed in other configurations of water molecules. 7. Analysis of the shielding constants and intramolecular SSCC for water molecules forming cages divided according to topological criteria shows that water molecules of different topologies are characterized by distinct NMR parameters. Therefore, it should also be possible to observe experimentally all topologically different water molecules present in crystalline clathrate hydrates (56, 5561 and 5462 in sI structure and 56, 5561 in sII structure). Summarizing, despite the approximations mentioned in the Introduction the results presented here provide important insight into influence of structure on NMR parameters of molecules forming CHs and are close to available experimental data. Although more detailed studies, addressing dynamics of the system would be desirable, our results provide a good starting point for the interpretation of NMR properties of clathrate hydrates of carbon dioxide. Acknowledgments Project operated within the Foundation for Polish Science MPD Programme cofinanced by the EU European Regional Development Fund. Computations conducted using supercomputing facilities of University of Warsaw, Poland (ICM and Department of Chemistry facilities), and University of Tromso, Norway. Profs. Kenneth Ruud and Michał Jaszun´ski are acknowledged for helpful discussions. This work was supported also by FRSE (FSS/2011/V/D3/W/0107). Appendix A. Supplementary data Supplementary data associated with this article can be found, in the online version, at http://dx.doi.org/10.1016/j.chemphys.2014. 01.007. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16]

T. Helgaker, M. Jaszun´ski, K. Ruud, Chem. Rev. 99 (1) (1999) 293. P. Siuda, J. Sadlej, J. Phys. Chem. A 115 (5) (2011) 612. J.A. Ripmeester, Ch.I. Ratcliffe, Energy Fuels 12 (2) (1998) 197. S. Alavi, P. Dornan, T.K. Woo, ChemPhysChem: Eur. J. Chem. Phys. Phys. Chem. 9 (6) (2008) 911. V. Chihaia, S. Adams, W.F. Kuhs, Theor. Study Chem. Phys. 297 (1–3) (2004) 271. E.D. Sloan, Nature 426 (6964) (2003) 353. E.D. Sloan, C.A. Koh, Clathrate Hydrates of Natural Gases, third ed., CRC Press, New York, 2008. J. Ripmeester, J. Tse, Ch. Ratcliffe, B. Powell, Nature 325 (6100) (1987) 135. L. Rufine, J.P.M. Trusler, J. Chem. Thermodyn. 42 (5) (2010) 605. Y.-T. Seo, H. Lee, J. Phys. Chem. B 108 (2) (2004) 530. R.K. McMullan, G.A. Jeffrey, J. Chem. Phys. 42 (8) (1965) 2725. O.S. Subbotin, T.P. Adamova, R.V. Belosludov, H. Mizuseki, Y. Kawazoe, J.-I. Kudoh, P.M. Rodger, V.R. Belosludov, J. Chem. Phys. 131 (11) (2009) 114507. Y.A. Dyadin, E.G. Larionov, A.Y. Manakov, F.V. Zhurko, E.Y. Aladko, T.V. Mikina, V.Y. Komarov, Mendeleev Commun. 5 (1999) 209. E. Dartois, B. Schmitt, Astron. Astrophys. 504 (3) (2009) 869. H. Mohammadi-Manesh, S. Alavi, T.K. Woo, M. Ashraffzaadeh, B. Najafi, Phys. Chem. Chem. Phys. 11 (39) (2009) 8821. E.K. Lee, J.D. Lee, H.J. Lee, B.R. Lee, Y.S. Lee, S.M. Kim, H.O. Park, Y.S. Kim, Y.-D. Park, Y.D. Kim, Environ. Sci. Technol. 43 (20) (2009) 7723.

41

[17] R. Kumar, P. Linga, I. Moudrakovski, J. Ripmeester, P. Englezos, AIChE J. 54 (8) (2008) 2132. [18] M. Kida, A. Hori, H. Sakagami, S. Takeya, Y. Kamata, N. Takahashi, T. Ebinuma, H. Narita, J. Phys. Chem. A 115 (5) (2011) 643. [19] V. Buch, J.P. Devlin, I.A. Monreal, B. Jagoda-C´wiklik, N. Uras-Aytemiz, L. C´wiklik, Phys. Chem. Chem. Phys. 11 (44) (2009) 10245. [20] J.A. Ripmeester, L. Ding, D.D. Klug, J. Phys. Chem. 100 (32) (1996) 13330. [21] I.A. Monreal, L. C´wiklik, B. Jagoda-C´wiklik, J.P. Devlin, J. Phys. Chem. Lett. 1 (1) (2010) 290. [22] T. Nakayama, S. Tomura, M. Ozaki, R. Ohmura, Y.H. Mori, Energy Fuels 24 (4) (2010) 2576. [23] C.A. Koh, A.K. Sum, E.D. Sloan, J. Appl. Phys. 106 (6) (2009) 061101. [24] R. Boswell, Science 325 (2009) 957. [25] S. Alavi, T.K. Woo, J. Chem. Phys. 126 (4) (2007) 044703. [26] P. Dornan, S. Alavi, T.K. Woo, J. Chem. Phys. 127 (12) (2007) 124510. [27] G. Tegze, T. Pusztai, G. Tóth, L. Gránásy, A. Svandal, T. Buanes, T. Kuznetsova, B. Kvamme, J. Chem. Phys. 124 (23) (2006) 234710. [28] P.S.R. Prasad, K.S. Prasad, N.K. Thakur, Curr. Sci. 90 (11) (2006) 1544. [29] S. Alavi, R. Susilo, J.A. Ripmeester, J. Chem. Phys. 130 (2009) 174501. [30] S. Alavi, K. Udachin, J.A. Ripmeester, Chem. – Eur. J. 16 (3) (2010) 1017. [31] J.A. Ripmeester, Ch.I. Ratcliffe, J. Phys. Chem. 90 (1) (1986) 1259. [32] U. Buck, I. Ettischer, M. Melzer, V. Buch, J. Sadlej, Phys. Rev. Lett. 80 (12) (1998) 2578. [33] P. Terleczky, L. Nyulaszi, Chem. Phys. Lett. 488 (4–6) (2010) 168. [34] S. Takeya, K.A. Udachin, I.L. Moudrakovski, R. Susilo, J.A. Ripmeester, J. Am. Chem. Soc. 132 (2) (2010) 524. [35] J.D. Bernal, R.H. Fowler, J. Chem. Phys. 1 (8) (1933) 515. [36] A.D. Becke, J. Chem. Phys. 98 (7) (1993) 5648. [37] C. Lee, W. Yang, R.G. Parr, Phys. Rev. B 37 (2) (1988) 785. [38] T.H. Dunning, J. Chem. Phys. 90 (2) (1989) 1007. [39] M.J. Frisch, G.W. Trucks, H.B. Schlegel, G.E. Scuseria, M.A. Robb, J.R. Cheeseman, J.A. Montgomery, Jr., T. Vreven, K.N. Kudin, J.C. Burant, J.M. Millam, S.S. Iyengar, J. Tomasi, V. Barone, B. Mennucci, M. Cossi, G. Scalmani, N. Rega, G.A.P. Sson, H. Nakatsuji, M. Hada, M. Ehara, K. Toyota, R. Fukuda, J. Hasegawa, M. Ishida, T. Nakajima, Y. Honda, O. Kitao, H. Nakai, M. Klene, X. Li, J.E. Knox, H.P. Hratchian, J.B. Cross, V. Bakken, C. Adamo, J. Jaramillo, R. Gomperts, R.E. Stratmann, O. Yazyev, A.J. Austin, R. Cammi, C. Pomelli, J.W. Ochterski, P.Y. Ayala, K. Morokuma, G.A. Voth, P. Salvador, J.J. Dannenberg, V.G. Zakrzewski, S. Dapprich, A.D. Daniels, M.C. Strain, O. Farkas, D.K. Malick, A.D. Rabuck, K. Raghavachari, J.B. Foresman, J.V. Ortiz, Q. Cui, A.G. Baboul, S. Clifford, J. Cioslowski, B.B. Stefanov, G. Liu, A. Liashenko, P. Piskorz, I. Komaromi, R.L. Martin, D.J. Fox, T. Keith, M.A. Al-Laham, C.Y. Peng, A. Nanayakkara, M. Challacombe, P.M.W. Gill, B. Johnson, W. Chen, M.W. Wong, C. Gonzalez, J.A. Pople. Gaussian 03, Revision C.02. Gaussian Inc, Wallingford, CT, 2004. [40] H. Cybulski, J. Sadlej, Chem. Phys. 323 (2–3) (2006) 218. [41] S. Huzinaga, J. Chem. Phys. 42 (4) (1965) 1293. [42] T. Helgaker, P. Wilson, R. Amos, N. Handy, J. Chem. Phys. 113 (8) (2000) 2983. [43] C. Angeli, K.L. Bak, V. Bakken, O. Christiansen, R. Cimiraglia, S. Coriani, P. Dahle, E.K. Dalskov, T. Enevoldsen, B. Fernandez, C. Hättig, K. Hald, A. Halkier, H. Heiberg, T. Helgaker, H. Hettema, H.J.Aa. Jensen, D. Jonsson, P. Jørgensen, S. Kirpekar, W. Klopper, R. Kobayashi, H. Koch, A. Ligabue, O.B. Lutnæs, K.V. Mikkelsen, P. Norman, J. Olsen, M. J. Packer, T.B. Pedersen, Z. Rinkevicius, E. Rudberg, T.A. Ruden, K. Ruud, P. Salek, A. Sanchez de Meras, T. Saue, S.P.A. Sauer, B. Schimmelpfennig, K.O. Sylvester-Hvid, P.R. Taylor, O. Vahtras, D.J. Wilson, H. Ågren, A. Dalton, A molecular electronic structure program, Release 2.0 (2005), see http://www.kjemi.uio.no/software/dalton/dalton.html. [44] D.J. Anick, J. Mol. Struct. (Thoechem) 587 (1–3) (2002) 87. [45] N. Bjerrum, Matematisk-Fysiske Meddelelser kongelige danske videnskabernes selskab 27 (1) (1951) 1. [46] A. Stirling, J. Phys. Chem. B 115 (49) (2011) 14683. [47] H. Omi, T. Ueda, K. Miyakubo, T. Eguchi, Appl. Surf. Sci. 252 (3) (2005) 660. [48] K. Jackowski, M. Jaszun´ski, W. Makulski, J. Vaara, J. Magn. Reson. 135 (2) (1998) 444. [49] F.G. Morin, M.S. Solum, J.D. Withers, D.M. Grant, D.K. Dalling, J. Magn. Reson. 48 (1982) 138. [50] K. Jackowski, W. Makulski, Chem. Phys. Lett. 341 (2001) 369. [51] A.A. Auer, J. Gauss, J.F. Stanton, J. Chem. Phys. 118 (23) (2003) 10407. [52] A. Antusek, K. Jackowski, M. Jaszun´ski, W. Makulski, M. Wilczek, Chem. Phys. Lett. 411 (1–3) (2005) 111. [53] T. Kupka, M. Nieradka, M. Stachów, T. Pluta, P. Nowak, H. Kjær, J. Kongsted, J. Kaminsky, J. Phys. Chem. A 116 (14) (2012) 3728. [54] G. Graner, C. Rossetti, D. Bailly, Mol. Phys. 58 (3) (1986) 627. [55] J. San Fabián, E. Díez, J.M. García de la Vega, R. Suardíaz, J. Chem. Phys. 128 (8) (2008) 084108. [56] R.K. Harris, E.D. Becker, S.M. Cabral de Menezes, P. Granger, R.E. Hoffman, K.W. Zilm, Pure Appl. Chem. 80 (1) (2008) 59. [57] A.A. Auer, J. Chem. Phys. 131 (2) (July 2009) 024116. [58] R.E. Wasylishen, J.O. Friedrich, S. Mooibroek, J.B. Macdonald, J. Chem. Phys. 83 (2) (1985) 548. [59] A.A. Auer, J. Gauss, J. Chem. Phys. 115 (4) (2001) 1619. [60] N.M. Sergeyev, N.D. Sergeyeva, Yu.A. Strelenko, W.T. Raynes, Chem. Phys. Lett. 277 (1–3) (1997) 142. [61] W. Makulski, J. Mol. Struct. 839 (1–3) (2007) 90.