© D.A. Rochman et al., Published by EDP Sciences, 2023

1. Introduction

High-level radioactive nuclear waste is one of the most important public concerns regarding the safety and viability of the nuclear industry and its acceptance. Questions are regularly raised about their safe storage (intermediate- and long-term disposal), handling, transport, or reprocessing. Regardless of its rational danger to the environment, it is therefore our responsibility to address this apprehension and to demonstrate that the nuclear community is making sensible efforts to keep such waste under control. Additionally, the economical aspect of waste treatment can be considered for both the industry and the consumer, leading to a virtual limit on practical solutions.

An additional aspect of the majority of the high-level radioactive waste, contained in Spent Nuclear Fuel (or SNF), is that their characteristics cannot be systematically measured. It is currently estimated that at the beginning of the 2020s, there are about 285 000 SNF assemblies in the United States and almost as much in Europe, leading to a total of more than 500 000 SNF assemblies, without accounting for Russia and Asia [1, 2]. With such a large increasing number of SNF assemblies worldwide, it is not realistic to base our knowledge on systematic measurements, but one should rather rely on our prediction capabilities through modelling and simulations, coupled with a limited number of measurements. The importance of adequate calculations, with a reliable estimation of uncertainties and biases, is therefore of paramount importance.

When characterizing the SNF, the origin of several observables is the SNF nuclide concentrations (also called nuclide inventory): short as well as long-lived actinides and fission products such as ²³⁵U, ²³⁹Pu or ¹³⁷Cs. A synonym expression is “source terms”, and their knowledge affects the criticality aspect of a storage facility, their design (e.g. with the SNF decay heat), and the protection of the environment. Such nuclide concentrations vary as a function of the assembly design, its irradiation conditions (spatially and timely), and cooling time. They are nevertheless not easily observed, and it is preferred to use different quantities for the classification of SNF, being closer to macroscopic or integral quantities, such as initial fissileenrichment (being a single value not changing as a function of irradiation time), average assembly burnup at the end of life, assembly average decay heat, reactor type, mass, and cooling time.

The knowledge of the SNF has also an economical aspect. For instance, countries having selected long-term solutions such as deep geological repositories need to optimize the cost of such facilities with their safety aspect [3]. In simple terms, the amount of SNF to be stored directly influences the cost of such a solution. And the number of SNF assemblies as well as their characteristics defines the number of required canisters to be stored underground (a canister can contain for instance 8 or 12 SNF; dual-purpose canisters can contain more SNF). Decay heat is one such characteristic, each canister having a hard limit for the total decay heat it (as well as the hosting rock or backfilling) can handle [4]. In a recent Nuclear Energy Agency Joint Workshop, it was recognized that decay heat suffers from a lack of realistic uncertainties and open questions related to conservatism, biases and required margins still remain [5]. It becomes then clear that the knowledge of the SNF decay heat directly influences the number of required canisters: it is expected that an average variation of 1% for the decay heat of a total SNF stockpile modifies the total amount of canisters also by ≈1%. Therefore, a “suspected” bias in the SNF decay heat predicted by best-estimate codes of a few percent, possibly leading to penalty factors, can have an important financial impact.

The present study aims at helping characterize the content of SNF (nuclide concentrations) as well as one of its integral observable parameters, the decay heat. By estimating biases and uncertainties, this study provides a sound approach to help to make informed decisions for the safe and economical handling of SNF. By doing so, this work is fulfilling one of the goals of the European EURAD Work Package 8 project, dedicated to defining the best practice for SNF characterization [6].

In the following, the term “decay heat” is used to express the decay power (as a function of cooling time) for SNF assemblies or samples. It is expressed in Watt (or Watt per ton) and should therefore be called a “power”, and not a “heat” (describing a transfer of energy, and expressed in calories). The term “heat” was historically used, as it was noticed that “When a nuclear reactor is shut down, following some period of operation, there are various nuclear species and processes that remain which are capable of generating heat” [7]. It is still much in use today, and we will therefore continue with the term “decay heat”, understanding that it means “decay power”.

1.1. Overview of SNF characteristics

As mentioned, the present work fits within the intention of demonstrating our capabilities in addressing part of the SNF issues, related to its knowledge (physical content of the used fuel) through modelling, prediction and validation. Although many physical characteristics of the SNF are relevant (e.g. thermo-mechanical and chemical properties), we focus on the neutronic aspect of the SNF characterization, such as its nuclide inventory (also called nuclide concentration or source terms), and quantities helping to model and classify various SNFs: irradiation history, enrichment, fuel type, cooling time, fuel utilization (or burnup) and release energy from radioactive decay (or decay heat). Examples of such quantities are presented in Figures 1 and 2. These figures present histograms of specific quantities, extracted from realistic nuclear power plants, operated over a few decades. These quantities are calculated and are therefore not absolutely certain. They are nevertheless well representative of existing SNF that the civil society will have to deal with. One can already observe a classification of SNF in terms of plants: Boiling and Pressurized Water Reactors (or BWR and PWR, respectively) and fuel types (uranium oxide or UO₂, and mixed oxide, or MOX). Such a list is not exhaustive, as there are other types of fuels and plants, but the present research will focus on this set of cases.

Fig. 1. Left: example of assembly burnup at the end of life, based on realistic irradiation history for a number of existing PWR and BWR power plants. Right: same but for average fissile enrichments. Vertical coloured segments indicate the burnup of available PIE samples (long segments: considered in this work; short segments: extracted from SFCOMPO [8]).

Fig. 2. Right: example of assembly decay heat at 10 years cooling time, based on realistic irradiation history for a number of existing PWR and BWR plants. Left: same for the assembly average 235U and 239Pu concentrations.

In Figure 1 left, the assembly-averaged burnup, in terms of MWd/kg¹, at the end of life (EOL) (i.e. when the assembly is finally discharged, after a number of reactor cycles) is presented for both PWR and BWR. Although this quantity is not the only key criterion to characterize SNFs, it is one of the most used (local peak burnup can also be of relevance, and is often linked to the average value). This figure already shows the diversity in the SNF cases: the burnup values span from a few MWd/kg (for the partially burned assemblies) up to 65 MWd/kg (for the fully burned assemblies). These differences in burnup rates will strongly influence other SNF quantities, such as nuclide concentrations and decay heat and lead to specific requirements for prediction capabilities. Such prediction capabilities can be tested by simple comparison between calculated and measured quantities. It is consequently of prime importance to access experimental data, as proposed in the SFCOMPO database for nuclide concentrations [8]. Such database includes a number of cases (called PIE data, for Post Irradiation Examination) and their spread as a function of the sample burnup (in general being a small part of a considered SNF, but there is a limited number of cases where a full assembly was also dissolved and nuclide concentrations were measured) is also presented in this figure with small vertical segments (differentiated in colours for PWR and BWR). As observed, available PIE samples from the SFCOMPO database cover well the spread of SNF in terms of burnup (although one can question the representativity of a local sample with respect to a full assembly).

The right part of Figure 1 presents the average fissile concentrations for the same realistic (fresh) assemblies. For UO₂ cases, the fissile content is the concentration in ²³⁵U, and for MOX, it is the sum of ²³⁵U and ^239, 241Pu. This figure also indicates that this characteristic of SNF is relatively diverse, spreading from low concentrations (less than 1%, often corresponding to assemblies used in the first reactor cycle), up to 5 or 6%. Also indicated in the figure is the occurrence of the SFCOMPO samples, separated for UO₂ and MOX fuel types (short vertical segments). If the UO₂ fuel case seems to overlap well between existing SNF and measured samples, it is not the case for the MOX fuel. MOX assemblies are by design more complex than their UO₂ counterpart: instead of a single enrichment for the assembly regions under high neutron flux (i.e except for the bottom and top of assemblies), the MOX assemblies often have three distinct enrichment zones, for instance including low, medium and high ²³⁹Pu contents. It becomes then difficult, and in practice not realistic, to find SFCOMPO MOX samples representing the three enrichment zones, and therefore being representative for a specific MOX assembly type. A similar comment can be made for BWR cases, being axially and radially more heterogeneous than their PWR counterparts. Another difficulty in the case of MOX, as shown in the figure, is the small number of SFCOMPO samples². There is therefore a mismatch between realistic MOX SNF characteristics and the available measurements from the SFCOMPO database, used to perform validation steps.

Another view on the existing SNF is presented in Figure 2. Two other quantities of interest for the SNF characterization are given in the form of histograms: calculated actinide concentrations (²³⁹Pu and ²³⁵U), as well as fuel assembly decay heat for a cooling time of 10 years after the EOL. Nuclide concentrations are also an important key factor for the characterization of SNF, as they directly influence decay heat values, as well as k_eff or reactivity values [9]. As observed in the figure, there is an adequate overlap between the SFCOMPO samples and the realistic assemblies, for both nuclides. Similarly, the cases considered in this study also conveniently overlap with the concentrations from the SNF. One can observe that the spread of the average ²³⁹Pu and ²³⁵U concentrations from the SNF is relatively large for both nuclides, resulting from the fuel type, initial enrichment and assembly burnup.

If the nuclide concentrations generally are not part of the “integral” quantities used to characterize the SNF, the decay heat is. An example of assembly average decay heat is presented in Figure 2 right. Values at the cooling time of 10 years after EOL are presented. Such values are naturally calculated based on realistic cases. Again, a large range of decay heat rates can be obtained; such variety will also be present at the encapsulation time (when SNF are packed in their final disposal canisters), rendering the estimation of low and high decay heat values relevant. The vertical bars in Figure 2 right indicate the measured values from the limited existing dataset (see Sect. 4.4.2 for further details): it is apparent that the measurements do not cover high decay heat values, therefore limiting the possibility of direct validation.

In conclusion, the aspects presented in Figures 1 and 2 indicate the complexity of SNF characterization. Some perspectives were not touched on, which are still strongly linked to source terms: criticality and radiological protection. A comprehensive study for SNF characterization is out of the present goal: as presented in the next section, the present work is limited to the estimation of source terms, with a direct relationship to the decay heat.

1.2. Main goal on the present work

The main goal of the work is to provide recommendations for nuclide inventory and decay heat predictions: for best practices in simulation options, as well as for expected uncertainties and biases.

As the decay heat is a key quantity for the transport, storage of SNF in dedicated casks and canisters, it is crucial to estimate the bias and uncertainty of a specific calculation scheme. For the nuclide inventory, as it is at the origin of all integral (or macroscopic) SNF quantities, such as neutron and gamma-ray emission, as well as reactivity calculations, it is at least as important as the decay heat itself. The analysis will be divided into two parts:

1. presentation of C/E − 1 values and determination of biases and uncertainties for nuclide concentrations, see Sections 4.1–4.3 (with C being the calculated quantity and E the measured one).

2. Similar study for decay heat, see Section 4.4.

An estimation of biases for nuclide concentrations and decay heat is a necessity for quantifying the confidence in calculations. It is therefore the first quantity presented. Then simulation options can play a very important role in the reachable accuracy of the calculated quantities. A few examples can be mentioned: Monte Carlo versus deterministic neutron transport, method of solving the burnup equations (or Bateman equations), multi-group or (quasi-) continuous energy nuclear data, nuclear data uncertainty, details of the power history during irradiation, 2D versus 3D modelling, single- or multi-assembly simulations, etc. In simple terms, it means that if a decay heat value is calculated for a specific assembly, at a specific cooling time, and obtained with a specific calculation tool and input (including a specific level of detail), we are able to justify an uncertainty (being equivalent to one standard deviation) and a bias (defined as the ratio of C/E).

In this paper, the term “bias” is used to express the difference between measured and calculated values. Such a definition assumes that the measurement leads to the true value, which might not be correct. We nevertheless use this term, understanding that it erroneously mixes measured and true values.

In the following, a number of cases are considered (PIE samples and decay heat calorimetric measurements, as presented in Figs. 1 and 2) and analyses are carried out for important cases. Both uncertainties and biases will be proposed in specific cases, together with justifications from the present analysis and other published work. As much as possible, we will follow the definitions provided in reference [10] for the uncertainties and related quantities; the usual definition of uncertainty in this work corresponds to one standard deviation.

1.3. Link with other current and past activities

An important activity is the compilation and analysis efforts within the Working Party on Nuclear Criticality Safety (WPNCS) of isotopic compositions from SNF within the SFCOMPO database and its technical review group of the Nuclear Energy Agency (NEA) [8, 11]. The present work would not have been possible without this international database, especially the results presented in Figure 4. The SFCOMPO database, continually improved with the addition of new PIE samples, is at the start of almost all code validations based on publicly available data.

Fig. 3. Histograms of C/E − 1 for all nuclides (left), and for a selection of nuclides (center and right).

Fig. 4. Interquartile ranges for the C/E − 1 nuclide concentrations, considering a total of more than 12 000 measured concentrations. The blue colour is given to important nuclides. See Tables 3 and 4 for numerical values. The limits of plotted boxes are equal to the median ± IQR/2, as presented in Tables 3 and 4.

Among the past activities, the WPNCS published 2011 a state-of-the-art report concerning “Spent Nuclear Fuel Assay Data for Isotopic Validation” [12]. It recognized the importance of measured spent nuclear fuel nuclide compositions for a number of SNF activities and highlights the importance of the SFCOMPO database. If it acknowledges that the validation of codes and nuclear data highly depends on such data, it also points out the necessity of additional effort in assessing biases and uncertainties for applications such as burnup credit. On this aspect, the present work takes into account such recommendations and extends them to another linked activity: the estimation of SNF decay heat.

The Electric Power Research Institute (EPRI) recently published a dedicated study on the state of the art regarding the estimation of SNF decay heat [13]. This report will be commented on in Section 4.4.2, and similar conclusions and recommendations will be presented in the present work.

Besides the EURAD project, three other international activities are connected to the present work. The first one is the WPNCS Subgroup 10 (SG10), on “Nuclear Data Uncertainties Quantification on Spent Fuel Inventory”, dedicated to the study of the ARIANE GU3 sample (also studied in the present work). Even if the SG10 has not yet concluded the work at the time of the present work, it indicates the importance of nuclear data and their uncertainties in the validation of PIE samples, and consequently for decay heat estimation. The second activity is in direct link with the present work, being the WPNCS Subgroup 12 (SG12), dedicated to “Spent nuclear fuel decay heat: assessing the confidence level in experimental and computational estimations (SNF-DH)”. The goal of this working group is to “evaluate” the decay heat of existing SNF, a goal which goes beyond the current study, but clearly presents overlaps. The last international activity to be mentioned is the Coordinated Research Project (CRP) of the International Atomic Energy Agency, on “Spent Fuel Characterisation”. This project tackles different aspects of the SNF characterizations, including decay heat and nuclide concentrations.

Finally, many national projects provide enormous efforts for the characterization of SNF, in relation to the neutronic aspects, or not. Needless to mention the activities in Sweden (SKB), Finland (Posiva), Switzerland (Nagra), France (ANDRA), Spain (ENRESA) and Belgium (SCK-CEN) cite only some of the European organizations.

2. Studied cases

A selection of SNF samples and assemblies was performed by the different partners of this work. It concerns PIE samples for nuclide inventory and results of decay heat measurements. The list can be summarized as follows:

• 17 PIE samples were used to calculate nuclide compositions, uncertainties, and biases, for 4 ARIANE samples (GU1, GU3, BM1, BM3), 8 ENRESA samples (not available in SFCOMPO), the U1 PROTEUS sample (also not available in SFCOMPO), 2 Takahama samples (SF95-4 and SF95-5), and 2 Gundremmingen samples (B23-A1-I2680 and C5-B23-K2680), with the addition of 2 computational cases (S1.PWR and krko.PWR [14, 15]),

• 271 calorimetric measurements from CLAB (labelled herein as CLAB-2006), GE-Morris, HEDL facilities and the “SKB-Vattenfall” blind test were analyzed for the assembly average decay heat [16–19],

• and the decay heat from the 17 (not measured) PIE samples.

Detailed studies were independently published for a number of cases, see references [20–32, 50]; summaries are presented in Tables 1 and 2. Regarding sensitivity studies presented in Table 1, considerable efforts have been made to address the effect of the irradiation history, as-built manufacturing data and effect of various calculation model approximation on the SNF characterization. This analysis is not yet completed and conclusions are therefore not be reported in the present paper. The results, together with other data in the literature, are used to perform a statistical analysis presented in the following sections. Based on comparison between calculated and measured values, different quantities will be recommended, such as simulation methods, biases and uncertainties for SNF decay heat and nuclide concentrations.

3. Analysis method

The analysis methods are based on individual simulations performed by the various institutes. Each partner decided on the cases to be studied, the source of information, and the simulation scheme. On one side such an approach makes the comparison of results and understanding of possible differences complicated, but on the other side, it reflects what currently happens in various institutes.

For the PIE samples, identical cases were selected by different institutes, such as the BWR ENRESA samples, SF95-5 from the Takahama-3 PWR reactor or the I2680 (A1-B23) sample from the BWR Gudremmingen reactor. In addition to these explicit choices by the EURAD partners, a number of external institutes also provided additional results, such as for the GU3 and ENRESA samples [23, 33]. As mentioned, in the case of the ARIANE GU3 sample, the current activity from the OECD WPNCS Subgroup 10 will bring additional calculation cases for this sample, based on different irradiation history assumptions [34].

Two samples (BM1 and A1-B23) were analyzed based on different assumptions, such as modelling a single pincell, a two-dimensional assembly (with reflective boundaries), and a full three-dimensional model, taking into account the heterogeneous environment [24, 50]. For uncertainties due to nuclear data, different geometrical assumptions were also tested in the case of the GU1 sample [22].

3.1. Codes and libraries

As mentioned, a variety of codes and simulation methods were applied for the PIE samples and decay heat cases. A short description can be found below.

• CASMO5, SIMULATE and SNF. The Studsvik codes were used to model a number of samples (see Tab. 1), and details can be found in references [20, 22, 23, 35, 50]. When possible, the assumption of reflective boundaries was avoided, and a number of different nuclear data libraries were also tested. Concerning the propagation of nuclear data uncertainties, the SHARK-X in-house tool was used with CASMO5 and three different libraries were considered: ENDF/B-VIII.0, JENDL-4.0 and JEFF-3.3. Details on the uncertainty propagation method can be found in the mentioned references and references [36, 37].

• Serpent 2 was used to model the Gudremminden sample as indicated in Table 1. The serpent is a three-dimensional continuous-energy Monte Carlo particle transport code developed at VTT in 2004 [38]. Serpent uses a matrix exponential solution based on the Chebyshev Rational Approximation Method (CRAM) for solving the Bateman equations [39]. For particle tracking, a combination of conventional ray-tracing-based surface tracking and the rejection sampling-based delta-tracking method is used. Neutron interactions are read from continuous-energy ACE format cross-section libraries similar to MCNP.

• The EVOLCODE system [40] developed by CIEMAT was used to simulate the SF95-5 pellet of the Takahama-3 reactor by means of a single pincell model with reflecting boundaries. The propagation of nuclear data uncertainties was done following the methodology described in reference [41]. The JEFF-3.3 library was used as a reference [42]. Sensitivity analyses have been done to explore the impact of the nuclear data library (ENDF/B-VIII.0 [43]), solver (MCNP/CINDER), resolution model (spectrum description and EVOLCODE irradiation steps) and geometry model (cylinder pincell instead of squared pincell); these results are presented in reference [31].

• The SCALE/Polaris code from the SCALE Code System, version 6.2.4, was used for modelling the nuclide inventiry of the BWR samples (details can be found in Ref. [33]). A two-dimensional single assembly was modelled with the exact irradiation history including the downtimes between cycles. The SCALE v7-252 multigroup cross-section library based on ENDF-B/VII.1 library was used. The moderator density was derived from the segment void fraction. Nuclide concentrations at the measurement time were obtained from an ORIGEN decay calculation using the calculated Polaris inventory at the end of irradiation. The sample power at each irradiation step was calculated from the nodal burnup, adjusted to the sample burnup (from the weighted average of Nd values) and divided by the irradiation time of each burnup step. Some of the six samples were analyzed together (ENRESA-1 and 2, ENRESA-3 and 7) since they are very close and share the same irradiation history (the nodal burnup was adjusted to the average burnup of both samples).

• In addition to the previous use of Polaris, the Polaris and Sampler super-sequence was used for both the modelling of selected PIE samples (for the validation with nuclide inventory), fuel assemblies (for decay heat validation), and uncertainty propagation. The generic modelling assumptions of the application of Polaris and Sampler in this study are described in reference [44]. A sampler is used for the uncertainty propagation, generating and running hundreds of input files of Polaris and analyzing the outputs. The cross-section covariances are primarily based on the ENDF/B-VII.1 nuclear data library. The nuclear data samples are available in SCALE, and the design and irradiation uncertainties are assumed, similar to the assumptions from the previous reference. The outputs are uncertainties of the calculated values.

• For descriptions of ALEPH, MURE and TRITON, see references [45–47].

The present study is complemented by results from the literature, allowing it to cover a wide range of C/E values.

3.2. Uncertainty methods and variables

Uncertainties coming from calculations can have various origins. They contribute to the understanding of our simulation capabilities and can eventually explain biases. For instance, if the calculated uncertainties for specific nuclide concentrations are larger than the observed bias, one can eventually assume that the differences between C and E come from limited knowledge of some model inputs or assumptions. On the contrary, small calculated uncertainties leave part of the biases unexplained.

Different origins of calculated uncertainties were considered, and they can basically be categorized into the following origins:

• irradiation history. Quantities such as fuel and coolant temperatures, void fraction, segment or assembly burnup and power history over the irradiation period(s) are considered assumptions during the sample irradiation. They are seldom measured and are usually provided as part of the sample irradiation conditions. They can for instance be calculated from a core simulator. These quantities are therefore not perfectly known and uncertainties are usually assumed for each of them.

• Manufacturing tolerances. Like the irradiation history, some quantities such as the sample enrichment, impurity concentrations, dimensions or densities are also assumed. Although uncertainties can be difficult to assess (due to their proprietary aspects), a number of assumptions (based on expert’s judgment) can be made.

• Model assumptions. As the description of the irradiation condition is often limited to the segment or assembly of interest, a number of assumptions concerning the sample environment: position of control rods or surrounding assemblies. In some cases, a limited amount of information is known (such as average burnup of neighbouring assemblies, and position in the core in different cycles), but in other cases no information is available. Assumptions such as reflective boundaries are therefore necessary and can influence the operating temperatures and the neutron spectrum and hence the calculated nuclide concentrations. Their impact is more difficult to assess but is certainly not negligible.

• Simulation methods. Choices for solving the neutron transport equation, depletion equations, possible rod relocation from one assembly to another, or burnup-induced changes (swelling, diffusion, etc.) can impact calculated quantities. Their impact is also difficult to quantify and cannot be neglected.

• Nuclear data. They correspond mainly to cross-sections, fission yields (as well as energy and angular distributions and multiplicities of outgoing particles) and decay data (branching ratios, decay modes, half-life and particle emissions). They are considered separately from other “input quantities” due to their impact and their well-defined physical origins which are independent to the considered samples. The covariance information is usually included in the nuclear data libraries used by the simulation codes (as for SCALE modules Polaris and TRITON), or are directly obtained from the original ENDF-6 nuclear data libraries (as for CASMO5, SNF, Serpent, CINDER, MURE and EVOLCODE).

Different solutions exist to perform so-called “uncertainty quantification”, or “uncertainty propagation”. One method is to first establish a trusted simulation model, and then repeat calculations based on specific inputs or quantities variations (for instance varying the ²³⁵U enrichment). The spread of the calculated quantities can then be quantified by a standard deviation. Another possibility is to repeat the same simulation a certain number of times, each time changing different simulation methods or model assumptions. A different method is based on the use of sensitivity vectors (relative variation of an output quantity with respect to the variation of individual specific input parameters), with the assumption of linear behaviour, in combination with covariance matrices (so-called sandwich method; see for instance Ref. [48] for the definition and application of the sandwich method with nuclear data). Finally, the simulation of the same sample irradiation can be performed by different institutes, based on different tools and understanding of the irradiation.

None of these solutions is comprehensive; however, they can contribute to the assessment of the total calculated uncertainties. For the samples studied here, all the above solutions were applied, based on the preferences of each participant (it is assumed that the different methods lead to results which can be compared; if extreme quantities are found, such as for uncertainties, the origin of the discrepancies will be indicated).

4. Results

As presented in the introduction, one of the main goals of the present work is to perform a dedicated statistical analysis to provide reliable estimates of the uncertainties and biases for both nuclide concentrations and decay heat. To reach this goal, we are taking advantage of the results for the samples presented in Table 1, as well as of a number of results from the literature, namely the C/E − 1 values from references [44, 49, 51, 52]. In these references, a large number of C/E − 1 values are given for a variety of PWR and BWR samples, including UO₂ fuel, MOX fuel and a large range of sample burnup. In total, more than 12 000 C/E − 1 values were extracted, including the inventory of many nuclides from actinides to fission products. Results are presented in a global manner in Figure 3 (left) as a histogram for all studied nuclides and all samples. This is a rather simple way of presenting the results, implying that all C/E − 1 values can be compared together. A more discriminant way is presented in Figure 3 (centre and right) by selecting specific nuclides. Independently of its relevance, such a histogram allows visualization that the vast majority of the values are close to zero. The median is +0.50% and the median absolute deviation is 5.4%. The mean is +2.9% and the standard deviation is 29.1%.

These quantities indicate that without using the information on the sample, i.e. its irradiation conditions and the nuclide of interest, most C/E − 1 are within a range between −26.1 and 32.0% when considering the mean and standard deviation (91% of the cases are within this range). Given that all calculated values are correlated (such as irradiation conditions, nuclear data, and burnup normalization to ¹⁴⁸Nd), these statistical values still quantify our capability to reproduce measured nuclide concentrations.

The next step in the analysis is to add information on the type of nuclide that is measured. This is partially presented in Figure 3 (right) for a limited number of nuclides of importance. In such cases, the median or mean values previously presented can be separated into single quantities. This is presented in detail in Figure 4 and Tables 3 and 4. Discussions are presented in the following sections. Evidently, further information on the samples (e.g. fuel type, reactor type, burnup) can be used to disentangle even more the C/E − 1 values.

4.1. Nuclide concentrations

As presented in the introduction of this section, more than 12 000 C/E − 1 values are taken into account. Figure 3 (right) presents the expected type of distributions one can obtain for a number of nuclides. The generalization for more nuclides is presented in Figure 4 and in Tables 3 and 4 (the data are represented in interquartile ranges: median and quarters). In these tables, n is the number of samples having measured values for the considered nuclide, the mean and median follow the usual definitions, “std” stands for standard deviation (1σ), and “mad” for median average deviation. Additionally, the interquartile quantities are also provided (IQR, Q1 and Q3). To indicate how many values are included within ±1σ or ±2σ, the two last columns of the tables provide the percentage of values included within these ranges. Finally, the column ΔC_max provides the maximum calculated uncertainties obtained in this work and presented in the references previously mentioned. Such uncertainties come from nuclear data, irradiation parameters and manufacturing tolerances (see dedicated references for details).

As presented in these tables, many quantities can be extracted to characterize the C/E − 1 distributions. Depending on the nuclide of interest, users can propose different methods to estimate biases; for instance, in the case of ²⁴⁴Cm, the average deviation between measured and calculated concentrations is 2%, with a higher calculated value (if one considers the mean of the distribution). But specific cases can strongly vary from the average, as quantified for instance by the standard deviation, which is 24%. To cover more than 90% of the calculated cases, one can say that C/E − 1 for ²⁴⁴Cm vary between −46% and +50% (by considering 2σ). Based on these values, and depending on the application, the user can define a confidence interval: for decay heat calculations, the decay heat at a specific cooling time can be underestimated by as much as 22% (times the ²⁴⁴Cm contribution to decay heat). One can for instance impose a conservative penalty factor of 1.22 on the calculated concentration, in order to account for a possible underestimated ²⁴⁴Cm amount. The same approach would hold for limits on neutron emission for transport cask. The inverse consideration can hold in the case of a strong neutron absorber (in burnup credit application), such as ¹⁴⁹Sm: the overestimation can reach more than 15%, by considering the mean and 2σ of the distribution. One could therefore justify imposing a penalty factor of 0.85 for the calculated ¹⁴⁹Sm inventory. Cases can also be more complicated, e.g. for fertile nuclides (leading to fissile nuclides after neutron capture, such as ²⁴⁰Pu), which can remove neutrons, but also lead to a fissile nuclide. As a side remark, sensitivity analysis for particular concentrations can help in better defining correction factors.

Due to the variety of possibilities (and interpretations of conservatism in considered applications), no recommendation is put forward regarding the definition of bias and penalty factor, as one needs to consider specific applications. But the biases will nevertheless be used for the estimation of decay heat, as presented in the following sections. Additional studies can be performed, for instance by analyzing the C/E − 1 values as a function of specific parameters, such as initial enrichment, reactor type, sample burnup, etc. This will imply smaller statistics per group, and extensive studies for various parameters, as partly performed in many publications (see for instance Ref. [44] for such type of study for decay heat and nuclide concentrations). Such effort was not performed in the present case.

4.2. Sample burnup

The burnup is a calculated quantity typically assessed from core calculations and it cannot be directly experimentally observed. It can nevertheless be derived from measurements, such as from the ¹⁴⁸Nd and ¹³⁷Cs concentrations. The latter can be obtained from gamma-ray spectroscopic measurements. A comparison of the burnup from these two origins – i.e. core calculations vs nuclide concentrations – can point to a mismatch due to various origins. In order to eliminate such differences, the calculated burnup (from either origin) is often adjusted so that it matches the one derived from one observable, for instance from the ¹⁴⁸Nd concentration. In PIE calculations, such an approach is almost always applied. This implies that the sample burnup cannot be known beyond the calculated and measured uncertainties of the observable, i.e. the ¹⁴⁸Nd concentration.

The measured uncertainties reported for the ¹⁴⁸Nd concentration are often comparable to the calculated ones, although they can reach 5% for specific samples. For the samples studied here the calculated uncertainties on ¹⁴⁸Nd range between 0.8 and 2.2%. They include the contribution of nuclear data as well as irradiation parameters and manufacturing tolerances and they represent 1σ of the calculated distributions. In addition, Table 4 indicates that for the 299 considered cases the ¹⁴⁸Nd concentration has an average C/E bias of −0.1% (due to normalization) and a C/E spread of 2.9% (1σ). This last value is certainly due to the experimental uncertainty on the ¹⁴⁸Nd concentration: when adjusting the calculated sample burnup to reproduce the measured ¹⁴⁸Nd concentration, one tries to be within the experimental ¹⁴⁸Nd uncertainty.

Two previous studies were compared with the present results. Reference [36] presents uncertainties due to nuclear data for the burnup of realistic LWR assemblies: a maximum of 2.3% is predicted for low burnup values (about 10 MWd/kg), down to uncertainties lower than 1% at higher burnup (up to 60 MWd/kg). These uncertainties are for assembly average values, which cannot be directly compared to the previous values (being for a specific sample, corresponding to a unique rod and vertical segment): they correspond to integral burnup values, which are likely to be better estimated than local values. Reference [37] reports burnup uncertainties of nuclear data origin for BWR assemblies only, with indications as a function of segment burnup. For segments located at the bottom or top of the assemblies, the local burnup was as low as 5 MWd/kg, and uncertainties (due to nuclear data) can reach 3%.

In the case of fuel samples with burnup between 60 and 100 MWd/kg, the calculation of the sample burnup is challenging due to the rod relocation. From an experimental aspect, the characterization of high-burnup samples is also more complicated due to the fuel-cladding interaction and the migration of specific nuclides. In this case, it is assumed that the sample characteristics are less known than for samples with lower burnup values.

Finally, it was observed in the present work that various simulation hypotheses (such as reflective boundaries, and two or three-dimensional models) or measurements in different laboratories for the same samples, can impact the calculated burnup by 2–3%.

Following these observations, it is recommended that the minimum estimated uncertainty on sample burnup is 5%, independent of burnup values.

4.3. Calculated nuclide uncertainties

Calculated uncertainties are presented in Tables 3 and 4 as ranges from the minimum to the maximum uncertainties obtained for the samples studied here.Differentiation as a function of sample characteristics is not performed here, and uncertainties are analysed in a general manner. It was nevertheless observed that the major factor in the change of uncertainties is not a sample characteristic, but rather prior assumptions for uncertainties, either for nuclear data or for engineering quantities (see for instance Refs. [36, 53]). Based on the uncertainties and sensitivities calculated in this work, the following remarks can be made.

• The first remark is that generally speaking, nuclear data is the largest contributor to uncertainties on nuclide concentrations (cross-sections and fission yields; the impact of thermal scattering data was shown to be negligible [21]). The only noticeable exceptions are ²³⁹Pu, for which the moderator temperature can be almost as important as nuclear data in the case of PWR samples, and ²³⁵U, ²³⁹Pu and ²⁴¹Am for which the void coefficient is also an important source of uncertainties for BWR samples. In the case of actinides, cross-sections are the major sources of uncertainties, and for fission products, it is mainly fission yields (with some noticeable exceptions such as ¹³⁴Cs and ¹⁵⁴Eu).

• Large variations of calculated uncertainties can be obtained, sometimes up to a factor of 10, especially for specific fission products. Fission products are generally very sensitive to fission yields, and less to thermal capture cross-sections (with some exceptions, such as some Cs, Nd, Gd and Eu isotopes). In addition, the information for fission yield covariances in nuclear data libraries is incomplete due to the lack of correlation matrices. This implies that each user can define its own correlation matrix, which largely influences some fission product uncertainties (e.g. Cs, Ag, Sb).

• For actinides, the spread of uncertainties is more limited than for fission products, except for ²³⁷Np and some neutron-rich curium isotopes. In the case of Cm, strong differences of uncertainties are observed in nuclear data libraries for nuclides in the chain of the Cm built-up (as for instance ²⁴²Pu). For ²³⁷Np, the main component of uncertainty change is the type of considered sample: UO₂ or MOX. It was observed that for MOX samples, the calculated uncertainty for the ²³⁷Np concentration is much higher than for the UO₂ sample, due to the difference in production routes.

Secondary remarks can be made. The effect of thermal scattering data, more especially H in H₂O (which is a special type of nuclear data) is minimal on nuclide concentrations and decay heat. Details can be found in reference [20]. Finally, nuclide uncertainties can differ if one considers a single assembly model, or a full core model (e.g. on ¹⁴⁸Nd) [50]. These variations are certainly linked to the calculation method and normalization, which are intrinsically different between single assembly and full core models.

From these tables and remarks, one can observe that calculated uncertainties of nuclide concentrations relevant for spent fuel applications can be relatively large. Some cases are analysed in the following. For other nuclides, maximum uncertainties presented in Tables 3 and 4 are recommended.

4.3.1. ¹⁴⁸Nd

The case of ¹⁴⁸Nd is important, as the calculated nuclide concentration is often used as a reference to normalize sample burnup values. As presented in reference [31], the calculated uncertainty can be as high as the extreme value of 7%. Such high value comes from the analysis of the SF95-5 sample, due to the impact of fission yields based on the JEFF-3.1.1 library (the impact of cross-sections is about 1.8% in this study). This high impact of fission yields reflects the (unfortunate) lack of correlations in the nuclear data library, with possible inconsistency between independent and cumulative yields and uncertainties, potentially leading to high uncertainties on calculated concentrations.

If one removes this high value, the uncertainties from other sources can reach 2.5%. The majority of uncertainties are in the range of 1–2%, mainly due to nuclear data (fission yields), but also burnup variations [25, 44].

Other variations tend to show important changes in ¹⁴⁸Nd concentrations. In the case of the Gundreminden samples [25], simple changes of nuclear data libraries (from JEFF-3.2 to ENDF/B-VI.8) indicate a change of 1.2% in the ¹⁴⁸Nd concentrations. In reference [31], a change from the JEFF-3.3 library to the ENDF/B-VIII.0 library indicates a modification of 2.2% for the ¹⁴⁸Nd concentration. The impact of the sample mass balance for the GU1 sample induces a change of about 2% [22]; in the case of the GU3 sample, differences in measurements between the two independent laboratories lead to C/E differences of 5%; finally, the considered geometry can change uncertainties from 0.5 to 1.1%, only by performing a pincell or full core simulation [50].

As already mentioned, the impact of the (calculated) sample burnup value directly modifies the ¹⁴⁸Nd content. The sample burnup is often varied by scaling the cycle-wise power levels, such that the EOL discharge yields the target ¹⁴⁸Nd concentration obtained from measurements. Such variations are based on mapping the ¹⁴⁸Nd spread due to the experimental uncertainty: for instance, one can vary the sample calculated power so that the calculated ¹⁴⁸Nd concentrations cover the measured one, plus or minus one experimental standard deviation. This method bonds the calculated uncertainty due to the sample burnup to the experimental ¹⁴⁸Nd values, which in principle are two independent quantities. Additionally, this method can only be applied for measured concentrations, and not for samples (or assemblies) where the ¹⁴⁸Nd content is not experimentally known. For this reason, it is believed that the impact of the (assumed) sample burnup can be underestimated, especially for SNF without nuclide concentration measurements.

Based on these observations, it is recommended to accept a calculated uncertainty for the ¹⁴⁸Nd concentration not smaller than 2% for measured cases, and 4% for other cases.

4.3.2. ¹³⁷Cs

¹³⁷Cs is also an important nuclide as the measurement of its gamma-ray during its decay can lead to an estimation of the SNF burnup value. The calculated uncertainties vary from 2 to 8%, and the main sources of such values are nuclear data, more specifically fission yields. The value of 8% comes from reference [31] and is due to problems related to the covariance matrix of independent fission yields. Other studies indicate values ranging from 2 to 6%: the JEFF-3.3 library tends to lead to lower uncertainties than ENDF/B-VIII.0. Other parameters can influence the ¹³⁷Cs concentration, such as the sample enrichment (between 1.5 and 2% [25, 31, 44]). Additionally, differences between laboratory measurements for the GU3 sample indicate a change in C/E of about 1%, and model assumptions (pincell, assembly or full core) can change the C/E by up to 5%.

Based on these observations, it is justified to accept a calculated uncertainty for the ¹³⁷Cs concentration not smaller than 5%.

4.3.3. ⁹⁰Sr

⁹⁰Sr and its decay product ⁹⁰Y are important nuclides as they significantly contribute to the decay heat for PWR SNF assembly from 2 to 50–100 years. ⁹⁰Sr has a half-life of 28.6 years and decays to ⁹⁰Y (with a half-life of 64 hours). Only the ⁹⁰Sr concentration is measured and if a bias exists, the same bias will be observed for ⁹⁰Y. Uncertainties due to nuclear data range from 0.7 to 1.5% (with higher values obtained with the JEFF-3.3 library) if one excludes the high value of almost 7% calculated in reference [31]. ⁹⁰Sr is mostly sensitive to fission yields, as in the case of ¹³⁷Cs. Other sources of uncertainties, such as operating conditions or manufacturing tolerances can account for 1–1.5% [23, 35, 50]. Additionally, it was also observed that differences in ⁹⁰Sr concentrations obtained when changing libraries (JEFF-3.2 and ENDF/B-VII.1) can reach 1.7% [25].

Based on these observations, it is justified to accept a calculated uncertainty for the ⁹⁰Sr concentration not smaller than 3%.

4.3.4. ²³⁵U

The concentration of main actinides such as ²³⁵U is of prime importance for a large range of SNF characterization. Table 3 and the relevant references indicate calculated uncertainties from 1 to 7%. The maximum value of 7% comes from a modelling effect and the discretization of the neutron flux (for one-group cross-section convolution) in a number of energy bins [31]. Apart from this high value, the impact of the void coefficient for the ENRESA study leads to an uncertainty of almost 6% for the ²³⁵U concentration. Such value might be on the conservative side, as the considered variation of the void coefficient was large [33]. Another important factors are the sample enrichment (in the case of UO₂ samples), as well as the moderator temperature (both leading to 2–3%).

The impact of nuclear data varies for different libraries: about 2% for ENDF/B-VIII.0 and 4% for JENDL-4.0 (almost exclusively due to cross-sections). There is a difference between UO₂ and MOX samples, with lower uncertainties in the case of MOX samples. Finally, the repetition of calculations by different institutes, in the case of the GU3 sample, has led to a variation of about 1.5% for the ²³⁵U concentration.

Based on these observations, it is recommended to accept uncertainties not smaller than 4% for the ²³⁵Uconcentration.

4.3.5. ²³⁹Pu

The uncertainty for the ²³⁹Pu concentration is similar to ²³⁵U. There is a difference between MOX and UO₂ fuel (as for ²³⁵U): about 3% for the BM1 sample, where the moderator temperature plays an important role. It was observed that for MOX samples, the impact of the modelling can be a major source of changes for the ²³⁹Pu concentration. It was observed that simulating a single-pin cell, a full assembly, or a full core can affect the amount of ²³⁹Pu by more than a factor of 2 [50]. In the case of UO₂ fuel, calculated uncertainties can reach 6% (due to the void coefficient for the ENRESA samples; again, this value might be on the conservative side), 4.7% due to the coolant density, or 3.3% for the Gundremingen sample.

We, therefore, recommend a minimum calculated uncertainty of 4%, similar to the case of ²³⁵U.

4.3.6. ²⁴⁴Cm

With a half-life of 18.1 years, ²⁴⁴Cm is an important isotope for both the SNF neutron emission and MOX SNF decay heat up to 30–50 years of cooling time. The uncertainty due to nuclear data on the ²⁴⁴Cm nuclide concentration is about 10% for UO₂ fuel, and about 7% for MOX fuel. Other sources of uncertainties can contribute up to 10% [35], but were also reported to impact the ²⁴⁴Cm concentration by 3–4% [23, 25, 50]. Reference [25] also indicates that the impact of changing the nuclear data library, from JEFF-3.2 to ENDF/B-VII.1, modifies the ²⁴⁴Cm nuclide concentration by 11%.

We, therefore, recommend a minimum calculated uncertainty of 10%, similar to the case of ²³⁵U.

4.3.7. Link to decay heat

Many of the nuclides presented in Tables 3 and 4 are source terms for the calculation of the decay heat. For instance, the decay of fission products such as ¹⁰⁶Rh, ¹⁴⁴Pr or ¹³⁴Cs will be an important contribution to the SNF decay heat for a relatively short cooling time (less than a few years), decays of ⁹⁰Sr/Y and ^137(m)Cs/Ba will be major contributors for the first tens of years of cooling time, whereas actinides such as ²³⁸Pu or ²⁴¹Am represent the main sources of heat for 100–200 years of cooling time. These contributions might be different in the case of MOX fuel, with a stronger emphasis on actinides, but the same nuclides are at the top of the contribution list.

One can easily foresee that biases in nuclide concentrations might lead to biases for the SNF decay heat. For instance, an underestimation of the ¹³⁴Cs concentration (as shown in Tab. 4 with an average bias of −8.5%) will lead to lower decay heat in the vicinity of 2–3 years of cooling time. The amplitude of the decay heat bias will then be proportional to the contribution of the nuclides to the total decay heat. Naturally, compensation between nuclide contributions can appear.

There is an analogy with the calculation of a k_eff value for a criticality benchmark: a bias in a cross-section can lead to a bias in k_eff. And compensations are known to exist, even if the correct k_eff is calculated.

Naturally, the link between calculated nuclide compositions and calculated decay heat might also not exist in practice. This happens in the case of decay heat calculations using standards such as the ANS or DIN standards [54–56]. This can also occur when a simulation is targeted to assess the decay heat and to compare it with a measured value, without particular regards for the nuclide concentrations. Additional details are provided in the next section.

4.4. Decay heat

The decay heat is one of the most important observables for SNF in the context of short-, intermediate- and long-term storage. If the amplitude of such decay heat varies as a function of the SNF characteristics, the global aspect is similar for all SNF; three examples are presented in Figure 5 for different cases. These three SNF correspond to realistic cases: existing fuel assemblies, under real irradiation conditions; they are well representative of the large variety of actual SNF with a wide range of enrichment and burnup values. As observed, their amplitude and, to some extent, their shapes vary as a function of these two quantities and cooling time. For this study, the cooling period between 1 and 1000 years is considered. An estimation of the decay heat for these types of SNF, i.e. the code predictability (estimated biases, eventually with uncertainties) is touched on in this section.

Fig. 5. Example of calculated decay heat for a UO2 and MOX assemblies. Despite some differences, similar nuclides are the most important contributors, as presented later in the text and in Figure 8.

We will approach it with three estimation angles: based on standards, based on measurements, and based on the nuclide inventory as developed in the previous sections. This will allow covering pragmatic studies as usually performed by the different actors of the SNF management. As presented in the next sections, biases are not straightforward to define in the case of SNF decay heat, as measured values are seldom defined. In this respect, the use of measured nuclide concentrations can be a substitute for measured SNF decay heat over larger cooling periods.

4.4.1. Decay heat estimation based on standards

This is possibly the most convenient way to calculate decay heat for a specific SNF. Three standards can be considered: American [54], German [55, 56]. A minimum amount of information is needed to obtain a calculated value of the decay heat: SNF mass, average burnup, average enrichment, irradiation power and time (including inter-cycle cooling times), and “fission splits” (amount of fission from ^235, 238U and ^239, 241Pu). The calculation methods will not be described here, but a few remarks need to be taken into account. Some standards are conservative by nature; they do not systematically take into account the decay of fission products and actinides (with some noticeable exceptions, as for the decay of ¹³⁴Cs); and they have their own range of applicability. Examples of decay heat curves obtained from the best estimate code (in this case, the SNF code [58]), from the American standard ANS-5.1-2014 and the German DIN standard, are presented in Figure 6. Both SNFs correspond to the CLAB measurement campaign reported in reference [18]. The calculations based on the ANS standard are more conservative than the other calculated values. In the case presented here, the DIN standard is very close to the calculations with the SNF code. A comparison of the results obtained by ANS-5.1-2014 and best estimate SCALE calculations for realistic PWR fuel assembly cases has shown more than 20% over-prediction by the standard [15]. It might not be relevant to define bias values when true values (usually equal to measurements) are not defined. As observed here, measured values are defined for specific cooling times, which renders the bias estimation not applicable in a general case (for both standards and best estimate calculations).

Fig. 6. Examples of calculated decay heat from best estimate codes and from the two standards.

The use of standards will not be further detailed here, but users need to keep in mind that they present a convenient way to obtain an estimation of the decay heat from spent nuclear fuel. They might still differ from measured values and best estimates, with various degrees of deviation. Due to their conservative and approximate nature, they are of limited help for code validation.

4.4.2. Decay heat and bias estimation based on measured values

As presented earlier, measurements of decay heat values from SNF assemblies exist in a limited number of cases. In the present work, the cases reported from the CLAB, GE-Morris and HEDL facilities were calculated and compared to measurements (see Tab. 1 and Ref. [30]). Details of the measurements can be found in references [18,59,60]. Recently, a comprehensive study has been published by EPRI, analyzing the current state of knowledge, as well as identifying gaps, for decay heat estimation [13]. Given the relevance of this recent work for the present analysis, we will refer to this report.

A summary of studies performed in this work (see Refs. [30, 44] for details) as well as some analysis from the open literature (see Refs. [13, 61]) is presented in Figure 7. In this figure, one C/E value represents the average of all cases presented in the specific reference. For instance, the value of C/E = 1.010 ± 0.026 for the Polaris calculations of the CLAB BWR cases (first point on the left) is the average of 81 decay heat values from 49 different BWR SNF, and the uncertainty is the standard deviation of the 81 values, assuming each measurement independent from each other. The points for a given type of SNF are not based on the same number of experiments.

Fig. 7. Plots of the average C/E values for the decay heat from various references.

Fig. 8. Contributions to the total decay heat for three different SNF. Note that the same nuclides are the main contributors to the decay heat, only the fraction is changing (e.g. 244Cm).

A simple average of the values presented in this figure leads to an average of 1.002 ± 0.015 (1.002 for BWR, and 1.002 for PWR). Such values can be misleading; a few remarks can add a different perspective:

1. apart from the 4 SNF assemblies from the HEDL campaign (representing 6 values), all other measurements come from two calorimeters: at CLAB and GE-Morris. It is likely that values obtained at the same facilities suffer from similar experimental bias effects. Such biases can be averaged to zero when enough different sources of measurement exist, which is not the case here.

2. Even if individual C/E values presented in Figure 7 are spread within approximately ±5%, the standard deviations are relatively large. This indicates that individual cases still present a non-negligible disagreement between calculated and measured values. In the case of the CLAB data, which are believed to be the most accurate, two standard deviations for the C/E values represent almost 10%.

3. Apart from the 5 SNF decay heat measurements from reference [62], all other calculations were not blindly performed, i.e. measurement values were available at the time of the simulations. In the case of the blind simulations (performed by 31 combinations of codes and institutions), a vast majority of results indicated an underestimation of the calculation of about −2%. This value is inconsistent with the average C/E of 1.002 previously mentioned.

The recommendation from the EPRI report also indicates that code predictions are “within a few percent” for decay heat.

Nonetheless, higher predictive power is for the time being

challenging, and the present study tends to validate this observation.

Regarding simulation uncertainties, studies performed in this project and presented in references [22, 23, 31, 35, 50]

indicate that uncertainties due to nuclear data, manufacturing tolerances and irradiation history are not lower than 2.5%, and can

reach 6–7%, depending on the cooling time. Nuclear data uncertainties are usually the largest contributor.

Based on the above remarks, the comparison between (best estimate) code simulations and measurements indicates that the code prediction cannot be expected to be better than 5% (one standard deviation).

4.4.3. Decay heat and bias estimation based on nuclide compositions

The previous comparisons were based on direct observations of SNF decay heat. Naturally, decay heat is related to the nuclide compositions at the end of irradiation, and their decay scheme. Considering that the decay data for the nuclides of relevance for decay heat (within the cooling time of interest, between 1 and 1000 years) are mostly well known, which is generally the case, a correct estimation of the nuclide concentration leads to a correct estimation of the decay heat.

As presented in Sections 4.1–4.3, the estimation of nuclide concentrations can be performed based on PIE validation, and Tables 3 and 4 can be used as a base of such estimation. One approximation performed in this argumentation is that the C/E for the nuclides of the PIE samples are the same for the full assemblies. This is evidently an approximation, given that PIE samples were originally located at different heights and radial positions within a number of assemblies.

Another precondition is that the measured nuclide concentrations are of relevance for the decay heat estimation. For instance, the measurement of ⁹⁹Tc is of no importance for the decay heat, whereas ⁹⁰Sr is. It is therefore of importance to estimate the contributing nuclides to the decay heat. Such a list of contributing nuclides can be found in many reports and publications, and examples are presented in Figure 8. This selection of SNF is representative of realistic cases, even if other concentrations or burnup values could have been presented. It is nevertheless important to notice that between 1 and 1000 years, a short list of similar nuclides is enough to represent more than 90% of the total decay heat for most of the SNF cases. The order of such nuclides can vary from one SNF characteristic to another, but nuclides are not changing, even when considering MOX or UO₂ fuel.

The most important contributors are: ⁹⁰Y+⁹⁰Sr, ¹³⁷Cs+^137m Ba, ¹³⁴Cs, ¹⁴⁴Pr, ¹⁰⁶Rh, ^242, 244Cm, ^{238, 239, 240}Pu and ²⁴¹Am. With these 13 nuclides, 90% (or more) of the decay heat contribution is covered, and it is convenient to recognize that the same nuclides are the main decay heat contributors for a large variety of SNF characteristics. Among these 13 nuclides, only 4 are not directly measured during the PIE analysis: ⁹⁰Y, ^137m Ba, ¹⁴⁴Pr and ¹⁰⁶Rh; their parents or daughters are nevertheless measured: ⁹⁰Sr, ¹³⁷Cs, ¹⁴⁴Ce, and ¹⁰⁶Ru. Given that the decay data are correct, it is then possible to consider that a bias in the estimation of the concentration for these 13 nuclides, their daughters or parents can be used for the estimation of the decay heat bias, proportionally to their contribution to the decay heat at a specific cooling time.

The definition of the bias for the nuclide concentration is left to the appreciation of the user, but a convenient quantity is the mean C/E (or C/E − 1) presented in Tables 3 and 4. By simply using the mean biases for the nuclide concentrations, and the contributions of the nuclide to the decay heat (such contribution is obtained with a best-estimate code), one can obtain the decay heat bias as a function of cooling time. To illustrate this approach, the case of a PWR UO₂ SNF, enriched at 3.95% and with a burnup of 53 MWd/kg is presented in Table 5, for cooling times between 1 and 1000 years. This simple example allows an understanding of the impact of each nuclide and possible compensation effects. Each individual nuclide bias in Table 5 is based on the mean biases presented in Tables 3 and 4. The advantage of considering mean values is that the effect of individual samples can be averaged out, assuming that the collected sample is well representative of the population of different design and irradiation conditions: the only bias left is the one from the simulation code (with its nuclear data library), without approximations due to modelling (e.g. irradiation history, sample location). These results might be of interest with respect to characterizing a code performance in the context of its licensing. If one additionally considers the standard deviations presented in the mentioned tables, the bias on the decay heat will be significantly larger. For instance, at 10 years of cooling time, instead of having a mean bias of −0.7% in the above example, the use of the mean nuclide bias plus one standard deviation will lead to a decay heat bias of +11% (similarly, the mean nuclide bias minus one standard deviation leads to −12%).

Such decay heat biases are well representative of the values for the vast majority of existing SNF. Even in the case of MOX fuel, the decay heat biases are not substantially different: Figure 9 presents the decay heat biases obtained for two SNF types: one UO₂, 3.95% enrichment and 53 MWd/kg, and one MOX, 4.0% enrichment and 45 MWd/kg. As mentioned, limited changes can be observed from one SNF type to another, but the trends are very similar. Such results (i.e. the nuclide contributions) are obtained based on the best-estimate code SNF, and differences can potentially be observed with different codes. The shape of the decay heat biases can be explained from the mean nuclide biases: below 5 years of cooling time, the positive bias is mainly due to ¹⁰⁶Rh (and ¹⁰⁶Ru); above 100 years, the negative bias is due to ²³⁸Pu and ²⁴¹Am. For a more exhaustive study, other simulation codes need to be used, for instance, based on the summation method (such as Serpent). Differences can appear for short cooling time (less than 1 year), where a large number of fission products can contribute to the total decay heat.

Fig. 9. Decay heat biases obtained from the mean biases on nuclide inventory.

Compared to the estimation of the decay heat biases using direct measurements (presented in Sect. 4.4.2), the present bias is smaller than the 5% previously recommended. Such difference comes from the difference of approach, and possibly from the fact that a limited number of direct measurements exists. As mentioned previously, results presented in Figure 9 are based on the mean bias for each nuclide; if one considers standard deviations (to cover all nuclide concentrations), the derived decay heat biases strongly change.

4.4.4. Decay heat uncertainty

The uncertainties for the SNF decay heat as a function of cooling time were calculated and presented in many references. Different values can be obtained, depending on the assumed sources of uncertainties, such as nuclear data or operation variables. In such cases, uncertainties, usually defined as one standard deviation, can be obtained by repeating the same calculation a large number of times, each time randomly changing a specific set of input variables (or nuclear data) according to a covariance matrix. Another possibility is to perform a limited number of similar simulations, for instance changing the entire nuclear data library at once, the simulation code, or specific modelling assumptions. These different approaches were applied in this work, and a summary is presented below. Results are provided for a cooling period between 1 and 1000 years.

• In references [14, 26], the case of a PWR UO₂ assembly is considered, with 4.95% enrichment and a burnup value of 60 MWd/kg. Changes in simulation codes indicate a maximum effect on the assembly decay heat of 4%. The modification of the nuclear data library indicates a change smaller than 2.5%. Readers interested in the results of the sensitivity studies should directly access the two mentioned references.

• In reference [32], a similar case of a PWR UO₂ assembly (4.0% enrichment, 51 MWd/kg) is studied, and changes in simulation codes modify the SNF decay heat by a maximum of 2%. For the same assembly, results in reference [28] indicate that changes in boundary conditions (from over-reflective to full-vacuum cases) moderately affect the assembly decay heat below 20 years of cooling time (less than 3%), but induce strong changes for longer cooling time (up to 25%). Such extreme changes in boundary conditions are likely to be unrealistic, and more reasonable changes (such as realistic surrounding assemblies vs. reflective boundaries) affect the SNF decay heat by less than 5%.

• The impact of nuclear data were quantified in the case of the SF95-5 PIE sample, considering the covariance matrices from the JEFF-3.3 and EAF-2010 libraries.A maximum of 5% impact is found, close to 2 years of cooling time.

• For one of the BWR Gundremingen samples studied in this work, the major impact on the sample decay heat comes from the knowledge of the sample burnup, changing the decay heat by a maximum of 2% [25].

• In the case of measured decay heat values as studied in reference [35], the impact of nuclear data as well as operational variables ranges between 2 and 3%. The average SNF burnup was identified as one of the major parameters affecting the decay heat.

• For the PWR samples GU1, GU3 and BM1 [22, 23, 50], similar variations were applied, in terms of nuclear data, operating conditions and manufacturing tolerances. It was shown that decay heat uncertainties were globally higher than 2% for the UO₂ samples and higher than 1.5% for the MOX sample. Nuclear data is one of the main sources of uncertainties, together with the temperature of the moderator. Regarding nuclear data, significant differences were found based on different libraries, mainly due to differences in uncertainties for the ¹³⁴Cs fission yield (a maximum decay heat uncertainty of 5–6% was found close to 2 years of cooling time with the ENDF/B-VIII.0 and JENDL-4.0 libraries).

• For the ENRESA BWR samples [33], similar conclusions as for the PWR samples were found, with additionally a strong contribution from the void coefficient. One of these samples was also modelled with different assumptions (related to the interpretation of the description of the sample irradiation), leading to changes as high as 12% for the decay heat.

• Finally, it is worth looking at references [36, 37], presenting decay heat uncertainties due to nuclear data for thousands of realistic SNF. Dependance with the SNF burnup and enrichment values can be observed, but the trends are globally the same as for the GU1, GU3, BM1 and ENRESA samples: a peak close to 2 years of cooling time between 6 and 7%, followed by a decrease generally close to 2%.

Based on these remarks, SNF decay heat uncertainties due to nuclear data, manufacturing tolerances, modelling effects and irradiation history are recommended to be higher than 4%, regardless of the SNF characteristics, for the period between 1 and 1000 years. Note that different values are recommended in reference [13].

5. Recommendations

In this section, a list of simple recommendations, supported by the present study, is proposed. It aims at guiding the modelling of SNF and at proposing adequate biases and uncertainties for SNF nuclide concentrations and SNF decay heat. A number of these recommendations is in agreement with the PIRT report [13] and the NEA WPNCS report [12].

5.1. Best practice

A number of best practices for SNF modelling can be followed in order to minimize the simulation biases:

1. use the latest nuclear data libraries. Some historical libraries might not contain all the necessary information for decay heat calculations, such as complete decay data schemes, or up-to-date fission yields.

2. Avoid single pincell simulations (i.e. a single rod calculation in two dimensions). In such a case, the actinide concentrations can be significantly modified compared to a full two-dimensional assembly simulation (especially in the case of heterogeneous fuel lattices, such as typical BWRs or MOX designs). The user needs also to be aware of the impact of the boundary conditions (neighbour assemblies or control rods), as well as the crucial impact of the burnup normalization through the ¹⁴⁸Nd content. In the latter case, the normalization of the sample power to match the ¹⁴⁸Nd content can be useful for the validation of a specific code (and its nuclear data library) using a PIE sample, as it adjusts the irradiation conditions to observables (i.e. ¹⁴⁸Nd content). For a full assembly, conclusions based on the previous validation effort are not identically applicable due to the impact of the irradiation history.

3. If possible, use the maximum information from a full core simulator. The main advantage is that the assembly average burnup is correlated with the core power, a quantity believed to be well known (in addition, for a specific core loading, all calculated assembly burnup values are correlated between themselves due to the core power normalization).

5.2. Uncertainties and biases

A simple summary of the previous section is presented in Table 6. For other nuclide concentrations, details are presented in Tables 3 and 4. A detailed comparison of these recommendations with other published analyses is out of the scope of the present paper, but it can be noticed that the values in Table 6 are more conservative than the ones reported in reference [13]. Additional work is necessary to understand the origins of these differences. The present values for the decay heat uncertainties and biases are nevertheless in agreement with the conclusions of reference [62] (blind calculations for five assemblies).

5.3. Needs

Based on the existing international efforts, the availability of measurements and databases, the following needs are identified in order to improve the estimation of the important SNF quantities (nuclide concentrations and decay heat):

1. the number of SNF decay heat measurements is extremely limited and based on a small number of experimental facilities (only one still exists nowadays). This is limiting the validation capabilities, and therefore the confidence in the predictive power of existing simulation codes. New SNF decay heat measurements are therefore recommended, overlapping with the characteristics of existing SNF (high burnup, high enrichment, large range of cooling time). If possible, measurements of the decay heat and of nuclide concentrations should be performed for the same assemblies. Only these sources of independent experimental information can validate at the same time calculated compositions and calculated decay heat.

2. Experimental uncertainties for nuclide concentrations are often not covering the full aspect of the experimental knowledge (often leading to too small reported uncertainties). The consequence is that such uncertainties are left aside during the simulation process, leading to a possible underestimation of recommended uncertainties. It is therefore suggested to perform a more complete uncertainty assessment during the measurement procedure.

3. Calculated quantities such as SNF nuclide concentrations and decay heat can be affected by the normalization procedure (e.g. to ¹⁴⁸Nd) and by the boundary conditions. It is therefore recommended to perform simulations based on the maximum information regarding the irradiation history, and, as much as possible, to avoid power normalization through the ¹⁴⁸Nd content.

4. The normalization of the sample power to the measured ¹⁴⁸Nd content presents advantages, but it does not disentangle between fission contributors (such as ²³⁵U or ²³⁹Pu). The use of a ratio such as ⁹⁰Sr over ¹⁵⁰Nd, with different fission yields for the two main fission contributors, can help to better assess if individual contributions are correctly captured during the simulation.

5. Finally, specific needs for nuclear data can also be expressed: the necessity of providing fission yield correlations in the evaluated libraries, better knowledge of the ²⁴²Pu(n,γ) cross-section (for the estimation of ²⁴⁴Cm), of the independent fission yields of isotopes of mass 134 (for the estimation of ¹³⁴Cs), of the ¹³⁴Cs(n,γ) cross-section (also for the estimation of ¹³⁴Cs), and of the production and disappearance chain of ²⁴⁰Pu. Last but not least, the precise knowledge of the cumulative fission yields leading to isotopes such as ¹³⁷Cs, ¹⁴⁸Nd or ⁹⁰Sr is of high importance of the characterization of the spent fuel.

6. Conclusions

In this work, a number of PIE samples were analyzed in detail, leading to estimations of C/E ratios for nuclide concentrations. In addition, existing measurements for SNF decay heat were compared to calculated values. The combination of the analysis of these measured quantities has led to recommendations regarding best practices for simulations, best-estimate values, uncertainties and biases for calculated SNF nuclide concentrations and decay heat.

The current simulation practices and existing experimental data are believed to limit our prediction capabilities for SNF characterization, for instance to a combination of uncertainties and biases not better than 4% for SNF decay heat between 1 and 1000 years of cooling time. The need for additional measurements for SNF decay heat is clearly expressed, in order to reduce the future cost of SNF storage and long-term disposal.

Conflict of interests

The authors declare that they have no competing interests to report.

Funding

It was funded by the European Union’s Horizon 2020 Research and Innovation Programme under grant agreement No. 847593 (EURAD, Work Package 8 on “Spent Fuel Characterisation and evolution until disposal”). The work of one author (D. Rochman) was partly supported by swissnuclear, the association of the Swiss nuclear power station operators, with the COLOSS-2 project.

Data availability statement

All used data are publicly available.

Author contribution statement

All the co-authors contributed equally to this work (modelling, simulations, analysis and writing).

Appendix A

Calculated nuclides

References

All Tables

All Figures

	Fig. 1. Left: example of assembly burnup at the end of life, based on realistic irradiation history for a number of existing PWR and BWR power plants. Right: same but for average fissile enrichments. Vertical coloured segments indicate the burnup of available PIE samples (long segments: considered in this work; short segments: extracted from SFCOMPO [8]).
In the text

	Fig. 2. Right: example of assembly decay heat at 10 years cooling time, based on realistic irradiation history for a number of existing PWR and BWR plants. Left: same for the assembly average 235U and 239Pu concentrations.
In the text

	Fig. 3. Histograms of C/E − 1 for all nuclides (left), and for a selection of nuclides (center and right).
In the text

	Fig. 4. Interquartile ranges for the C/E − 1 nuclide concentrations, considering a total of more than 12 000 measured concentrations. The blue colour is given to important nuclides. See Tables 3 and 4 for numerical values. The limits of plotted boxes are equal to the median ± IQR/2, as presented in Tables 3 and 4.
In the text

	Fig. 5. Example of calculated decay heat for a UO2 and MOX assemblies. Despite some differences, similar nuclides are the most important contributors, as presented later in the text and in Figure 8.
In the text

	Fig. 6. Examples of calculated decay heat from best estimate codes and from the two standards.
In the text

	Fig. 7. Plots of the average C/E values for the decay heat from various references.
In the text

	Fig. 8. Contributions to the total decay heat for three different SNF. Note that the same nuclides are the main contributors to the decay heat, only the fraction is changing (e.g. 244Cm).
In the text

	Fig. 9. Decay heat biases obtained from the mean biases on nuclide inventory.
In the text