Advanced Modeling and Simulation Methods for Evaluation of Thermal Neutron Scattering Materials ∗

.


Introduction
Several recent efforts have been made to validate thermal neutron scattering libraries (e.g., [1]). Although this work is certainly important, efforts must be made to ensure accurate modeling of double differential scattering cross sections (DDXSs) so they can be used in the evaluation procedure. Previous validation efforts using DDXS data [2] relied on processing Evaluated Nuclear Data File (ENDF) files with NJOY [3] and then simulating the experiment using MCNP [4]. This makes for a sound first-order comparison, but the ultimate goal is to remove as many approx- * Notice: This manuscript has been authored by UT-Battelle LLC, under contract DE-AC05-00OR22725 with the US Department of Energy (DOE). The US government retains and the publisher, by accepting the article for publication, acknowledges that the US government retains a nonexclusive, paid-up, irrevocable, worldwide license to publish or reproduce the published form of this manuscript, or allow others to do so, for US government purposes. DOE will provide public access to these results of federally sponsored research in accordance with the DOE Public Access Plan (http://energy.gov/downloads/doe-public-access-plan).
* * e-mail: chapmancw@ornl.gov imations in validation as possible. The LEAPR module of NJOY uses several approximations to generate the thermal scattering file (e.g., incoherent approximation, cubic approximation, atom-site approximation), and the MCNP input is only a simplified model of the beamline.
Several code packages created by the instrument scientists at the Spallation Neutron Source (SNS) have been developed to remove these approximations. Two codes used in this analysis are OCLIMAX [5], which provides more accurate processing of thermal neutron scattering properties, and MCViNE [6], which provides more accurate simulations of neutron scattering experiments. These codes are discussed in detail in Section 2; comparisons of these codes against both the previous code packages and experimental data are given in Section 3, and some concluding thoughts and recommendations for future work are provided in Section 4.

SNS codes
OCLIMAX is used to calculate neutron scattering properties, similar to the LEAPR module of NJOY or the FLASSH code [7]. OCLIMAX calculates the dynamic structure factor (DSF) from phonon frequencies and polarization vectors, which are calculated using density functional theory (DFT) codes. OCLIMAX is also capable of reading in atomic trajectories from molecular dynamics (MD) simulations to calculate the phonon density of states (PDOS), which is then used to calculate the DSF in the incoherent approximation. Options are available to calculate both the elastic and inelastic options, as well as coherent and incoherent scattering. Output formats include the full DSF, the PDOS (both total and subdivided by atom type), the structure factor S(Q), and the mean square displacement.
MCViNE is a Monte Carlo neutron ray-tracing program (similar to McStas [8]) used to perform computational simulations of neutron scattering experiments. A beam profile is created from a user-provided incident energy and from chopper settings of the experiment. Then the sample environment is modeled using the geometry of the sample in the experiment, as well as a full DSF file, or a PDOS if the incoherent approximation is allowed. The data are then stored in a binary histogram format for easy on-the-fly sampling. Then the experiment is simulated, and a resulting experimental DSF is calculated.

Comparison to Previous Methodology
Using these SNS codes provide two major benefits when compared to the previous codes, NJOY and MCNP. They are both freely available (OCLIMAX: https://sites.google.com/site/ornliceman/download; MCViNE: https://mcvine.ornl.gov/), and both use fundamental physics and use fewer approximations, leading to more accurate results. OCLIMAX calculates the scattering properties from first principles, and MCViNE has high-fidelity models of the beamline geometry, beam profiles, and instrument resolutions. The downsides to these codes are that, compared to NJOY and MCNP, they are computationally slow, and a uniform (Q,E) grid is required for input into MCViNE. This last point is particularly troublesome, as most ENDF thermal scattering files are not uniformly spaced in (Q,E) space, so these files must be regenerated to be used as MCViNE inputs.
Comparatively, NJOY and MCNP are both computationally fast, and they have much smaller data footprints compared to OCLIMAX and MCViNE. However, as mentioned in the introduction, NJOY has several built-in approximations. Additionally, MCNP is not freely available, and extensive modeling would be required to accurately recreate the SNS beamline geometry and instrument resolutions.

Workflow
To see the improvements between the two methodologies, a series of DFT simulations were carried out using Cam-bridge Serial Total Energy Package (CASTEP) [9] on ultra high molecular weight polyethylene (UHMWPE) and yttrium hydride. To compare with the old methodologies, these phonon calculations were run in OCLIMAX to generate a 1D PDOS, which is a necessary input to run NJOY. The PDOS was then processed using NJOY, and the resulting ACE file was used to simulate the experiment using the simplified beamline mock-up in MCNP. The incident beam profile in the MCNP input was the same profile simulated using MCViNE. An approximate resolution function (generated from https://rez.mcvine.ornl.gov/) was used to post-process the MCNP results.
To compare the new code packages, the CASTEP phonon information was read directly into OCLIMAX, which generated the full 2D S(Q,E). This was then used in the MCViNE simulation of the beamline. These two simulation methodologies were compared against recently measured data of UHMWPE at the Wide Angular-Range Chopper (ARCS) beamline and yttrium hydride at the SE-QUOIA beamline. Specifically, the measured spectra were divided into slices of DDXS over several angles for easy comparison to MCNP results. This is because the MCNP input only tallied angles of the scattered neutrons instead of the full detector array. Both samples were measured at 5 K at a range of incident energies, although it was necessary to generate the ACE files at 20 K because of an issue with how NJOY handles very low-temperature processing. This approach is not expected to show any significant discrepancies.

Results
When these simulations were run, there was an issue with how MCViNE handled the OCLIMAX-generated DSF file. This meant the MCViNE simulations had to be run without any multiple neutron scattering. The issue stems from how MCViNE bins the S(Q,E) data near the elastic peak. While this issue has since been resolved, the analysis here shows only single-scattering simulations. In both sets of experiments, the experimental DDXSs and simulated results are normalized so that the maximum value is 1.

Polyethylene
The comparison of measurements of UHMWPE at ARCS against the MCNP and MCViNE simulations are shown in Fig. 1.
As shown in these plots, the impact of using singleneutron scattering in the MCViNE simulations on the lower incident energy spectra (80, 175 meV) is clear. While the shape of the elastic peak in general is closely aligned with the experimental data, the inelastic spectra is underpredicted. Conversely, at higher energies the inelastic spectra at exit energies between 200-300 meV agree with the data, while the up-scattering quasi-elastic spectra between 400-450 meV is overpredicted.
It is also clear that the elastic peak for the MCNP simulations does not align up with either the experimental data or the MCViNE simulations. This can most easily be seen in the 80 and 250 meV plots, where the peak of the MCNP simulations is offset from the data, and the shape of the quasi-elastic peak does not align with the data. This is most likely the result of the approximate instrument resolution that was applied.

Yttrium hydride
The comparison of measurements of yttrium hydride (specifically YH 1.86 ) at SEQUOIA against the MCNP and MCViNE simulations are shown in Fig. 2. Here, the effects of running single-neutron scattering are not quite as pronounced as the UHMWPE plots, although the effects still exist at small scattering angles. This is largely because the YH 1.86 sample was much thinner than the UHMWPE sample, so the probability of multiple neutron scattering is significantly reduced. Here, the impact of the imperfect resolution function is not quite as pronounced until 600 meV, but curiously, for the first time, the MCViNE shows a shift in the inelastic spectra that does not align with the measured data.
In both the UHMWPE and yttrium hydride plots, the MCViNE spectra appear to be more jagged, particularly in the upscattering regime. This is most likely a result of how the data were spliced: the (Q,E) grid calculated by MCViNE does not necessarily align directly along a given scattering angle. As such, the (Q,E) grid was transformed to an (angle, E f ) grid for each incident energy for cases in which the closest angle to the MCNP simulated angle was used.

Conclusion
When results from the two methodologies are compared, it is clear that while OCLIMAX and MCViNE are better able to accurately model the elastic peak shape, NJOY and MCNP have a more accurate spectra between the phonon peak modes. Also, while NJOY and OCLIMAX took about the same amount of time to process the data, MCNP was able to simulate the experiment significantly faster than MCViNE.
Because of the way that MCViNE handles the input S(Q,E) data, it may not be possible to allow for a nonuniform (Q,E) grid, but there may be possible workarounds, ranging from how the data are handled (interpolation of nonuniform grid onto a uniform grid) to how the code handles the histogram binning (allowing for nonuniform data). Finally, the issue concerning using multiple neutron scattering has been resolved and will be investigated in a future analysis.