2. METHODOLGY
2.1. Structure generation
The NO, N2O and NH3 are sufficiently small molecules to form hydrates in the sI structure under certain pressure–temperature conditions.5 The primitive lattice unit consists of two 512 cages with 12 pentagonal faces, and six 51262cages with 12 pentagonal and 2 hexagonal faces, totally 46 water molecules.38 In the initial guesses, the structures for these three hydrates were based on high-resolution neutron diffraction data for deuterated methane hydrates, the hydrogen atoms are disordered just like those in the ice Ih structure and their positions were assigned randomly but in accordance with the Bernal–Fowler ice rules.5 The cage occupancy was assumed to be 100% with a guest gas molecule at the center of every cage.
2.2. Full elasticity tensor determination
In this work, we first determined the equilibrium cell parameters for the stress-free states. Then full sets of second order elastic constants (SOECs) were calculated using the homogeneous finite strain method incorporated with the first-principle total energy method. The detailed description of this method can be found in the early publication.39 In short, elastic constants are defined by expanding the Gibbs free energy as a Taylor series versus Lagrangian strains. The second-order Taylor expansion coefficients as combinations of SOECs of the crystal are obtained from a polynomial fit to the calculated energy-strain relation. Using SOECs, the bulk modulus (B ) and shear modulus (G ) of the clathrate hydrate crystals are calculated via the Voigt-Reuss-Hill averaging scheme.40 Afterwards, (E ) and Poisson’s ratio (v ) can be calculated from B and G .
The sI hydrate structure is cubic, thus, three independent elastic constants c 11, c 12, andc 44 are required which can be determined from three Lagrangian strains e 1,e 2, and e 3. Though a stress strain analysis could have been applied, an energy ­strain analysis is believed to be preferable.41 The total energy variation of the system can be expanded as a Taylor series of the elastic strain.42, 43 For a system at zero pressure, one has
\(\Delta{E}=\frac{v}{2}\sum_{i=1}^6\sum_{j=1}^6c_{ij}e_ie_j+O(e_i^3)\) (1)
Here V is the volume of the undistorted lattice cell at which the ionic positions were fully optimized, ΔE is the energy change resulting from the strain vector e = (e 1,e 2, e 3,e 4, e 5,e 6) in Voigt notation, and cij is the matrix elements of elastic constants. The primitive lattice vectors ai (i = 1…3) of a crystal are transformed to the new vectorsai ’ under the strain by
\(\par \begin{pmatrix}a_{1}^{{}^{\prime}}\\ a_{2}^{{}^{\prime}}\\ a_{3}^{{}^{\prime}}\\ \end{pmatrix}=\par \begin{pmatrix}a_{1}\\ a_{2}\\ a_{3}\\ \end{pmatrix}\bullet(+\varepsilon)\) (2)
Where\(\text{\ ε}\ \)is the strain tensor. It relates to the strain vector e by
\(\varepsilon=\par \begin{pmatrix}e_{1}&\frac{e_{6}}{2}&\frac{e_{5}}{2}\\ \frac{e_{6}}{2}&e_{2}&\frac{e_{4}}{2}\\ \frac{e_{5}}{2}&\frac{e_{4}}{2}&e_{3}\\ \end{pmatrix}\) (3)
Three modes corresponding to three different sets of small strains were used. A volume-conserving tetragonal strain e = (δ ,-δ ,δ 2/(1-δ 2),0,0,0) was applied according to Equation (3) yields:
\(\Delta \text{E} =V(c_{11}-c_{12})\delta^{2}+O(\delta^{4})\)  (4)
Then the [100] and [010] strain e = (δ,δ,0,0,0,0) was applied which leads to:
\({\Delta}\text{E}=V(c_{11}+c_{12})\delta^2+O(\delta^2)\) (5)
Last, c 44 was calculated from the [111] shear strain e = (0,0,0,δ,δ,δ)
\(\Delta \text{E} = \frac{3V}{2}c_{44}+O(\delta^{2})\) (6)
Combining (4), (5), and (6) the relevant elastic constantsc 11, c 12, andc 44 are obtained by fitting to a polynomial ofδ varying from -0.03 to 0.03. Equations of elastic properties are listed in the following.44
The bulk modulus (B ) can be calculated from:
\(B=\frac{(c_{11}+c_{12})}{3}\) (7)
Shear modulus:
\(G_{R\text{euss}}=\frac{5(c_{11}-c_{12})c_{44}}{4c_{44}+3(c_{11}-c_{12})}\text{\ \ \ \ \ }\)(8)
\(G_{\text{Voigt}}=\frac{(c_{11}-c_{12}+3c_{44})}{5}\text{\ \ \ \ \ \ }\)(9)
\(G=\frac{G_{\text{Voigt}}-G_{\text{Reuss}}}{2}\text{\ \ \ \ \ \ \ \ \ \ \ }\)(10)
Poisson’s ratio:
\(\upsilon=\frac{(\frac{3}{2})B-G}{G+3B}\text{\ \ \ \ \ \ \ \ \ \ \ }\text{\ \ \ }\)(11)
Young’s modulus:
\(=2G(1+\upsilon)\) (12)
Longitudinal wave speed:
\(V_{p}=\left(\frac{B+(\frac{4}{3})G}{\rho}\right)^{\frac{1}{2}}\text{\ \ \ \ \ \ \ \ \ }\)(13)
Transverse wave speed:
\(V_{s}=\left(\frac{G}{\rho}\right)^{\frac{1}{2}}\text{\ \ \ \ \ \ \ \ \ \ \ \ \ \ }\ \)(14)
Zener anisotropy ratio:
\(A_{Z}=\frac{2c_{44}}{c_{11}-c_{12}}\text{\ \ \ \ \ \ \ \ }\text{\ \ \ \ \ \ \ }\)(15)