Windowed multipole sensitivity to target accuracy of the optimization procedure

Abstract

This paper compares the accuracy of the windowed multipole direct Doppler broadening method to that of the ENDF-B/VII.1 libraries that come with MCNP6. Various windowed multipole libraries were generated with different maximum allowed relative errors. Then, the libraries were compared to the MCNP6 data via resonance integral and through single assembly Monte Carlo analysis. Since the windowed multipole uses resonance parameters, resonance integrals are only affected by the number of resonances included in the library and not by the order of the background fitting function. The relative performance of each library with varying maximum allowed error was evaluated. It was found that setting a maximum target relative error of 0.1% in the library provided highly accurate data that closely matches the MCNP6 data for all temperatures of interest, while still having suitable computational performance. Additionally, a library with a maximum relative error of 1% also provided reasonable accuracy on eigenvalue and reaction rates with a noticeable improvement on performance, but with a few statistically significant differences with the MCNP6 data.

Keywords:

• windowed multipole
• cross section
• nuclear data
• accuracy
• performance

1. Introduction

The windowed multipole (WMP) method is a new technique to perform Doppler broadening of resolved resonance cross section data in a highly efficient manner, based on the multipole formalism. Reich–Moore and Multi-Level Breit Wigner resolved resonance data can be converted into the multipole formalism via a process of partial fraction expansions [1,2]. The resulting pole and residue form can be analytically Doppler broadened. The primary drawback of the original multiple formalism is the computational cost, since each evaluation requires a Faddeeva function (a scaled complex complimentary error function) evaluation for each pole. This is a time consuming process if there are many poles. A brief overview of the mathematics is described in Section 2.

A simplification was found in that, for each energy point, only some poles contribute appreciably to the solution and fewer still fluctuate significantly. Therefore, instead of evaluating each pole exactly, some can be replaced with a polynomial [3]. This simplification, called the WMP formalism, improves computational performance by a significant margin at the slight cost of accuracy [4]. Furthermore, this format uses very little memory relative to pointwise data. This technique is explained in detail in Section 3.

There are three primary tunable parameters inside of the WMP method: inner window average size, curve fit order, and target accuracy. The first two components only affect the computational efficiency and memory requirements of a library. The last component, the target accuracy, is what maximum error the optimization routine will admit in the data. Reducing the accuracy will allow for faster but more approximate cross section reconstruction. This study will analyze how sensitive realistic reactor problems are to the library target accuracy.

To provide a suitable testbed to evaluate the target accuracy, three problems were studied. The first was a comparison of the resonance integrals between several libraries of differing target accuracies with that of the ENDF-B/VII.1 data from MCNP6. The results are shown in Section 4.2. The next problem evaluated the computational performance in an infinite hydrogen scattering medium. The time to compute a cross section was evaluated and the WMP libraries are compared in Section 4.3. Finally, the WMP method was implemented into OpenMC. A single assembly was run with the MCNP6 data and WMP libraries of varying accuracy. The impact on eigenvalue and energy-dependent tallies are presented in Section 4.5.

2. Multipole formalism

The multipole formalism is a mathematically exact alternate representation of Reich–Moore and Multi-Level Breit Wigner data [1,2]. In this form, the cross sections are represented by sums of poles and residues. For example, the general form of the 0 K cross section for reaction x is shown in Equation (1(1) $$\begin{array}{ccc}\hfill {\sigma }_x\left(E\right)& =& \frac{1}{E}\sum _j\Re \left[\frac{r_{j,x}}{p_j-\sqrt{E}}\right].\hfill \end{array}$$(1) ). In this equation, p_j are the poles and r_{j, x} are the residues corresponding to reaction x. For each set of quantum numbers, there corresponds a set of resonance levels. Each one of these resonance levels can be decomposed into a sum of poles and residues. The subindex j represents this decomposition, and spans all quantum numbers, all levels, and all poles necessary to represent the entirety of the cross section data (1) $$\begin{array}{ccc}\hfill {\sigma }_x\left(E\right)& =& \frac{1}{E}\sum _j\Re \left[\frac{r_{j,x}}{p_j-\sqrt{E}}\right].\hfill \end{array}$$(1) The key advantage to the multipole formalism is when Equation (1(1) $$\begin{array}{ccc}\hfill {\sigma }_x\left(E\right)& =& \frac{1}{E}\sum _j\Re \left[\frac{r_{j,x}}{p_j-\sqrt{E}}\right].\hfill \end{array}$$(1) ) is Doppler broadened. This integration is performed analytically and results in Equation (2(2) $$\begin{array}{ccc}& & {\sigma }_x(E,T)=\hfill \\ & & \frac{1}{2E\sqrt{\xi }}\sum _j\Re \left[i{r}_{j,x}\sqrt{\pi }W\left(z_0\right)-\frac{r_{j,x}}{\sqrt{\pi }}C\left(\frac{p_j}{\sqrt{\xi }},\frac{u}{2\sqrt{\xi }}\right)\right],\hfill \end{array}$$(2) ) (2) $$\begin{array}{ccc}& & {\sigma }_x(E,T)=\hfill \\ & & \frac{1}{2E\sqrt{\xi }}\sum _j\Re \left[i{r}_{j,x}\sqrt{\pi }W\left(z_0\right)-\frac{r_{j,x}}{\sqrt{\pi }}C\left(\frac{p_j}{\sqrt{\xi }},\frac{u}{2\sqrt{\xi }}\right)\right],\hfill \end{array}$$(2) where $$\begin{array}{ccc}\hfill u& =& \sqrt{E},\hfill \\ \hfill 2\sqrt{\xi }& =& \sqrt{\frac{kT}{A}},\hfill \\ \hfill z_0& =& \frac{u-p_j}{2\sqrt{\xi }},\hfill \\ \hfill C\left(\frac{p_j}{\sqrt{\xi }},\frac{u}{2\sqrt{\xi }}\right)& =& 2p_j{\int }_0^{\infty }d{u}^{\text{'}}\frac{exp\left[{\displaystyle -\frac{{(u+u^{\text{'}})}^2}{4\xi }}\right]}{p_j^2-u^{\text{'}2}}.\hfill \end{array}$$There are two important functions in Equation (2(2) $$\begin{array}{ccc}& & {\sigma }_x(E,T)=\hfill \\ & & \frac{1}{2E\sqrt{\xi }}\sum _j\Re \left[i{r}_{j,x}\sqrt{\pi }W\left(z_0\right)-\frac{r_{j,x}}{\sqrt{\pi }}C\left(\frac{p_j}{\sqrt{\xi }},\frac{u}{2\sqrt{\xi }}\right)\right],\hfill \end{array}$$(2) ). The first, W, is the Faddeeva function. There are several algorithms that can evaluate this function quickly using various approximations [5,6]. The second one, C, has been found to be negligible except at very low energy [1,7].

Unfortunately, performing Faddeeva function evaluations for all poles is not efficient. Some isotopes have very large numbers of poles. Each resonance contributes 2(l + 1) poles, where l is the neutron orbital angular momentum. In the case of ²³⁸U in the ENDF-B/VII.1 library [8], this results in 11,520 poles. To make this technique sufficiently fast to be useful, approximations must be made which led to the WMP formalism.

3. Windowed multipole

One key feature of a singular pole and residue is that the region of high fluctuation, and thus the region most difficult to approximate, is constricted to a narrow band near $\sqrt{E}\approx p_j$ [7]. By evaluating some poles exactly and approximating the rest as a curve fit, the number of poles that must be summed substantially decreases.

To efficiently determine which poles are important at run time, the entire energy range is chopped into windows. In prior studies, it was found that having the windows equally sized in momentum proved superior to equally spaced in lethargy or energy [4]. During preprocessing, an optimization routine is run on each window to identify the index of the first and last poles required to maintain a specific accuracy criterion, forming the outer window. These indices are stored inside of the library for use at runtime. The remaining components of the library (poles and pointwise data) are curve fit with a polynomial. These polynomials are relatively smooth, since they do not contain the fluctuating portions of resonances. This makes them mostly temperature insensitive except in the thermal region. The general process is depicted in Figure 1. These indices and curve fits are valid for the entire domain of the inner window.

Figure 1. Library generation mechanism for windowed multipole.

The curve fit used took the form of Equation (3(3) $${\sigma }_{cf}\left(E\right)=\sum _{i=0}^N{C}_i{E}^{i/2-1}.$$(3) ). This curve fit was chosen as this series contains the 1/E term from Equation (1(1) $$\begin{array}{ccc}\hfill {\sigma }_x\left(E\right)& =& \frac{1}{E}\sum _j\Re \left[\frac{r_{j,x}}{p_j-\sqrt{E}}\right].\hfill \end{array}$$(1) ), the 1/v term common to absorption cross sections at low energies, and a constant term (3) $${\sigma }_{cf}\left(E\right)=\sum _{i=0}^N{C}_i{E}^{i/2-1}.$$(3)

4. Results

The general goal was to analyze how sensitive the integral accuracy and the performance of the WMP formalism was to the target accuracy of the optimization procedure. First, several libraries were generated. Their properties are described in Section 4.1. Then, the libraries were compared to the ENDF-B/VII.1 data sourced from the MCNP6 distribution. The resonance integrals were compared as a function of initial library accuracy and temperature. This is described in Section 4.2. The relative performance of the libraries are studied in Section 4.3 and the memory requirements are quantified in Section 4.4. Finally, the libraries are loaded into a Monte Carlo code and used to simulate a single assembly. Tallies and eigenvalues are compared in Section 4.5.

4.1. Library

In order to test how the target accuracy actually affected results, four libraries of varying target accuracy were generated from the ENDF-B/VII.1 data. The two isotopes processed were ²³⁵U and ²³⁸U. The accuracy settings used were 0.01%, 0.1%, 1%, and 10% maximum allowed relative error over the temperature range from 300 to 3000 K. A parametric search over curve fit order and inner window average size was performed, and the library which took the shortest average time to evaluate a cross section was used. These are listed in Table 1.

Table 1. Library window and curve fit order.

	Max relative	Inner window	Curve fit
Isotope	error	(eV)	order
^235U	0.01 %	1	6
	0.1 %	1	2
	1 %	1	2
	10 %	1	2
^238U	0.01 %	10	3
	0.1 %	10	2
	1 %	10	2
	10 %	10	2

An interesting thing to note is that the optimization process selected the same inner window size for each target accuracy and, except for the 0.01% case, the same curve fit order. This indicates that the major difference between the four libraries is the number of resonances to be treated explicitly.

4.2. Resonance integral comparison

For each isotope, the resonance integral of the total cross section was calculated. Specifically, the calculation performed is shown in Equation (4(4) $${\text{RI}}_t={\int }_{E_1}^{E_2}\frac{{\sigma }_t\left(E\right)}{E}dE.$$(4) ). In the case of ²³⁵U, the limits of the integration were 0.1 to 2250 eV. For ²³⁸U, the limits were 0.1 eV to 20 keV. For the MCNP6-sourced data, integration was performed exactly for each line segment. The WMP data were evaluated at each energy point from the MCNP6-sourced data and integrated identically to provide a suitable comparison. The same process was run for all of the temperatures available: 0.1, 250, 293.6, 600, 900, 1200, and 2500 K (4) $${\text{RI}}_t={\int }_{E_1}^{E_2}\frac{{\sigma }_t\left(E\right)}{E}dE.$$(4)

4.2.1. Resonance integral results for ²³⁵U

As shown in Figure 2, all the WMP libraries become increasingly inaccurate with increasing temperature. However, the inaccuracies are still small considering the actual target accuracy. The 10% library yields a maximum error below 0.5%, which was expected since the WMP treats the local resonances exactly. The maximum errors in the library always occur at low cross section values, and the impact on resonance integral increases slightly with temperature as neglected resonances start contributing to the background term. In the case of the 0.01% library, the error at 2500 K of the resonance integrals is −0.003%. The value is −0.011% for the 0.1% library. The reduction in resonance integral going from the MCNP6-sourced data to the WMP library was found to be consistent for all temperatures, all reactions (σ_t, σ_s, σ_a, σ_f, and σ_c), and all library configurations tested.

Figure 2. Multipole library resonance integral compared to MCNP6-sourced data, ²³⁵U.

4.2.2. Resonance integral results for ²³⁸U

The exact same process was repeated for ²³⁸U, yielding Figure 3. Most notably, the errors do not begin at 0% as is the case for ²³⁵U, but at −0.004%. This likely indicates a slight discrepancy between NJOY and our processing code. The ²³⁸U library is far less sensitive to the target accuracy, with the 10% library being 50 times more accurate than the equivalent ²³⁵U library. Similar to the ²³⁵U results, the trend is negative and this holds for all reactions, temperatures, and libraries.

Figure 3. Multipole library resonance integral compared to MCNP6-sourced data, ²³⁸U.

4.3. Performance

These libraries were loaded into a simple hydrogen scattering slowing down benchmark. Neutrons are started at the end of the resolved resonance range and terminated once they reach 0.1 eV. The time to compute a full set of cross sections at each collision is calculated, and presented in Figure 4 for ²³⁵U and in Figure 5 for ²³⁸U. All runs were done on an Intel^® i7-970 @ 3.20GHz.

Figure 4. Cross section evaluation time in microseconds, ²³⁵U.

Figure 5. Cross section evaluation time in microseconds, ²³⁸U.

There are two notable points in these graphs. First, the more accurate the library, the more time it takes to run. The reason for this is that the higher accuracy libraries include more poles treated exactly inside of each window, requiring more Faddeeva evaluations per cross section evaluation. The second significant point is that the high-accuracy libraries have a positive slope as a function of temperature. The Faddeeva function used is most efficient when the value of |z₀| in W(z₀) (in Equation (2(2) $$\begin{array}{ccc}& & {\sigma }_x(E,T)=\hfill \\ & & \frac{1}{2E\sqrt{\xi }}\sum _j\Re \left[i{r}_{j,x}\sqrt{\pi }W\left(z_0\right)-\frac{r_{j,x}}{\sqrt{\pi }}C\left(\frac{p_j}{\sqrt{\xi }},\frac{u}{2\sqrt{\xi }}\right)\right],\hfill \end{array}$$(2) )) is large since it can be approximated with fewer terms. Since $\left|z_0\right|\propto 1/\sqrt{T}$, evaluation of the cross section is slower at high temperatures. Furthermore, the high-accuracy libraries have more poles and so the relative fraction of computational time for the Faddeeva function is higher.

4.4. Memory requirements

One advantage of the WMP method is the reduction in memory requirements as compared to pointwise. There are two reasons for this. The first is that the entire resolved energy range (several hundreds of thousands of points for some isotopes) has been replaced with a few thousand complex numbers and curve fits. The second is that since the resolved energy range has been eliminated, the other cross sections (inelastic, etc.) do not need to be unionized to those points. The memory requirements for each isotope was calculated. For this calculation, the memory listed is the sum of the sizes of all arrays for the resolved resonance data and for all secondary distributions. This does not include angular distributions and other components shared by both libraries. The memory requirements for the WMP library are shown in Table 2.

Table 2. Memory requirements for windowed multipole.

Library	^235U Mem. kB	^238U Mem. kB
10%	543	542
1%	543	542
0.1%	543	542
0.01%	841	589

The memory requirements for a single-temperature pointwise data-set from the MCNP6-sourced libraries varied with temperature. The memory requirements are shown in Table 3.

Table 3. Memory requirements for the MCNP6-sourced data.

Temperature	^235U Mem. kB	^238U Mem. kB
0.1 K	12,993	21,014
250 K	4478	9014
293.6 K	4240	8693
600 K	3370	7391
900 K	3005	6769
1200 K	2795	6374
2500 K	2429	5514

Overall, the WMP library footprint was mostly insensitive to the target accuracy except when the curve fit polynomial order increased in the 0.01% case. Furthermore, a WMP library, which is valid for all temperatures from 300 to 3000 K, took significantly less memory than any single-temperature ACE library.

4.5. Assembly OpenMC comparison

Lastly, the WMP libraries for ²³⁵U and ²³⁸U were loaded into the OpenMC Monte Carlo code [9] and a single assembly from the benchmark for evaluation and validation of reactor simulations (BEAVRS) was run [10]. The remaining components of the library, such as the unresolved range and angular distributions, were sourced from the MCNP6 data. The assembly chosen was enriched to 3.2%, and had no burnable absorber. The simulation was performed with the settings shown in Table 4. For this particular test, the 0.01% library was not used since that level of accuracy could not be maintained below 10⁻² eV. Research is ongoing as to the reason why.

Table 4. Assembly OpenMC run settings.

Parameter	Value
Inactive batches	200
Active batches	500
Neutrons per batch	4 million

For each tally, the eigenvalue is listed in Table 5. Notably, the WMP library results appear to be converging towards the MCNP6 data results. In the 0.1% case, the eigenvalues are only 5 pcm different, making the results identical within two standard deviations. This indicates a high degree of accuracy for the 0.1% library.

Table 5. Assembly OpenMC eigenvalue.

Data	Eigenvalue
MCNP6	1.21201 ± 0.00002
10% multipole	1.21216 ± 0.00002
1% multipole	1.21209 ± 0.00002
0.1% multipole	1.21206 ± 0.00002

The assembly averaged fission neutron production rate was tallied in 16 energy bins. The tally relative error is plotted in Figure 6 for the 0.1% library, Figure 7 for the 1% library, and Figure 8 for the 10% library. All errors were plotted with the same y-axis to allow for ease of comparison.

Figure 6. Fission production rate relative error, 0.1% target library.

Figure 7. Fission production rate relative error, 1% target library.

Figure 8. Fission production rate relative error, 10% target library.

In this particular case, we find that the 10% library has significant fluctuations across the spectrum, and is as much as 0.35% off in some of the bins. Only the last three of the 16 bins were within two standard deviations. The 1% library fairs much better with all tallies within 0.1%. Furthermore, only 5 of 16 bins were outside of two standard deviations. In the case of the 0.1% accuracy library, only four of 16 bins were outside of two standard deviations.

5. Conclusion

With the new WMP formalism, there had been relatively little testing of the accuracy and performance of the algorithm. With the target maximum relative error as a tunable parameter, it becomes imperative to see just how sensitive problems are to tuning. This paper has measured on a few sample problems the impact of this parameter.

It was found that for resonance integrals, the target accuracy has relatively little impact for ²³⁸U. No library exceeded 0.025% error at any temperature. For ²³⁵U, however, the impact was larger, around 0.5% error for the worst library at the highest tested temperature. This can be reduced to below 0.1% or less by using a more accurate library. For both isotopes, increasing the temperature decreased the accuracy of the library since approximated resonances start contributing at greater distances.

In performance, the higher accuracy libraries required far longer to evaluate than the low accuracy libraries. This was due to the increased number of Faddeeva function evaluations required. As temperature increases, the Faddeeva function becomes more costly timewise, further increasing the computational cost. For all libraries, time never exceeded 10 µ s per evaluation. In the case of the 0.1% library, each evaluation cost between 2 and 4 µ s for ²³⁸U and ²³⁵U, respectively. Both libraries used substantially less memory than any single-temperature pointwise library compared against.

Lastly, in the assembly test, the 0.1% tolerance library matched fission neutron production rate within two standard deviations for 12 of 16 energy bins. The 1% library was within two standard deviations in 11 of 16 energy bins. In the case of the 10% library, observed differences on reaction rates were as large as 0.35%.

With all of this combined together, a library prepared to 0.1% tolerance should provide sufficient accuracy while at the same time be suitably fast. However, that level of accuracy might not be needed for all isotopes since the 1% library also showed suitable accuracy. Further studies are required on a full suite of problems with varying spatial size and isotopic concentrations to evaluate the full impact of varying target accuracy as a way to improve performance.

Acknowledgements

This work was supported by the Consortium for Advanced Simulation of Light Water Reactors (CASL) under US Department of Energy [contract No. DE-AC05-00OR22725]; by the DOE’s Center for Exascale Simulation of Advanced Reactors (CESAR) under US Department of Energy [contract No. DE-AC02-06CH11357]; by a US Department of Energy Nuclear Energy University Programs Graduate Fellowship [contract No. DE-NE-0000102].

Disclosure statement

No potential conflict of interest was reported by the authors.