On the estimation of nuclide inventory and decay heat: a review from the EURAD European project

. In this work, a study dedicated to the characterization of the neutronics aspect of the Spent Nuclear Fuel (SNF), as part of the European project EURAD (Work Package 8), is presented. Both measured nuclide concentrations from Post Irradiation Examination samples and decay heat from calorimetric measurements are compared to simulations performed by diﬀerent partners of the project. Based on these detailed studies and data from the published literature, recommendations are proposed with respect to best practices for SNF modelling, as well as biases and uncertainties for a number of important nuclides and the SNF decay heat for a cooling period from 1 to 1000 years. Finally, speciﬁc needs are presented for the improvement of current code prediction capabilities.


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 (intermediateand 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, lead- * e-mail: dimitri-alexandre.rochman@psi.ch ing 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 235 U, 239 Pu or 137 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 fissile enrichment (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; dualpurpose 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".

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 2 , 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.
In Figure 1 left, the assembly-averaged burnup, in terms of MWd/kg 1 , 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 2 cases, the fissile content is the concentration in 235 U, and for MOX, it is the sum of 235 U and 239,241 Pu. 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 2 and MOX fuel types (short vertical segments). If the UO 2 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 2 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 239 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 2 . 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 ( 239 Pu and 235 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 239 Pu and 235 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.

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: 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.

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.
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.

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 krsko.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][17][18][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][21][22][23][24][25][26][27][28][29][30][31][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.

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 ARI-ANE 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].

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 inhouse 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 threedimensional 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-tracingbased 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 ORI-GEN decay calculation using the calculated Polaris Quantities in italics indicate comparisons with experiments; "PIE" means nuclide inventory from Post Irradiation Examination, "DH" is for decay heat, and "UQ" is for uncertainty propagation. The term TSL corresponds to the thermal scattering data, and the term γ, n corresponds to the study of gamma and neutron emission.
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][46][47].
The present study is complemented by results from the literature, allowing it to cover a wide range of C/E values.

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 crosssections, 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 EVOL-CODE).
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 235 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).

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 2 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 148 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 pre-sented 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.

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 244 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 244 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 244 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 244 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 149 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 149 Sm inventory. Cases can also be more complicated, e.g. for fertile nuclides (leading to fissile nuclides after neutron capture, such as 240 Pu), which can remove neutrons, but also lead   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. 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.

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 148 Nd and 137 Cs concentrations. The latter can be obtained from gammaray 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 148 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 148 Nd concentration. The measured uncertainties reported for the 148 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 148 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 148 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 148 Nd concentration: when adjusting the calculated sample burnup to reproduce the measured 148 Nd concentration, one tries to be within the experimental 148 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.

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 239 Pu, for which the moderator temperature can be almost as important as nuclear data in the case of PWR samples, and 235 U, 239 Pu and 241 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 134 Cs and 154 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 237 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 242 Pu). For 237 Np, the main component of uncertainty change is the type of considered sample: UO 2 or MOX. It was observed that for MOX samples, the calculated uncertainty for the 237 Np concentration is much higher than for the UO 2 sample, due to the difference in production routes.
Secondary remarks can be made. The effect of thermal scattering data, more especially H in H 2 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 148 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.

148 Nd
The case of 148 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 148 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 148 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 148 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 148 Nd content. The sample burnup is often varied by scaling the cyclewise power levels, such that the EOL discharge yields the target 148 Nd concentration obtained from measurements. Such variations are based on mapping the 148 Nd spread due to the experimental uncertainty: for instance, one can vary the sample calculated power so that the calculated 148 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 148 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 148 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 148 Nd concentration not smaller than 2% for measured cases, and 4% for other cases.

137 Cs
137 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 137 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 137 Cs concentration not smaller than 5%.

90 Sr
90 Sr and its decay product 90 Y are important nuclides as they significantly contribute to the decay heat for PWR SNF assembly from 2 to 50-100 years. 90 Sr has a halflife of 28.6 years and decays to 90 Y (with a half-life of 64 hours). Only the 90 Sr concentration is measured and if a bias exists, the same bias will be observed for 90 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]. 90 Sr is mostly sensitive to fission yields, as in the case of 137 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 90 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 90 Sr concentration not smaller than 3%.

235 U
The concentration of main actinides such as 235 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 235 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 2 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 2 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 235 U concentration.
Based on these observations, it is recommended to accept uncertainties not smaller than 4% for the 235 U concentration.

239 Pu
The uncertainty for the 239 Pu concentration is similar to 235 U. There is a difference between MOX and UO 2 fuel (as for 235 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 239 Pu concentration. It was observed that simulating a single-pin cell, a full assembly, or a full core can affect the amount of 239 Pu by more than a factor of 2 [50]. In the case of UO 2 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 235 U.

244 Cm
With a half-life of 18.1 years, 244 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 244 Cm nuclide concen-tration is about 10% for UO 2 fuel, and about 7% for MOX fuel. Other sources of uncertainties can contribute up to 10% [35], but were also reported to impact the 244 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 244 Cm nuclide concentration by 11%.
We, therefore, recommend a minimum calculated uncertainty of 10%, similar to the case of 235 U.

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 106 Rh, 144 Pr or 134 Cs will be an important contribution to the SNF decay heat for a relatively short cooling time (less than a few years), decays of 90 Sr/Y and 137(m) Cs/Ba will be major contributors for the first tens of years of cooling time, whereas actinides such as 238 Pu or 241 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 134 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][55][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.

Decay heat
The decay heat is one of the most important observables for SNF in the context of short-, intermediate-and longterm 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.
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.

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] and Japanese [57]. 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,238 U and 239,241 Pu). 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 134 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 stan-dard [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).
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.

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.   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%.  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).

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 99 Tc is of no importance for the decay heat, whereas 90 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 2 fuel.
The most important contributors are: 90 Y+ 90 Sr, 137 Cs+ 137m Ba, 134 Cs, 144 Pr, 106 Rh, 242,244 Cm, 238,239,240 Pu and 241 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: 90 Y, 137m Ba, 144 Pr and 106 Rh; their parents or daughters are nevertheless measured: 90 Sr, 137 Cs, 144 Ce, and 106 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 2 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 2 , 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 106 Rh (and 106 Ru); above 100 years, the negative bias is due to 238 Pu and 241 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.
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.

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 2 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 2 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 fullvacuum 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 2 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 134 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].

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].

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 148 Nd content. In the latter case, the normalization of the sample power to match the 148 Nd content can  Figure 9 The uncertainty represents one standard deviation (1σ).
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. 148 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).

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).

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 148 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 148 Nd content. 4. The normalization of the sample power to the measured 148 Nd content presents advantages, but it does not disentangle between fission contributors (such as 235 U or 239 Pu). The use of a ratio such as 90 Sr over 150 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 242 Pu(n,γ) cross-section (for the estimation of 244 Cm), of the independent fission yields of isotopes of mass 134 (for the estimation of 134 Cs), of the 134 Cs(n,γ) cross-section (also for the estimation of 134 Cs), and of the production and disappearance chain of 240 Pu. Last but not least, the precise knowledge of the cumulative fission yields leading to isotopes such as 137 Cs, 148 Nd or 90 Sr is of high importance of the characterization of the spent fuel.

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).