Nuclear data effect on BWR stability parameter uncertainty at stable core conditions

The objective of this study is to quantify the influence of nuclear data uncertainties on both steady-state and dynamic core parameters in BWRs, with a specific focus on the uncertainty quantification of stability parameters such as the decay ratio and resonance frequency at stable core conditions. This investigation is carried out within the framework of the OECD/NEA benchmark on the Oskarshamn-2 stability event that occurred in February 1999. The nuclear data uncertainties are propagated, using the in-house SHARK-X platform, through the 2-D lattice calculations using CASMO5 and followed by downstream static and dynamic calculations using and SIMULATE3 and SIMULATE-3 K, respectively. For steady-state, the uncertainty is estimated in terms of k-eff, radial and axial power peaking factors, void fraction, and average fuel temperature, while for the transient, the estimated uncertainties are obtained, for the first time, for the decay ratio and resonance frequency. The results indicate that nuclear data uncertainties have a non-negligible impact on most of the analyzed parameters. More specifically, the uncertainty in k-eff, radial and axial power peaking factors, nodal void fraction peak, and nodal fuel temperature peak for steady-state could reach about 0.5 %, 0.4 %, 2.8 %, 0.1 %, and 1.5 %, respectively, For transient, the uncertainty due to nuclear data for the two stability parameters could reach about 9 % and 2 % for the decay ratio and resonance frequency, respectively.


Introduction
Over the past few years, an OECD/NEA International benchmark has been initiated with the focus on the stability event that occurred at the Swedish Oskarshamn-2 nuclear power plant (Kozlowski, et al., 1999).This benchmark has been divided into three phases.The first phase involved analyzing the stability event that occurred on February 25th, 1999.The aim of this phase was to gain insights into the root cause of the stability event and to evaluate the effectiveness of different stability analysis methods/codes.The second phase focused on analyzing extreme scenarios that included nonlinear behavior with limit-cycle oscillations.This phase aimed to assess the safety margins of the reactor during extreme operating conditions and to identify potential improvements in safety features and system design.The third and final phase of the benchmark involved validating the developed models against stability measurements that were carried out a few weeks before and after the event, as well as performing uncertainty analyses.This phase aimed to evaluate the accuracy of the developed model and to quantify the uncertainties associated with the stability analysis.
Numerous studies have been conducted within the framework of this benchmark, as documented in references (Gajev et al., 2013;Balestra et al., 2013;Kozlowski et al., 2014;Dokhane et al., 2016;Dokhane et al., 2017;Dokhane et al., 2022).These studies have contributed significantly to our understanding of the mechanism behind the triggering of the instability event and the explanation of different observed trends.
The Paul Scherrer Institute has been an active participant in this benchmark and has made dedicated efforts to study this interesting and challenging event (Dokhane et al., 2016;Dokhane et al., 2017).In this context, a stability methodology was developed and validated, based on the three-dimensional reactor kinetics state-of-the-art code SIMULATE-3 K (S3K) (Dokhane et al., 2016).The results of these studies demonstrated the capability of S3K to not only simulate such complex behavior of the O2 stability event, but also to investigate possible scenarios that could have occurred, whether with or without operator intervention.For instance, it was also found that the core could have evolved into limitcycle oscillations if a stabilization of the feedwater flow and temperature had occurred just before the scram signal.Subsequently, an additional study was conducted, which showed the high sensitivity of the core dynamics to small changes in the core design.More specifically, it was found that the introduction of minor changes in the core through the reshuffling of the two hottest channels, and even reshuffling of less limiting channels, has small effects on the steady state core behavior.However, a significant change in the core dynamic behavior is observed through a drastic reduction of the power oscillation amplitude and the simulation demonstrates that, if such minor core design changes had been introduced in the beginning of cycle 24, the stability event of February 25, 1999 might have been potentially avoided (Dokhane et al., 2022).
Despite these efforts, no study so far, had been dedicated to uncertainty analysis, which is part of phase 3 of the benchmark.Therefore, the goal of this paper is to extend the PSI effort to cover uncertainty analysis, with a main focus on quantifying the effect of nuclear data uncertainties on the stability parameters, namely the decay ratio and resonance frequency, which to the authors' knowledge, is the first study of its kind.Such study of high relevance since it illustrates the impact of the propagation of uncertainties through the multiphysics coupling of neutronics, two-phase flow, and fuel heat transfer.The uncertainty propagation is carried out using the PSI methodology, called SHARK-X (Wieselquist et al., 2013), to propagate nuclear data uncertainties first in the assembly calculations using the CASMO5 code then followed by downstream static and dynamic calculations using SIMULATE and SUMULATE-3 K, respectively.This methodology has been successfully applied for both steady-state and transient calculations (Rochman et al., 2020;Dokhane et al., 2018;Dokhane et al., 2021).

Oskarshamn-2 stability tests
The Swedish Oskarshamn-2 reactor was designed by ABB ATOM, with an external-loop comprising four recirculation loops and a core with 444 fuel assemblies, belonging to four different bundle designs (Kozlowski, et al., 1999).It began operation in 1975 with a 1700 MW thermal power, then was uprated to 1800 MW in the early 1980 (Kozlowski, et al., 1999).As mentioned earlier, the OECD/NEA benchmark aims to verify the accuracy of the developed models by comparing their results with those of stability measurements obtained during two measurement campaigns.The first campaign took place on December 12, 1998, ten weeks before the stability event, while the second campaign was conducted on March 13, 1999, two weeks after the event.Both campaigns consisted of five tests conducted under different operating conditions.Table 1 summarizes the operating conditions and measured decay ratio (DR) values for the ten stability tests (Balestra et al., 2013); (Magedanz, 2010)).It is worth noting that, in (Dokhane et al., 2016) it was reported that some uncertainties in the operating conditions for REC.3 of the second campaign should be expected.Additionally, it is important to mention that only the DR was measured during the campaigns, and the values were obtained from an online stability monitor using Auto-regressive (AR) models of Average Power Range Monitor (APRM) signals.The reported DR values were determined based on a simple read by selecting the maximum value during the minutes of stable power and flow conditions (Kozlowski, 2011).However, simulations conducted by the plant after the measurements indicated that some of the measured DR values (Table 1) were likely higher than the reported values.As a result, the uncertainties associated with the measured DR values are unknown, but they are likely to be significant.

PSI analysis methodologies
This section provides a brief overview of the two methodologies used in this study.The first methodology focuses on stability analysis, while the second methodology deals with the propagation of nuclear data uncertainties in CASMO5 and downstream to SIMULATE3 and SIMU-LATE-3 K calculations.

BWR stability analysis methodology
The stability analysis methodology, developed at PSI in recent years, is based on the SIMULATE-3 K, a best-estimate LWR core dynamics code, which incorporates a full core coupling between neutronics and thermal hydraulics (Dokhane et al., 2013).The neutronics component uses a 3D time-dependent two-group nodal diffusion solver that integrates the frequency-transform method for time integration and the transverseintegrated method for spatial integration during the transient.The thermal hydraulics component uses a six-equation solver that considers mass, energy, and momentum conservation equations for both vapor and liquid (Grandi and Belblidia, 2011).
The coupling between thermal-hydraulics and neutronics is performed via fuel pin heat generation rate, which is directly determined from fission power.The thermal hydraulic module provides the neutronics with the appropriate hydraulic data, such as the moderator density and the average fuel temperature, that enter in the cross-section evaluation and thus permit calculations of the local feedback.The coupling is asynchronous, in the sense that the results from one module at the previous time step are used as input to the other module for the current time step (Grandi et al., 2011).
The stability analysis methodology follows several steps.First, lattice calculations are conducted using CASMO5, which is a two dimensional transport physics code based on the method of characteristics, to prepare the nuclear library (Rhodes et al., 2006).Subsequently, core static calculations are performed under the desired operating conditions to obtain initial core parameters for S3K.Finally, the S3K core model is obtained directly from the S3 restart file and supplemented by a hydraulic vessel model, which includes all peripheral system components.The S3K calculation involves two steps: a static calculation that initializes the core parameters and vessel hydraulic conditions at each investigated point (Dokhane et al., 2016), and a dynamic transient simulation that triggers a pressure pulse perturbation.Upon completion of the dynamic simulation, the decay ratio (DR) and resonance frequency (RF), the two stability parameters, are usually estimated from the total power time series (Dokhane et al., 2013).

Nuclear data uncertainty propagation methodology
The nuclear data uncertainty propagation is accomplished through the utilization of the in-house SHARK-X platform, which comprises a series of modules designed specifically to propagate ND uncertainties, provided in the form of covariance matrix libraries, in CASMO5 (C5) assembly calculations using a stochastic sampling or direct perturbation methods (Wieselquist et al., 2013).In the current study, the former approach is used and entails the random perturbation of nuclear data based on their joint probability distribution, derived from ENDF/BVII.1 19-group CM library and defined with multivariate normal distributions.This library consists of nuclear reaction pairs from 178 isotopes and six nuclear data, including fission, capture, elastic and inelastic scattering, fission spectrum, and nubar.It is noteworthy that, delayed neutron data are not perturbed in this study.It should be noted however that, as described in (Wieselquist et al., 2012), there is no guarantee samples drawn according to the multivariate normal distribution will be positive.Therefore, the final stage of the sampling algorithm in SHARK-X is to check for negative relative perturbations and sets any found to 0. This induces a small positive bias in the sample mean, but prevents code crashes, which could result from providing negative cross sections.Eventually, the SHARK-X methodology could be partly enhanced by adding the capability to use the multivariate lognormal distribution to reduce appearance of negative non-physical values.However, so far, no such updates were done for the SHARK-X tool, and such developments can be subject to future studies at PSI.For each sample, CASMO-5 calculations are performed for all assemblies in the core, using a randomly perturbed nuclear data, followed by downstream SIMULATE3 and SIMULATE-3 K calculations and analyses.In this study, CASMO5 calculations are performed using version 2.03 and ENDF-B/VII.0 (E7.0) 586-group library.Once all samples and associated calculations are completed, statistics are performed on the mean (first moment) and variance (second moment) of the SIMULATE 3 and S3K output distributions.

Results
The results presented here are obtained through the analysis of the ten stability tests, described in section 2 (Table 1), which were conducted before and after the stability event of February 1999 in Oskarshamn-2.The analysis was performed using SHARK-X, which employs a stochastic sampling technique to perturb the nuclear data set.For each sample, a corresponding C5 calculation was performed using ENDF/B-VII.0,586-group library, followed by downstream S3 and S3K calculations.In total, 300 samples were generated and the results are presented in the following sections, which include quantities of interest for both steady-state and transient conditions.

Steady-state
The steady-state results, calculated by S3K, are presented in terms of k-eff, axial and radial power peaking factors (Fax and Frad) at assembly level, maximum nodal (MN) void content (MNVoid), and maximum nodal fuel temperature (MNTFU), in Table 2 and Figs. 1 and 2. Table 2 summarizes the statistics for two representative stability tests of campaign 1 and 2 (Test No.1), while Figs. 1 and 2 represents the uncertainty or the Relative Standard Deviation (Rel.Std(%)) and the associated statistical error of the estimated uncertainty (SE Rel.Std.) in terms of keff, Frad, Fax, MN void content, and MN fuel temperature for all the stability tests of campaign 1 and 2, respectively.Fig. 3, represents the statistical uncertainty of both mean and standard deviation of k-eff and also the 95 % confidence interval as a function of sample size for Test No.1 of the first stability campaign (12.12.1998).Clearly, Fig. 3 illustrates the convergence of results using 300 samples.As can be observed in Table 2, the bias between the average value (over 300 samples) and the reference value (unperturbed case) is quite small for all parameters, e.g. a maximum of 68 pcm for k-eff.
Concerning uncertainty in k-eff, it is about 550 pcm (0.54 %) for all tests, which is of the same range in a study performed for a Swiss BWR (Rochman et al., 2020).The uncertainty (Rel.Std(%)) in terms of Frad, Fax, MN void content, and MN fuel temperature is found to have a maximum of 0.37 %, 2.79 %, 0.08 %, and 1.45 %, respectively (Fig. 1).It should be noted that, the uncertainties in the axial power peaking factor and maximum nodal fuel temperature of stability tests of the second campaign (13.03.99) are systematically lower than those of the first campaign (12.12.98).This trend could be related to burnup effect.
The statistical error of the estimated uncertainty is estimated using a bootstrap replacement approach, with 100 replicates (Efron and Tibshirani, 1986).As can be seen in Fig. 2, the statistical error for all cases is around 4 %, as expected by the chi-squared approximation (1/ ̅̅̅̅̅̅̅̅̅̅̅̅̅ ̅ 2*299 √ ≈ 4 %), where it is assumed that the target parameter follows the normal distribution.
In Fig. 4, the 2D map of steady-state results in terms of the average fuel assembly power (in MW) and the associated uncertainty, based on 300 samples, are presented for all the 444 fuel assemblies in the core.As can be seen, as expected, the fuel assembly average power for test 1 of first campaign (REC1-121298) is systematically slightly higher than that of test 1 of campaign 2 (REC1-130399) because the core total power at the former is slightly higher than that of the latter (see Table 1).Concerning the associated uncertainty, i.e. in the FA power, both tests show similar 2D uncertainty distributions.More specifically, the maximum power uncertainties are registered at the most peripheral FAs of the core with a maximum value about 1.7 %, then followed by the second peripheral inner ring with a maximum FA power uncertainty about 0.8 % and finally the inner region of the core shows a lower uncertainty, lower than 0.4 %.One possible reason behind the high power uncertainty at the peripheral FAs is the fact that the FA power in this region is very low and therefore become more sensitive to any change.It should be noted that, there is a clear increase of power uncertainty for FAs surrounding the inserted control rods in the core (five CRs inserted in REC1-121298 and one in REC1-130399).This is again could be by the lower power level at these FAs.

Transient
Using S3K, the transient calculations for 300 samples were performed with a duration of 15 s, after applying a pressure pulse perturbation, as explained in section 3.1.The transient results are presented in terms of time-dependent total power and the associated uncertainty.More importantly, the results of the stability characteristics of the core are presented in terms of the decay ratio and resonance frequency and their associated uncertainties.

Total power
Fig. 5 shows an example, corresponding to Test No. 1 (REC 1) of 12.12.98 and 13.03.99 stability campaigns, of the time evolution of the total reactor power for the 300 samples along with the associated uncertainty and its 95 % confidence interval estimated using a bootstrap replacement approach with 100 replicates.The red area in the figure (left-hand side figures) represent the spread of the results of the 300 samples due to the uncertainties of the nuclear data.The two dashed red and green curves represent the two extreme cases, among the 300 samples, corresponding to those with maximum and minimum power peak.The black curve corresponds to the unperturbed (reference) case.
As can be seen in Fig. 5, the reactor core is stable for all samples but with different peak levels and decaying rate of the oscillations.The uncertainty in time-dependent power is found to be oscillating too with a maximum value about 0.18 % and with peaks obtained around time frames when power is at maximum and minimum values (Fig. 5 (Right)).The reason behind is mostly due to the shift in amplitude and time of the peak creating larger spread in power these time frames.The uncertainty decreases with time and its value at the end of the simulation becomes almost similar to that at the beginning of the simulation.

Stability parameters
The uncertainty quantification results of the stability parameters due to nuclear data based on 300 samples are presented in Fig. 6 (average value), 7 (relative standard deviation, Rel.Std (%)) and 8 (statistical error of the estimated uncertainty (SE Rel.Std.(%)).As can be observed in Fig. 6, the discrepancy between the average value of both DR and RF of each stability test, based on 300 samples, and those of the reference case (unperturbed cases) is almost negligible.
Regarding the corresponding uncertainty, as can be seen in Fig. 7, the uncertainty in DR due to the perturbation of all nuclear data could reach a maximum of 8.50 %, while for RF the uncertainty is lower and does not exceed 2 %.It should be noted that, there is almost a systematic trend in which the uncertainty in DR is lower in tests of campaign 1 (12.12.98) compared to those of campaign 2 (13.03.99), which could be due to burnup effect.Therefore, further investigation is needed to unveil the main reasons behind such trend.However, this is beyond the current scope and will be performed in a future work.
The statistical error of the estimated uncertainty (SE Rel.Std.) in the stability parameters is presented in Fig. 8 and it is estimated using the bootstrap replacement approach, with 100 replicates.As can be seen, the statistical error for all cases is around 4 %, with maximum and minimum values of 5.28 % and 3.09 %, respectively.
The above transient results are presented in terms of global decay ratio and resonance frequency, based on total core power.In order to have more resolution on the distribution of the stability parameter         A. Dokhane et al. values and their associated uncertainty within the core, an in-house Matlab code has been used with the main purpose to collect first the time-dependent power of each fuel assembly for 300 samples, then to evaluate the corresponding stability parameters.Fig. 9 shows the 2D core map of the mean decay ratio value and the associated uncertainty, based on 300 samples, for the first tests of the two stability campaigns (REC1-121298 and REC130399).Concerning the DR distribution, for both tests, the radial variation of the DR is very small (between 0.42 and 0.44 for REC1-121298; and between 0.43 and 0.45 for REC1-130399) with a systematic slightly higher values at the center of the core.
Regarding the DR uncertainty results, it is clearly seen that and in agreement with Fig. 7, the DR uncertainty of the individual FAs in test of campaign 1 is lower (between 6.8 and 7.7 %) than that of campaign 2 (between 7.5 and 8.4 %).Again this could be related to burnup effect, however further investigations are needed to unveil the main reasons behind such trend.However, this is beyond the current scope and will be performed in a future work.
Fig. 10 shows the 2D core map of the mean resonance frequency value and the associated uncertainty, based on 300 samples, for the first tests of the two stability campaigns.Concerning the RF distribution, for both tests, the radial variation of the RF is even very small, compared to that of DR, between 0.535 and 0.540 Hz with a homogenous value distribution throughout the core.
Regarding the RF uncertainty results, it is clearly seen that the distribution is almost homogenous with an maximum uncertainty about 1.9 % and 1.8 % for REC1 of campaign 1 and 2, respectively.
The calculated DR results, including uncertainty (2σ) due to nuclear data, are presented along the measured values in Fig. 11.In addition, the minimum and maximum of the calculated value are presented also.As can be seen, for the majority of the five tests of 12.12.98stability campaign, the measured DR value is within the uncertainty of the calculated value due to nuclear data.However, for the 13.03.99 stability campaign, the measured value is within the uncertainty interval of the calculation for three tests (1, 4 and 5), while the remaining two tests (2 and 3) show a considerable discrepancy with lower calculated values compared to measurements.This might be due to the potential high uncertainty in the manner how the measurement of DR was performed and also due to uncertainties in the operating conditions for REC3, as explained in section 2.

Conclusions
The current research is part of the OECD/NEA benchmark on the Oskarshamn-2 stability event, which took place in February 1999 during cycle 24.The primary objective of the study is to evaluate the impact of uncertainties in nuclear data on both steady-state and transient core parameters during stability tests conducted before and after the stability event.To achieve this goal, the in-house SHARK-X platform is utilized, which propagates nuclear data uncertainties, provided in the form of covariance matrices, in 2-D lattice calculations using CASMO5.Downstream static and dynamic calculations using SIMULATE3 and SIMU-LATE-3 K are also performed.
In terms of steady-state, the uncertainty is assessed in terms of k-eff, radial and axial power peaking factors, void fraction, and average fuel temperature.For the transient case, the study provides for the first time an estimation of uncertainties in the decay ratio and resonance frequency.The findings reveal that nuclear data uncertainties have a nonnegligible impact on most of the analyzed parameters.More specifically, the uncertainty in k-eff, radial and axial power peaking factors, nodal void fraction peak, and nodal fuel temperature peak for steady-state could reach about 0.5 %, 0.4 %, 2.8 %, 0.1 %, and 1.5 %, respectively, On the other hand, the uncertainty in stability parameters such as the decay ratio and resonance frequency can be as high as 9 % and 2 %, respectively, with a tendency to have higher uncertainties for higher burnup.
In terms of comparing the calculation results with measurements, it A. Dokhane et al. is found that the calculated decay ratio falls within the uncertainty range for the first stability campaign.However, for the second campaign, there are significant discrepancies in two out of five tests, with the calculated value being lower than the measured value.This discrepancy could potentially be attributed to the high uncertainty in the measurement reporting of the DR values.Overall, the study highlights the impact of accounting for nuclear data uncertainties in the evaluation of parameters for steady-state and, for the first time, transient conditions for a realistic BWR at stable core conditions.It would be interesting to assess such impact when the core is in unstable conditions, as the case during which the stability event of February 1999 occurred.The impact could be significant both quantitatively and qualitatively since the core is found to be very sensitive to small changes during that event.This investigation is planned to be performed in a future work.This study highlights also the need for further research to explain certain observed trends.Specifically, there is a need to investigate the reasons behind the higher uncertainties in power obtained for the peripheral FAs in the core.Additionally, the systematic lower uncertainties in the power axial peaking factor and the maximum nodal fuel temperature of tests of the second campaign (13.03.99) compared to those of the first campaign (12.12.98).Furthermore, the systematic higher uncertainties in DR for tests of the second campaign compared to those of the first campaign.These trends could be related to burnup effect, which will also be explored in a future investigation.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.A. Dokhane et al.

Fig. 3 .
Fig. 3. k-eff and associated STD convergence and confidence Interval for test No.1 of 12.12.98.

Fig. 11 .
Fig. 11.Decay Ratio results with uncertainty and comparison to measured values.

Table 1
Operating Conditions and Measured DR of Stability Tests.

Table 2
Stability testsNo. 1 of 12.12.98 and 13.03.99Campaigns: UQ of steady-state parameters based on 300 samples.