Table of Contents Author Guidelines Submit a Manuscript
Science and Technology of Nuclear Installations
Volume 2018, Article ID 5706409, 9 pages
Research Article

Effect of Molten Corium Behavior Uncertainty on the Severe Accident Progress

1Department of Nuclear Engineering, Hanyang University, Seoul 04763, Republic of Korea
2Institute of Nano Science and Technology, Hanyang University, Seoul 04763, Republic of Korea

Correspondence should be addressed to Sung Joong Kim;

Received 22 March 2018; Revised 13 June 2018; Accepted 25 June 2018; Published 1 August 2018

Academic Editor: Jerzy Foit

Copyright © 2018 Wonjun Choi et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.


Uncertainty of a severe accident code output needs to be handled reliably considering its use in safety regulation of nuclear industry. In particular, severe accident codes are utilized for probabilistic safety assessment (PSA), where the uncertainty of severe accident progress should be considered carefully due to its influence on human reliability analysis. Therefore, in this study, the uncertainty analysis of severe accident progress was performed using MELCOR code, and a total of 200 data sets of in-vessel uncertainty parameters were generated by Latin hypercube sampling method. The rank regression analysis was also performed to investigate the effect of uncertainty parameters on the severe accident progress. Sensitivity coefficients (SCs) in MELCOR such as molten clad drainage rate and zircaloy melt breakout temperature showed significant influence on relocation time and dryout time of lower plenum. However, the influence of uncertainty parameter diminished as the accident progressed.

1. Introduction

Since Fukushima accident in 2011, numerous studies about the severe accident (SA) phenomena and their progress have been investigated experimentally and numerically. However, due to the large scale of nuclear power plants, not many experimental studies were conducted for the accident sequence. Therefore the researches about the SA progress were conducted mostly by the numerical analysis [15]. As far as the numerical analysis of SA progress is concerned, lumped parameter SA codes such as MELCOR [6], MAAP [7], and SCDAP/RELAP [8] have been utilized in general. However, it is widely accepted that the SA code results may generate large uncertainty because rather simplified models were adopted in the codes due to the incomplete knowledge about SA phenomena. Because the code results can be reutilized for calculating source terms in the frame of PSA level 3, it is very important to investigate the uncertainty range and to identify relative importance of input parameters [9, 10].

Among the code outputs, identifying timing of certain SA sequence is important in terms of human reliability analysis (HRA). As a practical example, the success probability of human related actions can be defined as a function of available time in HRA [11]. Thus, the timing of SA progress is important to define the available time for certain human related actions. Therefore, the uncertainty of the SA progress in terms of the major sequence time needs to be analysed for better understanding on the nuclear safety. For this reason, the major objective of this study is to perform the uncertainty analysis about occurrence time of SA progression phenomena. A representative SA scenario was selected as short term station black out (STSBO). STSBO is the scenario that all power, even direct current (DC) batteries, is unavailable. Because the safety feature using batteries is also unavailable, the accident progress can be faster than other accident scenarios.

2. Numerical Methodologies

2.1. Description of MELCOR Code [6]

MELCOR code is an integrated SA code for various kinds of nuclear power reactors. Since 1982, MELCOR code has been developed by Sandia National Laboratories (SNL) and used for plant risk assessment and source term analysis. To simulate various SA phenomena, a number of modular packages were coupled in a unified frame. By the interaction of each modular package, MELCOR code can simulate thermal-hydraulic response of the plant, heat-up, oxidation, degradation and relocation of the core, and fission product release and transport. Because of the lack of knowledge about the SA phenomena and the interest of quick code execution time, parametric equations were used to model complicated physical processes in the past. However most recent MELCOR models are mechanistic thank to the improved calculation speed of modern computers. In many of the mechanistic models, adjustable optional parameters can be also available for uncertainty analyses and sensitivity studies. For that reason, MELCOR version 2.1.6342 was used.

2.2. MELCOR Input Model of OPR1000

As a reference plant, Optimized Power Reactor 1000 MWe (OPR1000), which consists of majority of Korean operating NPP, was selected. The nodalization of OPR1000 input model is shown in Figure 1. The input model includes reactor coolant system (RCS) composed of 25 control volumes (CVs) and simple secondary side composed of two steam generators (SGs) and a turbine. RCS consists of two hot legs (CV310 and 410), four cold legs (CV380, 390, 480 and 490), pressurizer (CV500), and reactor pressure vessel (RPV) with 6 detailed volumes (CV130, 150, 180, 170, 190, and 260). At the top of a pressurizer, pressurizer safety relief valves (PRSVs) and safety depressurization system (SDS) were modeled to simulate release of RCS coolant in case of depressurization. The containment (CNMT) consists of 20 CVs and the total free volume is 77,440 m3. Figure 2 shows the nodalization of the core. Total 7 ~ 8 radial rings and 14 axial nodes were modeled for the detailed simulation of the reactor core. 1st ~ 3rd axial core nodes were coupled with CV150 (LP), and 4th ~ 14th axial core nodes were coupled with CV170 (core region). In the core region, only 4th ~ 13th axial core nodes contained the nuclear fuel. The supporting structures were located at the 1st, 2nd, 3rd, 6th, 9th, and 12th axial node, and support plate was modeled in 3rd axial node. Total initial mass of uranium fuel, zirconium cladding, and stainless steel in the core were 85,692.0 kg, 23,910.58 kg, and 20,024.4 kg, respectively.

Figure 1: The nodalization of OPR1000 MELCOR input model.
Figure 2: The core nodalization of OPR1000 MELCOR input model.
2.3. Description of Uncertainty Parameters and Distribution

In SOARCA project [9], the uncertainty analysis of source term was performed using MELCOR code. Based on PIRT process by SNL and nuclear regulatory commission (NRC), total 24 epistemic MELCOR parameters were selected for uncertain input parameters of Surry power plant. To study the effect of molten corium behavior on the severe accident progress, 5 of uncertain parameters in COR package were selected in this study. The detailed distribution of uncertain parameters referred to the uncertainty analysis of Surry plant in SOARCA project. Table 1 shows the selected uncertain parameters and distributions. The lower/upper bound of zircaloy melt breakout temperature (SC1131(2)) was defined by melting point of zircaloy (2100 K) and the value of hydrogen uncertainty study [12] based on Phebus experiments [13]. The distribution was updated from the distribution of Peach bottom plants with the S/Q simulation results [14]. Because the technical base for the boundary values of molten clad drainage rate (SC1141(2)) and radial solid/molten debris relocation time constant (SC1020) were insufficient, the boundaries of these uncertainty parameters were selected based on order of magnitude. The default value of molten clad drainage rate was determined with CORA-13 experiment. Therefore, log triangular distribution, that is, the mode value is this default value, was used for the distribution of molten clad drainage rate. The uniform distribution was used for radial solid/molten debris relocation time constant due to the lack of experimental data. The distribution of Effective temperature at which the eutectic formed from UO2 and ZrO2 melts was created based on the mean and standard deviation of 6 test results from VERCORS experiment [15].

Table 1: The information of PDF about selected uncertainty parameters.
2.4. Description of Sampling and Quantification of Results

Symbolic nuclear analysis package (SNAP) is a graphical user interface for simplifying the task of creating input files for the analytic codes and helping to visualize code results [16]. Input files of MELCOR code can also be created by SNAP. For uncertainty and sensitivity analysis, SNAP provides the additional plugin called “DAKOTA”. Using DAKOTA plugin, users can create the random samples of each sensitivity coefficient by LHS method and combine them into an input data set. The correlation between uncertain parameters was not considered in this study. Variance inflation factor and correlation coefficient between sampled input parameters were checked to prevent significant correlation between input parameters by random sampling. In this study, total 200 cases were generated and simulated until 259,200 s (72 hrs). For the quantification of the effect of the selected uncertain parameters on the accident progress, the time of initial melt relocation, LP dryout, RPV failure, and CNMT leak were investigated as major output parameters. In addition, the statistics of the major output parameter were calculated. To investigate the individual effect of uncertain input parameters, linear rank regression analysis was also performed.

3. Results and Discussion

3.1. Description of Base Case

To investigate the uncertainty of SA progress, the aforementioned STSBO was selected as a base case. All of safety features were assumed to fail by losing AC power, and steam driven auxiliary feed water (AFW) pump was also modeled to fail. Table 2 shows the SA sequence of the base case. The STSBO was initiated by losing AC power and subsequently the reactor was tripped. With the reactor tripped, reactor coolant pumps (RCPs) and main feed water (MFW) pumps stopped. Following the reactor trip, the decay heat from the core was removed by the natural circulation of primary coolant during early phase of the accident. Because the AFW pumps were unavailable by the assumption, all SGs were exhausted at 3,762 s by the heat from the RCS. The RCS pressure increased with the loss of secondary heat removal, and the PSRV was opened at 4,881 s as shown in Figure 3. The reactor coolant was released to the CNMT through the PSRV, and thus it caused the core uncovery at 5,359 s. The exposed core was heated up due to the lack of coolant, and additional oxidation heat generated between high temperature steam and high cladding temperature (>1,100 K) accelerated the core heat-up from 8,991 s. The coolant in the core region was dried out at 9,990 s and dryout of the LP region occurred at 12,866 s due to the relocated core debris as shown in Figure 4. After the first LP dryout, thru-wall yielding mode of RPV failure was initiated by the heat of relocated core debris.

Table 2: The accident sequence of base case.
Figure 3: RCS pressure in the base case (STSBO).
Figure 4: Core water level in the base case (STSBO).

When the RPV failure was initiated, the RCS pressure reached to 16.25 MPa (Figure 3). Due to the high pressure, the core debris was released into the CNMT via high pressure melt ejection (HPME) phenomena. Figure 5 shows the ejected debris mass from the RPV, debris mass in the cavity, and mass of debris ejected to CNMT atmosphere by the HPME. At first ejection with RPV failure (14,462 s), 6,341 kg of the 8,819 kg of debris ejected from the RPV was ejected to CNMT atmosphere. Because the RCS pressure decreased after the first ejection, the ejected debris mass following the first ejection only contributes to the increase of debris mass retained in the cavity. After about 19,200 s, 119,870 kg of core debris was ejected to the cavity due to the collapse of remained fuel structure. Figure 6 shows the CNMT pressure during the STSBO accident. Before the RPV failure, the CNMT pressure was repeated to increase and decrease by the opening and closure of the PSRV. After the RPV failure, 0.26 MPa of pressure peak was observed due to the direct containment heating (DCH) by the HPME. Because the debris within the cavity released the heat and noncondensable gas by molten corium concrete interaction (MCCI), the CNMT pressure increased continuously and reached to design pressure (0.44 MPa) at 60,100 s. The CNMT pressure at the end of calculation (259,200 s) was calculated as 1.17 MPa.

Figure 5: Ejected debris mass after RPV failure.
Figure 6: CNMT pressure in the base case (STSBO).
3.2. Statistics of Major Accident Progress

Among a total of 200 sample cases, only 162 cases were calculated until 259,000 s (72 hrs). Because of COR or CAV package errors, 37 cases were stopped as soon as the RPV failure was reported. A case was calculated until the initial melt relocation. However, it was terminated before the RPV failure. For that reason, the statistics of accident progress such as mean, coefficient of variance (C.V.), median, 5th percentile (P5), and 95th percentile (P95) were calculated through the available cases. Table 3 shows the statistical values of accident progress. The selected uncertainty parameters affected the accident progress after cladding melt. Therefore, the statistical values of initial melt relocation, first LP dryout, RPV failure, and CNMT overpressurization time were investigated. The mean values of initial melt relocation, first LP dryout, RPV failure, and CNMT overpressurization were calculated as 10,445 s, 13,421 s, 15,850 s, and 60,100 s, respectively. All mean values of accident progress were larger than those of the base case. The values of C.V. show that the variance of accident progress increased as the accident progressed. The median values of initial melt relocation, first LP dryout, RPV failure, and CNMT overpressurization were estimated as 10,491 s, 13,466 s, 15,731 s, and 58,932 s, respectively. Median values of initial melt relocation and first LP dryout time were larger than mean values of those. It means that the distributions were biased toward the high values. On the other hand, RPV failure and CNMT overpressurization time were biased toward the low values. Figure 7 shows these biased distributions of accident progress. Counts on Y-axis corresponded to the number of samples that was included in specific range about the time of accident progression phenomena on X-axis. The 5th/95th values of initial melt relocation, first LP dryout, RPV failure, and CNMT overpressurization were 10,244 s/10,601 s, 12,831 s/13,811 s, 14,911 s/17,252 s, and 55,882 s/69,272 s, respectively. Using bias-corrected and accelerated bootstrap method, 95% confidence interval of 5th/95th values was also calculated. The confidence intervals for 5th/95th values of initial melt relocation, first LP dryout, RPV failure, and CNMT overpressurization were 10,236 ~ 10,246 s/10,586 ~ 10,651 s, 12,676 ~ 12,931 s/13,781 ~ 13,841 s, 14,826 ~ 14,973 s/17,192 ~ 17,304 s, and 55,621 ~ 56,187 s/68,798 ~ 69.639 s, respectively.

Table 3: The statistics values of accident progress.
Figure 7: Histogram of major accident progress.
3.3. Regression Analysis of Uncertainty Parameters on Accident Progress

Although the ordinary linear regression can show the potential relations between input and output values, capturing any complex relationship is challenging in reality [9]. Because the accident progress is the cumulative result of various phenomena, the relationship between in-vessel uncertain parameters and accident progress may be nonlinear or complex. Thus, in this study, rank regression technique was used to investigate the potential relations between in-vessel uncertainty and accident progress. The rank regression is a kind of linear regression, which uses rank values rather than original values. Using the rank values, the rank regression can capture the potential relations including nonlinear influence [17, 18]. For the rank regression results, standardized rank regression coefficient (SRRC) and partial rank correlation coefficient (PRCC) were summarized in Figure 8.

Figure 8: Linear regression results of accident progress.

The coefficient of determination values shows how the regression model explains the relationship between input and output parameters. In this study, adjustable R square (Adj. R2) was used for the coefficient of determination. The Adj. R2 values of the regression model about initial melt relocations show that the in-vessel uncertainty parameters explain the potential relationships about initial melt relocation in a limited way. However, as the accident progressed, the influence of in-vessel uncertainty parameters on the regression model decreased. In cases of RPV failure and CNMT overpressurization, Adj. R2 values were significantly low (Adj. R2 < 0.3). Thus the regression model could not explain the relationship between input and output. These results might be caused by the complex behavior of molten corium relocation. In this case, although the uncertainty parameters affected the tendency of molten corium behavior, the effect on the results can be insignificant due to the cumulative random variation.

In case of initial melt relocation, SRRC and PRCC show strong negative correlation (r < -0.7) between SC1141(2) and the timing of initial melt relocation. Figure 9 shows the scatter plot of relocation time versus SC1141(2). The distribution was divided into two categories by 10,400 s of the relocation time. The strong negative correlation resulted from the upper relocation time category. More cases in the upper relocation time category were observed in lower value of SC1141(2) due to the blockage effect. During the candling process, the molten core materials can be refrozen by the candling heat transfer and form the blockage. Because the lower value of SC1141(2) indicates the slower downward flow of molten core materials, the possibility of blockage formation increased with the lower value of SC1141(2). Therefore, the increment of blockage could delay the relocation time.

Figure 9: The scatter plot of relocation time versus SC1141(2).

SRRC of first LP dryout shows that SC1141(2) and SC1131(2) were the most significant parameters affecting the first LP dryout. SRRC and PRCC indicated moderate negative relationship (-0.3 > r > -0.7) between these two uncertainty parameters and first LP dryout time. Figure 10 shows the scatter plot of first LP dryout time versus SC1131(2). Using the color scatter, the values of SC1141(2) were shown in the scatter plot. The linear fitting line shows the negative relationship between SC1131(2) and first LP dryout time. In case of SC1141(2), more blue color scatters were observed below the fitting line. However, this trend was not clear due to the number of samples below the fitting line. The hydrogen production mass and relocated mass before first LP dryout can describe the effect of SC1131(2) and SC1141(2) on first LP dryout time. Figure 11 shows the SRRC and PRCC values of in-vessel hydrogen production and relocated mass before first LP dryout. The high Adj. R2 value (Adj. R2 > 0.7) shows that the rank regression model can properly explain the relationship between uncertainty parameters and hydrogen production and relocated mass. PRCC indicated strong negative correlation between SC1141(2) and in-vessel hydrogen production and strong positive correlation between SC1141(2) and relocated mass before first LP dryout. While the lower hydrogen production reduced the total heat accumulated in the core, the higher relocated mass transferred more heat to water in LP. Therefore, the effect of SC1141(2) on first LP dryout can be reduced by the countering effect. In case of SC1131(2), only moderate positive correlation between SC1131(2) and in-vessel hydrogen production was observed. For that reason, the effect of SC1131(2) on in-vessel hydrogen production directly influenced the first LP dryout time.

Figure 10: The scatter plot of first LP dryout time versus SC1141(2).
Figure 11: Linear regression results of potential factor about first LP dryout time.

4. Conclusions

To investigate the uncertainty of accident progress induced by in-vessel uncertainty parameters, uncertainty quantification was performed with regard to the timing of initial melt relocation, first LP dryout, RPV failure, and CNMT overpressurization. A total of 200 data sets of sensitivity coefficients in COR package were generated by SNAP-DAKOTA plugin using the LHS method. 5th/95th range of initial melt relocation, first LP dryout, RPV failure, and CNMT overpressurization were 10,244 s/ 10,601 s, 12,831 s/13,811 s, 14,911 s/17,253 s, and 55,882 s/69,272 s, respectively. The rank regression analysis was also performed to investigate the effect of in-vessel uncertainty parameters. The regression results show that zircaloy melt breakout temperature (SC1131(2)) and molten clad drainage rate (SC1141(2)) exhibited significant influence on initial melt relocation time and first LP dryout time, respectively. Although only in-vessel uncertainty parameters among input parameters were modified, the effect on RPV failure and CNMT overpressurization time could not be described by the rank regression model. To clearly understand the effect of in-vessel uncertainty on accident progress, various regression techniques and the stochastic studies about the core behavior should be performed for the future work.


Adj. R2:Adjustable R square
AFW:Auxiliary feed water
C.V.:Coefficient of variance
CV:Control volume
DCH:Direct containment heating
HPME:High pressure melt ejection
HRA:Human reliability analysis
LHS:Latin hypercube sampling
LP:Lower plenum
MCCI:Molten corium concrete interaction
MFW:Main feed water
P5:5th percentile
P95:95th percentile
PDF:Probability distribution function
PRCC:Partial rank correlation coefficient
PSA:Probabilistic safety assessment
PSRV:Pressurizer safety relief valve
r:Value of coefficient
RCS:Reactor coolant system
RCP:Reactor coolant pump
RPV:Reactor pressure vessel
SA:Severe accident
SC:Sensitivity coefficient
SDS:Safety depressurization system
SG:Steam generator
SANP:Symbolic nuclear analysis package
SNL:Sandia national laboratories
SRRC:Standardized rank regression coefficient
STSBO:Short term station black out.

Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.


This work was supported by the National Research Foundation of Korea (NRF) funded by the Ministry of Science & ICT (MSIT, Korea) [Grant no. NRF-2017M2A8A4018213] and the Nuclear Safety Research Program through the Korea Foundation of Nuclear Safety (KoFONS) using the financial resource granted by the Nuclear Safety and Security Commission (NSSC) of the Republic of Korea (no. 1805001).


  1. T. Sevón, “A MELCOR model of Fukushima Daiichi Unit 1 accident,” Annals of Nuclear Energy, vol. 85, pp. 1–11, 2015. View at Publisher · View at Google Scholar · View at Scopus
  2. R.-J. Park, S.-B. Kim, and H.-D. Kim, “Evaluation of the RCS depressurization strategy for the high pressure sequences by using SCDAP/RELAP5,” Annals of Nuclear Energy, vol. 35, no. 2, pp. 150–157, 2008. View at Publisher · View at Google Scholar · View at Scopus
  3. S. Lee, K. S. Ha, H.-Y. Kim, and S. J. Kim, “Validation of RCS depressurization strategy and core coolability map for independent scenarios of SBLOCA, SBO, and TLOFW,” Journal of Nuclear Science and Technology, vol. 51, no. 2, pp. 181–195, 2014. View at Publisher · View at Google Scholar · View at Scopus
  4. S. Seo, Y. Lee, S. Lee, H.-Y. Kim, and S. J. Kim, “Effectiveness and adverse effects of reactor coolant system depressurization strategy with various severe accident management guidance entry conditions for OPR1000,” Journal of Nuclear Science and Technology, vol. 52, no. 5, pp. 695–708, 2015. View at Publisher · View at Google Scholar · View at Scopus
  5. Y. Lee, W. Choi, and S. J. Kim, “Efficacy assessment of independent severe accident mitigation measures in OPR1000 using MELCOR code,” Journal of Nuclear Science and Technology, vol. 54, no. 1, pp. 89–100, 2017. View at Publisher · View at Google Scholar · View at Scopus
  6. SNL, “MELCOR Computer Code Manuals, Vol. 1: Primer and Users’ Guide,” SAND2015-6691R, Sandia National Laboratories, 2015. View at Google Scholar
  7. EPRI, MAAP4 applications guidance desktop reference for using MAAP4 software, Electric Power Research Institute, 2010.
  8. INEEL, “SCDAP/RELAP5-3D code manual bolume3: user's guide and input manual,” Idaho National Engineering and Environmental Laboratory, 2003. View at Google Scholar
  9. U. S. NRC, State-of-the-Art Reactor Consequence Analyses Project: Uncertainty Analysis of the Unmitigated Long-Term Station Blackout of the Peach Bottom Atomic Power Station, NUREG/CR-7155, SAND2012-10702P, United States Nuclear Regulatory Commission, 2012.
  10. SNL, “Fukushima Daiichi Unit 1 Uncertainty Analysis-Exploration of Core Melt Progression Uncertain Parameters-Volume II,” Tech. Rep., SAND2015-6612, Sandia National Laboratories, 2015. View at Publisher · View at Google Scholar
  11. J. C. Joe, R. B. Shirley, D. Mandelli, R. L. Boring, and C. L. Smith, “The Development of Dynamic Human Reliability Analysis Simulations for Inclusion in Risk Informed Safety Margin Characterization Frameworks,” Procedia Manufacturing, vol. 3, pp. 1305–1311, 2015. View at Publisher · View at Google Scholar
  12. R. O. Gauntt, N. E. Bixler, and K. C. Wagner, “An uncertainty analysis of the hydrogen source term for a station blackout accident in Sequoyah using MELCOR 1.8.5,” Tech. Rep. SAND2014-2210, 2014. View at Publisher · View at Google Scholar
  13. N. Girault and F. Payot, “Insights into iodine behaviour and speciation in the Phébus primary circuit,” Annals of Nuclear Energy, vol. 61, pp. 143–156, 2013. View at Publisher · View at Google Scholar · View at Scopus
  14. P. Hofmann et al., “ZrO2 Dissolution by Molten Zircaloy and Cladding Oxide Shell Failure. New Experimental Results and Modeling,” in Wissenschaftliche Berechte, INVCIT(98)-P026, 1999. View at Google Scholar
  15. Y. Pontillon, P. P. Malgouyres, G. Ducros et al., “Lessons learnt from VERCORS tests.: Study of the active role played by UO2-ZrO2-FP interactions on irradiated fuel collapse temperature,” Journal of Nuclear Materials, vol. 344, no. 1-3, pp. 265–273, 2005. View at Publisher · View at Google Scholar · View at Scopus
  16. U. S. NRC, “Symbolic Nuclear Analysis Package (SNAP),” NUREG/CR-6974, 2009. View at Google Scholar
  17. K.-I. Ahn, S.-H. Park, H.-D. Kim, and H. S. Park, “The plant-specific uncertainty analysis for an ex-vessel steam explosion-induced pressure load using a TEXAS-SAUNA coupled system,” Nuclear Engineering and Design, vol. 249, pp. 400–412, 2012. View at Publisher · View at Google Scholar · View at Scopus
  18. S. Y. Park, K.-I. Ahn, and Y. M. Song, “Uncertainty Analysis of the Potential Hazard of MCCI during Severe Accidents for the CANDU6 Plant,” Science and Technology of Nuclear Installations, vol. 2015, pp. 1–6, 2015. View at Publisher · View at Google Scholar