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)