Next Article in Journal
Analysis of Voltage and Reactive Power Algorithms in Low Voltage Networks
Next Article in Special Issue
Early Warning Weather Hazard System for Power System Control
Previous Article in Journal
Quantitative Comparison of Infrared Thermography, Visual Inspection, and Electrical Analysis Techniques on Photovoltaic Modules: A Case Study
Previous Article in Special Issue
Risk Management Scenarios for Investment Program Delays in the Polish Power Industry
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Uncertainty Study of the In-Vessel Phase of a Severe Accident in a Pressurized Water Reactor

1
Faculty of Electrical Engineering and Computing, University of Zagreb, Unska 3, 10000 Zagreb, Croatia
2
Innovative Systems Software, 3585 Briar Creek Lane, Ammon, ID 83406, USA
3
Energy Software S.L, C/Catalunya 13, Taradell, 08552 Barcelona, Spain
*
Author to whom correspondence should be addressed.
Energies 2022, 15(5), 1842; https://doi.org/10.3390/en15051842
Submission received: 22 January 2022 / Revised: 18 February 2022 / Accepted: 28 February 2022 / Published: 2 March 2022
(This article belongs to the Special Issue Risk Management in the Energy Sector)

Abstract

:
A comprehensive uncertainty analysis in the event of a severe accident in a two-loop pressurized water reactor is conducted using an uncertainty package integrated in the ASYST code. The plant model is based on the nuclear power plant (NPP) Krško, a Westinghouse-type power plant. The station blackout scenario with a small break loss of coolant accident is analyzed, and all processes of the in-vessel phase are covered. A best estimate plus uncertainty (BEPU) methodology with probabilistic propagation of input uncertainty is used. The uncertain parameters are selected based on their impact on the safety criteria, the operation of the NPP safety systems and to describe uncertainties in the initial and boundary conditions. The number of required calculations is determined by the Wilks formula from the desired percentile and confidence level, and the values of the uncertain parameters are randomly sampled according to appropriate distribution functions. Results showing the thermal hydraulic behaviour of the primary system and the propagation of core degradation are presented for 124 successful calculations, which is the minimum number of required calculations to estimate a 95/95 tolerance limit at the 3rd order of the Wilks formula application. A statistical analysis of the dispersion of results is performed afterwards. Calculation of the influence measures shows a strong correlation between the decay heat and the representative output quantities, which are the mass of hydrogen produced during the oxidation and the height of molten material in the lower head. As the decay heat increases, an increase in the production of hydrogen and the amount of molten material is clearly observed. The correlation is weak for other input uncertain parameters representing physical phenomena, initial and boundary conditions. The influence of the order of the Wilks formula is investigated and it is found that increasing the number of calculations does not significantly change the bounding values or the distribution of results for this particular application.

1. Introduction

Uncertainties in assessing the nuclear power plant (NPP) behaviour according to predictions of severe accident codes affect the on-site and off-site risk assessment. The calculation of core damage and radioactive release depends on uncertainties associated with code models, the plant’s numerical representation and determination of the NPP operating parameters.
A significant experience exists in the application of best estimate thermal hydraulic system codes for quantification of uncertainty during different types of calculations: the analyses of a large break loss of coolant accident (LOCA) in experimental facilities [1] and NPPs [2,3], assessment of thermal hydraulic phenomena in boiling water reactors with a focus on plant systems and fuel performance, two-phase flow and sub channel analysis [4]. In the area of severe accidents there are fewer examples of the use of uncertainty quantification, and these include estimating the influence of uncertainty of input parameters and assumptions on fission product releases [5,6,7] and hydrogen production [8]. These analyses cover a station blackout sequence at the U.S. NPP Sequoyah, Westinghouse 4-loop pressurized water reactor (PWR) [5], 20 accident sequences covering a wide range of conditions in eight reactor designs that are operated in the European Union [6] and the analyses of experiment PHEBUS FPT-1 conducted in the PHEBUS facility used to investigate reactor accident scenarios in PWR reactors [7,8]. Severe accident codes have reached a high level of development and are used, among other applications, to support severe accident management programs [9,10] and to assess modifications in nuclear power plants [11]. Uncertainty analysis can be a tool to help improve these activities to make current and future power plants more resilient to the potentially severe consequences of unmitigated NPP accidents. Preliminary phenomena identification and ranking tables for severe accidents have been established [12], but they are still not as detailed as for design basis accident analyses [13,14,15].
Uncertainty methods are usually divided in two groups [16]. The first group consists of methods based on the propagation of input uncertainties, which can be further divided into probabilistic and deterministic methods. The second group includes methods that evaluate uncertainties based on the extrapolation of output uncertainties. The propagation of the input uncertainty approach involves the selection of input uncertain parameters with their uncertainty specification and performing multiple best estimate code calculations, so the term BEPU (best estimate plus uncertainty) is often used. The probabilistic methods [17,18,19,20], one of which was used in the presented analysis, use statistical principles to characterize input uncertainties, while the deterministic approach [16,21] does not use probability distributions, but reasonable bounding intervals of input parameters instead. In the extrapolation of the output uncertainty procedure [22,23], the uncertainty is determined based on a comparison of the calculation output results and the relevant experimental data.
Sources of uncertainty are, in general, uncertainties related to calculation models, insufficient knowledge of physical processes, especially during severe accidents, temporal and spatial (nodalization) discretization of the system, the plant’s initial and boundary conditions, etc. The code and the plant representation uncertainties are minimized by using a system code with best estimate thermal hydraulics and mechanistic treatment of a severe accident progression, and the plant nodalization qualified at the steady state level. The ASYST code [24], an advanced analysis software for nuclear safety applications [25,26,27,28,29], is used in the analysis presented in the paper. The code is developed as part of ASYST Development and Training Program (ADTP), an international nuclear technology project managed by Innovative Systems Software (ISS). It consists of three modules: the SAMPSON THA module for thermal hydraulic calculation that uses best estimate two-fluid, nonequilibrium models and correlations, the SCDAPSIM module for calculation of reactor core severe accident progression and the SCDAPSIM COUPLE module for the finite element thermal analysis of the reactor pressure vessel lower head. It also uses a variety of other member-developed computational packages.
The severe accident core behaviour is described to a certain extent by using correlations derived based on the separate effect tests and integral experimental results [30,31]. The uncertainties related to an insufficient knowledge of physical phenomena, such as fuel rod cladding oxidation, core fragmentation, metallic meltdown, molten pool formation, etc., are covered by varying user-defined parameters in these correlations. Initial and boundary conditions can vary up to 2% in a nuclear power plant. The major issue is the uncertainty in calorimetric error when measuring thermal power, but variations in other important NPP parameters such as reactor coolant pump (RCP), accumulator and break junction friction coefficients can also have significant consequences during a transient calculation.
The uncertainty package integrated in the ASYST code is based on the BEPU methodology with probabilistic propagation of input uncertainty [32]. It consists of the following steps:
  • Development of the nuclear power plant computation model for the steady state and the transient calculation;
  • Selection, identification and ranking of the relevant phenomena based on the safety criteria;
  • Selection of uncertain code parameters to represent those phenomena and determination of applicable probability density function (PDF);
  • Random sampling of uncertain parameters based on their PDFs and performing multiple code calculations determined by the percentile and confidence level using the Wilks formula;
  • Post-processing of results, determination of uncertainty bands, quantification of dispersion of output values and determining the relationship between input and output variables using influence measures.
A severe accident scenario in a two-loop PWR plant caused by a station blackout and coolant leakage through degraded RCP seals was selected for the uncertainty analysis. The integrity of the seal is compromised due to the unavailability of the charging pump, which normally provides its cooling. The in-vessel phase is analyzed thoroughly and the relevant figure of merits (FOM) are hydrogen production and the corium level in the lower plenum, which most accurately characterize the processes of oxidation of the fuel rod cladding and the core damage, respectively. However, other important output results are also presented and commented on.

2. Computation Model, Initial and Boundary Conditions

2.1. Nuclear Power Plant Nodalization

The plant model is based on the NPP Krško, a two-loop PWR plant of Westinghouse design [33,34]. The nodalization of the plant’s primary and secondary systems [35], including the radial cross-section of the core is shown in Figure 1.
The ASYST code, a six-equation one-dimensional, nonequilibrium and nonhomogeneous system code, uses detailed thermal hydraulic NPP representation, which is input-compatible with RELAP5/SCDAPSIM format [36]. The model consists of 554 thermal hydraulic control volumes (CV), 641 junctions, 361 heat structures with 1923 mesh points, 748 control variables, 197 variable and 221 logical trips. The reactor pressure vessel (RPV) is represented with control volumes 101 to 175. Primary sides of steam generators (SG) are modelled with control volumes 215–245 and 315–345 for SG1 and SG2, respectively. The secondary side of the first steam generator is modelled with CVs 411–429 and the second steam generator with CVs 511–529. The pressurizer is represented with CVs 61–69. The thermal hydraulic model also includes models of reactor coolant pumps (CVs 265 and 365), main feedwater and auxiliary feedwater lines, accumulators (CVs 701 and 801), high pressure and low pressure safety injection lines, steam lines, power operated and safety valves, etc., that is, all the relevant components for a thorough safety analysis. It is assumed that the seals of both pumps will degrade when left without cooling, and the rupture is modelled to be located on the suction side of the pump, on the control volume 259 for the first pump (RCP1) and on the CV 359 for the second pump (RCP2).
Fuel rods, control rods, the core baffle, grid spacers and structures in the upper plenum are modelled as SCDAPSIM components, which can break, melt or oxidize under severe reactor conditions. The hemispherical RPV lower head is modelled using the COUPLE code, a two-dimensional, finite element heat conduction code incorporated in the ASYST code.
Core fuel assemblies, 121 in total, are divided in five regions by grouping similarly powered fuel assemblies together (Figure 1). Each region is assigned a separate thermal hydraulic channel. These five channels are connected by transverse junctions (12 junctions in height between the two adjacent channels). The core will begin to melt in the central upper part, and the melting will progress radially and axially downward. The core needs to be divided into several regions in order to cover the spatial domain of the degradation and melting. Accumulation of impermeable debris could block the flow of the coolant in the thermal hydraulic channel. Thus, it is essential to model the crossflow junctions between the channels to allow the coolant flow to be redirected in the radial direction from the blocked channel to the adjacent channel with a free flow surface.

2.2. Selection of Uncertain Parameters

The uncertain parameters used in the analysis can be divided into three groups (Table 1). The first group contains parameters related to the calculation of the core behaviour during a severe accident. Although ASYST is a mechanistic code that automatically calculates the degradation of the core and is based on stand-alone models, user defined parameters are still necessary to describe certain processes in the core. These are processes that have not been fully investigated, i.e., for which there is a lack of experimental data or knowledge about their progression. In general, physical processes during severe accidents are described by mathematical correlations, determined on the basis of experimental research, which may contain constants (parameters) that the user can change or use the recommended values by the code developers. Changing the values of these parameters affects the course of the accident. Phenomena covered by ASYST, the calculation of which can be influenced by changing user-defined parameters, include oxide layer stability, metallic meltdown, molten pool growth, core fragmentation and cladding deformation behaviour. Table 1 provides a detailed description of uncertain parameters that describe these phenomena (parameters 1–12) and gives their reference values and probability density functions. Variations in parameters for which there are no strict constraints are described by a normal distribution, while for parameters that are limited by the core design, they are described by a uniform distribution. The mean value (relative to reference value) and the variance are required data inputs for the parameters described by the normal distribution. The minimum and the maximum values (relative values) are required for the parameters described by the uniform distribution. Accordingly, most parameter variations are represented with the normal distribution, except those that describe fuel rod cladding deformation. Their parameter uncertainties are represented with the uniform distribution so that the input values of the cladding strain do not exceed the given fuel assembly geometry limits.
The grid spacer uncertain parameters are their mass, height, plate thickness and radius of the contact area between the spacer and the fuel rod cladding. The grid spacers can affect the course of the accident in two ways. On the one hand, they can slow down the movement of corium and enhance the formation of the molten pool, and on the other hand, due to chemical reactions between zirconium (in the cladding) and nickel (in the grid spacer) occurring at relatively low temperatures of 1200–1400 K, they can lead to earlier core melt. The fuel rod uncertain parameters include fuel rod pellet, inner and outer cladding radii, fuel rod plenum length and the plenum void volume.
In addition, the fraction of theoretical fuel density and helium gas inventory are also varied. As for the dimensions of the fuel pellet and the cladding, there is a small tolerance specified by the manufacturer, so we wanted to examine the impact of variations in fuel rod dimensions on the processes of core damage and oxidation. As for the helium inventory, variations in its mass will affect the integrity of the cladding. The cracking of the cladding allows the penetration of steam into the gap of the fuel rod and oxidation of the cladding from the inside and thus increases the production of hydrogen.
The second group of uncertain parameters covers thermal hydraulic initial and boundary conditions. Break junction and steam generator safety valve friction coefficients vary in a wide range because the uncertainties related to the pressure drop at the break and in the valves are large. The pressure drops in the components that are active during the normal operation of the power plant are adjusted in such a way as to obtain a satisfactory steady state. For other components, such as breaks and safety valves, there is no particular way of determining friction coefficients, so in our case, variations in the coefficients are taken to be ±20%, which is still a rather conservative assumption because we did not want to exaggerate their influence on the calculation results.
Since the analyzed scenario is the loss of AC power, all major safety systems are out of operation except for passive hydro accumulators. Therefore, the variable initial conditions for the accumulators are modelled: the pressure, temperature and volume of water inside them. Accumulators inject significant amounts of water at short time intervals and can successfully quench the core, but that can also negatively affect the integrity of the core due to possible thermal shock and fuel rod shattering. The last components for which variations in input parameters are assumed are reactor coolant pumps. The RCP moment of inertia and friction torque coefficients are varied to examine the effect of increased or decreased primary coolant flow on the values of primary pressure and temperature, and the consequent behaviour of the core.
The third group of uncertain parameters consists of time-dependent decay heat data. The uniform distribution of decay heat is assumed within the range of 98% and 102% of the base value. The reasons for the power variations are the uncertainties related to the calorimetric error when measuring thermal power, and the standard deviations in fission product yields, half-lives and the amount of energy released during the radioactive decay [37]. The generation of heat has a key impact on the thermal response of the core, its damage and fuel rod cladding oxidation, which exponentially depends on temperature. Decay heat data are given in a table according to the ANS79-1 standard [38] and change uniformly depending on the weighting factor determined by the probability density function.

2.3. Determination of the Sample Size

The uncertain parameters are sampled randomly in accordance with their distribution function. They vary simultaneously in each calculation and the user only needs to specify uncertain parameters and their probability density function data.
The number of code runs (sample size) depends on three factors: the confidence interval γ, the confidence level β and the “order statistics” m. Using the well-known Wilks formula [39] approach, the sample size N for m = 1 can be determined by the equation
β = 1 γ N
The equation defined in this way gives the minimum number of calculations required to estimate the one-sided tolerance limit. Thus, for a 95th percentile (γ = 0.95) with 95% confidence level (β = 0.95), which are standard values, the number of code runs for this 95/95 tolerance limit is found to be 59.
If the order of application of the Wilks formula m is greater than 1, the sample size is determined by the equation given by Nutt and Wallis [40]
β = 1 i = N m + 1 N N ! i ! ( N i ) ! γ i ( 1 γ ) N i
Using this type of statistical formula, it is achieved that the number of calculations does not depend on the number of uncertain parameters but only on the percentile and confidence level of the desired uncertainty bound. For the 95/95 tolerance limit, the sample size depending on the order of the Equation (2) is as follows:
  • For the 1st order there are 59 code runs required;
  • For the 2nd order there are 93 code runs required;
  • For the 3rd order there are 124 code runs required;
  • For the 4th order there are 153 code runs required;
  • For the 5th order there are 181 code runs required, etc.
The interpretation of the order statistics (order of application of the Wilks formula) m is not clearly defined, as discussed by Frepoli [3]. It refers to the number of outcomes of the conducted analysis that are to be further evaluated. For a design basis, LOCA, the licensing criteria for emergency core cooling systems for light-water reactors (U.S. NRC regulations [41]), identify three criteria and their limiting values: peak cladding temperature, local maximum oxidation and core-wide oxidation. Thus, if we were to strictly adhere to the regulations as recommended by Guba et al. [42], the number of calculations for 95/95 tolerance limit should be 124. In the area of severe accidents, there is no precise recommendation for the selection of critical parameters, although the results of the analysis such as fission product releases, generated hydrogen mass and core condition during the degradation process are essential for planning accident management measures. On the other hand, in contrast to Guba et al., Wallis [43] concludes that there is only one outcome of interest from an accident analysis and that is a simple statement (pass or fail) whether prescribed safety limits (regulatory criterion) have been exceeded or not. In order for the regulatory criterion to be met, all output variables must be below the permitted values.
The sample size also depends on the confidence level and the confidence interval. For example, for the first order and the 97/97 tolerance limit the sample size is 116, and for the 99/99 limit the sample size is 459. Since the sampling process is completely random regardless of the tolerance criteria, the sample size can be regarded as an approximate value, and its higher value guarantees better accuracy in determining the bounding values of output results. We decided to use a sample size of 124 (γ = 0.95, β = 0.95, m = 3) for the presented analysis. Prior to that, a preliminary analysis was performed with the sample size of 59, with uncertain parameters covering only the modelling of the core damage (the first group of parameters in Table 1). A short discussion and a comparison of the main results of these two analyses is given in Section 4.3.

3. Severe Accident Analysis

3.1. Description of the SBO Accident

The analyzed transient is a station blackout, which includes a loss of electric power supply and the simultaneous leakage of coolant through the breaks at both intermediate legs, at control volumes representing the suction side of the reactor coolant pumps. The break size is 18 mm, a small break LOCA scenario. The transient lasts 15,000 s and is preceded by the 1000 s steady state calculation. The transient time is long enough to cover the entire sequence of events for an in-vessel phase of a severe accident, including the late phase of corium accumulation in the lower head.
After the loss of electric power and opening the breaks, the coolant is released from the reactor coolant system (RCS) and the primary pressure decreases. On the secondary side, the injection of main and auxiliary feedwater into the steam generators is interrupted, so the steam generators gradually dry out. During this period, the primary system is cooled by the natural circulation established between the core and steam generators, but after losing the heat sink on the secondary side (after the steam generators have been emptied), the natural circulation is terminated. The primary pressure then rises again, until the moment of the loop seal clearing when it finally begins to decrease. That process is standard during a small break LOCA accident when water accumulates in the intermediate leg, the lowest piping elevation, between the reactor coolant pump and the steam generator. The pressure upstream increases due to the heating of the primary fluid, and when it becomes higher than the back pressure of the water, there is a breach of the flow blockage through the intermediate leg, drainage of the pipe, and a decrease in pressure.
Since no emergency core cooling system other than the accumulators is available, water starts to boil causing a decrease in heat removal from the core. The heating of the core is additionally supported by the fuel rod cladding oxidation. As the temperature rises, chemical reactions begin in the core, iron-zirconium, nickel-zirconium, silver–zirconium, etc. They cause liquefaction and melting of the fuel assemblies and core supporting structures. The initial melting of a local character eventually leads to the subsequent formation of a molten pool, usually formed on the grid spacers. With the absence of cooling, the molten pool continues to grow in the radial and axial directions. In the final stage of the in-vessel accident progression there is a relocation of the molten material to the lower head. The corium level will increase over time and thermally and mechanically load the reactor vessel wall. The rupture of the wall and release of corium in the containment is not explicitly taken into account because the reactor pressure vessel wall can rupture as early as a few hundred seconds after the beginning of the relocation, but also much later. That depends on the conditions in the lower head, numerical models of interaction between the corium and the wall, wall deformation and the failure criteria. It is a complex analysis involving the interaction of the primary system and the containment, which goes beyond the scope of this study. Therefore, the idea of the paper was, among other things, to monitor the long-term outflow of corium from the core into the lower plenum, in order to gain a better insight into the dynamics of the in-vessel phase of a severe accident. In this way, an assessment of the uncertainty of the initial conditions for the analysis of the ex-vessel phase of the accident was achieved.

3.2. Calculation Results

The presented description of the transient coincides with the results of the analysis. Figure 2, Figure 3, Figure 4, Figure 5, Figure 6, Figure 7, Figure 8 and Figure 9 show the results of the base case and all 124 analyzed cases with the variation in input parameters.
The generation of a low flow signal in the cold legs, after losing power to the primary pumps at 0 s, will lead to the trip of the reactor and the turbine. Simultaneously with the loss of AC power, opening the breaks in the RCS causes a loss of coolant and decrease in primary pressure, Figure 2. The cold legs are dried out at 2000 s, but in the hot legs the void fraction is about 30%, which allows the establishment of the natural circulation between the reactor vessel and the steam generators. The drying of the hot legs follows quickly after drying of the steam generators and the termination of cooling of the core through natural circulation at approximately 4500 s. The pressure then increases until 6000 s when the loop seal clearing occurs. In the late phase of the accident there are intermittent increases in pressure caused, on the one hand, by the evaporation of water in the lower plenum after the pouring of hot material from the core and, on the other hand, by the evaporation of water freshly injected from the accumulators that comes into contact with the hot melt in the reactor vessel. The collapsed water level in the core and RCS coolant mass are shown in Figure 3 and Figure 4, respectively. The time of the core water level decrease roughly coincides with the decrease in the primary pressure. As the pressure decreases, so does the saturation temperature, the water in the core evaporates, which causes a decrease in the collapsed water level. The collapsed water level (COR LVL) is not the actual level but only an indication of the amount of water in the core. It is expressed as a percentage of the volume fraction of the liquid phase (voidf) and is calculated according to the equation
C O R   L V L ( % ) = i = 1 5 ( j = 1 12 1 12 N F A i 121 v o i d f i j ) 100
The core area is represented with five thermal hydraulic channels, each divided in 12 control volumes of the same height. The cross-sectional area of the channel is proportional to the number of fuel assemblies (NFA) located in the channel. In the first channel there are nine fuel assemblies, in the second 20 assemblies, in the third 52 assemblies, in the fourth 16 assemblies and in the fifth 24 fuel assemblies.
The mass of coolant in the RCS decreases continuously until 8500 s, but it is noticeable that this decrease is faster when the pressure in the primary system is higher, which is expected because the higher the pressure, the greater the force acting on the fluid at the break. After 8500 s, there is an increase in the amount of water in the primary system due to the injection from the accumulators. The accumulator pressure is shown in Figure 5. A pressure drop means that water has been ejected. This occurs when the pressure in the surge line decreases below the pressure of nitrogen, which is above the water in the accumulator and ensures a constant accumulator pressure before the activation. The initial pressures are different in different scenarios since the initial accumulator conditions are taken as varying parameters.
Figure 6 shows the production of hydrogen. Hydrogen generation is the result of the exothermic chemical reaction of zirconium in the fuel rod cladding and water vapor at high temperatures. The layer of zirconium dioxide is formed on the outside cladding surface during this oxidation process. The majority of hydrogen is produced in the interval between 7500 s and 8500 s (~70% of total production) when the collapsed core water level is close to or equal to zero (in steam-rich conditions), i.e., in the period between the core dryout and the start of accumulator operation. The average rate of hydrogen generation is then 0.25 kg/s. In general, an injection from accumulators reduces a further increase in oxidation. There are two groups of cases depending on the accumulator activation time and the decay heat value. When the accumulators start to inject water 50–100 s earlier, or the amount of decay heat is lower, the total hydrogen production is 15% lower than in the cases that have higher core heat production and a later start of the accumulators. This second group of cases is characterized by an extended period of high oxidation rate at a time of about 8500 s (Figure 6). The oxidation process can be divided in two phases. The first phase corresponds to the already mentioned period from 7500 s to 8500 s when the oxidation is more intense and predominantly the central and upper parts of the core are oxidized. After 8500 s the oxidation rate slows down and steam production in the core is reduced. Steam is produced by the evaporation of water injected from the accumulators, and the oxidation front moves to the lower parts of the core. Melting of the fuel rods and relocation of the corium have little effect on the oxidation rate. In some cases, however, increased oxidation was indeed observed at the time of the core melt and molten pool formation. This can be attributed to the cracking and disappearance of the protective zirconia layer on the fuel rod cladding and the exposure of the zirconium metal alloy to hot steam. Nevertheless, this process does not always have the same effect because in some cases, especially those with slower oxidation, no substantial oxidation was observed during the core damage.
The maximum core temperature is shown in Figure 7. Its sharp increase between 7500 s and 8000 s coincides with the processes of core drying and rapid cladding oxidation. The maximum value of the temperature is about 2850 K. It is determined by the solidus and liquidus temperatures of the UO2-ZrO2 quasi-binary system, which are 2830 K and 2873 K, respectively. The analysis of the late phase of the transient shows that the injection of cold water from the accumulators into the reactor vessel did not cause a decrease in the maximum temperature. Sudden temperature drops occur in some cases after 8000 s due to the relocation of the high-temperature corium from the core to the lower head. The maximum core temperature thus decreases after relocation and rises afterwards as the core is heated again by the decay heat.
The radius of the molten pool within the core is shown in Figure 8. By melting the core structures the radius increases until the moment of relocation to the lower plenum, after which it decreases, or drops to zero, and then begins to grow again as the new material melts. A key event in the molten pool formation is the transition of highly oxidized and partially damaged, or liquefied, fuel rods into low porosity debris. This occurs due to the thermal shock experienced by the fuel rods after the injection of cold water from the accumulators. Namely, zirconium dioxide, formed by the oxidation of zirconium on the outside of the fuel rod cladding, and uranium dioxide are low-strength ceramic materials that crack when exposed to sudden changes in temperature. Water, i.e., steam, hardly penetrates into the newly formed layer of debris, which makes it difficult to cool it down. The generation of decay heat is greater than the heat removal by water/steam flow, so the debris melts and a layer of melt is formed. Water injection from accumulators will, thus, not cool the core, but vice versa, it will accelerate the degradation process, which roughly corresponds to the scenario that occurred during the accident at the Three Mile Island power plant [44].
Figure 9 shows corium height in the lower head. The relocation of corium to the lower head represents a turning point in the further course of the accident because after that there is a significant scattering in the results. The interaction of corium and water in the lower plenum causes swift evaporation and an increase in the primary pressure. The pressure dynamics, but more generally the thermal hydraulic plant feedback, ultimately determine the operating time of the accumulators and other events caused by their activation. Emptying of the in-core molten pool affects both the thermohydraulic behaviour and the progression of a severe accident. Thus, after the corium relocation, there is no longer a clear pattern of events between different scenarios, but the key events take place at different times, leaving the impression of a stochastic distribution of results.
The molten pool will slump to the lower head if the crust at the periphery of the core surrounding the corium fails, which takes place when the stresses in the crust become higher than its ultimate strength. The initial relocation occurs between 9000 s and 10,300 s for different scenarios. The stability of the crust is a function of the crust thickness, differential pressure between the fluid inside the crust and the fluid outside of the crust and the weight of material that is supported by the crust above the molten pool. If the molten pool is surrounded by coolant that is considerably cooler than the corium, then a thick crust forms around the pool and holds it in place. Otherwise, if the coolant is at about the same temperature as the corium, then the crust will be thin and may fail. Variations in the input parameters affect the fuel and coolant temperatures, as well as the differential pressure around the crust, and therefore lead to different times of the crust failure.

4. Statistical Analysis and Discussion of Results

The performed calculations generated a large amount of data. Statistical analysis is needed to critically evaluate the results. The scatter of results is estimated using the relative standard deviation, and the relationship between the input parameters and the output variables using the Pearson correlation coefficient.

4.1. Dispersion of Output Values

The dispersion of a set of values around the mean, or the expected value, is quantified using the standard deviation. As its unit corresponds to the unit of a particular variable, it is not possible to directly compare the scatter of results of different output variables. Therefore, a relative value of standard deviation (RSD) equal to the ratio of the standard deviation to the mean was used instead in the analysis
R S D = 1 N i = 1 N ( x i x ¯ 1 ) 2
where N is the size of the population (sample size), xi i-th population value and x ¯ the population mean.
Relative standard deviations of thermal hydraulic variables: collapsed core water level (CORWAT), reactor coolant system mass (RCSMASS), integral break flow (BRKFLINT), pressurizer (PRZPRES) and accumulator (ACCPRES) pressures are shown in Figure 10, while for variables describing core damage: maximum core temperature (BGMCT), production of hydrogen (HYDR), corium height in the lower head (HGTCOR) and radius of the in-core molten pool (REPOOL), RSD values are shown in Figure 11. A small value of the relative standard deviation implies that the data are grouped close to the mean, while its large value implies that the data are scattered over a wide range around the mean, or some data deviate substantially from others.
The first thing that can be noticed in the figures are the sharp increases in the relative standard deviations of the collapsed core water level, hydrogen production, radius of the molten pool and the height of corium in the lower head at different times. They occur due to large differences in the results, at the same time points, between individual scenarios, where these differences are several orders of magnitude. There are three visible peaks on the core water level RSD curve, in 7650 s, in 8400 s and in 9100 s. The first peak refers to the initial drying of the core, the second to the first activation of the accumulators and the third peak to the second activation. These events do not occur simultaneously in all scenarios and are accompanied by sudden changes in water levels, which is why the values of the standard deviation are so high. Thus, we limited the maximum values of RSD to 2 because after reaching the peak value, they decrease fast afterwards. In this way it is then easier to monitor the behaviour of variables throughout the transient. The peak values are not representative (they are given in Table 2 though), but more important are the moments when they are occurring. The maximum value of the relative standard deviation of the hydrogen production occurs in 6500 s when the oxidation of the fuel rod cladding begins, which, again, does not occur simultaneously in all scenarios. The maximum value of the relative standard deviation for the molten pool radius occurs in 8400 s, and for the lower head corium height in 9000 s. The reason for the appearance of the first peak is the transition of the layer of debris into the molten pool shortly after the formation of core debris caused by injection of water from the accumulators. The peak value of the corium height in 9000 s corresponds to the moment of the first corium relocation to the lower head for the scenarios with the earliest cracking of the crust surrounding the corium within the core region. The relative standard deviation of the maximum temperature in the core is among the lowest compared to the other results (the lowest being for the coolant mass released through the breaks). The amount of water that enters the reactor vessel from the accumulators is not enough to quench the core, so this has very little effect on the maximum temperature transient behaviour.
The time profiles of the relative standard deviations of the four thermal hydraulic variables: reactor coolant system mass, integral break flow, pressurizer and accumulator pressures show similar behaviour. They are all affected by the activation of accumulators and injection of water, which changes the dynamics of mass and energy balance in the RCS, as discussed in the previous section. The values are small before the start of the accumulator injection, and after that they tend to increase steadily, except for the deviation of the coolant mass in the primary circuit, which decreases in the last 2000 s of the transient. The lowest value of the relative standard deviation is for the integral break flow, followed by the accumulator and the primary pressure, while the highest value of the RSD is for the already mentioned primary coolant mass. These variables do not have particularly pronounced sudden changes in their values, but an increase is visible after 10,000 s, especially for the values of primary pressure and reactor coolant mass RSD, which exceed relative standard deviations of collapsed core water level, hydrogen production, maximum core temperature and height of the molten material in the lower head.

4.2. Correlation between Input Uncertain Parameters and Output Data

The influence, or importance, measures are used to assess the impact of input parameter variation (Xi) on the variation in a selected output variable (Y). The dimensionless influence measure form is
Y X i Δ X i Δ Y
There are several ways to determine measures of influence, and the most commonly used are the Pearson correlation coefficient, Spearman’s rank correlation coefficient, the Kendall rank correlation coefficient, etc. All these approaches are equally valid [1], and the Pearson correlation coefficient was selected for the analysis because it uses the data values and not the ranks, so the information on the distance between the data is not lost. It is also relatively easy to implement it using an Excel spreadsheet. The Pearson correlation coefficient r is given by the following equation
r = i = 1 N ( x i x ¯ ) ( y i y ¯ ) i = 1 N ( x i x ¯ ) 2 i = 1 N ( y i y ¯ ) 2
where N is the sample size, xi, and yi are the individual sample points, and x ¯ and y ¯ the sample means.
The Pearson correlation coefficient is positive when both variables x and y increase or decrease at the same time, and is negative when one variable increases and the other one decreases. The lower and upper limits of the Pearson correlation coefficient are −1 and 1, respectively. According to Jaeger et al. [4], if the Pearson coefficient is less than 0.3 absolute, there is a weak correlation between the variables, and if it is between 0.3 and 0.5, there is a medium correlation and if it is larger than 0.5, the correlation is strong.
In order to unambiguously determine the relationship between the input parameters and output results, the Pearson correlation coefficient should be calculated for output variables that are either constant or increase over time. Hydrogen production and corium height in the lower head were thus selected as logical candidates. They are furthermore used to characterize two essential processes during a severe accident, oxidation of the fuel rod cladding and the core degradation. Time-dependent correlation coefficients between the selected input parameters and the two output parameters, accumulated hydrogen mass and the corium height, are shown in Figure 12, Figure 13, Figure 14 and Figure 15.
Pearson correlation coefficients describing the relationships between selected input uncertain parameters, listed in Table 1, used in modelling severe accident core behaviour (CORE parameters) and hydrogen mass, are shown in Figure 12. The values of the correlation coefficients indicate that there is a weak correlation between the input parameters and the output hydrogen production. Around 8500 s, shortly after the accumulators are activated, most of the coefficients change abruptly, with some changing the trend from negative to positive values, and vice versa. However, the results after 8500 s have no major significance because then the geometry of the upper half of the core is lost, while the input parameters refer to modelling of the initial, or only partially damaged core. An increase in the value of the CORE parameter 2 (fraction of oxidation of the fuel rod cladding for the stable ZrO2 oxide layer) leads to an increase in the hydrogen production. When the cladding temperature exceeds the limiting temperature (CORE parameter 1) and the oxidation is of low intensity (CORE parameter 2), the oxide layer will fail, metallic zirconium below that layer will be exposed to steam and the oxidation rate will increase. The influence of CORE parameters 10 (rupture strain for the fuel rod cladding) and 23 (helium inventory in the fuel rod) can be attributed to the same phenomenon. Lower cladding rupture strain and higher helium mass inside the fuel rod will cause earlier cladding rupture, shattering of the zirconium dioxide layer, increased exposure of the metal to steam and possible double-sided oxidation after the penetration of steam into the fuel rod. The impact of CORE parameters 16 (radius of the contact area between the grid spacer and the fuel rod cladding) and 21 (fuel rod outer cladding radius) is more straightforward. If the first parameter increases, the area of the cladding exposed to steam decreases and thus the oxidation rate decreases as well. The oxidation rate also decreases when the second parameter (cladding radius) increases because the flow area between the fuel rods is smaller and less steam is available for oxidation.
The relationships between selected thermal hydraulic input parameters (THA parameters), including decay heat and hydrogen generation, are shown in Figure 13. By far the greatest influence, a strong correlation, on hydrogen production is that of the amount of decay heat. The intensity of oxidation increases exponentially with increasing temperature, which in turn is directly related to heat generation. The increase in initial water volume in the accumulators (THA parameter 3) causes an increase in the mass of produced hydrogen because the evaporation of that water, once it enters the core, contributes to the oxidation process. The influence of the last displayed parameter (THA parameter 1—break junction friction coefficient) cannot be simply explained since this coefficient affects the complex mass and energy balance within the primary system [35]. Discharge of the coolant through the break directly affects the amount of fluid in the primary system, and indirectly the primary pressure, accumulator activation time, heat transfer in the core and in the steam generators. Therefore, the small value of the Pearson correlation coefficient is understandable.
The influence of input parameters on the corium height in the lower head is not as direct as for the hydrogen production because the requirement for relocation is first the formation of a molten pool in the core, followed by cracking of the peripheral crust supporting the pool when the necessary conditions are met. Figure 14 shows the influence of selected CORE parameters, and Figure 15 shows the influence of decay heat and two thermal hydraulic parameters on the corium height. Again, there is a weak correlation, although a little higher than for the hydrogen production, between the two sets of data, except for the decay heat, which has a medium to strong influence on the amount of melt in the lower head. The CORE parameter 5 (surface temperature for freezing of drops of liquefied fuel rod cladding) is characterized by a positive correlation, which means that when the temperature for freezing of the liquefied material is higher, then the accumulation of material in the melt layer will be larger. If the temperature at which fuel rods disintegrate during quench increases, then the core will stay intact longer and degradation and melting will be weaker, which is determined by the CORE parameter 8 (temperature above saturation at which rod fragmentation occurs during quench). The strain limit for the rod-to-rod contact (CORE parameter 12) has a quite substantial impact on the mass of corium in the lower head. As it increases, it blocks the downward movement of the melt along the fuel rods, more and more material can accumulate inside the core, and once the crust fails, there will be a larger relocation. The CORE parameter 13 is the mass of grid spacer per the fuel rod. Its correlation to the corium height is positive because a larger amount of grid spacer material leads to a stronger chemical reaction between the grid and the fuel rod cladding, which, in the end, results in more liquefied material. Finally, the fuel pellet radius (CORE parameter 19) does not have a serious influence on the core relocation and accumulation in the lower head. At the time of relocation, most of the core is completely damaged and the fuel assemblies no longer exist, except their lower part that is better cooled than the upper part, but have been turned into debris or melted. Therefore, this parameter should not have any major impact, which is confirmed by the calculation.
The heat generation inside the core has a strong influence on the radius of the in-core molten pool (r = 0.9) before the relocation, but for the corium height in the lower head, the Pearson correlation coefficient is twice as small. Unfortunately, since the molten pool radius decreases sharply every time when corium relocates to the lower head and increases afterwards as new material is added to the pool, the Pearson correlation coefficient for the relationship between the input parameters and the in-core molten pool is not applicable. Similarly, when accumulators inject more water into the RCS, there is more damage to the fuel assemblies due to thermal shock; however, this is not clearly visible in the diagram in Figure 15 (THA parameter 3). The correlation coefficient between the accumulator water inventory and the debris volume inside the core, before the corium relocation, again shows, although small, a clear enough relationship between these two parameters. Regarding the THA parameter 1 (break junction friction coefficient), the same conclusion can be drawn as for the hydrogen production, that its impact cannot be separated from all other impacts on the mass and energy balance of the system. The Pearson correlation coefficient increases in the late phase of the accident; however, the mass of the released fluid from the primary system is then very small and has no major influence on the process of material accumulation in the lower head of the reactor vessel.

4.3. Comparison of Scenarios with Smaller and Larger Sample Sizes

In addition to the presented analysis, which includes 124 calculations, a preliminary analysis was performed with 59 calculations and a reduced number of uncertain parameters. It contained only 23 uncertain parameters describing the behaviour of the core during a severe accident and listed in Table 1.
Figure 16 shows bounding values of selected output variables for these two calculations. Interestingly, there are no large differences, only of the order of a few percent, for the bounding values at each time step. In the early phase of the accident, before the start of the core damage, there is a small influence of the variation in the thermohydraulic parameters and the core power on the minimum and maximum values of the thermal hydraulic variables, such as primary pressure and the fluid inventory. In the period between 4500 s and 6000 s, i.e., between the termination of the natural circulation in the primary system and the loop seal clearing, the differences between the extrema are 4–5% in the case of the analysis with 124 calculations, and 1% in the analysis with 59 calculations. During the process of core melting, and especially during the relocation of corium, there is no clear trend in the progression of results due to large differences between individual calculations. The difference between the minimum and maximum values in the analysis with 124 calculations is a bit higher than in the analysis with 59 calculations, although at certain time intervals, the latter analysis yields higher maximum, or lower minimum values.
Dispersion of results is determined using the relative standard deviation and the coefficient of range. Mean and maximum values of relative standard deviations and coefficients of range of output variables are shown in Table 2. The numbers in parentheses indicate the number of performed calculations. The coefficient of range is defined as the ratio of difference between the highest and the lowest value of a variable to their sum, (xmaxxmin)/(xmax + xmin). Thus, this coefficient takes into account only minimum and maximum values, regardless of the dispersion of other results.
Table 2. Relative standard deviations and coefficients of range for output variables, mean and maximum values.
Table 2. Relative standard deviations and coefficients of range for output variables, mean and maximum values.
VariableRelative Standard DeviationCoefficient of Range
Mean Value (59)Mean Value (124)Maximum Value (59)Maximum Value (124)Mean Value (59)Mean Value (124)Maximum Value (59)Maximum Value (124)
Pressurizer pressure0.220.170.420.400.380.340.760.82
Collapsed core water level0.750.607.758.050.580.641.01.0
RCS fluid mass0.310.290.730.720.430.450.920.93
Integral break flow0.020.020.040.040.040.040.090.10
Accumulator pressure0.160.110.240.230.260.290.430.44
Production of hydrogen0.080.231.685.780.190.241.01.0
Maximum core temperature0.020.030.160.270.050.070.290.36
Radius of the in-core molten pool0.620.677.756.510.980.991.01.0
Corium height in the lower head0.650.717.7511.180.770.781.01.0
The mean values of relative standard deviations do not indicate any major differences in the results of the two analyses. The only exception is the hydrogen production, for which the mean relative standard deviation in the analysis with 124 runs is three times higher than the mean relative standard deviation in the analysis with 59 runs. This is due to differences in the calculations that occur in the initial phase of the oxidation. However, when looking at the total mass of produced hydrogen, the differences are small. In the case with 59 calculations the mean value is 297.9 kg and the standard deviation 12.6 kg. The confidence interval is (294.7 kg, 301.1 kg) for the 95% confidence limit. In the case with 124 calculations the mean value is 302.4 kg and the standard deviation 14.9 kg. The confidence interval is (299.8 kg, 305.0 kg) for the 95% confidence limit.
Mean standard deviations of the collapsed core water level, RCS fluid mass, pressurizer and accumulators’ pressures are larger for the case with fewer calculations. Regarding the peak values, RSDs for two variables, production of hydrogen and the corium height in the lower head, are clearly higher for the analysis with 124 calculations. These peak values appear at the onset of the oxidation and the relocation of corium, respectively, but the relative standard deviations decrease quickly afterwards, as shown in Figure 11. There are no such large differences between the peak values of the other variables, while the peak value for the radius of the in-core molten pool is higher in the analysis with 59 calculations.
The differences between the coefficients of range are far lower than for the relative standard deviations of the two analyses. This is because the coefficient of range cannot be greater than one. These coefficients have the highest values for variables related to the core damage, radius of the molten pool and the corium height in the lower head, for which the largest differences between individual calculations were determined. The peak values of coefficients of range are equal to one for variables that can be zero during the calculation, hydrogen production, collapsed core water level, corium height in the lower head and the radius of the in-core molten pool.

5. Summary and Conclusions

A comprehensive uncertainty analysis of a severe accident requires the selection of uncertain parameters that cover both thermohydraulic and core damage phenomena, their probability density functions and qualified BEPU methodology. Small variations in the parameters related to the fuel rod cladding deformation, oxide layer stability, core fragmentation, stability of the crust supporting the molten pool, fuel rod dimensions, decay heat, junction friction coefficients and initial accumulator conditions have different influence on the output results of performed calculations in specific phases of the accident.
The course of the accident can be divided into four phases. In the first phase, from the start of the accident until the core dryout, the scattering of results, primary pressure and temperature, is small. In the second phase, which lasts until the moment of accumulator activation, there is a substantial oxidation of the fuel rod cladding and hydrogen production, but the scattering of the results is still limited. The third phase of the accident begins at the time of water injection from the accumulators, which causes interruption of the rapid increase in the oxidation rate, but also leads to the fuel rods’ fragmentation, formation of debris and the in-core molten pool. During this phase, deviations in results are somewhat greater, especially of the molten pool dimensions. The fourth phase that begins with relocation of the corium, is characterized by the largest variations in results. Drainage of the molten pool and relocation of the corium into the lower head affect significantly the reactor coolant system behaviour, accumulator operation and the dynamics of the core damage progression. Even small differences in the timing of the drainage occurrence in different scenarios greatly change the further course of the accident.
Calculation of the Pearson correlation coefficient shows that there is no straightforward relationship between the output variables (hydrogen production and the corium height in the lower head) and the selected input uncertain parameters. Although for almost all input parameters there is a weak correlation between input values and output results, except for the decay heat where there is a strong correlation, the physical relationship between each parameter and hydrogen production, i.e., the amount of corium, can be clearly established. Pearson correlation coefficients related to hydrogen production change suddenly at the time of disintegration of fuel assemblies, but since most of the input parameters are associated with the modelling of the initial, or only partially damaged core, this effect is not important. The effect of thermohydraulic parameters is complex because they primarily affect the global balance of mass and energy in the primary system, and in this sense their influence is not particularly pronounced. The largest influence on hydrogen production, and partly on the process of the corium relocation, is had by the decay heat because the oxidation rate increases exponentially with temperature, which in turn is directly proportional to the heat released in the core. The relocation of the corium to the lower head is preceded by the molten pool formation, and it would be natural to look for a correlation between input parameters and the pool diameter, but as the pool diameter can increase and decrease over time, this dependence could not be expressed using the Pearson correlation coefficient over the full time interval. When looking only at the period before the first corium relocation, the correlation between the decay heat and the molten pool radius is very strong, while for the corium height in the lower head, this correlation is twice as small.
A comparison of the bounding values of the analyses results with 59 calculations and 124 calculations does not indicate large differences between the two data sets. The dispersion of results, quantified by the relative standard deviation and the coefficient of range, is also relatively small. The only exceptions are the hydrogen production and the corium height in the lower head for which deviations are larger in the case with 124 calculations. However, these differences are more pronounced in the initial phases of hydrogen production and corium relocation and are rapidly decreasing as the accident progresses. It can be concluded that increasing the number of calculations is not crucial to obtain a relevant sample if the choice of uncertain parameters and the sample size are consistent with qualified uncertainty analysis methodology.

Author Contributions

Conceptualization, S.Š. and D.G.; methodology, C.A. and M.P.-F.; software, M.P.-F.; validation, C.A.; formal analysis, S.Š.; investigation, D.G.; resources, C.A.; data curation, S.Š.; writing—original draft preparation, S.Š. and D.G.; writing—review and editing, D.G. and M.P.-F.; visualization, S.Š.; supervision, D.G. All authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Not applicable.

Acknowledgments

The authors express their gratitude to the Krško NPP for providing the plant operational and design data.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Bazin, P.; de Crécy, A.; Glaeser, H.; Skorek, T.; Joucla, J.; Probst, P.; Chung, B.; Oh, D.-Y.; Kyncl, M.; Pernica, R.; et al. BEMUSE Phase III Report–Uncertainty and Sensitivity Analysis of the LOFT L2-5 Test, NEA/CSNI/R(2007)4; OECD NEA: 2007. Available online: https://www.oecd-nea.org/jcms/pl_18442/bemuse-phase-iii-report-uncertainty-and-sensitivity-analysis-of-the-loft-l2-5-test?details=true (accessed on 25 February 2022).
  2. Reventós, F.; Batet, L.; Pérez, M. BEMUSE Phase V Report–Uncertainty and Sensitivity Analysis of a LB-LOCA in ZION Nuclear Power Plant. NEA/CSNI/R(2009)13; OECD NEA: 2009. Available online: https://www.oecd-nea.org/jcms/pl_18866/bemuse-phase-v-report-uncertainty-and-sensitivity-analysis-of-a-lb-loca-in-zion-nuclear-power-plant?details=true (accessed on 25 February 2022).
  3. Frepoli, C. An Overview of Westinghouse Realistic Large Break LOCA Evaluation Model. Sci. Technol. Nucl. Install. 2008, 498737, 15. [Google Scholar] [CrossRef]
  4. Jaeger, W.; Sanchez Espinoza, V.H.; Montero Mayorga, F.J.; Queral, C. Uncertainty and Sensitivity Investigations with TRACE-SUSA and TRACE-DAKOTA by Means of Post-Test Calculations of NUPEC BFBT Experiments. NUREG/IA-0462. U.S. NRC. 2017. Available online: https://www.nrc.gov/reading-rm/doc-collections/nuregs/agreement/ia0462/index.html (accessed on 25 February 2022).
  5. Gosh, T. State-of-the-Art Reactor Consequence Analyses (SOARCA) Project Sequoyah Integrated Deterministic and Uncertainty Analyses; Prepared by: Severe Accident Analysis Department Sandia National Laboratories: NUREG/CR-7245, U.S. NRC. 2019. Available online: https://www.nrc.gov/reading-rm/doc-collections/nuregs/contract/cr7245/index.html (accessed on 25 February 2022).
  6. Ang, M.L.; Grindon, E.; Dutton, L.M.C.; Garcia-Sedano, P.; Santamaria, C.S.; Centner, B.; Auglaire, M.; Routamo, T.; Outa, S.; Jokiniemi, J.; et al. A risk-based evaluation of the impact of key uncertainties on the prediction of severe accident source terms—STU. Nucl. Eng. Des. 2001, 209, 183–192. [Google Scholar] [CrossRef]
  7. Elsalamouny, N.; Kaliatka, T. Uncertainty Quantification of the PHEBUS FPT-1 Test Modelling Results. Energies 2021, 14, 7320. [Google Scholar] [CrossRef]
  8. Darnowski, P.; Mazgaj, P.; Włostowski, M. Uncertainty and Sensitivity Analysis of the In-Vessel Hydrogen Generation for Gen-III PWR and Phebus FPT-1 with MELCOR 2.2. Energies 2021, 14, 4884. [Google Scholar] [CrossRef]
  9. Sonnenkalb, M. Application of the Integral Code MELCOR for German NPPs and use within Accident Management and PSA Projects. In Proceedings of the Technical Meeting on Severe Accident and Accident Management, Tokyo, Japan, 14–16 March 2006. [Google Scholar]
  10. Nowack, H.; Chatelard, P.; Chailan, L.; Hermsmeyer, S.; Sanchez, V.; Herranz, L. CESAM–Code for European severe accident management, EURATOM project on ASTEC improvement. Ann. Nucl. Energy 2018, 116, 128–136. [Google Scholar] [CrossRef]
  11. Šadek, S.; Grgić, D.; Šimić, Z. Application of ASTEC, MELCOR, and MAAP Computer Codes for Thermal Hydraulic Analysis of a PWR Containment Equipped with the PCFV and PAR Systems. Sci. Technol. Nucl. Install. 2017, 2017, 8431934. [Google Scholar] [CrossRef] [Green Version]
  12. Magallon, D.; Mailliat, A.; Seiler, J.-M.; Atkhen, K.; Sjövall, H.; Dickinson, S.; Jakab, J.; Meyer, L.; Buerger, M.; Trambauer, K.; et al. European expert network for the reduction of uncertainties in severe accident safety issues (EURSAFE). Nucl. Eng. Des. 2005, 235, 309–346. [Google Scholar] [CrossRef]
  13. Wilson, G.E.; Boyack, B.E. The role of the PIRT process in experiments, code development and code applications associated with reactor safety analysis. Nucl. Eng. Des. 1998, 186, 23–37. [Google Scholar] [CrossRef]
  14. Shaw, R.A.; Larson, T.K.; Dimenna, R.K. Development of a Phenomena Identification and Ranking Table (PIRT) for Thermal-Hydraulic Phenomena during a PWR Large-Break LOCA; NUREG/CR-5074, EGG-2524; EG&G Idaho, Inc.: Gaithersburg, MD, USA, 1988. [Google Scholar]
  15. Boyack, B.E.; Motta, A.T.; Peddicord, K.L.; Alexander, C.A.; Andersen, J.G.M.; Blaisdell, J.A.; Dunn, B.M.; Ebeling-Koning, D.; Fuketa, T.; Hache, G.; et al. Phenomena Identification and Ranking Tables (PIRTs) for Loss of Coolant Accidents in Pressurized and Boiling Water Reactors Containing High Burnup Fuel; NUREG/CR-6744, LA-UR-00-5079, U.S. NRC. 2001. Available online: https://www.nrc.gov/reading-rm/doc-collections/nuregs/contract/cr6744/index.html (accessed on 25 February 2022).
  16. IAEA. Best Estimate Safety Analysis for Nuclear Power Plants: Uncertainty Evaluation; Reports Series No. 52. 2008. Available online: https://www.iaea.org/publications/7768/best-estimate-safety-analysis-for-nuclear-power-plants-uncertainty-evaluation (accessed on 25 February 2022).
  17. Boyack, B.; Duffey, R.; Griffith, P.; Lellouche, G.; Levy, S.; Rohatgi, U.; Wilson, G.; Wulff, W.; Zuber, N. Quantifying Reactor Safety Margins–Application of Code Scaling, Applicability and Uncertainty Evaluation Methodology to a Large-Break, Loss-of-Coolant Accident; NUREG/CR-5249, EGG-2552; U.S. NRC: Washington, DC, USA, 1989.
  18. Glaeser, H. GRS Method for Uncertainty and Sensitivity Evaluation of Code Results and Applications. Sci. Technol. Nucl. Install. 2008, 2008, 798901. [Google Scholar] [CrossRef] [Green Version]
  19. Wickett, T.; Sweet, D.; Neill, A.; D’ Auria, F.; Galassi, G.; Belsito, S.; Ingegneri, M.; Gatta, P.; Glaeser, H.; Skorek, T.; et al. Report of the Uncertainty Methods Study for Advanced Best Estimate Thermal Hydraulic Code Applications; NEA/CSNI/R(97)35, OECD NEA: 1998; Volumes 1–2. Available online: https://www.oecd-nea.org/upload/docs/application/pdf/2020-01/csni-r1997-35-vol2.pdf (accessed on 25 February 2022).
  20. Perez Ferragut, M. Integration of a Quantitative-Based Selection Procedure in an Uncertainty Analysis Methodology for NPP Safety Analysis. Ph.D. Thesis, UPC-Universitat Politècnica de Catalunya, Barcelona, Spain, September 2011. [Google Scholar]
  21. Ludmann, M.; Sauvage, J.-Y. LB LOCA analysis using the deterministic realistic methodology—Application to the 3-loop plant. In Proceedings of the 7th International Conference on Nuclear Engineering (ICONE-7), Tokyo, Japan, 19–23 April 1999; Japan Society of Mechanical Engineers: Tokyo, Japan, 1999. [Google Scholar]
  22. D’Auria, F.; Giannotti, W. Development of a Code with the Capability of Internal Assessment of Uncertainty. Nucl. Technol. 2000, 131, 159–196. [Google Scholar] [CrossRef]
  23. D’Auria, F.; Debrecin, N.; Galassi, G.M. Outline of the Uncertainty Methodology Based on Accuracy Extrapolation. Nucl. Technol. 1995, 109, 21–38. [Google Scholar] [CrossRef]
  24. ASYST VER 3 User Reference Manuals; SDTP/ADTP; Innovative Systems Software: Ammon, ID, USA, 2020.
  25. Allison, C.M.; Hohorst, J.K. Role of RELAP/SCDAPSIM in Nuclear Safety. Sci. Technol. Nucl. Install. 2010, 2010, 425658. [Google Scholar] [CrossRef]
  26. Nakata, A.; Naitoh, M.; Allison, C.M. Need of a Next Generation Severe Accident Code. J. Nucl. React. Technol. Tri Dasa Mega 2019, 21, 119–126. [Google Scholar] [CrossRef] [Green Version]
  27. Nakata, A.; Allison, C.M.; Hohorst, J.K.; Naitoh, M.; Matsubara, K.; Pericas, R. Development and Preliminary Assessment of the new ASYST Integral Analysis BEPU Code using the PBF SFD-ST Bundle Heating and Melting Experiment, a Typical BWR Under Fukushima-Daiichi-Accident-Like Thermal Hydraulic Conditions and PWR for a Steam Line Break in the Containment. In Proceedings of the ICAPP 2020, Abu Dhabi, United Arab Emirates, 15-19 March 2020. [Google Scholar]
  28. Trivedi, A.K.; Novog, D.R.; Allison, C.M. Implementation of Solar Salt as Fluid in ASYST4.1 and Validation for a Natural Circulation Loop. In Proceedings of the ICONE 28, Online, 4–6 August 2021. [Google Scholar]
  29. Trivedi, A.K.; Perez-Ferragut, M.; Hohorst, J.K.; Allison, C. Integrated uncertainty analysis of LBLOCA in AP1000 using RELAP/SCDAPSIM. In Proceedings of the 12th International Conference of the Croatian Nuclear Society, Zadar, Croatia, 3–6 June 2018. [Google Scholar]
  30. Trambauer, K.; Haste, T.J.; Adroguer, B.; Hozer, Z.; Magallon, D.; Zurita, A. In-vessel Core Degradation Code Validation Matrix, Update 1996–1999; NEA/CSNI/R(2000)21, OECD NEA: 2001. Available online: https://www.oecd-nea.org/upload/docs/application/pdf/2020-01/csni-r2000-21.pdf (accessed on 25 February 2022).
  31. Sehgal, B.R. (Ed.) Nuclear Safety in Light Water Reactors-Severe Accident Phenomenology, 1st ed.; Elsevier Inc.: Amsterdam, The Netherlands, 2012; ISBN 9780123884466. [Google Scholar]
  32. Perez, M.; Reventos, F.; Wagner, R.; Allison, C. Integrated Uncertainty Analysis using RELAP/SCDAPSIM/MOD4.0. In Proceedings of the 13th International Topical Meeting on Nuclear Reactor Thermal Hydraulics (NURETH-13), Kanazawa City, Japan, 27 September–2 October 2009. [Google Scholar]
  33. Benčik, V.; Grgić, D.; Šadek, S.; Vlahović, Š. NEK RELAP5/MOD3.3 Steady State Qualification Report for Cycle 29; NEK ESD-TR-13/16, FER-ZVNE/SA/DA-TR02/16-0; Faculty of Electrical Engineering and Computing: Zagreb, Croatia, 2016. [Google Scholar]
  34. Grgić, D.; Benčik, V.; Šadek, S. NEK RELAP5/MOD3.3 Post-RTDBE Nodalization Notebook; NEK ESD-TR-02/13, FER-ZVNE/SA/DA-TR03/13-0; Faculty of Electrical Engineering and Computing: Zagreb, Croatia, 2013. [Google Scholar]
  35. Šadek, S.; Grgić, D.; Benčik, V.; Vlahović, Š. Analysis of the upflow conversion modification and influence on selected LOCA accidents in a PWR plant. Nucl. Eng. Des. 2020, 369, 110854. [Google Scholar] [CrossRef]
  36. SCDAP/RELAP5/MOD3.2 Code Manual, Volume 1–5; Prepared by: SCDAP/RELAP5 Development Team. NUREG/CR-6150, INEL-96/0422, U.S. NRC. 1998. Available online: https://www.nrc.gov/docs/ML0103/ML010310397.pdf (accessed on 25 February 2022).
  37. Hermann, O.W.; Parks, C.V.; Renier, J.P. Technical Support for a Proposed Decay Heat Guide Using SAS2H/ORIGEN-S Data; NUREG/CR-5625, ORNL-6698. U.S. NRC. 1994. Available online: https://www.osti.gov/servlets/purl/43752 (accessed on 25 February 2022).
  38. Schrock, V.E. A Revised ANS Standard for Decay Heat from Fission Products. Nucl. Technol. 1979, 46, 323–331. [Google Scholar] [CrossRef]
  39. Wilks, S.S. Determination of sample sizes for setting tolerance limits. Ann. Math. Stat. 1941, 12, 91–96. [Google Scholar] [CrossRef]
  40. Nutt, W.T.; Wallis, G.B. Evaluation of nuclear safety from the outputs of computer codes in the presence of uncertainties. Reliab. Eng. Syst. Saf. 2004, 83, 57–77. [Google Scholar] [CrossRef]
  41. Code of Federal Regulations 10 CFR 50.46: Acceptance Criteria for Emergency Core Cooling Systems for Light-Water Nuclear Power Reactors. U.S. NRC. 1988. Available online: https://www.nrc.gov/reading-rm/doc-collections/cfr/part050/part050-0046.html (accessed on 25 February 2022).
  42. Guba, A.; Makai, M.; Pál, L. Statistical aspects of best estimate method—I. Reliab. Eng. Syst. Saf. 2003, 80, 217–232. [Google Scholar] [CrossRef]
  43. Wallis, G.B. Uncertainties and probabilities in nuclear reactor regulation. Nucl. Eng. Des. 2007, 237, 1586–1592. [Google Scholar] [CrossRef]
  44. Broughton, J.M.; Kuan, P.; Petti, D.A.; Tolman, E.L. A Scenario of the Three Mile Island Unit 2 Accident. Nucl. Technol. 1989, 87, 34–53. [Google Scholar] [CrossRef]
Figure 1. ASYST nodalization of NPP Krško, including a five-region core model.
Figure 1. ASYST nodalization of NPP Krško, including a five-region core model.
Energies 15 01842 g001
Figure 2. Pressurizer pressure.
Figure 2. Pressurizer pressure.
Energies 15 01842 g002
Figure 3. Collapsed water level in the core.
Figure 3. Collapsed water level in the core.
Energies 15 01842 g003
Figure 4. Primary system mass.
Figure 4. Primary system mass.
Energies 15 01842 g004
Figure 5. Accumulator pressure.
Figure 5. Accumulator pressure.
Energies 15 01842 g005
Figure 6. Production of hydrogen.
Figure 6. Production of hydrogen.
Energies 15 01842 g006
Figure 7. Maximum core temperature.
Figure 7. Maximum core temperature.
Energies 15 01842 g007
Figure 8. Radius of the in-core molten pool.
Figure 8. Radius of the in-core molten pool.
Energies 15 01842 g008
Figure 9. Corium height in the lower head.
Figure 9. Corium height in the lower head.
Energies 15 01842 g009
Figure 10. Relative standard deviations for thermal hydraulic variables.
Figure 10. Relative standard deviations for thermal hydraulic variables.
Energies 15 01842 g010
Figure 11. Relative standard deviations for variables related to core damage description.
Figure 11. Relative standard deviations for variables related to core damage description.
Energies 15 01842 g011
Figure 12. Pearson correlation coefficient between parameters related to core damage propagation and hydrogen production.
Figure 12. Pearson correlation coefficient between parameters related to core damage propagation and hydrogen production.
Energies 15 01842 g012
Figure 13. Pearson correlation coefficient between thermal hydraulic input parameters and hydrogen production.
Figure 13. Pearson correlation coefficient between thermal hydraulic input parameters and hydrogen production.
Energies 15 01842 g013
Figure 14. Pearson correlation coefficient between parameters related to core damage propagation and corium height in the lower head.
Figure 14. Pearson correlation coefficient between parameters related to core damage propagation and corium height in the lower head.
Energies 15 01842 g014
Figure 15. Pearson correlation coefficient between thermal hydraulic input parameters and corium height in the lower head.
Figure 15. Pearson correlation coefficient between thermal hydraulic input parameters and corium height in the lower head.
Energies 15 01842 g015
Figure 16. Bounding values of output variables for 59 and 124 calculations: (a) Pressurizer pressure; (b) collapsed water level in the core; (c) primary system mass; (d) integral mass flow rate at the break (loop 1); (e) accumulator pressure; (f) production of hydrogen; (g) maximum core temperature; (h) radius of the in-core molten pool; (i) corium height in the lower head.
Figure 16. Bounding values of output variables for 59 and 124 calculations: (a) Pressurizer pressure; (b) collapsed water level in the core; (c) primary system mass; (d) integral mass flow rate at the break (loop 1); (e) accumulator pressure; (f) production of hydrogen; (g) maximum core temperature; (h) radius of the in-core molten pool; (i) corium height in the lower head.
Energies 15 01842 g016
Table 1. Input uncertain parameters.
Table 1. Input uncertain parameters.
PhenomenonNo.DescriptionReference ValueProbability Density Function
Severe accident core behaviour (CORE)1Temperature for failure of the oxide layer on the outer cladding surface2500 KNormal (1.00, 10−4)
2Fraction of oxidation of the fuel rod cladding for the stable oxide layer0.2Normal (1.00, 10−4)
3Cladding hoop strain threshold for double-sided oxidation0.07Normal (1.00, 10−4)
4Fraction of cladding surface area covered with drops that results in the blockage that stops local oxidation0.2Normal (1.00, 10−4)
5Surface temperature for freezing of drops of liquefied fuel rod cladding1750 KNormal (1.00, 10−4)
6Velocity of drops of cladding material slumping down outside the surface of the fuel rod0.5 m/sNormal (1.00, 10−4)
7Multiplication factor on fuel pellet diameter that defines minimum thickness that the crust must have in order to support the molten pool1.0Normal (1.00, 10−4)
8Temperature above saturation at which rod fragmentation occurs during quench100 KNormal (1.00, 10−4)
9Gamma heating fraction0.026Normal (1.00, 10−4)
10Rupture strain for the fuel rod cladding0.18Uniform (0.99, 1.01)
11Strain for transition from the sausage type cladding deformation to the localized deformation0.2Uniform (0.99, 1.01)
12Strain limit for the rod-to-rod contact0.22Uniform (0.99, 1.01)
13Mass of the grid spacer per the fuel rod3.23·10−3 kgNormal (1.00, 10−4)
14Height of the grid spacer0.0336 mNormal (1.00, 10−4)
15Plate thickness of the grid spacer4·10−4 mNormal (1.00, 10−4)
16Radius of the contact area between the grid spacer and the fuel rod cladding3.39·10−3 mNormal (1.00, 10−4)
17Fuel rod plenum length0.186 mNormal (1.00, 10−4)
18Fuel rod plenum void volume9.56·10−6 m3Normal (1.00, 10−4)
19Fuel rod pellet radius4.096·10−3 mUniform (1.00, 1.005)
20Fuel rod inner cladding radius4.178·10−3 mUniform (1.00, 1.005)
21Fuel rod outer cladding radius4.75·10−3 mUniform (1.00, 1.005)
22Fraction of fuel theoretical density0.95Uniform (0.99, 1.01)
23Helium inventory in the fuel rod6.2·10−5 kgNormal (1.00, 10−4)
Thermal hydraulic system behaviour (THA)1Break junction friction coefficient1.0Uniform (0.80, 1.20)
2Steam generator safety valve friction coefficient1.0Uniform (0.80, 1.20)
3Accumulator initial water volume35.94 m3Uniform (0.98, 1.02)
4Accumulator initial pressure4.928 MPaUniform (0.98, 1.02)
5Accumulator initial temperature322 KUniform (0.98, 1.02)
6RCP moment of inertia2.697·103 kgm2Uniform (0.98, 1.02)
7, 8RCP friction torque (RFT) coefficients (TFC1, TFC2):
RFT = TFC1∙S + TFC2∙S2,
where S is the ratio of the current pump rotational velocity to the rated pump rotational velocity
3556.61 Nm,
−715.41 Nm
Uniform (0.98, 1.02),
Uniform (0.98, 1.02)
Decay heatPowerCore decay heat vs. timeANS79-1 standard dataUniform (0.98, 1.02)
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Šadek, S.; Grgić, D.; Allison, C.; Perez-Ferragut, M. Uncertainty Study of the In-Vessel Phase of a Severe Accident in a Pressurized Water Reactor. Energies 2022, 15, 1842. https://doi.org/10.3390/en15051842

AMA Style

Šadek S, Grgić D, Allison C, Perez-Ferragut M. Uncertainty Study of the In-Vessel Phase of a Severe Accident in a Pressurized Water Reactor. Energies. 2022; 15(5):1842. https://doi.org/10.3390/en15051842

Chicago/Turabian Style

Šadek, Siniša, Davor Grgić, Chris Allison, and Marina Perez-Ferragut. 2022. "Uncertainty Study of the In-Vessel Phase of a Severe Accident in a Pressurized Water Reactor" Energies 15, no. 5: 1842. https://doi.org/10.3390/en15051842

APA Style

Šadek, S., Grgić, D., Allison, C., & Perez-Ferragut, M. (2022). Uncertainty Study of the In-Vessel Phase of a Severe Accident in a Pressurized Water Reactor. Energies, 15(5), 1842. https://doi.org/10.3390/en15051842

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop