Issue 
EPJ Nuclear Sci. Technol.
Volume 8, 2022



Article Number  8  
Number of page(s)  21  
DOI  https://doi.org/10.1051/epjn/2022004  
Published online  15 June 2022 
https://doi.org/10.1051/epjn/2022004
Regular Article
Generation of thermal scattering files with the CINEL code
CEA, DES, IRESNE, DER, Cadarache, 13108 Saint Paul Les Durance, France
^{*} email: gilles.noguere@cea.fr
Received:
8
February
2022
Received in final form:
30
March
2022
Accepted:
14
April
2022
Published online: 15 June 2022
The CINEL code dedicated to generate the thermal neutron scattering files in ENDF6 format for solid crystalline, free gas materials and liquid water is presented. Compared to the LEAPR module of the NJOY code, CINEL is able to calculate the coherent and incoherent elastic scattering cross sections for any solid crystalline materials. Specific material properties such as anharmonicity and texture can be taken into account in CINEL. The calculation of the thermal scattering laws can be accelerated by using graphics processing unit (GPU), which enables to remove the short collision time approximation for large values of momentum transfer. CINEL is able to generate automatically the grids of dimensionless momentum and energy transfers. The Sampling the Velocity of the Target nucleus (SVT) algorithm capable of determining the scattered neutron distributions is implemented in CINEL. The obtained distributions for free target nuclei such as hydrogen and oxygen are in good agreement with analytical results and MonteCarlo simulations when incident neutron energies are above a few eV. The introduction of the effective temperature and the rejection step to the SVT algorithm shows improvements to the neutron upscattering treatment of hydrogen bound in liquid water.
© S. Xu et al., Published by EDP Sciences, 2022
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
1 Introduction
Neutron transport through materials is crucial for nuclear reactor core design, neutron detection and shielding. If the incident neutron energy is much larger than the atomic binding energy, the target materials can be approximated by a collection of unbounded atoms. In this case, the treatment of the neutron scattering with materials can be well approximated by the twobody collision kinetics. However, when neutron energies lie below a few eV, the atomic vibration behaviors and/or crystalline structure play a predominant role in the neutron scattering reactions. Lowenergy neutrons have been used to probe the crystalline structure of oxides of uranium [1–4] and other materials [5, 6] since the sixties.
The low energy (below a few eV) neutron scattering theory is well detailed in the literature [7–9]. The key quantity is the double differential neutron scattering cross section (DDXS) which represents the probability of finding a scattered neutron with a given energy at a specific solid angle. The DDXS has been expressed as a function of the dynamic structure factor $S\overrightarrow{Q},\omega $ by Van Hove [10] in 1950s, where $\hslash \overrightarrow{Q}$ and ħω are respectively the momentum transfer and energy transfer of neutrons. Under the incoherent and Gaussian approximations, the dynamic structure factor $S\overrightarrow{Q},\omega $ can be calculated from the phonon density of states (PDOS) of materials. The PDOS representing the structure dynamics of materials can be obtained from ab initio density functional theory (DFT) or molecular dynamics (MD) calculations. $S\overrightarrow{Q},\omega $ are later converted to dimensionless thermal scattering laws (TSL) or S(α, β). S(α, β) are stored in ASCII files by following the ENDF6 format requirements [11], to prepare the neutron scattering data libraries.
The TSL in ENDF format have been generated by using the LEAPR module [12] of the nuclear data processing code NJOY [13] which is an industry standard of the same class codes. Stemming from the incoherent and Gaussian approximations, LEAPR is able to calculate the TSL for various kinds of materials. For solid crystalline materials, S(α, β) are further divided in an elastic and inelastic parts, with coherent and incoherent terms. However, the calculations of the coherent elastic scattering cross sections in LEAPR are limited to crystalline materials with specific symmetries. In addition, to reduce the computational time and/or for coding convenience, LEAPR adopts the short collision time (SCT) approximation and the atom site approximation in which the meansquared displacement (MSD) are supposed the same for different types of atoms in the same material.
To overcome the limitations of the LEAPR module, new codes have been developed in the scientific community: NJOY+NCrystal [14], FLASSH [15] and OCLIMAX [16, 17]. They are able to calculate the coherent elastic scattering cross sections of any kind of solid crystalline materials. NJOY+NCrystal allows to generate the coherent and incoherent elastic cross sections in the same ENDF file (mixed elastic format [18]). FLASSH and OCLIMAX improve the calculations of TSL by introducing the coherent onephonon correction [8] to the incoherent approximation and eliminating the cubic approximation.
In this work we present the CINEL code which is dedicated to generate the thermal neutron scattering files in ENDF6 format for solid crystalline, free gas materials and liquid water. Compared to the LEAPR module of the NJOY code, CINEL is able to compute the coherent elastic scattering cross sections of any crystals. The generated ENDF files support the mixed elastic format. CINEL is able to generate automatically the α and β grids used for the calculations of S(α, β). Specific properties such as anharmonicity and texture can be taken into account in CINEL. The SCT approximation is removed thanks to the use of graphics processing unit (GPU) acceleration. CINEL is able to use the PDOS obtained by MD or DFT calculations.
The main expressions involved in the lowenergy neutron scattering formalism are reported in Section 2. The implementation of CINEL is briefly presented in Section 3. CINEL is composed of three modules named Cubic, INELastic and SVT, respectively. Validations of Cubic and INELastic by comparing with other codes and databases are presented in Section 4. The automatic generation of α and β grids, mixed elastic cross section, taking into account the anharmonicity and texture are illustrated with examples. Section 5 presents the validations of the SVT module by comparing with MonteCarlo calculations and analytical results. Conclusions and perspectives are given in the last section.
2 Theory
2.1 General formalisms
As mentioned above, the key quantity in lowenergy neutron scattering is the double differential cross section (DDXS). For unpolarised neutrons and samples, based on Born approximation or Fermi^{,}s Golden rule [7], the DDXS can be obtained by: $$\frac{{\text{d}}^{2}{\sigma}_{{\overrightarrow{k}}_{i}\Rightarrow {\overrightarrow{k}}_{f}}}{\text{d}\Omega \text{d}{E}_{f}}=\frac{{k}_{f}}{{k}_{i}}S\overrightarrow{Q},\omega ,$$(1)
where ${\overrightarrow{k}}_{i}$ and ${\overrightarrow{k}}_{f}$ are respectively the wave vectors of incident and scattered neutrons, k_{i} and k_{f} represent the wavenumbers with $k=\left\overrightarrow{k}\right,\hslash \overrightarrow{Q}$ is the momentum transfer of neutron defined by $\hslash \overrightarrow{Q}\equiv \hslash {\overrightarrow{k}}_{i}{\overrightarrow{k}}_{f}$ in which Η is the Planck constant divided by 2π, ħω is the energy transfer with ħω ≡ E_{i} − E_{f}^{1}, and $S\overrightarrow{Q},\omega $ is the dynamic structure factor or scattering function defined by: $$S\overrightarrow{Q},\omega \equiv \frac{1}{2\pi \hslash}{\displaystyle \sum _{j,j\text{'}=1}^{N}\overline{{b}_{j\text{'}}{b}_{j}}{\displaystyle {\int}_{\infty}^{\infty}\langle j\text{'},j\rangle \mathrm{exp}}}j\omega t\text{dt.}$$(2)
In equation (2), N is the number of particles in the scattering system under consideration, b_{j} represents the scattering length of the particle $j,\text{\hspace{0.17em}}\overline{{b}_{j}}$ is the average value, 〈A〉 represents the operator expectation value in the scattering system, and 〈j′, j〉 (notation taken from Ref. [9], which is referred as intermediate function in Ref. [8]), represents the correlation between the position of the particle j at time t and the position of the particle j′ at time 0: $$\langle j\text{'},j\rangle \equiv \langle \mathrm{exp}i\overrightarrow{Q}\cdot {\overrightarrow{R}}_{j\text{'}}0\mathrm{exp}i\overrightarrow{Q}\cdot {\overrightarrow{R}}_{j}t\rangle ,$$(3)
where $\overrightarrow{R}t$ is the timedependent Heisenberg operator.
It is assumed that there is no correlation between the scattering lengths and the positions of the particles in the scattering system [7], then $$\overline{{b}_{j\text{'}}{b}_{j}}=\{\begin{array}{cc}\overline{{b}_{j\text{'}}}\cdot \overline{{b}_{j}}& j\ne j\text{'}\\ {b}_{j}^{2}& j=j\text{'}\end{array}.$$(4)
Equation (4) enables to decompose the dynamic structure factor in equation (2) into two distinct parts: $$S\overrightarrow{Q},\omega ={S}_{\text{coh}}\overrightarrow{Q},\omega +{\mathrm{S}}_{\text{inc}}\overrightarrow{Q},\omega ,$$(5)
where $${S}_{\text{coh}}\overrightarrow{Q},\omega \equiv \frac{1}{2\pi \hslash}{\displaystyle \sum _{j,j\text{'}=1}^{N}\overline{{b}_{j\text{'}}}}\cdot \overline{{b}_{j}}{\displaystyle {\int}_{\infty}^{\infty}\langle j\text{'},j\rangle \mathrm{exp}i\omega t\text{d}t},$$(6)
and $${S}_{\text{inc}}\overrightarrow{Q},\omega \equiv \frac{1}{2\pi \hslash}{\displaystyle \sum _{j=1}^{N}\overline{{b}_{j}^{2}}{\overline{{b}_{j}}}^{2}}{\displaystyle {\int}_{\infty}^{\infty}\langle j,j\rangle \mathrm{exp}i\omega t\text{d}t},$$(7)
${S}_{\text{coh}}\overrightarrow{Q},\omega $ is the coherent scattering function, representing the correlation between the positions of the same particle at different time (when j = j′) and the correlation between the positions of different particles at different time (when j ≠ j′). Therefore, ${S}_{\text{coh}}\overrightarrow{Q},\omega $ gives interference effects [8]. ${S}_{\text{inc}}\overrightarrow{Q},\omega $ is the incoherent scattering function, which represents the correlation between the positions of the same particle at different time. The incoherent scattering does not give interference effects (there is no terms 〈j′, j〉 for j ≠ j′ in Eq. (7)).
Based on equation (2), the value of the scattering function $S\overrightarrow{Q},\omega $ increases with the number of particles in the scattering system. In practice, the normalization is performed with respect to the number of particles given by the chemical formula or the unit cell [7]. In addition, the principle of detailed balance must be fulfilled for the dynamic structure factor $S\overrightarrow{Q},\omega $ [7]: $$S\overrightarrow{Q},\omega =\mathrm{exp}\frac{\hslash \omega}{{k}_{B}T}S\overrightarrow{Q},\omega ,$$(8)
where k_{B} is the Boltzmann constant, T is the temperature of the material under consideration.
2.2 Scattering models for solid crystalline materials
2.2.1 Harmonic, incoherent and cubic approximations
For solid crystalline materials, the operator representing the position of atom $\overrightarrow{R}t$ in the intermediate function (Eq. (3)) can be represented as the displacement of the atom $\overrightarrow{u}t$ from its equilibrium position $\overrightarrow{d}$ at time t: $$\overrightarrow{R}(t)=\overrightarrow{d}+\overrightarrow{u}(t).$$(9)
In addition, the equilibrium position $\overrightarrow{d}$ is a simple vector, which commutes with all the operators. Hence, the intermediate function can be obtained by: $$\begin{array}{c}\langle j\text{'},j\rangle =\mathrm{exp}i\overrightarrow{Q}\cdot {\overrightarrow{d}}_{j\text{'}}{\overrightarrow{d}}_{j}\\ \times \langle \mathrm{exp}i\overrightarrow{Q}\cdot {\overrightarrow{u}}_{j\text{'}}0\mathrm{exp}i\overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}t\rangle .\end{array}$$(10)
For practical purposes, it is desired to permute the exponential and the operator expectation in $\langle \mathrm{exp}i\overrightarrow{Q}\cdot {\overrightarrow{u}}_{j\prime}0\mathrm{exp}i\overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}t\rangle $ in equation (10). To this end, the harmonic approximation is supposed, in which atoms are assumed to have harmonic vibrations, i.e., the interatomic forces are linear with respect to the displacements of atoms and all the higher terms related to the anharmonic vibrations are neglected. The intermediate function 〈j′, j〉 in equation (10) becomes: $$\begin{array}{c}\langle j\text{'},j\rangle =\mathrm{exp}i\overrightarrow{Q}\cdot {\overrightarrow{d}}_{j\text{'}}{\overrightarrow{d}}_{j}\mathrm{exp}{W}_{j\text{'}}\overrightarrow{Q}\mathrm{exp}{W}_{j}\overrightarrow{Q}\\ \times \mathrm{exp}\langle i\overrightarrow{Q}\cdot {\overrightarrow{u}}_{j\text{'}}0\mathrm{exp}i\overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}t\rangle ,\end{array}$$(11)
where ${W}_{j}\overrightarrow{Q}$ is the DebyeWaller function, which measures the meansquared displacement (MSD) of the atom j along the direction $\overrightarrow{Q}$: $${W}_{j}\overrightarrow{Q}=\frac{1}{2}\langle {\overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}0}^{2}\rangle .$$(12)
Note that corrections related to the anharmonic effects can be lately introduced to the intermediate function in equation (11), in the case of small atomic displacements compared to interatomic distances.
In equation (11), $\langle \overrightarrow{Q}\cdot {\overrightarrow{u}}_{j\prime}0\overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}t\rangle $ represents the correlation of linear displacement along $\overrightarrow{Q}$ of two atoms at different time. Though the correlation for distinct atoms j ≠ j′ is well formulated, the calculation is computationally cumbersome. Thus, in this work, a widely used approximation is utilized, called incoherent approximation, in which the correlation for distinct atoms is neglected for inelastic scattering, i.e., $$\langle \overrightarrow{Q}\cdot {\overrightarrow{u}}_{j\text{'}}0\overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}t\rangle =0\text{\hspace{1em}if\hspace{0.17em}\hspace{0.17em}}j\ne j\text{'}.$$(13)
For j = j′, the last term in equation (11) does not vanish. Its computation is quite difficult $\langle \overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}0\overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}t\rangle $. To circumvent this difficulty, the cubic approximation is used in this work. It assumes that for atoms of a solid crystalline material under consideration, the interatomic forces along all directions are isotropic. Under the cubic approximation, we have: $$\begin{array}{c}\langle \overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}0\overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}t\rangle =\frac{\hslash {Q}^{2}}{2{M}_{j}}{\displaystyle {\int}_{0}^{\infty}\frac{{\rho}_{j}(\omega )}{\omega}}\\ \times \mathrm{coth}\frac{\hslash \omega}{2{k}_{B}T}\mathrm{cos}\omega t+i\text{\hspace{0.17em}}\mathrm{sin}\omega t\text{d}\omega ,\end{array}$$(14)
where M_{j} is the mass of the atom j, ρ_{j}(ω) represents the phonon density of states (PDOS) of the atom j satisfying $${\int}_{0}^{\infty}{\rho}_{j}\omega \text{d}\omega}=1.$$(15)
The cubic approximation has the considerable merit that it enables to connect the PDOS of atoms to their dynamics. By utilizing the symmetry of the PDOS versus ω, i.e., ρ_{j} (−ω) = ρ_{j} (ω), equation (14) can be rewritten as: $$\langle \overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}0\overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}t\rangle =\frac{{\hslash}^{2}{Q}^{2}}{2{M}_{j}{k}_{B}T}{\displaystyle {\int}_{\infty}^{\infty}{P}_{j}\omega \mathrm{exp}i\omega t\text{d}\omega},$$(16)
where $${P}_{j}\omega =\frac{{\rho}_{j}\omega}{\frac{2\hslash \omega}{{k}_{B}T}\mathrm{sinh}\frac{\hslash \omega}{2{k}_{B}T}}\mathrm{exp}\frac{\hslash \omega}{2{k}_{B}T}.$$(17)
An important condition emerges from equations (16) and (17): ρ_{j} (ω) must vary as ω^{2} as ω goes to zero, to ensure the convergence of the temperaturedependent function P_{j} (ω). The integral of P_{j} (ω) provokes the DebyeWaller coefficient Λ_{j} (T) [13]: $${\Lambda}_{j}T={\displaystyle {\int}_{\infty}^{\infty}{P}_{j}\omega \text{d}\omega}.$$(18)
Based on equations (12) and (16), the DebyeWaller function ${W}_{j}\overrightarrow{Q}$ represents the dynamics of the atom j at t = 0. In addition to equation (18), ${W}_{j}\overrightarrow{Q}$ can be obtained from the corresponding DebyeWaller coefficient Λ_{j}(T) as: $${W}_{j}\overrightarrow{Q}=\frac{{\hslash}^{2}{Q}^{2}}{4{M}_{j}{k}_{B}T}{\Lambda}_{j}T.$$(19)
Thanks to the incoherent and cubic approximations, the dynamic term $\langle \overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}0\overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}t\rangle $ and the DebyeWaller function ${W}_{j}\overrightarrow{Q}$ in the intermediate function (Eq. (11)) can be obtained with equations (16) and (19), respectively. The following step consists of expanding $\langle \overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}0\overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}t\rangle $ in equation (11) in a Taylor series: $$\begin{array}{c}\mathrm{exp}\langle \overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}0\overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}t\rangle \\ {\displaystyle \sum _{n=0}^{\infty}\frac{1}{n!}{\langle \overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}0\overrightarrow{Q}\cdot {\overrightarrow{u}}_{j}t\rangle}^{n}.}\end{array}$$(20)
Equation (20) is referred to phonon expansion, originally introduced by Sjolander [19]. The first term n = 0 corresponds to elastic scattering (no energy exchange between the incident neutron and the system under consideration). The terms n ≥ 1 represent inelastic scattering contributions.
The cross sections of the coherent elastic scattering, incoherent elastic scattering and inelastic scattering are presented in the following sections.
2.2.2 Coherent elastic scattering
Important crystallographic quantities used for the calculations of the coherent elastic scattering cross section for solid crystalline materials are presented. Any crystal can be characterized by a periodic unit cell containing N_{uc} atoms. In the direct lattice, it is defined by a set of unit vectors $\overrightarrow{a},\overrightarrow{b},\overrightarrow{c}$ of lengths a, b, c and volume ${V}_{uc}=\overrightarrow{a}\cdot \overrightarrow{b}\times \overrightarrow{c}$. The angles between them are conventionally denoted α, β and γ. The position of the jth atom located at the point (x_{j}, y_{j}, z_{j}) is given by: $${\overrightarrow{p}}_{j}\text{\hspace{0.33em}}={x}_{j}\overrightarrow{a}+{y}_{j}\overrightarrow{b}+{z}_{j}\overrightarrow{c}\text{\hspace{0.33em}}{x}_{j},\text{\hspace{0.33em}}{y}_{j},\text{\hspace{0.33em}}{z}_{j}\text{\hspace{0.33em}}\in \text{\hspace{0.33em}}\left[0,\text{\hspace{0.33em}}1\right].$$(21)
In the reciprocal lattice, the unit vectors become: $${\overrightarrow{\tau}}_{a}=\frac{2\pi}{{V}_{uc}}(\overrightarrow{b}\times \overrightarrow{c})\text{\hspace{0.33em}},{\overrightarrow{\tau}}_{b}=\frac{2\pi}{{V}_{uc}}\overrightarrow{c}\times \overrightarrow{a},\text{\hspace{0.33em}}{\overrightarrow{\tau}}_{c}=\frac{2\pi}{{V}_{uc}}\overrightarrow{a}\times \overrightarrow{b}$$(22)
and indices h, k, l denote planes (hkl) orthogonal to the reciprocal lattice vector: $${\overrightarrow{\tau}}_{hkl}\text{\hspace{0.33em}}=h{\overrightarrow{\tau}}_{a}+k{\overrightarrow{\tau}}_{b}+l{\overrightarrow{\tau}}_{c}.$$(23)
Given that the Fourier transform of a constant function is a Dirac delta function: $${\int}_{\infty}^{\infty}\mathrm{exp}}(i\omega t)\text{d}t=2\pi \delta (\omega ),$$(24)
and thanks to equations (6), (11) and (20), the coherent elastic dynamic structure factor ${S}_{\text{coh}}^{\text{el}}\overrightarrow{Q},\omega $ is given by: $$\begin{array}{r}\hfill {S}_{\text{coh}}^{\text{e}1}(\overrightarrow{Q},\text{\hspace{0.33em}}\omega )=\frac{\delta (\omega )}{\hslash}{\displaystyle \sum _{j,j\text{'}=1}^{{N}_{uc}}\overline{{b}_{j}}}\cdot \overline{{b}_{{j}^{\prime}}}\mathrm{exp}i\overrightarrow{Q}\cdot {d}_{{j}^{\prime}}^{\to}{d}_{j}^{\to}\\ \hfill \times \mathrm{exp}{W}_{{j}^{\prime}}\overrightarrow{Q}\mathrm{exp}{W}_{j}\overrightarrow{Q}\text{\hspace{0.33em}}.\end{array}$$(25)
Considering a crystal with the previous notations, ${S}_{\text{coh}}^{\text{el}}\overrightarrow{Q},\omega $ can be rewritten as: $${S}_{\text{coh}}^{\text{e}1}\overrightarrow{Q},\text{\hspace{0.33em}}\omega =\frac{{2\pi}^{3}\delta \omega}{\hslash {V}_{uc}}{\displaystyle \sum _{hkl}\delta}\overrightarrow{Q}{\overrightarrow{\tau}}_{hkl}F{\overrightarrow{\tau}}_{hkl}{}^{2}$$(26)
The neutron diffraction will occur from the planes (hkl) that are oriented at the correct angle to fulfill the Bragg condition. The coherent elastic cross section per atom emerges from the sum of all the neutron scattering contributions over the N_{uc} atoms of the unit cell and plans (hkl). The analytical expression taking into account the anharmonicity and the texture properties for crystal powder is given by: $${\sigma}_{\text{coh}}^{\text{el}}E=\frac{{\pi}^{2}{\hslash}^{2}}{m{N}_{uc}{V}_{uc}E}{\displaystyle \sum _{hkl}^{E\ge {E}_{hkl}}{d}_{hkl}}F{\overrightarrow{\tau}}_{hkl}{}^{2}{\mathcal{P}}_{hkl}{\Theta}_{hkl},$$(27) $$F({\overrightarrow{\tau}}_{hkl})={\displaystyle \sum _{j=1}^{{N}_{uc}}{b}_{j}}\mathrm{exp}\frac{{\hslash}^{2}{\tau}_{hkl}^{2}}{4{M}_{j}{k}_{B}T}{\Lambda}_{j}(T){\text{e}}^{i\overrightarrow{\tau}{\overrightarrow{p}}_{j}i{c}_{123}^{j}{\tau}_{hkl}^{3}},$$(28)
in which m is the neutron mass, ${E}_{hkl}={\hslash}^{2}{\tau}_{hkl}^{2}/8m$ represents the Bragg edges, d_{hkl} = 2π/τ_{hkl} stands for the distance between adjacent planes (hkl) and b_{j} is the bound coherent scattering length.
Compared to the reported equation in the literature [9, 14], in our work, the texture or preferred orientation correction is taken into account by introducing the cylindrically symmetric PoleDensity Distribution Function (PDDF) ${\mathcal{P}}_{hkl}{\Theta}_{hkl}$ which depends on the orientation angle Θ_{kkl}, i.e., angle between the preferred orientation vector $\text{P}\overrightarrow{\text{O}}$ and the plan vectors ${\overrightarrow{\tau}}_{hkl}$: $${\Theta}_{hkl}\text{\hspace{0.33em}}=\mathrm{arccos}\frac{\text{P}\overrightarrow{\text{O}}\cdot {\overrightarrow{\tau}}_{hkl}}{\Vert \text{P}\overrightarrow{\text{O}}\Vert \Vert {\overrightarrow{\tau}}_{hkl}\Vert}.$$(29)
Various types of PDDFs reported in the literature [20] are investigated, such as the models proposed by MarchDollase [21, 22], Altomare et al. [23], Cerny, Valvoda, and Chladek [24]. The impact of introducing the PDDFs to the calculations of ${\sigma}_{\text{coh}}^{\text{el}}E$ will be illustrated in Section 4.4.
The thirdcumulant coefficient ${c}_{123}^{j}$ is introduced in the form factor to account for the anharmonicity of atoms [25]. The influence of the anharmonicity of oxygen atoms in uranium dioxide on ${\sigma}_{\text{coh}}^{\text{el}}E$ will be investigated in Section 4.3.
2.2.3 Incoherent elastic scattering
${S}_{\text{inc,}j}^{\text{el}}Q,\omega $ representing the incoherent elastic dynamic structure factor of the atom j is obtained thanks to equations (7), (11), (19), (20) and (24), which is given by: $${S}_{\text{inc},j}^{\text{el}}Q,\text{\hspace{0.33em}}\omega =\frac{\delta (\omega )}{\hslash}\overline{{b}_{j}^{2}}{(\overline{{b}_{j}})}^{2}\mathrm{exp}\frac{{\hslash}^{2}{Q}^{2}}{2{M}_{j}{k}_{B}T}{\Lambda}_{j}(T).$$(30)
The incoherent elastic scattering cross section ${\sigma}_{\text{inc,}j}^{\text{el}}E$ is obtained by integrating the dynamic structure factor ${S}_{\text{inc,}j}^{\text{el}}Q,\omega $ in energy and scattering angle. The analytical result is given by: $${\sigma}_{\text{inc},j}^{\text{el}}(E)=\frac{{A}_{j}{k}_{B}T}{4E{\Lambda}_{j}(T)}1\mathrm{exp}\frac{4E{\Lambda}_{j}(T)}{{A}_{j}{k}_{B}T}{\sigma}_{\text{inc},j},$$(31)
where A_{j} = M_{j} /m is ratio of the mass of the atom j to the neutron mass and ${\sigma}_{\text{inc,}j}=4\pi \overline{{b}_{j}^{2}}{\overline{{b}_{j}}}^{2}$ is the bound incoherent scattering cross section of the atom j. The incoherent elastic cross section can thus be determined by the DebyeWaller factor Λ_{j}(T).
2.2.4 Inelastic scattering
Under the incoherent and cubic approximations, the inelastic dynamic structure factor S^{inel}( Q, ω) can be obtained by summing the coherent and incoherent inelastic components. Thanks to equations (6), (7), (11) and (20), ${S}_{j}^{\text{inel}}Q,\omega $ representing the inelastic dynamic structure factor of the atom j is given by: $$\begin{array}{l}{S}_{j}^{\text{inel}}(Q,\text{\hspace{0.33em}}\omega )\hfill \\ =\overline{\frac{{b}_{j}^{2}}{\hslash}}{\text{e}}^{\frac{{\hslash}^{2}{Q}^{2}}{2{M}_{j}{k}_{B}T}{\Lambda}_{j}(T)}{\displaystyle \sum _{n=1}^{\infty}\frac{1}{n!}}{\frac{{\hslash}^{2}{Q}^{2}}{2{M}_{j}{k}_{B}T}{\Lambda}_{j}(T)}^{n}{\mathcal{T}}_{j}{}_{,n}(\omega ),\hfill \end{array}$$(32)
where ${\mathcal{T}}_{j,n}\omega $ has the generic form (when n ≥ 2): $${\mathcal{T}}_{j}{}_{,n}(\omega )={\displaystyle {\int}_{\infty}^{\infty}{\mathcal{T}}_{j,1}}({\omega}^{\prime}){\mathcal{T}}_{j}{}_{,n1}(\omega {\omega}^{\prime})\text{d}{\omega}^{\prime}.$$(33)
The onephonon term is given by: $${\mathcal{T}}_{j}{}_{,\text{\hspace{0.33em}}1}\text{\hspace{0.33em}}(\omega )=\text{\hspace{0.33em}}\frac{1}{{\Lambda}_{j}(T)}{P}_{j}(\omega ).$$(34)
2.3 Scattering models for free gas materials
The momentum transfer $\hslash \overrightarrow{Q}$ and the energy transfer ħω are large when neutrons scatter with free gas particles. Hence, the interaction time is short and the time dependent Heisenberg operator $\overrightarrow{R}t$ which represents the position of the particle in the intermediate function (Eq. (3)) can be approximated by its first order term [26]: $$\overrightarrow{R}(t)=\overrightarrow{R}(0)+\overrightarrow{v}t,$$(35)
where $\overrightarrow{v}$ is the speed of the particle.
Since there is no interparticle interaction for a free gas material, the incoherent approximation is automatically valid. Thus, the intermediate function 〈j, j〉 is given by: $$\langle j,\text{\hspace{0.33em}}j\rangle \text{\hspace{0.33em}}=\langle \mathrm{exp}i\overrightarrow{Q}\cdot \overrightarrow{v}t\rangle .$$(36)
The free gas material is assumed to be monatomic. From the statistical point of view, the speeds of the free gas atoms $\overrightarrow{v}$ follow the MaxwellBoltzmann (MB) distribution with: $$\overrightarrow{\upsilon}\propto \mathrm{exp}{\upsilon}^{2}/{\upsilon}_{p}^{2},$$(37)
where ${v}_{p}^{2}=2{k}_{B}T/M$ is the square of the most probable speed, T represents the temperature of the gas and M is the mass of the gas atom.
By replacing the MB distribution (Eq. (37)) in equation 36), thanks to the Gaussian integral property and the Fourier transform of the Gaussian function, the intermediate function for monatomic free gas is given by: $$\langle j,j\rangle =\mathrm{exp}{Q}^{2}{w}^{2}t/2,$$(38)
where w(t) is the width function of the free gas: $$wt=t\sqrt{\frac{{k}_{B}T}{M}}$$(39)
w(t) in equation (39) should be corrected to account for quantum effects [27, 28]: $$\tilde{w}t=\sqrt{{w}^{2}t\frac{i\hslash}{M}t}.$$(40)
Based on equations (2), (38), (39) and (40), the dynamic structure factor for monatomic free gas can be obtained by: $$SQ,\omega =\frac{\overline{{b}^{2}}}{\hslash}\sqrt{\frac{M}{2\pi {k}_{B}T{Q}^{2}}}\mathrm{exp}\text{\hspace{0.17em}}\frac{M}{2{k}_{B}T{Q}^{2}}{\omega \frac{\hslash {Q}^{2}}{2M}}^{2}.$$(41)
2.4 Scattering models for liquid water
The intermediate function 〈j′, j〉 (Eq. (3)) of liquid materials obtained by a Fourier transform of the generalized pair distribution function $G\overrightarrow{r},t$ is introduced by Van Hove [10].
By distinguishing the contribution from the same particle to different particles, Van Hove split $G\overrightarrow{r},t$ into a self part ${G}_{s}\overrightarrow{r},t$ and a distinct part ${G}_{d}\overrightarrow{r},t$: $$G\overrightarrow{r},\text{\hspace{0.17em}}t={G}_{s}\overrightarrow{r},\text{\hspace{0.17em}}t+{G}_{d}\overrightarrow{r},\text{\hspace{0.17em}}t.$$(42)
For a given particle at position $\overrightarrow{r}$ and time t, ${G}_{s}\overrightarrow{r},t\text{d}\overrightarrow{r}\text{d}t$ represents the probability of finding the same particle between $\overrightarrow{r}$ and $\overrightarrow{r}+\text{d}\overrightarrow{r}$ in a time interval t and t + dt. Similarly, ${G}_{d}\overrightarrow{r},t\text{d}\overrightarrow{r}\text{d}t$ is the probability of finding a distinct particle between $\overrightarrow{r}$ and $\overrightarrow{r}+\text{d}\overrightarrow{r}$ in a time interval t and t + dt.
Since it is difficult to predict theoretically the contribution from the distinct part, we use the incoherent approximation, as for solid crystalline materials: $${G}_{d}\overrightarrow{r},\text{\hspace{0.17em}}t\cong 0.$$(43)
To obtain an analytical form of ${G}_{s}\overrightarrow{r},t$, Vineyard [29] introduces the Gaussian approximation that assumes a weak timedependent atomic position coupling: $${G}_{s}\overrightarrow{r},\text{\hspace{0.17em}}t\cong {2\pi}^{3/2}{w}^{3}t\mathrm{exp}{r}^{2}/2{w}^{2}t,$$(44)
where 3w^{2}(t) represents the mean square displacement of the particle after time t, in which w(t) is called the width function.
Under the incoherent and Gaussian approximations, the dynamic structure factor $S\overrightarrow{Q},\omega $ for liquid materials is given by: $$S\overrightarrow{Q},w=\frac{1}{2\pi \hslash}\overline{{b}^{2}}{\displaystyle {\int}_{\infty}^{\infty}\mathrm{exp}{Q}^{2}{w}^{2}t/2\mathrm{exp}iwt\text{d}t}.$$(45)
The Q^{2} dependence of equation (45) indicates that $S\overrightarrow{Q},\omega $ is isotropic with respect to the momentum transfer, leading to: $$S\overrightarrow{Q},w\equiv SQ,w.$$(46)
The width function w(t) plays a central role in the calculation of the dynamic structure factor. In this work, we focus on liquid water which is widely used as moderator in pressurized water reactors. An analytical expression of w(t) was proposed in the sixties [30] for liquid water. However, introducing such an expression in equation (45) leads to numerical issues at zero energy transfer. These mathematical problems can be solved by decoupling the mean square displacement into a diffusion part w_{diff}(t) and a solidlike vibrational part w_{vib}(t): $${w}^{2}t={w}_{\text{diff}}^{2}t+{w}_{\text{vib}}^{2}t.$$(47)
The corresponding diffusion S_{diff}(Q, ω) and solidlike S_{vib}(Q, ω) parts of the dynamic structure factor are calculated with equation (45) by introducing different models for w_{diff}(t) and w_{vib}(t). According to the convolution theorem for Fourier transforms, the composite dynamic structure factor is calculated as follows: $$SQ,\omega =\frac{\hslash}{\overline{{b}^{2}}}{\displaystyle {\int}_{\infty}^{\infty}{S}_{\text{diff}}Q,{\omega}^{\prime}{S}_{\text{vib}}Q,\omega {\omega}^{\prime}\text{d}{\omega}^{\prime}},$$(48)
where $\hslash /\overrightarrow{{b}^{2}}$ is a scale coefficient, the values of S_{diff}(Q, ω) and S_{vib}(Q, ω) depend on the density of states ρ(ω) of liquid water, which is partitioned in two parts: $$\rho \omega ={w}_{d}{\rho}_{\text{diff}}\omega +{w}_{\upsilon}{\rho}_{\text{vib}}\omega .$$(49)
The diffusion w_{d} and solidlike vibrational w_{υ} weights satisfy the condition: $${w}_{d}+{w}_{\upsilon}=1.$$(50)
The diffusion weight w_{d} can be calculated with the model involving the fluidicity [31, 32], proposed by Marquez Damian (private communications, 2021). The rototranslational diffusion model with random jump diffusion correction for S_{diff}(Q, ω) has been investigated in reference [33]. The formalisms are recalled here.
In view of using a descriptive analytical model, the diffusion term in equation (48) related to liquid water is approximated by a translational diffusion motion corrected for rotational diffusion contributions: $${S}_{\text{diff}}Q,\omega ={\displaystyle {\int}_{\infty}^{\infty}{S}_{\text{trans}}Q,{\omega}^{\prime}{S}_{\text{rot}}Q,\omega {\omega}^{\prime}\text{d}{\omega}^{\prime}}.$$(51)
The translational diffusion behavior of the water molecule was interpreted with the EgelstaffSchofield diffusion model [34], in which the width function is given by: $${w}_{\text{trans}}^{2}t=\frac{2\text{\hspace{0.17em}}c\text{\hspace{0.17em}}{w}_{d}\text{\hspace{0.17em}}{\hslash}^{2}}{M{k}_{B}T}\sqrt{{\frac{{k}_{B}T}{\hslash}}^{2}+{c}^{2}+\frac{1}{4}}c$$(52)
where M is the mass of the scattering atom and c is a dimensionless diffusion parameter. Its meaning has been revisited in reference [35] by introducing the SingwiSjolander residence time το [36] to account for the noncontinuous motions of the water molecules: $$c=\frac{M}{\hslash {w}_{d}}\mathcal{D}Q,$$(53)
with $$\mathcal{D}Q=\frac{D}{1+D{Q}^{2}{\tau}_{0}},$$(54)
in which D represents the selfdiffusion coefficient of the liquid water. The dynamic structure factor S_{trans}(Q, ω) can be derived analytically by introducing equation (52) in equation (45) with correction to fulfill the principle of detailed balance: $${S}_{\text{trans}}(Q,\omega )=\frac{\overline{{b}^{2}}}{{k}_{B}T}\frac{\mathcal{D}{Q}^{2}}{\pi}{\text{e}}^{\frac{M{\mathcal{D}}^{2}{Q}^{2}}{w{d}^{k}{B}^{T}}+\frac{\hslash \omega}{2{k}_{B}T}}\frac{\sqrt{{c}^{2}+1/4}}{\sqrt{\omega +(\mathcal{D}{Q}^{2})}}{K}_{1}(x),$$(55)
in which K_{1}(x) is the modified Bessel function of the second kind, with: $$x=\frac{\hslash}{{k}_{B}T}\sqrt{{c}^{2}+\frac{1}{4}}\sqrt{{\omega}^{2}+{(\mathcal{D}{Q}^{2})}^{2}.}$$(56)
The rotational diffusion correction in equation (51) is related to the dynamics of hydrogen bonds. Their contributions is approximated with the Sears expansion [37] expressed in terms of spherical Bessel functions j_{l}(QR_{cm}) of order l: $$\begin{array}{c}{S}_{\text{rot}}(Q,\omega )={j}_{0}^{2}(Q{R}_{\text{cm}})\delta (\omega )\\ +\frac{1}{\pi}{\displaystyle \sum _{l=1}^{\infty}(2l+1)}\\ \times {j}_{l}^{2}(Q{R}_{\text{cm}})\frac{l(l+1)(\hslash /6{\tau}_{R})}{{\omega}^{2}+[l(l+1)(\hslash /6{\tau}_{R})]2}.\end{array}$$(57)
This model assumes a continuous and isotropic rotation of the hydrogen atoms arround the center of mass of the molecule with a rotational relaxation time τ_{R}. The rotation of the water molecule is confined to a spherical surface of radius R_{cm}, whose value is often taken equal to the OH intramolecular distance.
The computation of the solidlike S_{vib}(Q, ω) part of the dynamic structure factor is equivalent to the inelastic part of dynamic structure factor for the solid crystalline materials.
2.5 Converting S(Q, ω) into thermal scattering laws
The dynamic structure factor for nuclear materials (e.g., beryllium, iron, liquid water) is routinely used to simulate the neutron behaviors in light water reactors. The neutron transport formalism implemented in dedicated MonteCarlo or deterministic code systems relies on the symmetric form of the thermal scattering laws S(α, β) [12, 38], which is defined as a function of the dimensionless parameters α: $$\alpha =\frac{{\hslash}^{2}{Q}^{2}}{2M{k}_{B}T}=\frac{{E}^{\prime}+E2\mu \sqrt{{E}^{\prime}E}}{A{k}_{B}T},$$(58)
and β: $$\beta =\frac{\hslash \omega}{{k}_{B}T}=\frac{E{E}^{\prime}}{{k}_{B}T}$$(59)
where µ = cos(θ) is the cosine of the scattering angle θ in the laboratory system and A is the ratio of the mass M of the scattering atom to the neutron mass. The change of sign in equation (59) is introduced to keep neutron energy gains positive in neutron transport calculations. The relationship between S(α, β) and S(Q, ω) is given by the following expression: $$\mathcal{S}(\alpha ,\beta )=\frac{{k}_{B}T}{{b}^{2}}{\text{e}}^{\frac{\hslash \omega}{2{k}_{B}T}}S(Q,\omega ).$$(60)
Note that the thermal scattering laws S(α, β) are symmetric in energy.
3 Implementation of CINEL
CINEL is developed in Python by using JupyterLab [39]. JupyterLab is a new generation of user interface [40] which is an opensource tool enabling to mix live code, texts, mathematical equations and interactive graphics. An example of JupyterLab interface is given in Figure 1 for illustration. This interactive development environment facilitates the verification and visualization of calculations.
CINEL consists of three modules: Cubic, INELastic and SVT. The data flow presenting the generations of the thermal scattering laws (TSL) in ENDF6 format by using the Cubic and INELastic modules, and the scattered neutron energy distributions via the SVT module are shown in Figure 2. These modules will be briefly presented in the following.
Fig. 1 JupyterLab interface. Live code, texts, mathematical equations and interactive graphics are mixed in Jupyter Notebooks which are integrated in JupyterLab together with blocks like terminal and text editor. 
Fig. 2 Data flow presenting the generation of the thermal scattering laws in ENDF6 format and scattered neutron energy distributions by using CINEL. 
3.1 Cubic and INELastic modules
Various physical quantities are necessary for the generation of the TSL such as scattering length, positions of atom in the unit cell, general settings and temperaturedependent parameters including the lattice parameters, PDOS and the phonon expansion order. Optional parameters are preferred orientation of crystal PO, cylindrically symmetric PDDF ${\mathcal{P}}_{hkl}{\Theta}_{hkl}$ and the thirdcumulant coefficient ${c}_{123}^{\text{O}}$.
The Cubic module calculates the DebyeWaller coefficients Λ(T) and the coherent elastic scattering cross sections ${\sigma}_{\text{coh}}^{\text{el}}$. Λ(T) is computed from the PDOS (Eq. (18)) since the cubic approximation is assumed. ${\sigma}_{\text{coh}}^{\text{el}}$ is calculated by directly looping through the indexes hkl until the interplanar distance d_{hkl} of the corresponding plan is smaller than 0.1 Å. The default threshold value is chosen to be 0.1 Å because the contribution from plans with d_{hkl} < 0.1 Å to ${\sigma}_{\text{coh}}^{\text{el}}E$ can be neglected, as discussed in reference [9]. Nonetheless, the threshold value can always be modified by the users. During the loop over hkl, plans with identical interplanar distance d_{hkl} and form factor F(τ_{hkl}) are grouped together since they contribute equally to ${\sigma}_{\text{coh}}^{\text{el}}E$ (Eq. (27)). The number of identical plans (hkl) is referred to multiplicity M_{hkl}. These physical quantities are calculated by CINEL and stored in ASCII file, as illustrated in Figure 3.
If the incident neutron energy E is provided, the diffraction angle θ_{hkl} corresponding to the Bragg edge E_{hkl} is given by: $${\theta}_{hkl}=\mathrm{arccos}(12{E}_{hkl}/E).$$(61)
The orientation angle Θ_{hkl} and ${\mathcal{P}}_{hkl}{\Theta}_{hkl}$ can be obtained if the preferred orientation of crystal and the PDDF model are given in input.
The INELastic module calculates the S(α, β) tables with the phonon expansion method. The EgelstaffSchofield diffusion model and the recombination of the diffusive and vibrational parts are also implemented in INELastic. Physical constants such as Boltzmann constant and common functions used by Cubic and INELastic are available via the library CINELLib.
To improve the performance of Cubic and INELastic, a JustinTime (JIT) compiler named Numba [41] is used. Numba allows to reduce the computational time by directly adding simple Python syntaxes to the initial functions without rewriting them in lowlevel languages. In addition, Numba supports Compute Unified Device Architecture (CUDA) graphics processing unit (GPU) programming, which enables these two modules to benefit from the powerful computing capability of GPU.
GPU programming, introduced in the early 2000s, has been developed tremendously during the last two decades. In our work, the use of GPU enables to accelerate the phonon expansion calculation. To test the speedup with GPU, the TSL of hydrogen bound in water molecule (H in H_{2}O) at room temperature were calculated with a phonon expansion order N_{phonon} = 2000. Two kinds of GPU are used: NVIDIA® QUADRO® K620 with 2 gigabyte memory and 384 cores; NVIDIA® QUADRO® K6000 with 12 gigabyte memory and 2880 cores. The comparisons of the computational time are presented in Table 1. High performing GPU allows to significantly reduce the computational time down to 2 minutes.
Approximations adopted in codes are summarized in Table 2. For strongly coherent crystalline materials such as pyrolytic graphite, the introduction of the coherent onephonon correction [8] to the incoherent approximation enables to better reproduce the experimental inelastic scattering cross sections as shown in reference [42]. The calculations of onephonon coherent inelastic scattering function necessitates the dispersion relations of materials. This correction will be implemented in the future version of CINEL. The cubic approximation will also be eliminated since it only necessitates the partial PDOS which can be calculated by DFT codes.
Fig. 3 Screenshot of cross section file generated by the Cubic module of CINEL for uranium dioxide at 296 K. Columns 1 to 11 represent respectively: index h, k, l, interplanar distance d_{hkl} (in Å), diffraction angle θ_{hkl} (in degree), orientation angle Θ_{hkl} (in degree), PoleDensity Distribution Function value ${\mathcal{P}}_{hkl}{\text{\Theta}}_{hkl}$, form factor square F(τ_{hkl})^{2} (in barn), multiplicity M_{hkl}, Bragg edge E_{hki} (in eV) and ${E}_{hkl}{\sigma}_{\text{coh}}^{\text{el}}{E}_{hkl}$ (in eV.barn), which serve as calculations of the coherent elastic cross sections. 
Comparison of the computational time of the TSLs of H in H_{2}O at room temperature with a phonon expansion order N_{phonon} = 2000.
Comparison of codes or platforms which enable to calculate the TSL.
3.2 SVT module
The SVT module is implemented with the Sampling the Velocity of the Target nucleus (SVT) algorithm [43, 44], in which the neutronnucleus scattering is approximated by the twobody kinetics as illustrated in Figure 4. The neutron mass and the mass of target nucleus are respectively m and M with M = Am. The speeds of neutron and target nucleus are respectively $\overrightarrow{v}$ and $\overrightarrow{V}$. Angles between two speeds are represented by θ. The subscripts of $\overrightarrow{\upsilon}$ or θ indicate their belonging frame and the superscript represents the physical quantities after scattering. For example, ${\overrightarrow{\upsilon}}_{\text{TR}}$ and ${{\overrightarrow{\upsilon}}^{\prime}}_{\text{TR}}$ are respectively the speeds of neutron in the targetatrest (TR) frame before and after scattering. For the sake of clarity, the subscripts of quantities in the laboratory (LAB) frame are omitted.
As presented in the data flow of CINEL (Fig. 2), the incident neutron energy E and the velocities of the target nuclei V are necessary for the calculation of the scattered neutron energy distributions. V are either sampled from trajectory files obtained by molecular dynamics (MD) calculations, or from a MaxwellBoltzmann (MB) distribution ${\mathcal{M}}_{T}(V)$ which is determined by the temperature or effective temperature of the nuclei. In the latter case, the cosine of the angle between the incident neutron and the target nucleus µ = cos(θ) before collision is sampled uniformly on [−1; 1].
The equations implemented in the SVT module to calculate the scattered neutron energy E′ and the cosine of the scattering angle µ′ = cos(θ′) are summarized in the following. The first step consists of calculating the relative velocity between the incident neutron and the target nucleus: $${\upsilon}_{\text{TR}}=\Vert \overrightarrow{\upsilon}\overrightarrow{V}\Vert =\sqrt{{\upsilon}^{2}2\mu \upsilon V+{V}^{2}},$$(62)
where $\upsilon =\sqrt{2E/m}$ represents the incident neutron velocity.
This step enables to transform the neutrontarget system from the LAB frame to the TR frame, as shown in Figure 4, from the top plot to the bottom left plot. Next step is to transform the studied system from the TR frame to the centerofmass (CM) frame in which the neutrontarget collision is assumed to be isotropic. The speed of the CM frame in the TR frame ${\overrightarrow{c}}_{\text{TR}}$ is computed by: $${\overrightarrow{c}}_{\text{TR}}=\frac{{\overrightarrow{\upsilon}}_{\text{TR}}}{1+A}.$$(63)
Then the speeds of neutron and target in the CM frame ${\overrightarrow{\upsilon}}_{\text{CM}}$ and ${\overrightarrow{V}}_{\text{CM}}$ can be obtained respectively by subtracting the speed of the CM frame: $${\overrightarrow{\upsilon}}_{\text{CM}}=\frac{A{\overrightarrow{\upsilon}}_{\text{TR}}}{1+A},$$(64)
and $${\overrightarrow{V}}_{\text{CM}}=\frac{{\overrightarrow{\upsilon}}_{\text{TR}}}{1+A}.$$(65)
Since the neutrontarget collision is assumed to be isotropic in the CM frame, the cosine of the scattering angle μ_{CM} = cos(θ_{CM}) is sampled uniformly on the interval [−1; 1]. The total momentum and the kinetic energy of the neutrontarget system are conserved before and after the collision. It can be easily shown that the scalar velocities of the neutron and the target remain invariant in the CM frame, i.e., ${\upsilon}_{\text{CM}}={{\upsilon}^{\prime}}_{\text{CM}}$ and ${V}_{\text{CM}}={{V}^{\prime}}_{\text{CM}}$, as illustrated in the bottom right plot of Figure 4.
The following step is to transform back the neutrontarget system from the CM frame to the TR frame. Based on geometric relations, speeds in the TR and CM frames lie on the same plan (P′). The scalar velocity of the postcollision neutron ${{\upsilon}^{\prime}}_{\text{TR}}$ and the cosine of the scattering angle ${{\mu}^{\prime}}_{\text{TR}}=\mathrm{cos}({{\theta}^{\prime}}_{\text{TR}})$ in the TR frame can be respectively calculated by: $${{\upsilon}^{\prime}}_{\text{TR}}=\frac{R(1+2{\mu}_{\text{CM}}A+{A}^{2})({\upsilon}^{2}2\mu \upsilon V+{V}^{2})}{1+A}$$(66)
and $${{\mu}^{\prime}}_{\text{TR}}=\frac{1+{\mu}_{\text{CM}}A}{\sqrt{1+2{\mu}_{\text{CM}}A+{A}^{2}}}.$$(67)
The final step consist of transforming the neutrontarget system back to the LAB frame to calculate υ′ and μ′ = cos(θ′). The target speed $\overrightarrow{V}$ lying on the plan (P) which does not coincide with the plan (P′) except for μ = ±1. Thus, an azimuthal angle ϕ_{CM} must be sampled. ϕ_{CM} is sampled uniformly on the interval ] − π; π] based on the assumption that the neutron scattering is isotropic in the CM frame. In summary, there are four variables sampled in the SVT algorithm [45]: V and µ are related to the motion of target nuclei, and µ_{CM} and ϕ_{CM} are used in the collision kinematics.
µ_{r} = cos(θ_{r}) is the cosine of the angle between the relative speed $\overrightarrow{\upsilon}\overrightarrow{V}$ and the speed of the target nucleus $\overrightarrow{V}$: $${\mu}_{r}=\frac{\mu \upsilon V}{\sqrt{{\upsilon}^{2}2\mu \upsilon V+{V}^{2}}}.$$(68)
The speed of the scattered neutron in the TR frame ${{\overrightarrow{\upsilon}}^{\prime}}_{\text{TR}}$ is decomposed into two components with one parallel to the relative speed $\overrightarrow{\upsilon}\overrightarrow{V}$ and the other perpendicular to it, i.e., ${{\overrightarrow{\upsilon}}^{\prime}}_{\text{TR}}={{\overrightarrow{\upsilon}}^{\prime}}_{\text{TR,}\parallel}+{{\overrightarrow{\upsilon}}^{\prime}}_{\text{TR,}\perp}$. Equally, the speed of the target nucleus is decomposed into $\overrightarrow{V}={\overrightarrow{V}}_{\parallel}+{\overrightarrow{V}}_{\perp}$. Hence, the speed of the scattered neutron in the LAB frame $\overrightarrow{{\upsilon}^{\prime}}$ can be obtained by: $${{\overline{\upsilon}}^{\prime}}_{\text{TR}}={{\overline{\upsilon}}^{\prime}}_{\text{TR}}+\overrightarrow{V}={{\overline{\upsilon}}^{\prime}}_{\text{TR,}\parallel}+{{\overline{\upsilon}}^{\prime}}_{\text{TR,}\perp}+{\overrightarrow{V}}_{\parallel}+{\overline{V}}_{\perp}.$$(69)
Based on geometric relations, the velocity of the scattered neutron in the LAB frame υ′ can be calculated by: $${\upsilon}^{\prime}={({{\upsilon}^{\prime}}_{\text{TR}})}^{2}+2{\mu}_{V}{{\upsilon}^{\prime}}_{\text{TR}}V+{V}^{2},$$(70)
where µ_{V} = cos(θ_{V}) represents the cosine of the angle between ${{\upsilon}^{\prime}}_{\text{TR}}$ and $\overrightarrow{V}$, which is obtained thanks to the decompositions of ${{\upsilon}^{\prime}}_{\text{TR}}$ and $\overrightarrow{V}$ in equation (69) and $${\mu}_{V}={{\mu}^{\prime}}_{\text{TR}}{\mu}_{r}+\sqrt{1{({{\mu}^{\prime}}_{\text{TR}})}^{2}}\sqrt{1{\mu}_{r}^{2}}\mathrm{cos}({\varphi}_{\text{CM}}).$$(71)
µ_{V} can be computed thanks to the calculations of ${{\mu}^{\prime}}_{\text{TR}}$ and µ_{r} from equations (67) and (68). Since ${{\upsilon}^{\prime}}_{\text{TR}}$ can be obtained from equation (66), the velocity of scattered neutron in the LAB frame υ′ in equation (70) can be determined. Finally, the cosine of the scattering angle µ′ = cos(θ′) can be calculated by: $$\begin{array}{c}{\mu}^{\prime}=\frac{({\overrightarrow{\upsilon}}_{\text{TR}}+\overrightarrow{V})({{\overrightarrow{\upsilon}}^{\prime}}_{\text{TR}}+\overline{V})}{\upsilon {\upsilon}^{\prime}}\\ =\frac{{\upsilon}_{\text{TR}}{{\upsilon}^{\prime}}_{\text{TR}}{{\mu}^{\prime}}_{TR}+{\upsilon}_{\text{TR}}{V}_{{\mu}_{r}}+{{\upsilon}^{\prime}}_{\text{TR}}{V}_{\mu V}+{V}^{2}}{\upsilon {\upsilon}^{\prime}}.\end{array}$$(72)
Fig. 4 Schematic representation of neutronnucleus scattering in LAB (top plot), TR (bottom left plot) and CM (bottom right plot) frames. Detail descriptions can be found in the text. 
4 Numerical validations of the cubic and INELastic modules
As presented in the data flow of CINEL in Figure 2, the DebyeWaller factor Λ(T) and the coherent elastic scattering cross section ${\sigma}_{\text{coh}}^{\text{el}}E.\text{\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}}\Lambda T$ computed from the Cubic module, and the symmetric TSL (S(α, β)) calculated from the INELastic module are merged into a file in ENDF6 format. The ENDF files are later processed by using the THERMR module of NJOY to obtain the inelastic scattering cross sections. The numerical validations of the Cubic and INELastic modules of CINEL were performed by comparing the calculated scattering cross sections of solid crystalline materials with the ENDF/BVIII.0 database and the NJOY+NCrystal library.
4.1 Automatic generation of α and β grids: UO_{2}
Uranium dioxide (UO_{2}) is the major component of the UOX fuel. In the ENDF/BVIII.0 database, UO_{2} as ideal fluorite structure with Fm3m symmetry is illustrated in Figure 5a. The atom positions in the unit cell are presented in Table 3.
The PDOS of ^{238}U in UO_{2} and ^{16}O in UO_{2} for the latest ENDF/BVIII.0 database are generated by using the ab initio lattice dynamics methods [48]. The two peaks of the PDOS of ^{238}U in UO_{2} in Figure 6 correspond to the acoustic modes of the uranium atoms. The peak at nearly 10 meV corresponds to the transverse acoustic (TA) mode while the peak at nearly 20 meV corresponds to the longitudinal acoustic (LA) mode. The peaks at nearly 30 meV, 55 meV and 70 meV are dominated by the optical modes of the oxygen atoms. The peak at roughly 30 meV stems from the transverse and longitudinal optical (TO1, LO1) modes. The peak at nearly 55 meV corresponds to the transverse optical mode (TO2) and the last peak at nearly 70 meV corresponds to the longitudinal optical mode (LO2).
CINEL is able to generate automatically the α and β grids for the calculations of the S(α, β) tables. Important steps are summarized in the following. The first step consists of determining the minimum and the maximum values of the α and β grids. According to equation (58), the minimum momentum transfer is obtained when µ = 1: $${\alpha}_{\text{min}}=\frac{{\sqrt{E\prime}+\sqrt{E}}^{2}}{A{k}_{B}T}.$$(73)
The scattered neutron energy is regarded as having an energy gain δE compared to the incident energy, i.e., E′ = E + δE, then equation (73) can be rewritten as: $${\alpha}_{\text{min}}=\frac{\delta E}{4A{k}_{B}T}.$$(74)
The momentum transfer achieves its maximum value when µ = −1: $${\alpha}_{\text{max}}=\frac{{\sqrt{E\prime}+\sqrt{E}}^{2}}{A{k}_{B}T}.$$(75)
Since the incident neutron energy has a threshold E_{thr}, α_{max} in equation (75) is obtained: $${\alpha}_{\text{max}}=\frac{4{E}_{\text{thr}}}{A{k}_{B}T}.$$(76)
In the case of the energy transfer β, based on equation 59), the minimum value β_{max} and the maximum value β_{max} are respectively: $${\beta}_{\text{min}}=0\text{\hspace{1em}when\hspace{0.17em}}E\prime =E,$$(77)
and $${\beta}_{\text{max}}=\frac{{E}_{\text{thr}}}{{k}_{B}T}\text{\hspace{1em}when\hspace{0.17em}}E={E}_{\text{thr}}\text{\hspace{0.17em}and\hspace{0.17em}}E\prime =0.$$(78)
The second step is to fix the number of grids N_{α} and N_{β} and the distribution modes: logarithmic, linear or mixed.
In the case of UO_{2}, the value of δE is set to be 2.8 meV. The energy threshold E_{thr} = 5 eV, N_{α} = 150 and N_{β} = 200. Note that these parameters can be modified by the users.
Results are presented in Figure 7. The total scattering cross sections of UO_{2} for a wide temperature range (from room temperature to 1200 K) obtained by using the Cubic and INELastic modules are in very good agreement with the ENDF/BVIII.0 database.
Fig. 5 Left hand plot presents the cubic unit cell of UO_{2}: a = b = c and α = β = γ = 90°. Blue and red spheres represent respectively uranium and oxygen atoms. Right hand plot presents the tetragonal unit cell of MgH_{2} or MgD_{2}: b = c and α = β = γ = 90°. Orange spheres represent magnesium atoms. Pink spheres represent hydrogen or deuterium atoms. The drawings of these two figures are with the aid of the Crystal Toolkit in the Materials Project [46]. 
Fig. 6 PDOS of ^{238}U in UO2 and ^{16}O in UO2, obtained from ab initio lattice dynamic methods [48]. 
Fig. 7 Comparison of the total scattering cross sections of UO2 for temperatures from 296 K to 1200 K. The resonances of uranium are not included in the calculated cross sections. 
4.2 Mixed elastic format: MgH_{2} and MgD_{2}
Magnetism hydride (MgH_{2}) is investigated as a potential candidate for cold neutron reflectors [51]. Magnetism deuteride (MgD_{2}) shows improvements on neutron slowing down since the neutron capture cross section of deuterium is lower than that of hydrogen. MgH_{2} and MgD_{2} share the same crystal structure which is illustrated in Figure 5b. The PDOS of components in MgH_{2} and MgD_{2} generated by Campi and Bernasconi from University of MilanoBicocca within the density functional perturbation theory [49,50], are presented in Figure 8.
To overcome the shortage of the ENDF6 format in which the coherent and incoherent elastic components cannot be encompassed in the same file, the scientific community has proposed recently a mixed elastic format. This new option is supported in CINEL when generating the ENDF files. The mixed elastic cross sections are illustrated with MgD_{2} in Figure 9. The comparisons of the CINEL calculated total cross sections of MgH_{2} and MgD_{2} with the NJOY+NCrystal library are presented in Figure 10. An excellent agreement with the NJOY+NCrystal library is obtained.
The total cross sections (elastic and inelastic) are in excellent agreement with the ENDF/BVIII.0 database and the NJOY+NCrystal library, which validates the Cubic and the INELastic modules. The symmetric $\mathcal{S}\alpha ,\beta $ tables generated by using the INELastic module are in excellent agreement with the ENDF/BVIII.0 database, as illustrated in Xu’s thesis [52].
Fig. 8 PDOS of ^{24}Mg in ΜgΗ_{2}, ^{1}H in ΜgΗ_{2}, ^{24}Mg in ΜgD_{2} and D in MgD_{2}, generated by using the density functional perturbation theory [49,50]. 
Fig. 9 Total scattering cross section (solid line) and the coherent/incoherent elastic/inelastic components (dot lines) of MgD_{2} at 293.6 K. 
Fig. 10 Comparison of the total scattering cross sections of MgH_{2} and MgD_{2} for temperatures from 20 K to 600 K. 
4.3 Anharmonicity of oxygen
The anharmonic motion of the oxygen atoms in UO_{2} with fluorite symmetry (illustrated in Fig. 5a) is taken into account by using the thirdcumulant coefficients ${c}_{123}^{j}$[25]. ${c}_{123}^{\text{U}}$ are zero due to the symmetry of UO_{2}. The impact of introducing the nonzero coefficient ${c}_{123}^{\text{O}}$ on the coherent elastic cross sections is illustrated in Figure 11. The positions of the Bragg diffraction peaks remain unchanged. However, the impact on the coherent elastic cross section ${\sigma}_{\text{coh}}^{\text{el}}E$ increases with the incident neutron energy E until ${\sigma}_{\text{coh}}^{\text{el}}E$ reaches the asymptotic part (when E is larger than about 0.1 eV). The impact of ${c}_{123}^{\text{O}}$ on the MonteCarlo simulations of uranium dioxide is discussed in reference [53].
Fig. 11 Comparison of coherent elastic scattering cross sections of UO_{2} with different thirdcumulant coefficient values: ${c}_{123}^{\text{O}}=0,{10}^{3},5\times {10}^{3}$. 
4.4 Texture: zirconium
Texture or preferred orientation presented in solid crystalline materials impacts the neutron diffraction compared to materials with random orientation. In this work, the texture is taken into account by the cylindrically symmetric PoleDensity Distribution Function (PDDF) P_{hkl}(Θ_{hkl}). The PDDFs proposed by MarchDollase [21, 22] and Altomare et al. [23] are respectively given by: $${\mathcal{P}}_{hkl}^{\text{MD}}{\text{\Theta}}_{hkl}={{P}_{1}^{2}{\text{cos}}^{2}{\text{\Theta}}_{hkl}+{P}_{1}^{1}{\text{sin}}^{\text{2}}{\text{\Theta}}_{hkl}}^{3/2},$$(79)
and $${\mathcal{P}}_{hkl}^{\text{Al}}{\text{\Theta}}_{hkl}=\mathrm{exp}{P}_{1}\text{cos}2\text{\Theta},$$(80)
where P_{1} is given by the users.
The comparison of the PDDFs ${\mathcal{P}}_{hkl}^{\text{MD}}{\text{\Theta}}_{hkl}$ and ${\mathcal{P}}_{hkl}^{\text{Al}}{\text{\Theta}}_{hkl}$ is illustrated in Figure 12. These PDDFs are applied to the zirconium crystalline with preferred orientation $\overrightarrow{\text{PO}}=\left[0,0,1\right]$ at room temperature. The obtained coherent elastic scattering cross sections of zirconium are compared to the calculations with random orientation in Figure 13. The positions of Bragg diffraction peaks are invariant. The modifications of the peak intensities correspond to the introduced PDDFs presented in Figure 12.
Fig. 12 Comparison of cylindrically symmetric PDDF proposed by MarchDollase with parameter P_{1}= 0.6 and Altomare et al. with P_{1} = 1.5. 
Fig. 13 Comparison of coherent elastic scattering cross section ${\sigma}_{\text{coh}}^{\text{el}}E$ of zirconium with random orientation and with preferred orientation $\overrightarrow{\text{PO}}=\left[0,0,1\right]$. ${\sigma}_{\text{coh}}^{\text{el}}E$ obtained by using the cylindrically symmetric PDDFs proposed by MarchDollase and Altomare et al. are respectively presented. 
5 Numerical validations of the SVT module
5.1 Comparison with the Tripoli4® calculations
In this part, the algorithm implemented in the SVT module of CINEL (cf. Sect. 3.2) is validated by comparing the scattered neutron energy distributions of ^{1}H and ^{16}O in form of free gas materials with the results obtained by using the MonteCarlo neutron transport code Tripoli4^{®} [54].
The physical model employed in the Tripoli4^{®} simulations is briefly presented. A capillary sample with a radius equal to 5pm and a height of 10 pm is placed in the center of the simulation model. The dimensions are chosen to minimize the impact of multiple scattering. The incident neutron energy is monoenergetic which is set to be 5 eV for free ^{1}H and 6.67 eV for free ^{16}O. The outgoing neutrons are detected at θ’ = 10°, 45° and 90°. The Tripoli4^{®} simulations are performed at room temperature with T = 294 K for ^{1}H and T = 298 K for ^{16}O, respectively.
The incident neutron energy E is set to be 5 eV and 6.67 eV for ^{1}H and ^{16}O in the SVT module, respectively. The velocities of the target nuclei are sampled from a MB distribution ℳ_{T} (V) with T = 294 K for ^{1}H and T = 298 K for ^{16}O: $${\mathcal{M}}_{T}V=\frac{4}{\sqrt{\pi}}{B}^{3}{V}^{2}\mathrm{exp}{B}^{2}{V}^{2},$$(81)
where $$B=\sqrt{\frac{M}{2{k}_{B}T}.}$$(82)
The scattered neutron energy distributions obtained by using the SVT module of CINEL and Tripoli4^{®} are compared in Figure 14. Overall good agreement is obtained for all the angles, which validates the twobody kinematic equations implemented in the SVT module of CINEL.
Fig. 14 Comparison of the scattered neutron energy distributions of ^{1}H and ^{16}O (free gas materials) obtained by using the SVT module of CINEL and Tripoli4 for incident neutron energies E = 5 eV and 6.67 eV, respectively. 
5.2 Comparison with analytical results of the free gas model
The DDXS for free gas materials is obtained by converting the dynamic structure factor $S\overrightarrow{Q},\omega $ (Eq. (41)) to dimensionless S(α, β): $$\frac{{\text{d}}^{2}{\sigma}^{\text{FGM}}}{\text{d}\mu \text{\hspace{0.17em}d}E\prime}=\frac{{\sigma}_{b}}{2{k}_{B}T}\sqrt{\frac{E\prime}{E}}\frac{1}{\sqrt{4\pi \alpha}}\text{exp}\frac{\alpha +\beta}{4\alpha}.$$(83)
The analytical results for free ^{1}H with incident neutron energies E = 1 and 5 eV are compared with the scattered neutron energy distributions obtained by using the SVT module in Figure 15. Excellent agreement is obtained for E = 5 eV, because the SVT algorithm satisfies the assumptions made in the analytical free gas model: the momentum and energy transfers are large, and the velocities of target nuclei follow a MB distribution, as presented in Section 2.3. Discrepancies are observed for large scattering angles when E = 1 eV. Since in this case, the first assumption is not satisfied. Similar results are found for free ^{16}O, as presented in Figure 16. These results enable to validate the sampling of the velocities from a M–B distribution in the SVT module.
The comparison with Tripoli4^{®} calculations enables to validate the twobody kinematic equations implemented in the SVT module (as presented in Sect. 3.2). We remark in Section 2.3 that a simplified $\mathcal{S}\alpha ,\beta $ formula for free gas model is obtained thanks to the assumption of MB distribution. Thus, the comparison with the analytical free gas models allows to validate the sampling of the velocities from a MB distribution in the SVT module. Therefore, we think that these two comparisons are necessary.
Fig. 15 Comparison of the scattered neutron energy distributions of free ^{1}H obtained by using the SVT module of CINEL and analytical results of the free gas model for incident neutron energies E = 1 and 5 eV. 
Fig. 16 Comparison of the scattered neutron energy distributions of free ^{16}O obtained by using the SVT module of CINEL and analytical results of the free gas model for incident neutron energies E = 0.1 and 1 eV. 
5.3 Comparison with velocities sampled from MD trajectory files
As presented in the data flow of CINEL (Fig. 2), the SVT module is able to sample the velocities from the MD calculated trajectory files. In this work, the trajectory files of liquid water generated by using the MD code GROMACS [55] are used to validate the velocity sampling.
The generation of the trajectory files for liquid water by using GROMACS is described in Scotta’s thesis [56]. The physical model used in the MD simulation consists of in total 512 water molecules within a cubic box with side length 2.48 nm. We calculated the number volume density of water at 294 K by ${\rho}_{V}=\frac{{N}_{A}}{{M}_{{\text{H}}_{2}\text{O}}}{\rho}_{{\text{H}}_{2}\text{O}}=33.5\text{\hspace{0.17em}}{\text{nm}}^{3}$ where Na is the Avogadro constant, ${M}_{{\text{H}}_{2}\text{O}}$ and ${\rho}_{{\text{H}}_{2}\text{O}}$ are respectively the molar mass and mass density of water. Then the size of the box is determined by (512/33.5)^{1/3} = 2.48 nm. TIP4P/2005f [57], a flexible water potential model taking into account the intermolecular and intramolecular interactions, is adopted to quantify the forces between atoms. Thanks to this potential, the positions and velocities of the hydrogen and oxygen atoms in the water molecules are obtained for a time step of 0.6 fs and they are stored in the trajectory files with a time duration of 100 ps.
Each water molecule is composed of two hydrogen atoms and one oxygen atom. For each time step, the partial velocities along the x, y and z directions are respectively V_{j,x}, V_{j,y} and V_{j,z}, with j representing the hydrogen atom H or oxygen atom O. The velocities V_{j} is thus computed via the partial components as follows: $${V}_{j}=\sqrt{{V}_{j,x}^{2}+{V}_{j,y}^{2}+{V}_{j,z}^{2}}.$$(84)
The distributions of velocities of H and O are presented in the left hand plot of Figure 17. The velocity distributions are fitted with a MB model as presented in equation (81). The fitted parameter is the thermodynamic temperature of material, which is 294 K for both H and O.
µ_{j} representing the cosine of the angle between the z components V_{j,z} and the velocities V_{j} for oxygen and hydrogen atoms, is calculated by: $${\mu}_{j}=\frac{{V}_{j,z}}{{V}_{j}}.$$(85)
µj is shown in the right hand plot of Figure 17, which is nearly an uniform distribution in [—1;1]
The velocities V_{j} and the corresponding cosine of angle µ_{j} are sampled from the MD trajectory files by using the SVT module. The scattered neutron energy distributions for H and O are compared with the results obtained with the velocities sampled from a MB distribution at 294 K, as presented in Figure 18. Excellent agreement is obtained, which validates the sampling of the velocities from the MD trajectory files implemented in the SVT module.
Compared to a MB distribution, velocities in the MD trajectories are strongly correlated between nearby time steps. These correlations can be illustrated with the velocity autocorrelation function (VACF), as shown in reference [56]. While the correlations between the velocities in the trajectories cannot be directly taken into account by the SVT algorithm. We have tried to incorporate these correlations into the SVT algorithm but our attempts did not succeed.
Fig. 17 Left hand plot presents the distributions of velocities of H and O and their MB fitting curves. The fitted temperatures for H and O are both 294 K. Right hand plot is the distribution of µ_{j}. 
Fig. 18 Comparison of the scattered neutron energy distributions of ^{16}O obtained by sampling the velocities of MD trajectory files and MB distribution with incident neutron energies E = 5 and 6.67 eV for H and O, respectively. 
5.4 Discussions on rejection option to improve the neutron upscattering
The free gas model or SVT treatment for neutronnucleus scattering at thermodynamic temperature T are not able to take into account the vibration dynamics or atom binding effects of materials, but the effective temperature T_{eff} must be involved, as discussed in the literature [60–63]. In this work, special attentions are given to the neutron upscattering: when neutron energy lies below a few eV, the neutronnucleus scattering may lead to an energy gain for the incident neutron. Based on the short collision approximation (SCT), after introducing the effective temperature T_{e}ff to the DDXS of free gas materials (Eq. (83)), the upscattering part of the DDXS must be further corrected to fulfill the principle of detailed balance [62]. The obtained DDXS involving T and T_{eff} is given by: $$\begin{array}{lll}\frac{{\text{d}}^{2}{\sigma}^{\text{SCT}}}{\text{d}\mu \text{\hspace{0.17em}d}E\prime}\hfill & =\hfill & \frac{{\sigma}_{b}}{2{k}_{B}T}\sqrt{\frac{E\prime}{E}}\frac{1}{\sqrt{4\pi \alpha {T}_{\text{eff}}/T}}\hfill \\ \hfill & \hfill & \times \mathrm{exp}\frac{{\left\beta \right\alpha}^{2}}{4\alpha {T}_{\text{eff}}/T}\frac{\beta +\left\beta \right}{2}\hfill \end{array}$$(86)
where T_{eff} is calculated from the PDOS [13,61]: $${T}_{\text{eff}}={w}_{d}T+\frac{\hslash}{2{k}_{B}{w}_{\upsilon}}{\displaystyle {\int}_{0}^{\infty}\omega {\rho}_{\text{vib}}\omega \text{coth}}\frac{\hslash w}{2{k}_{B}T}d\omega .$$(87)
Compared to the DDXS obtained by using the free gas model at T_{eff}, the DDXS based on the SCT can be regarded as multiplying a factor for the upscattering part (when E < E’): $$\frac{{\text{d}}^{2}{\sigma}^{\text{SCT}}}{\text{d}\mu \text{\hspace{0.17em}d}E\prime}=\frac{{\text{d}}^{2}{\sigma}^{\text{FGM}}}{\text{d}\mu \text{\hspace{0.17em}d}E\prime}{T}_{\text{eff}}\times \mathrm{exp}\frac{EE\prime}{{k}_{B}}\frac{1}{T}\frac{1}{{T}_{\text{eff}}}.$$(88)
It is confirmed in Section 5.2 that the SVT algorithm is equivalent to the analytical free gas model for large incident neutron energies. Based on the relation in equation (88), the SCT can be introduced in the SVT module by implementing an additional rejection step in case of neutron upscattering:
when E < E’, a variable R is sampled uniformly in [0;1]; ‘
by comparing ℛ and exp $\frac{E{E}^{\prime}}{{k}_{B}}\frac{1}{T}\frac{1}{{T}_{\text{eff}}}$, if ℛ is larger, the scattered neutron will be rejected.
To validate the rejection step implemented in the SVT module, neutron scattering with ^{1}H in ^{1}H_{2}O and ^{16}O in U^{16}O_{2} are investigated in this work. The first step consists of calculating their effective temperatures from the PDOS of the CAB model [58] and the PDOS reported in reference [59], respectively. The second step is to employ these effective temperatures in the SVT module and use the rejection option for the upscattering part. The final step is to perform MonteCarlo simulations with the Tripoli4^{®} code by using the S( α, β) tables. This step is able to provide accurate scattered neutron energy distributions for comparison. Tested models are summarized in Table 4.
The comparison of results obtained by using the SVT module and the Tripoli4 code for ^{1}H in ^{1}H_{2}O is illustrated in the top plot of Figure 19. At small scattering angle (θ’ = 10°), the SVT at T_{eff} with/without the rejection both fail to reproduce the quasielastic peak and the associated structures. However, the SVT treatment with rejection (green curve) reproduces better the scattered neutron distribution related to the upscattering part. At large θ’, the neutron upscattering can be neglected. Therefore, these two SVT treatments produce the same distributions.
Same comparisons are performed in the case of ^{16}O in U^{16}O_{2}. The obtained results are presented in the bottom plot of Figure 19. At small θ^{’}, both SVT treatments are incapable of reproducing the structures related to the atomic motions of the oxygen atoms bound in UO_{2}. This result is in consistency with the case of ^{1}H in ^{1}H_{2}O. The effective temperature of ^{16}O in U^{16}O_{2} (T_{eff} = 381 K) is close to T = 298 K. Therefore, the upscattering part produced by using the SVT algorithm with rejection does not show significant improvements compared to the SVT treatment at T_{eff}.
These results show interest to add a new rejection option in the Tripoli4^{®} code that allows to correct the upscattering by taking into account T and T_{eff} in the simulation. This new option will be of great interest to test the impact of different scattering models on integral benchmarks.
Simulations performed by using the MonteCarlo neutron transport code Tripoli4^{®} and the SVT module of CINEL for ^{1}H in ^{1}H_{2}O and ^{16}O in UO_{2}. The effective temperature T_{eff} of ^{1}H in ^{1}H_{2}O is computed via the PDOS of the CAB model [58]. T_{eff} of ^{16}O in UO_{2} is obtained via the PDOS of [59].
Fig. 19 Comparison of scattered neutron energy distributions of ^{1}H in ^{1}H_{2}O and ^{16}O in U^{16}O_{2}, obtained by using the SVT module of CINEL at T_{eff} with rejection, and Tripoli4® with S(a, β) for incident neutron energy E = 5 eV and 6.67 eV, respectively. 
6 Conclusions and prospectives
In this work we present a new code CINEL dedicated to generate the thermal scattering files in ENDF6 format for solid crystalline, free gas materials and liquid water. The lowenergy neutron scattering theory is reviewed by presenting the involved approximations and the equations to calculate the elastic (coherent and incoherent) scattering cross section and S(α, β). The implementation of the three modules of CINEL (Cubic, INELastic and SVT) is presented. The Cubic and INELastic modules are validated by comparing the CINEL results with the ENDF/BVIII.0 database and the NJOY+NCrystal library. Examples are given to illustrate the automatic generation of α and β grids, mixed elastic format, anharmonicity and texture in CINEL. The SVT module is validated by comparing the scattered neutron distributions with analytical and MonteCarlo simulations for free target nuclei such as hydrogen and oxygen.
Though CINEL presents improvements compared to the LEAPR module of the NJOY code, the incoherent approximation adopted in the CINEL code is reported to be inadequate for strongly coherent crystalline materials such as pyrolytic graphite. Therefore, the coherent onephonon correction will be implemented in the future version of CINEL to improve the incoherent approximation. The cubic approximation used in the current version of CINEL will also be eliminated. CINEL will also be extended to calculate the thermal scattering files for liquid materials other than water.
Conflict of interests
The authors declare that they have no competing interests to report.
Funding
This work is carried out in the framework of S. Xu’s thesis funded by the CEA Cadarache center.
Data availability statement
This article has no associated data generated and/or analyzed.
Author contribution statement
S. Xu: Conceptualization, methodology, software, validation, writing — review & editing, visualization. G. Noguere: Conceptualization, methodology, software, writing – review & editing, supervision.
Acknowledgements.
The authors would like to express their gratitude to their colleague P. Tamagno from the same laboratory for the GPUs used in this work. The authors also wish to express their appreciation to J.I. Marquez Damian and K. Ramić from European Spallation Source for many helpful discussions. Tripoli4 is a registered trademark of CEA. The authors would like to thank Electricite de France (EDF) for partial financial support.
References
 B.T.M. Willis, Neutron diffraction studies of the actinide oxides i. uranium dioxide and thorium dioxide at room temperature, Proc. R. Soc. Lond. A 274, 122–133 (1963) [CrossRef] [Google Scholar]
 B.T.M. Willis, Neutron diffraction studies of the actinide oxides ii. thermal motions of the atoms in uranium dioxide and thorium dioxide between room temperature and 1100 °C, Proc. R. Soc. Lond. A 274, 134–144 (1963) [CrossRef] [Google Scholar]
 B.O. Loopstra, Neutron diffraction investigation of U_{3}O_{8}, Acta Crystallogr. 17, 651–654 (1964) [CrossRef] [Google Scholar]
 L. Desgranges, G. Baldinozzi, G. Rousseau, J.C. Niepce, G. Calvarin, Neutron diffraction study of the in situ oxidation of UO2, Inorg. Chem. 48, 7585–7592 (2009) [CrossRef] [Google Scholar]
 H.A. Levy, P.A. Agron, M.A. Bredig, M.D. Danford, Xray and neutron diffraction studies of molten alkali halides, Ann. New York Acad. Sci. 79, 762–780 (1960) [Google Scholar]
 M. Atoji, Neutron diffraction studies of CaC_{2}, YC_{2}, LaC_{2}, CeC2, TbC2, YbC2, LuC2, and UC2, J. Chem. Phys. 35, 1950–1960 (1961) [Google Scholar]
 H. Schober, An introduction to the theory of nuclear neutron scattering in condensed matter, J. Neutr. Res. 17, 109–357 (2014) [CrossRef] [Google Scholar]
 G.L. Squires, Introduction to the Theory of Thermal Neutron Scattering, 3rd edn. (Cambridge University Press, 2012) [Google Scholar]
 X.X. Cai, T. Kittelmann, Ncrystal: a library for thermal neutron transport, Comput. Phys. Commun. 246, 106851 (2020) [CrossRef] [Google Scholar]
 L. Van Hove, Correlations in space and time and born approximation scattering in systems of interacting particles, Phys. Rev. 95, 249–262 (1954) [CrossRef] [Google Scholar]
 A. Trkov, D.A. Brown, Endf6 formats manual: Data formats and procedures for the evaluated nuclear data files, Tech. rep., Brookhaven National Lab. (BNL), Upton, NY (United States), 2018 [Google Scholar]
 R.E. MacFarlane, New thermal neutron scattering files for ENDF/BVI release 2. United States, 1994. Available at https://www.osti.gov/biblio/10192168 [Google Scholar]
 R. Macfarlane, D.W. Muir, R. Boicourt, A.C. Kahler III, J.L. Conlin, The enjoy nuclear data processing system, version 2016, Tech. rep., Los Alamos National Lab. (LANL), Los Alamos, NM (United States), 2017 [Google Scholar]
 K. Ramic, J.I. Marquez Damian, T. Kittelmann, D.D. Di Julio, D. Campi, M. Bernasconi, G. Gorini, V. Santoro, NJOY + NCrystal: an opensource tool for creating thermal neutron scattering libraries with mixed elastic support, Nucl. Instr. Methods Phys. Res. A 1027, 166227 (2022) [CrossRef] [Google Scholar]
 Y. Zhu, A.I. Hawari, Full law analysis scattering system hub (flassh), in Proc. PHYSOR (2018), pp. 22–26 [Google Scholar]
 Y.Q. Cheng, L.L. Daemen, A.I. Kolesnikov, A.J. RamirezCuesta, Simulation of inelastic neutron scattering spectra using OCLIMAX, J. Chem. Theory Comput. 15, 1974–1982 (2019) [CrossRef] [Google Scholar]
 Y.Q. Cheng, A.J. RamirezCuesta, Calculation of the thermal neutron scattering crosssection of solids using OCLIMAX, J. Chem. Theory Comput. 16, 5212–5217 (2020) [CrossRef] [Google Scholar]
 M.L. Zerkle, Mixed elastic scattering format proposal (2020). Available at https://indico.bnl.gov/event/7233/contributions/43822/attachments/31592/49906/Mixed_Elastic_Scattering_Format.pdf [Google Scholar]
 A. Sjolander, Multiphonon processes in slow neutron scattering by crystals, Arkiv Fysik 14, 330–333. [Google Scholar]
 H. Sitepu, Marchtype models for the description of texture in granular materials., Ph.D. thesis, Curtin University, 1998 [Google Scholar]
 A. March, Mathematische theorie der regelung nach der korngestah bei affiner deformation, Zeitsch. Kristallogr. 81, 285–297 (1932) [Google Scholar]
 W.A. Dollase, Correction of intensities for preferred orientation in powder diffractometry: application of the March model, J. Appl. Crystallogr. 19, 267–272 (1986) [CrossRef] [Google Scholar]
 A. Altomare, G. Cascarano, C. Giacovazzo, A. Guagliardi, Early finding of preferred orientation: a new method, J. Appl. Crystallogr. 27, 1045–1050 (1994) [CrossRef] [Google Scholar]
 R. Cerny, V. Valvoda, M. Chladek, Empirical texture corrections for asymmetric diffraction and inclined textures, J. Appl. Crystallogr. 28, 247–253 (1995) [CrossRef] [Google Scholar]
 B.T.M. Willis, R.G. Hazell, Reanalysis of singlecrystal neutrondiffraction data on UO2 using third cumulants, Acta Crystallogr. A 36, 582–584 (1980) [CrossRef] [Google Scholar]
 E.M. Tammar et al., Contribution a l'etude du facteur de structure dynamique des liquides simples, Ph.D. thesis, Metz (2000) [Google Scholar]
 P. Egelstaff, Neutron scattering studies of liquid diffusion, Adv. Phys. 11, 203–232 (1962) [CrossRef] [Google Scholar]
 R. Turner, The quasiclassical approximation for neutron scattering, Physica 27, 260–264 (1961) [CrossRef] [Google Scholar]
 G.H. Vineyard, Scattering of slow neutrons by a liquid, Phys. Rev. 110, 999–1010 (1958) [CrossRef] [Google Scholar]
 A. Rahman, K.S. Singwi, A. Sjolander, Theory of slow neutron scattering by liquids. i, Phys. Rev. 126, 986–996 (1962) [CrossRef] [Google Scholar]
 S.T. Lin, M. Blanco, W.A. Goddard, The twophase model for calculating thermodynamic properties of liquids from molecular dynamics: Validation for the phase diagram of lennardjones fluids, J. Chem. Phys. 119, 11792–11805 (2003) [CrossRef] [Google Scholar]
 T.A. Pascal, S.T. Lin, W.A. Goddard III, Thermodynamics of liquids: standard molar entropies and heat capacities of common solvents from 2PT molecular dynamics, Phys. Chem. Chem. Phys. 13, 169–181 (2011) [CrossRef] [Google Scholar]
 G. Noguere, J.P. Scotta, S. Xu, E. Farhi, J. Ollivier, Y. Calzavarra, S. Rols, M. Koza, J.I. Marquez Damian, Temperaturedependent dynamic structure factors for liquid water inferred from inelastic neutron scattering measurements, J. Chem. Phys. 155, 024502 (2021) [CrossRef] [Google Scholar]
 P.A. Egelstaff, P. Schofield, On the evaluation of the thermal neutron scattering law, Nucl. Sci. Eng. 12, 260–270 (1962) [CrossRef] [Google Scholar]
 J. Marquez Damian, J. Granada, F. Cantargi, J. Dawidowski, Generation of thermal scattering libraries for liquids beyond the gaussian approximation using molecular dynamics and njoy/leapr, Ann. Nucl. Energy 92, 107–112 (2016) [CrossRef] [Google Scholar]
 K.S. Singwi, A. Sjolander, Diffusive motions in water and cold neutron scattering, Phys. Rev. 119, 863–871 (1960) [CrossRef] [Google Scholar]
 V.F. Sears, Theory of cold neutron scattering by homonuclear diatomic liquids: II. Hindered rotation, Can. J. Phys. 44, 1299–1311 (1966) [CrossRef] [Google Scholar]
 M. Mattes, J. Keinert, Status of thermal neutron scattering data for graphite, Tech. rep., International Atomic Energy Agency (2005) [Google Scholar]
 B. Granger, J. Grout, Jupyterlab: Building blocks for interactive computing, Slides of presentation made at SciPy. [Google Scholar]
 T. Kluyver, B. RaganKelley, F. Pérez, B. Granger, M. Bussonnier, J. Frederic, K. Kelley, J. Hamrick, J. Grout, S. Corlay, P. Ivanov, D. Avila, S. Abdalla, C. Willing, J. development team, Jupyter notebooks—a publishing format for reproducible computational workflows, in Positioning and Power in Academic Publishing: Players, Agents and Agendas, edited by F. Loizides, B. Scmidt (IOS Press, 2016), pp. 87–90 [Google Scholar]
 S.K. Lam, A. Pitrou, S. Seibert, Numba: a llvmbased python jit compiler, in Proceedings of the Second Workshop on the LLVM Compiler Infrastructure in HPC (2015), pp. 1–6 [Google Scholar]
 I.I. AlQasir, Thermal neutron scattering in graphite, Ph.D. thesis, North Carolina State University (2008) [Google Scholar]
 R. Coveyou, R. Bate, R. Osborn, Effect of moderator temperature upon neutron flux in infinite, capturing medium, J. Nucl. Energy 2, 153–167 (1956) [Google Scholar]
 I. Lux, Monte Carlo Particle Transport Methods: Neutron and Photon Calculations (2018). Available at 10.1201/9781351074834. [Google Scholar]
 W. Rothenstein, Neutron scattering kernels in pronounced resonances for stochastic doppler effect calculations, Ann. Nucl. Energy 23, 441–458 (1996) [CrossRef] [Google Scholar]
 S.P. Ong, W.D. Richards, A. Jain, G. Hautier, M. Kocher, S. Cholia, D. Gunter, V.L. Chevrier, K.A. Persson, G. Ceder, Python Materials Genomics (pymatgen): a robust, opensource python library for materials analysis, Comput. Mater. Sci. 68, 314–319 (2013) [Google Scholar]
 S.P. Ong, W.D. Richards, A. Jain, G. Hautier, M. Kocher, S. Cholia, D. Gunter, V.L. Chevrier, K.A. Persson, G. Ceder, UO_{2} unit cell with fm3m symmetry. Available at https://materialsproject.org/materials/mp1597/ (asscessed 20210407) [Google Scholar]
 J.L. Wormald, N.C. Fleming, A.I. Hawari, M.L. Zerkle, Generation of the thermal scattering law of uranium dioxide with ab initio lattice dynamics to capture crystal binding effects on neutron interactions, Nucl. Sci. Eng. 195, 227–238 (2021) [CrossRef] [Google Scholar]
 S. Baroni, S. de Gironcoli, A. Dal Corso, P. Giannozzi, Phonons and related crystal properties from densityfunctional perturbation theory, Rev. Mod. Phys. 73, 515–562 (2001) [CrossRef] [Google Scholar]
 P. Giannozzi, S. Baroni, N. Bonini, M. Calandra, R. Car, C. Cavazzoni, D. Ceresoli, G.L. Chiarotti, M. Cococcioni, I. Dabo, A.D. Corso, S. de Gironcoli, S. Fabris, G. Fratesi, R. Gebauer, U. Gerstmann, C. Gougoussis, A. Kokalj, M. Lazzeri, L. MartinSamos, N. Marzari, F. Mauri, R. Mazzarello, S. Paolini, A. Pasquarello, L. Paulatto, C. Sbraccia, S. Scandolo, G. Sclauzero, A.P. Seitsonen, A. Smogunov, P. Umari, R.M. Wentzcovitch, QUANTUM ESPRESSO: a modular and opensource software project for quantum simulations of materials, J. Phys.: Condens. Matter 21, 395502 (2009) [CrossRef] [PubMed] [Google Scholar]
 J.R. Granada, J.I. Marquez Damian, J. Dawidowski, J.I. Robledo, C. Helman, G. Romanelli, G. Skoro, Development of neutron scattering kernels for cold neutron reflector materials. Available at http://arxiv.org/abs/2103.09145 [Google Scholar]
 S. Xu, Measurement of the double differential neutron cross section of UO_{2} and determination of the thermal scattering law as a function of temperature, Ph.D. thesis, AixMarseille Universitae (2021) [Google Scholar]
 S. Xu, G. Noguere, L. Desgranges, J.M. Damian, Atomic scale montecarlo simulations of neutron diffraction experiments on stoichiometric uranium dioxide up to 1664 k, Nucl. Instr. Methods Phys. Res. Sect. A 1002, 165251 (2021) [CrossRef] [Google Scholar]
 E. Brun, F. Damian, C. Diop, E. Dumonteil, F. Hugot, C. Jouanne, Y. Lee, F. Malvagi, A. Mazzolo, O. Petit, J. Trama, T. Visonneau, A. Zoia, Tripoli4, cea, edf and areva reference monte carlo code, Ann. Nucl. Energy 82, 151–160 (2015) [Google Scholar]
 M.J. Abraham, T. Murtola, R. Schulz, S. Paall, J.C. Smith, B. Hess, E. Lindahl, Gromacs: high performance molecular simulations through multilevel parallelism from laptops to supercomputers, SoftwareX 12, 19–25 (2015) [CrossRef] [Google Scholar]
 J. Scotta, Improvement of the thermal and epithermal neutron scattering data for the interpretation of integral experiments, Ph.D. thesis, AixMarseille Universitae (2017). [Google Scholar]
 M.A. Gonzaalez, J.L.F. Abascal, A flexible model for water based on TIP4P/2005, J. Chem. Phys. 135, 224516 (2011) [CrossRef] [Google Scholar]
 J. Maarquez Damiaan, J. Granada, D. Malaspina, Cab models for water: A new evaluation of the thermal neutron scattering laws for light and heavy water in endf6 format, Ann. Nucl. Energy 65, 280–289 (2014) [CrossRef] [Google Scholar]
 G. Noguere, J.P. Scotta, S. Xu, A. Filhol, J. Ollivier, E. Farhi, Y. Calzavarra, S. Rols, B. Fak, J.M. Zanotti, Q. Berrod, Combining density functional theory and monte carlo neutron transport calculations to study the phonon density of states of Uo_{2} up to 1675 k by inelastic neutron scattering, Phys. Rev. B 102, 134312 (2020) [CrossRef] [Google Scholar]
 W.E. Lamb, Capture of neutrons by atoms in a crystal, Phys. Rev. 55, 190–197 (1939) [CrossRef] [Google Scholar]
 G. Borgonovi, Neutron scattering kernels calculations at epithermal energies, Tech. rep., Gulf General Atomic, Inc., San Diego, Calif. (1970) [Google Scholar]
 T. Sutton, T. Trumbull, C. Lubitz, Comparison of some monte carlo models for bound hydrogen scattering, in: International Conference on Mathematics, Computational Methods & Reactor Physics, Saratoga Springs, New York USA (2009) [Google Scholar]
 J. Dawidowski, L. Rodriguez Palomino, J. Marquez Damian, J. Blostein, G. Cuello, Effective temperatures and scattering cross sections in water mixtures determined by deep inelastic neutron scattering, Ann. Nucl. Energy 90, 247–255 (2016) [CrossRef] [Google Scholar]
Cite this article as: Shuqi Xu and Gilles Noguere. Generation of thermal scattering files with the CINEL code, EPJ Nuclear Sci. Technol. 8, 8 (2022)
All Tables
Comparison of the computational time of the TSLs of H in H_{2}O at room temperature with a phonon expansion order N_{phonon} = 2000.
Simulations performed by using the MonteCarlo neutron transport code Tripoli4^{®} and the SVT module of CINEL for ^{1}H in ^{1}H_{2}O and ^{16}O in UO_{2}. The effective temperature T_{eff} of ^{1}H in ^{1}H_{2}O is computed via the PDOS of the CAB model [58]. T_{eff} of ^{16}O in UO_{2} is obtained via the PDOS of [59].
All Figures
Fig. 1 JupyterLab interface. Live code, texts, mathematical equations and interactive graphics are mixed in Jupyter Notebooks which are integrated in JupyterLab together with blocks like terminal and text editor. 

In the text 
Fig. 2 Data flow presenting the generation of the thermal scattering laws in ENDF6 format and scattered neutron energy distributions by using CINEL. 

In the text 
Fig. 3 Screenshot of cross section file generated by the Cubic module of CINEL for uranium dioxide at 296 K. Columns 1 to 11 represent respectively: index h, k, l, interplanar distance d_{hkl} (in Å), diffraction angle θ_{hkl} (in degree), orientation angle Θ_{hkl} (in degree), PoleDensity Distribution Function value ${\mathcal{P}}_{hkl}{\text{\Theta}}_{hkl}$, form factor square F(τ_{hkl})^{2} (in barn), multiplicity M_{hkl}, Bragg edge E_{hki} (in eV) and ${E}_{hkl}{\sigma}_{\text{coh}}^{\text{el}}{E}_{hkl}$ (in eV.barn), which serve as calculations of the coherent elastic cross sections. 

In the text 
Fig. 4 Schematic representation of neutronnucleus scattering in LAB (top plot), TR (bottom left plot) and CM (bottom right plot) frames. Detail descriptions can be found in the text. 

In the text 
Fig. 5 Left hand plot presents the cubic unit cell of UO_{2}: a = b = c and α = β = γ = 90°. Blue and red spheres represent respectively uranium and oxygen atoms. Right hand plot presents the tetragonal unit cell of MgH_{2} or MgD_{2}: b = c and α = β = γ = 90°. Orange spheres represent magnesium atoms. Pink spheres represent hydrogen or deuterium atoms. The drawings of these two figures are with the aid of the Crystal Toolkit in the Materials Project [46]. 

In the text 
Fig. 6 PDOS of ^{238}U in UO2 and ^{16}O in UO2, obtained from ab initio lattice dynamic methods [48]. 

In the text 
Fig. 7 Comparison of the total scattering cross sections of UO2 for temperatures from 296 K to 1200 K. The resonances of uranium are not included in the calculated cross sections. 

In the text 
Fig. 8 PDOS of ^{24}Mg in ΜgΗ_{2}, ^{1}H in ΜgΗ_{2}, ^{24}Mg in ΜgD_{2} and D in MgD_{2}, generated by using the density functional perturbation theory [49,50]. 

In the text 
Fig. 9 Total scattering cross section (solid line) and the coherent/incoherent elastic/inelastic components (dot lines) of MgD_{2} at 293.6 K. 

In the text 
Fig. 10 Comparison of the total scattering cross sections of MgH_{2} and MgD_{2} for temperatures from 20 K to 600 K. 

In the text 
Fig. 11 Comparison of coherent elastic scattering cross sections of UO_{2} with different thirdcumulant coefficient values: ${c}_{123}^{\text{O}}=0,{10}^{3},5\times {10}^{3}$. 

In the text 
Fig. 12 Comparison of cylindrically symmetric PDDF proposed by MarchDollase with parameter P_{1}= 0.6 and Altomare et al. with P_{1} = 1.5. 

In the text 
Fig. 13 Comparison of coherent elastic scattering cross section ${\sigma}_{\text{coh}}^{\text{el}}E$ of zirconium with random orientation and with preferred orientation $\overrightarrow{\text{PO}}=\left[0,0,1\right]$. ${\sigma}_{\text{coh}}^{\text{el}}E$ obtained by using the cylindrically symmetric PDDFs proposed by MarchDollase and Altomare et al. are respectively presented. 

In the text 
Fig. 14 Comparison of the scattered neutron energy distributions of ^{1}H and ^{16}O (free gas materials) obtained by using the SVT module of CINEL and Tripoli4 for incident neutron energies E = 5 eV and 6.67 eV, respectively. 

In the text 
Fig. 15 Comparison of the scattered neutron energy distributions of free ^{1}H obtained by using the SVT module of CINEL and analytical results of the free gas model for incident neutron energies E = 1 and 5 eV. 

In the text 
Fig. 16 Comparison of the scattered neutron energy distributions of free ^{16}O obtained by using the SVT module of CINEL and analytical results of the free gas model for incident neutron energies E = 0.1 and 1 eV. 

In the text 
Fig. 17 Left hand plot presents the distributions of velocities of H and O and their MB fitting curves. The fitted temperatures for H and O are both 294 K. Right hand plot is the distribution of µ_{j}. 

In the text 
Fig. 18 Comparison of the scattered neutron energy distributions of ^{16}O obtained by sampling the velocities of MD trajectory files and MB distribution with incident neutron energies E = 5 and 6.67 eV for H and O, respectively. 

In the text 
Fig. 19 Comparison of scattered neutron energy distributions of ^{1}H in ^{1}H_{2}O and ^{16}O in U^{16}O_{2}, obtained by using the SVT module of CINEL at T_{eff} with rejection, and Tripoli4® with S(a, β) for incident neutron energy E = 5 eV and 6.67 eV, respectively. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.