Theoretical calculations of neutron scattering cross sections for tetrahydrofuran-containing clathrate hydrates at low temperature

. Tetrahydrofuran-containing clathrate hydrate is considered as a potential cold and very cold neutron moderator material. The fully deuterated form is more promising because of the lower neutron absorption cross section of deuterium. In this work we present theoretical calculations of neutron scattering cross sections for hydrogenated and deuterated clathrates at low temperature. The scattering cross sections are generated by using the crystalline structure and phonon density of states obtained from density functional theory calculations. The theoretical neutron scattering cross sections serve to compare against existing or future experimental data.


Introduction
The European Spallation Source (ESS), currently under construction in Lund, Sweden, will be the most powerful pulsed neutron source in the world [1].The current source is located above the spallation target and designed for high cold and thermal brightness [2].A design study project named HighNESS is dedicated to develop a second source below the target which will provide cold, very cold and ultra cold neutrons with high intensity [3,4].This study includes exploration of promising novel reflector and moderator materials.For the purpose of cold neutron reflectors, nano-diaomond particles and magnesium hydride have been investigated [5][6][7].A promising cold neutron moderator class of materials also included in this study are clathrate hydrates, which are the focus of this work.
Clathrate hydrates are ice-like compounds having a cage structure.Small molecules, such as methane, can be enclathrated in the cage, leading to stabilisation of the structure.Clathrate hydrate hosting oxygen is considered as a promising candidate due to neutron inelastic magnetic scattering with oxygen [8].Tetrahydrofuran (THF)containing clathrate hydrate has also been considered as an advanced cold neutron moderator [9].The material can be produced directly from a stoichiometric THF/water mixture at normal pressure which makes it suitable for the large volume production needed for neutron moderators/reflectors.Rattling modes of THF molecules mostly overlap with the low energy translational modes of the cage molecules at about 7 and 10.5 meV, which are common to all type II clathrates [10].However, THF also gives rise with a lower energy peak at about 5 meV in * e-mail: shuqi.xu@ess.eu the phonon density of states [10] that might also contribute to the inelastic neutron cross section.Both THF and fully deuterated tetrahydrofuran (TDF)-clathrates [11] have been investigated by inelastic neutron scattering measurements [9,11].The lower neutron absorption cross section of deuterium makes the TDF-clathrate a more promising cold moderator candidate.
The objective of this work is to provide theoretical neutron scattering cross sections and scattering kernels for THF-and TDF-clathrates to be used for simulating measurements and design studies.Therefore, we performed density functional theory (DFT) calculations using CP2k [12] and phonopy [13] in order to obtain the crystalline structure and the phonon density of states (PDOS), which captures the dynamics of the atoms around their equilibrium positions.The DFT results were used to generate scattering kernels and cross-sections using the open source software package NCrystal [14][15][16].The calculated scattering kernels and cross sections can be used to compare with the neutron transmission, diffraction and inelastic neutron scattering measurements performed in another work package of this project [17].

Density functional theory calculations
The initial configuration of the THF-containing clathrate hydrate in structure II was taken from Ref. [18].The unit cell of the face centered cubic lattice contains 34 H 2 O molecules and two THF molecules inserted in the large cages.The unit cell was replicated to build the conventional cubic cell with 136 H 2 O and 8 THF molecules (see Fig. 1).The initial position of the THF molecules was chosen by performing classical molecular dynamics (MD) simulations with GROMACS [19], the SPC/E potential for  water [20], and the general Amber force field for THF and THF-water interactions [21].The initial structure was optimized by means of DFT calculations with CP2k [12].We employed the Perdew-Burke-Ernzerhof approximation for the exchange and correlation functional [22] and norm conserving pseudopotentials [23].Kohn-Sham orbitals were expanded in a triple-zeta-valence plus polarization basis set while the electronic density was expanded in plane waves up to a kinetic cutoff of 1000 Ry as implemented in CP2k.Brillouin Zone (BZ) integration for electrons was restricted to the Γ-point.Semiempirical van der Waals interactions were included according to Grimme (D3) [24].
The atomic positions were initially optimized at the experimental lattice parameter a = 17.13Å [25] and then at different volumes.The resulting energy versus volume points were fitted by a Birch-Murnaghan equation of state which yielded a theoretical lattice parameter of 16.7 Å.As the contraction of the lattice parameter with respect to ex- periments leads to a stiffening of phonon frequencies, we preferred to compute the PDOS at the experimental lattice parameter.
The phonons were computed using the CP2k code combined with phonopy [13]; the force constant matrix was computed from forces at finite atomic displacements in the conventional cubic cell (136 water molecules).The dynamical matrix, on a fine q-mesh in the BZ, was then built by phonopy.The coulombic long range contribution was computed analytically from effective charges and the dielectric constant obtained in turn from ionic forces and polarization at finite electric field within the Berry phase approach [26] implemented in the CP2k code.
The PDOS projected on different atomic species for the THF-containing clathrate hydrate are shown in Fig. 2. The upper panel of Fig. 2 shows the PDOS of hydrogen and oxygen bound in H 2 O.The dynamics of the water molecules inherit the behaviours observed in ice [27].At low energy, the translational modes of the H 2 O molecules dominate.In this spectral region, the experimental neutron weighted PDOS of type II clathrate hydrates has two characteristic peaks at about 7 and 10.5 meV [11,28].This double peak feature is less evident in the theoretical PDOS projected on the water molecules.The position of the two peaks is also blue-shifted to about 10 and 14.5 meV.This double peak feature is instead well reproduced by classical simulations with rigid water molecules as reported for instance in Refs.[11,29].A possible improvement of the theoretical PDOS would then result from the use of the DFT-PDOS for the intramolecular and librational modes and the classical result for the PDOS in the translational low frequency region.The band from 70 meV to 150 meV is due to molecular librations.At higher energy, intramolecular modes are dominant with O-H bending and O-H stretching at around 200 meV and 400 meV, respectively.
The PDOS of atoms bound in the THF molecules are  illustrated in the lower panel of Fig. 2. Several sharp peaks are observed in the energy range from 70 meV to 200 meV due to a weak coupling with the water molecules of the cage.Hydrogen bound in THF has excitations in a broad band from about 4 to 10 meV with contributions below the peak of the PDOS of hydrogen bound in H 2 O, also considering the blueshift of the theoretical peaks of water compared to experiments that we mentioned above.This property makes THF-clathrate a better cold neutron moderator compared to ice.The PDOS of fully deuterated TDF-clathrate projected on different atomic species are shown in Fig. 3.The PDOS are shifted to lower energies when substituting hydrogen by deuterium due to the increase of mass and moment of inertia, as in the case of ice [30].The PDOS of the cage can be transferred to other clathrates, e.g., oxygenclathrate, based on weak guest-host coupling assumption as discussed in Ref. [28].

Neutron scattering cross sections
The lattice parameters, atomic positions of the crystalline structure, and the PDOS were converted into input files for the NCrystal toolkit.Based on the Gaussian and incoherent approximations [31,32], NCrystal is able to calculate the neutron scattering cross sections and scattering kernels including the physics of both coherent and incoherent elastic, and inelastic scattering in a wide range of materials.Fig. 4 shows the total scattering cross section of the THF-clathrate, which is the sum of coherent elastic, incoherent elastic and inelastic cross sections.Incoherent elastic scattering is dominant for neutron energy below a few meV because of the large incoherent cross section of hydrogen [33].The coherent cross section of deuterium is more important than its incoherent part.Therefore, the coherent elastic scattering dominates in the case of TDFclathrate as illustrated in Fig. 5.The large Bragg cutoff (around 0.2 meV or 20 Å) is an advantage for cold neutron moderation because neutrons can be reflected within the cage structure thus increasing the interactions with the guest molecules.

Conclusions
In this work we present current advance of the calculation of the neutron scattering cross sections of THF-containing clathrate hydrate and its fully deuterated form.The crystalline structure and phonon density of states obtained from DFT calculations have been used to compute the neutron scattering cross sections by means of the open source software package NCrystal.The resulting cross sections at different wavelength can be used to compare against existing or future measurements on THF-and TDF-clathrates.Clathrate hydrates can also be considered as potential reflector materials because of the large Bragg cutoff.

Figure 1 .
Figure 1.Crystalline structure of THF-containing clathrate hydrate (structure II) obtained from DFT calculations.THF molecules are inserted in eight large cages.TDF-clathrate shares the same crystalline structure.

Figure 2 .
Figure 2. Normalized PDOS of THF-containing clathrate hydrate projected on different atomic species.The global PDOS can be obtained by summing the projected PDOS multiplied by the number of corresponding atoms in the unit cell.

Figure 3 .
Figure 3. Normalized PDOS of fully deuterated TDF-containing clathrate hydrate projected on different atomic species.The global PDOS can be obtained by summing the projected PDOS multiplied by the number of corresponding atoms in the unit cell.The shift of the excitations to low energies is due to the increase of mass and moment of inertia of deuterium.

Figure 5 .
Figure 5. Neutron scattering cross sections of type II fully deuterated TDF-containing clathrate hydrate.