<figure><img src="image_1.png" /></figure>

HAL
open science

# Computation of Sensitivity in Monte-Carlo Codes: Investigation of Convergence for Reactor Cores with High Dominance Ratio

M Rababah, A Bidaud

► To cite this version:

M Rababah, A Bidaud. Computation of Sensitivity in Monte-Carlo Codes: Investigation of Convergence for Reactor Cores with High Dominance Ratio. International Conference on Physics of Reactors, Apr 2024, San Francisco, United States. pp.739-747, 10.13182/physor24-43648 . hal-04702404

HAL Id: hal-04702404

https://hal.science/hal-04702404v1

Submitted on 29 Oct 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Computation of Sensitivity in Monte-Carlo codes: Investigation of Convergence for Reactor Cores with High Dominance Ratio.

Rababah M¹,*, Bidaud A¹

¹LPSC, Université Grenoble-Alpes, CNRS/IN2P3
53, rue des Martyrs, 38026 Grenoble, France

*[leave space for DOI, which will be inserted by ANS]*

## ABSTRACT

The utilization of the Monte Carlo method in reactor physics has expanded notably, and its use for sensitivity calculations has expanded with the integration of perturbation theory in numerous Monte Carlo-based codes, enhancing the sensitivity calculations significantly. Despite the method's growing application, research on the statistical convergence of sensitivities remains relatively limited.

This paper provides a review of various methods employed for sensitivity calculations and statistical convergence of these sensitivities. Our literature review indicates a general consensus among researchers, suggesting that reasonable calculation time is sufficient for obtaining converged sensitivities. Nevertheless, we identified a case in the literature, in high-dominance ratio reactor where sensitivities failed to converge, even with very large computational cost.

We present preliminary findings from simulating two reactor configurations: the Dummy-Core (which is a fictious SMR core with plate fuel assemblies) and the UAM-PWR Benchmark, with dominance ratios of 0.94 and 0.99, respectively. Our research agrees with literature findings in relatively low dominance ratio, but with high dominance ratios, neither the effective multiplication factor nor the power distribution sensitivities converge. This non-convergence persists despite augmenting the total neutron history. Notably, such computations demand substantial computational resources.

In conclusion, we outline our prospective research focusing on the use of propagation of sensitivities on higher modes to accelerate the convergence and enhance computational efficiency in Monte Carlo simulations for reactor physics.

*Keywords:* Monte Carlo, Neutron Transport, Sensitivity Convergence, Dominance Ratio, Latent Generation

## 1. INTRODUCTION

Despite being the reference neutron transport method, Monte Carlo used to be limited in its capacity to calculate sensitivities, or adjoint-based outputs, such as kinetics parameters. This has completely changed with the very fast development of perturbation methods and generalization in most of the codes in the last decade. This field is still growing and few researchers have focused on the question of the statistical

<footnote>*mohammad.rababah@lpsc.in2p3.fr</footnote>
convergence of those new objects, such as eigenvalues, in particular for problems with high dominance ratios.

The convergence of local variables such as local fluxes in Monte Carlo may be much slower than global ones, e.g. keff. The convergence of their sensitivities with the number of histories simulated is not well known and their actual convergence may be doubtful as it depends on « new » parameters of the Monte Carlo method used for the calculations of sensitivities such as the number of « latent » generations. Latent generation is the terminology used in SERPENT defined by the number of propagation generations through which events (perturbations) are weighted to compute sensitivities [1] [2].

In the first part of the paper, we provide a literature review of current research related to sensitivity and uncertainty analysis - used in Monte-Carlo codes - and their convergence. There is limited research on statistical convergences of sensitivities, especially in reactors with high dominance ratio. Reactors with a high dominance ratio pose challenges in achieving convergence, particularly for local variables sensitivities, such as power distribution sensitivity. The literature review suggest that 20 latent generations is reasonably enough for obtaining converged sensitivities. Nevertheless, we identified a case in the literature [3], in high-dominance ratio reactors where sensitivities failed to converge, even after 100 latent generations.

In the second part, we present preliminary results for sensitivity convergence in two different reactor configurations with distinct dominance ratio. We performed sensitivity calculations using SERPENT2.2 for two reactor configurations, the Dummy-Core and the UAM UOXPWR benchmark. The Dummy-Core, is a fictious PWR core designed by TechnicAtome [4]. The core is composed of 8 by 8 fuel assemblies and thus, relatively small. Two fuel assemblies are used, standard fuel assemblies and hafnium-plated controlled assemblies. The fuel pattern consists of a succession of standard Uranium silicate slabs composed of planar JHR-type fuel and contained in aluminum cladding, the core has similar thermal neutron energy spectrum as in typical PWRs. The UAM UOX-EPR benchmark is a well-known benchmark set by NEA-OECD, aims to unravel uncertainties in key reactor core parameters [5].

Such computations demand tremendous computational resources. Therefore, it is crucial to explore less computationally expensive methods. In the last part, we present our future work that aims to tackle the statistical convergence issue, by testing the use of higher modes projections of perturbations to calculate sensitivities [6].

## 2. LITERATURE REVIEW

The Monte Carlo method is widely used for reactor physics calculations and has seen significant improvements in sensitivity calculations in recent years. However, certain limitations have been identified in Monte Carlo criticality calculations, such as convergence issues of local variables (e.g., fission source) in large, loosely-coupled systems with a high dominance ratio. Additionally, a phenomenon known as the "clustering effect" has been observed, where neutrons tend to gather together, leading to undesirable fluctuation [7] [8]. Nevertheless, a lot of research has been already conducted to improve the accuracy of Monte-Carlo criticality calculations.

In contrast, limited research is available on the statistical convergence of sensitivity calculations. Kiedrowski et al [9] provides an overview of various methods and techniques used in Monte Carlo perturbation and sensitivity calculations, with a focus on techniques for K-eigenvalue problems. Perturbing the fission source in these calculations requires propagating the effect across iterations due to the iterative solution process in k-eigenvalue problems. This propagation introduces complexities and potential inaccuracies, particularly when the scattering term significantly impacts the shape of the neutron flux. Modern approaches address this challenge through explicit or implicit adjoint weighting.
Several techniques have been developed to estimate the effect of perturbed fission sources, including fission source correction using the correlated sampling and differential operator sampling methods [10] [11], adjoint-based perturbations for the k-eigenvalue using the iterated fission probability method [12], an extension to reaction rate ratios using generalized perturbation theory [13], and a recent development using a collision history approach allowing for the calculation of sensitivity coefficients of bilinear ratios and generalized responses [1] [9].

Both correlated and differential operator sampling techniques involve discretization and are prone to large variances. Nagaya and Mori [10] have improved these methods by splitting the change in k-eigenvalue into two terms expressing direct and indirect effects. The indirect effect is estimated by propagating weights between iterations, with the concept of "propagation batches" utilized to limit variance. Although increasing the number of propagation iterations improves the accuracy of the correction for the fission source perturbation, it also increases the statistical variance. Empirical studies suggest that 5 to 10 iterations yield acceptable correction before the variance becomes significant. Most methods adopt the concept of "propagation batches" with around 10 generations deemed sufficient for convergence.

Another article by Kiedrowski et al [14], reveals that accurate importance weighting factors for localized quantities require numerous generations to settle. Problems with a high dominance ratio may necessitate even more generations. However, global quantities, such as kinetics parameters, converge with far fewer iterations. It also indicates that the convergence in the adjoint flux is observed after around twenty latent generations, although the statistical noise increases substantially as the number of latent generations increases. According to the author, this is due to less frequent scoring in tallies, requiring substantially longer runs for statistically meaningful results. However, the difference between two latent generations and twenty generations is negligible (~1%). Thus, 10 to 20 generations are "said" to be likely appropriate for typical power reactors.

Jinaphanh et al [15] also observed an increase in variance with respect to increase number of latent generations. The author observed that the sensitivity coefficients are relatively stable with respect of latent number generation (up to 10 generations), yet the variance increases.

Concerning high dominance ratio configuration, in a paper by Lopez et Bidaud [3], sensitivity calculation was performed on UAM GENIII Benchmark. They observed that the sensitivities of multiplication factor converges faster than those of power distribution which require at least 100 latent generation to converge.

### 3. METHODOLOGY

In this study, we explore the statistical convergence of sensitivities in two distinct reactor configurations: the Dummy-Core and the UAM benchmark which have different dominance ratios. The dominance ratio was calculated using the fission matrix method deployed in Serpent2.2 and found to be 0.94 and 0.99 for the dummy core and for the UAM benchmark respectively.

We used the continuous energy Monte-Carlo code SERPENT2.2. to calculate the sensitivities. The simulations were executed over 100 latent generations with different neutron history to test the convergence speed of sensitivities and their uncertainties. The sensitivities were calculated in different detector locations representing diverse fuel assemblies to ascertain whether location influences convergence. Additionally, the effect of different perturbations on convergence speed was analyzed, focusing on total cross-section, fission cross-section, and others.

The simulations were run with 50000 neutrons per batch and for 25000 batches, with 100 latent generations, and 56g energy grid from the ENDF/B-VII.1 evaluated nuclear data library [16] that we use for uncertainty analysis, and the perturbation were run for all cross sections. Notably, such computations demand substantial computational resources; a single simulation with 100 latent generations, 50,000 neutrons per
batch, and 10,000 cycles required approximately 21 days of runtime on 72 CPUs allocated 78 GB of
memory.

Only sensitivities of energy-integrated sensitivities of both the effective multiplication factor and the power distribution to total cross-section are presented because of the similarity in the behavior of sensitivity to other cross sections.

4. RESULTS

For each subsection, several figures are presented to provide a comparative analysis between the convergence of multiplication factor and power distribution sensitivities with number of latent generations, to observe the change in sensitivity convergence with neutron history, and to shed the light on how sensitivities converge differently with different dominance ratio.

The sensitivity of power distribution is obtained via detector ratio, one detector representing a fuel assembly (in the case presented, the central fuel assembly) to the detector representing all the core.

**4.1. Dummy Core**

Figure 1 shows two graphs, the one on the left is the variation of sensitivity of the multiplication factor to the change in the total cross section as a function of latent generations from 1 to 100 and number of batches. We see that it requires around at least 20 latent generations for the sensitivity to converge. Moreover, there is small difference when changing the number of batches, knowing that number of neutrons per batch is 50000 neutrons for all cases. While sensitivity is flattening after 20 latent generation, the relative variance is always increasing with number of latent generations, Figure 1 right. But the increase in variance is relatively small and does not dramatically change the expected sensitivities.

Figure 2 left shows the sensitivity of the power distribution to the change in total cross section with different number of latent generations, with 50000 neutrons and 25000 batches. As we couldn't clearly see convergence of sensitivities, more neutron histories are probably needed. Figure 2 on the right, shows the change in relative variance, which is also increasing with number of latent generations, and is much higher than the one for multiplication factor sensitivities. Yet, despite the high computational cost of such calculation, it is not clear whether the sensitivity converges or not.

<figure><img src="image_3.png" /></figure>

<figure><img src="image_4.png" /></figure>
**Figure 1.** left: Sensitivity of keff in the dummy core to total cross section, right: relative variance of the sensitivtiy. The legend on both graphs represents the number of batches.

<figure><img src="image_5.png" /></figure>

**Figure 2.** left: Sensitivity of power distribution in the dummy core to total cross section, right: relative variance of the sensitivtiy. The legend on both graphs represents the number of batches.

## 4.2. UAM Benchmark

In Figure 3, we see the change in effective multiplication factor sensitivity to the total cross section varies quite largely with number of latent generations, and with number of neutron history, unlike the case of the dummy-core. Even with 100 latent generations the convergence speed is much slower considering the increase in relative variance with number of latent generations, Figure 3, right.

<figure><img src="image_6.png" /></figure>

**Figure 3** left: Sensitivity of multiplication factor in the UAM core to total cross section, right: relative variance of the sensitivtiy. The legend on both graphs represents the number of batches.

In Figure 4, even though sensitivities calculated at 20000 and 25000 batches seems to match, we clearly see the divergence of the power distribution sensitivity with number of latent generations. moreover, the relative variance behaves rather strangely. One explanation for the drop in variance in the first 20 latent generations, is the change in sign in sensitivity. But after 20 latent generations, the variance also does not stabilize despite a false plateau by 50 latent generations.
We also want to highlight that the relative variance behavior is neither consistent with the literature nor with results of dummy core, we expected to have a similar trend of increasing variance with number of latent generations, but here, it is not the case. This could be explained either by the fact that sensitivities predictions and their uncertainties are not accurate, or that even more neutron histories are required. It is also important to note that it took around 55 days to reach 25000 batches with 50000 neutrons per batch, which in the end, even at 100 latent generations did not lead to converged value, which is around -2.3, and we obtained it using direct sensitivity calculation. Therefore, it is crucial to explore more accurate and less computationally expensive method, which we describe in the following section.

<figure><img src="image_7.png" /></figure>

<figure><img src="image_8.png" /><figcaption>Figure 4. left: Sensitivity of power distribution in the UAM core central assembly to total cross section, right: relative variance of the sensitivity. The legend on both graphs represents the number of batches.</figcaption></figure>

## 5. DISCUSSION AND FUTURE WORK

It is clear that in a reactor core with high dominance ratio, the statistical convergence of sensitives, in particular of local variables, is an issue, beside the high computational cost in term of time and memory, we are neither sure whether these results are accurate, i.e., converge toward the “good” value, nor if they are precise, i.e., known with a small enough statistical uncertainty.

One way to accelerate convergence may be via using decomposition of power distribution perturbation on modes [6]. The effect of perturbations on fundamental fission source distribution can be represented by Equation (1)

$$ \frac{dS_0}{dP} = \sum_{i=0}^{\infty} S_i \cdot \frac{S_i^{\dagger} \frac{dF}{dP} S_0}{k_0 - k_i} \quad (1) $$

Where $\frac{dS_0}{dP}$ represent the derivatives of the fundamental fission source distribution $S_0$ with respect to the perturbation $P$. $\frac{dF}{dP}$ represents the effect to the perturbation in $P$ on the fission kernel. $F \cdot S_i \cdot S_i^{\dagger}$ and $k_i$ from the $i^{th}$ eigentriplet. It is assumed that the forward and adjoint eigenmodes respect the bi-orthogonality condition,
$$ \langle S_i^\dagger | S_j \rangle = \begin{cases} 1, & i = j \\ 0, & i \neq j \end{cases} $$

Aufiero et al [6]. have shown that the effect of perturbations on the fission kernel can be estimated via Monte Carlo Perturbation theory and scores of advanced tallies that are defined globally which may allow for a faster statistical convergence than local tallies.

As a first step, we try to check whether the convergence speed of ratio of a larger detector to the core detector is better, as mode sensitivities will be defined as ratios of full geometry tallies

Figure 5 shows the results of power distribution sensitivity to total cross section in the UAM core, the left and center figures shows the sensitivity and the relative variance when detector is scored on a group of fuel assemblies instead of one, the right figure shows the sensitivity when the detector represents the central fuel assembly.

We calculate the reference value by changing the density of 10% and directly calculating the sensitivity of power distributions. We see that there is an improvement in convergence speed when calculating sensitivity of a group of fuel assemblies rather than on a fuel assembly level. The sensitivity here converged faster with number of total neutrons, i.e. have smaller relative uncertainties, and diverged less with number of latent generations when using the larger size detector. We see that in the large detector, the sensitivity converges toward the reference value, where in the small detector sensitivity does not converge even at 100 latent generations.

Furthermore, we see that by 50 generations in the case of small detector, the sensitivities tend to plateau with variations that tends to be lower that the claimed uncertainties, which may lead the user to believe that sensitivities as a function of latent generation are converging. But in fact, they are not as the sensitivity eventually continues to grow with higher number of generations.

In the case of the bigger detector, the sensitivity reaches the reference value with a lower number of latent generations. With the highest number of neutron batches, the obtained value is compatible with the reference value given the uncertainty bars, moreover, the relative variance is one order of magnitude less compared to the relative variance of small detector. But when a smaller number of neutron histories are used, the user may be misled as the value once again tend to reach a plateau, that can be considered flat given the given uncertainties. Nevertheless, the value is not compatible within a 2 sigmas interval with the reference value.

Utilizing the sensitivities of modes to calculate sensitivities of power distribution will need some development. We plan to make those development and to test whether the uncertainties on the sensitivities are accurate and converge faster than those on local detectors sensitivities.

<figure><img src="image_9.png" /><figcaption>Figure 5. Sensitivity of power distribution (detector ratio) in UAM-PWR to change in total cross section. left: large detector Sensitivity, center: Large detector relative variance on sensitivtiy , right: small detector sensitivity. The legend on the graphs represents the number of batches.</figcaption></figure>
## CONCLUSIONS

There have been vast advancements in Monte Carlo methods capabilities to calculate sensitivities, which was facilitated by the rapid growth of computing efficiency and the innovative methods of to calculate sensitivities. Nevertheless, the literature review shows some gaps that needs to be improved, such as the statistical convergence of sensitivities in high dominance ratio configurations. Preliminary results have shown alignment with existing literature for configuration with relatively low dominance ratio. In particular the rather strange but widely observed increase of variance with the number latent generations has been reproduced. Nevertheless, the results on statistical convergence of sensitivities in high dominance ratio cores is not yet fully understood, but it is clear that both effective multiplication factor and power distribution do not converge easily. These simulations are very demanding in term of time and memory. Therefore, it is crucial to explore and test alternative methods to accelerate convergence, reduce uncertainty and improve computing efficiency.

## ACKNOWLEDGMENTS

We express our gratitude to the CNRS for funding this research through the joint project NEEDS, in which CNRS, CEA and IRSN work together to study the propagation of uncertainty in fuel burn-up calculations.

We also would like to express our warm regards and gratitude to Manuele Aufiero from Milano Multiphysics for his support and suggestions in this paper.

## REFERENCES

[1] M. Aufiero, A. Bidaud and M. Hursin, "A collision history-based approach to sensitivity/perturbation calculations in the continuous energy Monte Carlo code SERPENT," *Annals of Nuclear Energy*, 2015.

[2] M. Aufiero, M. Martin and M. Frat, "XGPT: Extending Monte Carlo Generalized Perturbation Theory capabilities to continuous-energy sensitivity functions," *Annals of Nuclear Energy*, 2016.

[3] B. A. Lopez P, "SENSITIVITY ANALYSIS AND ITS CONVERGENCE THROUGH MONTE CARLO CALCULATIONS FOR THE UAM GEN-III BENCHMARK: APPLICATION TO POWER DISTRIBUTIONS," in *The International Conference on Mathematics and Computational Methods Applied to Nuclear Science and Engineering*, North Carloina, 2021.

[4] S. Nicolas, "A dummy core for V&V and education and training purposes at TechniAtome: in and ex-core calculations," in *18th conference of the International Group on Research Reactors and IAEA workshop on safety reassessment of research reactors in light of the lessons learned from the Fukushima Daiichi accident*, 2017.

[5] k. Ivanova and M. Avramova, "Benchmarks for Uncertainty Analysis in Modelling (UAM) for the Design, Operation and Safety Analysis of LWRs - Volume I: Specification and Support Data for Neutronics Cases (Phase I)," *Organisation for Economic Co-Operation and Development*, 2013.

[6] M. Aufiero and M. Fratoni, "A new approach to the stabilization and convergence acceleration in coupled Monte Carlo-CFD calculations: The Newton method via Monte Carlo perturbation theory," 2017.
[7] I. Mickus, "Does neutron clustering affect tally errors in Monte Carlo criticality calculations?," *Annals of Nuclear Energy*, 2021.

[8] E. Dumonteil, "Patchy nuclear chain reactions," *Commun Phys* 4, 2021.

[9] B. C. Kiedrowski, "Review of Early 21st-Century Monte Carlo Perturbation and Sensitivity Techniques for k-Eigenvalue Radiation Transport Calculations," *NUCLEAR SCIENCE AND ENGINEERING*, 2017.

[10] Y. N. &. T. MORI, "Impact of Perturbed Fission Source on the Effective Multiplication Factor in Monte Carlo Perturbation Calculations," *Journal of Nuclear Science and Technology*, 2005.

[11] Y. Nagaya and F. Brown, "IMPLEMENTATION OF A METHOD TO ESTIMATE CHANGE IN EIGENVALUE DUE TO PERTURBED FISSION SOURCE DISTRIBUTION INTO MCNP," 2003.

[12] Raskach Fedorovich Kirill, "A Technique for Accounting for Multigroup Cross-Section Interdependence in Sensitivity Calculations," 2017.

[13] C. Perfetti and B. Rearden, "A new method for calculating generalized response sensitivities in continuous-energy Monte Carlo applications in SCALE," 2013.

[14] B. Kiedrowski and F. Brown, "AdjointWeighted Tallies for k-Eigenvalue Calculations with Continuous-Energy Monte Carlo," *Nuclear Science and Engineering*, 2011.

[15] A. Jinaphanh, C. Carabajal and A. Zoia, "Implementation and Testing of the GPT Sensitivities in TRIPOLI-4," in *International conference on physics of reactors*, Pittsburg (United States), 2022.

[16] M. B. Chadwick and et al., "ENDF/B-VII.1 Nuclear Data for Science and Technology: Cross Sections, Covariances, Fission Product Yields and Decay Data,", *Nuclear Data Sheets*, vol. 112, no. 12, pp. 2887-2996, 2011.