Effects of Hydrogen Bonding on Nuclear Data Development of Liquid Anhydrous HF

,


Introduction
Hydrogen fluoride (HF) is a material widely applied in nuclear fuel processing and manufacturing. As HF is commonly used directly with uranium, an accurate representation of the thermal scattering law (i.e. TSL, ( , )) and cross sections are vital to criticality safety. Due to its aggressive chemical character, the structure and dynamics of liquid HF are of scientific interest for experimental research. As a consequence, many experiments have been conducted to investigate the structure and dynamics of liquid HF [1][2][3][4]. Because the fluorine atom in the liquid HF molecule has a larger electronegativity than hydrogen, a negative massless charge is positioned near the fluorine atom along the covalent bond of hydrogen and fluorine. Therefore, this molecule has a significant tendency to form long chains of hydrogen bonds in the liquid phase [2][3][5][6][7]. Due to the difficulty of the simulations, the existing polarizable HF models have not been used much in recent years [8][9][10]. Therefore, the non-polarizable HF model is the most straightforward representation of liquid HF's characteristics.
A 2-site model has been investigated by the LAMMPS classical molecular dynamics (CMD) simulation package [11] for nuclear data development [12]. Cournoyer and Jorgensen (TIPS) [5], Jedlovszky and Vallauri (JV-NP)] [9], and Orabi and Faraldo-Gómez (OFG) [13] have presented 3-site nonpolarizable models for liquid HF. Kreitmeir et al. [4] adapted the JV-NP model for liquid HF in 2005, * Corresponding author: tahmed4@ncsu.edu resulting in the HF-Kr model. Mazack and Gao have performed a comparative investigation of the polarizable and non-polarizable models [10]. The 2-site flexible model [12] adopted in LAMMPS, is unable to capture the H-bond. Therefore, a 3-site flexible model was developed (referred to as NCSU HF) based on the OFG model [13] in GROMACS [14] to accurately capture the H-bond effects.

Potential Development
The total energy ( ( )) of liquid HF is the sum of intermolecular and intramolecular interactions. In GROMACS, the NCSU HF model was developed using a newly reparametrized potential with the addition of bond stretching, LJ, and Coulombic interaction among the hydrogen (H), massless charge (M), and fluorine atoms (F). The total energy is described as where and are the instantaneous and equilibrium H-F distances respectively, and is the optimized bond force constant. The others parts of equation (1) describe nonbonded interactions of indicated atoms. The geometric mean resultant combining rules are applied at = √ and = √ for the interactions. We used the experimental and simulation results [2][3][15][16][17] as a benchmark to set the parameters of the new NCSU HF model. The potential parameters are given in Table  1 below. The equilibrium F-H distance is set to 0.917 Å, which is close to the simulation and experimental value (0.93±0.02 Å [2][3], and 0.9169 Å [18]) in the liquid phase. The fractional charges are placed at all atoms in a molecule by setting the permanent dipole ( ) to 2.149 D [13] which is larger than the gas phase permanent dipole 1.83 D [19]. When the bond length = 0.917 Å is used as an input parameter, the average bond length over time for a simulation of liquid HF is 0.92 Å which is in good agreement with the simulation and experimental values (0.93±0.02 Å [2][3], and 0.9169 Å [18]). So, the current parameterization gives the best results. The 2-site HF model is described in the ref. [12].

Computational Protocol
To capture the H-bond, the NCSU HF model was developed for high temperature and pressure and implemented in the GROMACS [14]. For this model, the initial structure of 1000 molecules was produced using fftool [20] and PACMOL [21]. The Nose-Hoover thermostat and the Parrinello-Rahman barostat were applied. At 0.90 nm, the LJ interactions were terminated, and the long-range electrostatic forces were calculated using the particle mesh Ewald technique. Initially, the system was equilibrated for 1000 ps using a 0.01 fs timestep. After the equilibration, the system was simulated for a total of 2500 ps. The 2500 ps was split into 500 loops (5 ps each) to find the exact timeaveraged values. The simulations have been performed at six thermodynamic state points are summarized in Table 2.

CMD RESULTS
The densities of the above mentioned six state points are shown from CMD simulations and compared with the experimental data [15] in Fig. 1. It is important to note that the experimental densities fluctuate between 0.5% to 0.7%, and the temperatures vary ±1 K [15]. When the experimental densities of liquid HF are compared to this simulation, they match up well and show the expected linear trend for the liquid phase over the temperature and pressure ranges. Although the 2-site flexible model [12] gives a better agreement with the benchmarked experiment, the NCSU HF model also gives a reasonable agreement with them. The simulated potential energy at 373 K and 10 atm is -23.141±0.241 KJ/mol which is consistent with -23.08 KJ/mol [16][17]. In general, the radial distribution function (i.e. RDF, ( )) represents how the density of atoms in a system changes as a function of the distance from a reference particle. It is the ratio of the average local number density of particles at a distance, and the bulk density of particles. The expression for the RDF is as where ( ) and ( ) are the position and velocity of the atom at time respectively. Fig. 3 shows the MSD of hydrogen, fluorine, and HF molecule. The slope of all MSDs increases due to the increase in temperature. Due to the bonded system, the MSD of hydrogen, fluorine, and HF molecules confirms the diffusion coefficients are identical. The density of state (DOS) describes the system's excitation states. Hydrogen DOS is calculated using the normalized velocity autocorrelation function (VACF) of liquid HF. In general, the DOS (in equ. 7) is the Fourier transform of the normalized VACF. In CMD simulations, the VACF is calculated from the atomic trajectories. The normalized velocity autocorrelation function ( ), is calculated as Fig. 4 illustrates the variation in the normalized velocity autocorrelation function of the hydrogen atom in temperature for liquid HF. Fig. 5 depicts the self-diffusion coefficient of liquid HF with a comparison of simulation [12][13] (one of the MD models in ref. [13] based on the model in Ref. [5]) and experimental [23] results and the current CMD simulation results agree well with the literature.

TSL of Hydrogen in HF
Based on the CMD simulations and DOS data calculated from the normalized VACF (see Fig. 4), the TSL for H in HF and associated cross sections were evaluated using the Full Law Analysis Scattering System Hub (FLASSH) code [24]. Hydrogen diffusion in HF was modeled using Schofield model parameterized by DOS and diffusivity data. The diffusive and solid components were convolved, generating the incoherent TSL ( , ), where which is the fundamental input to the double differential cross section given by where is the bound cross section, and ′ are the incoming and outgoing neutron energies, respectively, is temperature, and and are dimensionless momentum and energy exchange [25][26].   In order to quantify the impact of these cross sections, the NCSU HF model TSLs were implemented in the ICSBEP HEU-SOL-THERM-039 benchmark where HF is the primary moderating material [15]. In this benchmark, six different configurations were arranged. Using the ENDF/B-VIII.0 libraries to calculate these critical configurations, the evaluations significantly overpredict [27]. However, with the addition of the TSL, the is greatly improved with an average reduction of 1170 pcm.

CONCLUSIONS
In this study, a new classical molecular dynamics simulation based NCSU HF model was established. The densities, potential energy, and bond length were compared to the experimental data and found to be in good agreement. The model demonstrates that the NCSU HF model is just as accurate as the other models and effective for nuclear data development. The DOS, TSL, and related thermal neutron scattering cross sections were produced based on this model. The impact from the HF TSL significantly improved the ICSBEP HEU-SOL-THERM-039 benchmark results.