Validation of Thermal Neutron Scattering Cross Sections for Heavy Water based on Molecular Dynamics Simulation

. Recently, the thermal scattering libraries of ENDF/B-VIII.0 and JEFF-3.3 for light and heavy water were released with a new water model (CAB model) proposed by Damian. For the CAB model, the molecular dynamics code GROMACS was used to more accurately describe the realistic motions of water molecules. In this paper, to consider the coherent component we also generated the thermal scattering cross section of the deuterium and oxygen bound in the heavy water molecules using the GROMACS code and EPSR code. In addition, the frequency spectrum was also calculated using the GROMACS code. Thermal scattering cross sections based on the newly calculated Sköld correction factor and the frequency spectrum were generated by NJOY2016 code. Finally, the performance of the generated thermal scattering cross sections were validated by performing an ICSBEP benchmark simulation using MCNPX code.


Introduction
Light and heavy water are two of the most important materials used as a moderator in thermal-neutron reactors. In the thermal energy region (~5 eV), to perform neutronics simulations based on nuclear data, it is important to accurately describe changes in the reaction cross section and the energy-angular distribution of secondary neutrons, using a thermal scattering library. Recently, the newly evaluated nuclear data files ENDF/B-VIII.0 and JEFF-3.3 have been released, which includes measurable improvements with respect to the thermal scattering libraries of light and heavy water. A new water model of the CAB model, which was developed at the Neutron Physics Department at Centro Atómico Bariloche (CAB model) [1] has been applied in the thermal scattering libraries of ENDF/B-VIII.0 and JEFF-3.3. The significant difference between the CAB model and existing models is that the CAB model was produced using a molecular dynamics (MD) simulation with GROMACS code [2] and NJOY code [3] which are comprehensive computer code packages for producing nuclear cross sections in the ENDF format. In particular, in the CAB model for heavy water, the Sköld correction factor [4] and frequency spectrum of the deuterium and oxygen bound in the heavy water molecules were calculated in the MD simulation using GROMACS code, taking the coherency and vibration of the molecules of heavy water into consideration.
In this paper, we also generated the Sköld correction factors and frequency spectra of deuterium and oxygen to generate a thermal scattering library for heavy water using GROMACS v.5.1.4. and EPSR (Empirical Potential Structure Refinement) code [5] and compared the results against an ICSBEP [6] benchmark simulation using the MCNPX 2.7.0 code [7].

Thermal scattering law for heavy water
In the ENDF/B notation, the thermal incoherent scattering cross section is given by Eq. (1).
where E and   are the incident and secondary neutron energies,   is the characteristic bound cross section, k is the Boltzmann constant and T is the Kelvin temperature of the material. Also, the thermal scattering law depends on parameters of  and . Here  is the momentum transfer parameter defined as and  is the energy transfer parameter defined as where A is the mass ratio of the scattering nuclide to the neutron and cos  is the scattering angle in the laboratory system.

Sk̈ld approximation for heavy water
In comparison with existing libraries, the thermal scattering libraries of ENDF/B-VIII.0 and JEFF-3.3 for heavy water for the CAB model have two significant improvements. Firstly, it has considered the coherent scattering cross section of heavy water by introducing the Sköld approximation with Sköld correction factors. Secondly, the thermal scattering cross section dealing with the oxygen atom in the heavy water molecule (i.e. the bound O in the D 2 O library) is newly generated to take into account the interference effects of oxygen atoms.
To generate the scattering cross section which is considered the coherent component of the deuterium and oxygen, the thermal scattering cross section for heavy water uses the Sköld approximation: where    and    are the Sköld correction factors for deuterium and oxygen, respectively.
The   (),   () and   () are the static structure factors for each atom in the heavy water molecules.

GROMACS and EPSP code simulations
To apply the Sköld approximation, the CAB model introduced the Sköld correction factors by calculating the static structure factors using the MD simulation with the GROMACS code. In this paper, we also calculated the Sköld correction factors using not only GROMACS but also EPSR code. The GROMACS code is a versatile package for performing the MD simulation, which simulates the motions of a system composed of hundreds to millions of particles using Newtonian equations. In sequence, EPSR uses a Monte Carlo method which evolved from the Reverse Monte Carlo, which was developed to build a structural model of disordered materials such as a glass or liquid.
As a result, Fig.1 shows the Sköld correction factors. When the Sköld correction factor from the GROMACS and EPSR codes were compared with the ENDF/B-VIII.0 data, both Sköld correction factors for deuterium showed high consistency. In case of the oxygen, however, the Sköld correction factors from the GROMACS and EPSR code showed lower peaks than the ENDF/B-VIII.0 data.

Frequency spectrum
The thermal scattering law which contains dynamic and structural information about the target system and determines the energy and angular distribution of the secondary neutrons can be written by Eq. (9).
where ̂ is time measured in the unit of ℏ  ⁄ seconds. The intermediate scattering function γ() is given by Eq. (10).
where () is the frequency spectrum. In the case of heavy water, it consists of a diffusion spectrum, a continuous spectrum and discrete oscillators, which is used as an input of the LEAPR module in the NJOY code for generating the thermal scattering cross section. Figure 2 shows the difference in the continuous spectra of the deuterium and oxygen in heavy water. The solid line represents the continuous spectrum used in the thermal scattering cross section in ENDF/B-VIII.0, and the dotted line is the spectrum calculated by GROMACS in this paper. The frequency spectrum was calculated by using just the GROMACS code because the EPSR code doesn't have a function for calculating the frequency spectrum.
As shown in Fig. 2, although both spectra are based on the GROMACS code, they show some differences. The differences might be from differences in the simulation conditions, such as simulation time, applied ensemble, conditions of equilibration and so on. Nevertheless, the frequency spectrums in this work have a shape analogous to frequency spectrum used in the thermal scattering library in ENDF/B-VIII.0. In particular, the peaks of both spectra show substantial accordance. The peaks at around 0.6 meV and 2~3 meV represent the motions of intermolecular bending and stretching of the deuterium and oxygen, respectively. In contrast with the oxygen, the deuterium has a third peak at around 4~6 meV, which indicates the motion of the librations.  Finally, we generated the thermal scattering cross sections for heavy water using the NJOY2016 code with the new Sköld correction factors and the frequency spectrum calculated by MD simulations. As shown in Fig.3, three kinds of thermal scattering cross sections are compared with the ENDF/B-VIII.0 data, which is used as a reference. The red lined scattering cross sections (FS+ENDF/B-VIII.0) are generated by using the same Sköld correction factor and the newly calculated frequency spectrum (FS) in this work. The green lined thermal scattering cross sections (FS+EPSR) were generated using the Sköld correction factor newly calculated with the EPSR code and the calculated frequency spectrum. Lastly, the blue lined scattering cross sections (FS+GROMACS) were generated using the Sköld correction factor newly calculated using the GROMACS code and the calculated frequency spectrum.
As a result, when comparing the reference and FS+ENDF/B-VIII.0, both scattering cross sections agreed well above the 1 meV range. At energy ranges lower than 1 meV, the FS+ENDF/B-VIII.0 scattering cross sections of the deuterium and oxygen had a slightly higher value than the reference. However, the differences were not outstanding.
In case of the scattering cross sections of the FS+EPSR and FS+GROMACS, the differences from the reference are more remarkable, which means the scattering cross sections were more significantly affected by the Sköld correction factors. Also, Fig. 3 shows the tendency that the scattering cross sections of the deuterium generated using the GROMACS and EPSR Sköld correction factors have larger values than the reference below 3meV.
However, for the scattering cross section for oxygen from 2 meV to 9 meV, values lower than the reference were observed. Also, all scattering cross sections had a common feature, in that two dips were observed around 1 meV and from 10 to 30 meV. The first dip comes from the coherence of each atom and the second dip is caused by the effect of oxygen interference in the heavy water molecules. In contrast with the second dip, each scattering cross section showed considerable differences in the first dip range. However, because the neutron flux in thermal systems is much higher in the energy range of the second dip, the description of the second dip is more important.

Criticality benchmark problems
To estimate the effects of the thermal scattering libraries, 59 heavy water moderated/reflected experiments were taken from the International Handbook of Evaluated Criticality Safety Benchmark Problems (ICSBEP). All of benchmark calculations were performed using Monte Carlo Transport code MCNPX 2.7.0. The ENDF/B-VII.1-based KNE71 library [8] was used for all nuclides except for the thermal scattering cross sections of D in D 2 O and/or O in D 2 O. All the MCNP benchmark simulations were carried out at 293.6K. Also, the MCNP runs were terminated after statistical uncertainty was reduced to below 20 pcm. Fig. 4 shows the benchmark uncertainties and differences between the k eff calculated from the benchmark k eff with the thermal scattering libraries of ENDF/B-VIII.0, and generated thermal scattering libraries in the categories of low enriched uranium (LEU) and highly enriched uranium (HEU). The performance of the generated thermal scattering cross sections are similar to the thermal scattering cross section of ENDF/B-VIII.0.
As mentioned in Section 3.1, because the thermal neutron flux spectrum is biased to the energy range of the second dip, the benchmark results were comparable although the generated scattering cross sections showed relatively large differences in the low energy range. Also, the results of the root mean square (RMS) errors for heavy water problems are described in Table 1. The total RMS errors of ENDF/B-VIII.0 and the generated libraries were 0.801~0.808%, which also shows the generated libraries have performance similar to the ENDF/B-VIII.0 library in aggregate.

Conclusions
In this work, we generated thermal scattering libraries for heavy water using the MD simulation and NJOY code. To take into account the coherent scattering cross section which is included in the scattering cross section of heavy water, we calculated the Sköld correction factors of the deuterium and oxygen bound in heavy water molecules, using the MD simulation with the GROMACS and EPSR code. In sequence, the frequency spectra of each atom were also calculated using GROMACS code to consider the intermolecular vibration of heavy water molecules. To generate the thermal scattering cross section for heavy water, finally, the Sköld correction factor and frequency spectrum were used as inputs in the NJOY code. In comparison with the ENDF/B-VIII.0 data, the generated thermal cross sections showed some discrepancies at energy ranges lower than 3 meV. Nevertheless, we confirmed that the performance of the generated libraries was similar to the thermal scattering libraries of ENDF/B-VIII.0 for heavy water by comparison with the results of the criticality benchmark simulation of ICSBEP.