Science and Technology of Nuclear Installations
Volume 2017 (2017), Article ID 8431934, 16 pages
https://doi.org/10.1155/2017/8431934

Research Article

Application of ASTEC, MELCOR, and MAAP Computer Codes for Thermal Hydraulic Analysis of a PWR Containment Equipped with the PCFV and PAR Systems

Siniša Šadek,¹ Davor Grgić,¹ and Zdenko Šimić²

¹University of Zagreb, Faculty of Electrical Engineering and Computing, 10000 Zagreb, Croatia
²European Commission Joint Research Centre, Directorate G – Nuclear Safety and Security, Westerduinweg 3, 1755 LE Petten, Netherlands

Correspondence should be addressed to Zdenko Šimić; zdenko.simic@ec.europa.eu

Received 20 December 2016; Accepted 12 March 2017; Published 14 May 2017

Academic Editor: Tim Haste

Copyright © 2017 Siniša Šadek 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.

Abstract.

The integrity of the containment will be challenged during a severe accident due to pressurization caused by the accumulation of steam and other gases and possible ignition of hydrogen and carbon monoxide. Installation of a passive filtered venting system and passive autocatalytic recombiners allows control of the pressure, radioactive releases, and concentration of flammable gases. Thermal hydraulic analysis of the containment equipped with dedicated passive safety systems after a hypothetical station blackout event is performed for a two-loop pressurized water reactor NPP with three integral severe accident codes: ASTEC, MELCOR, and MAAP. MELCOR and MAAP are two major US codes for severe accident analyses, and the ASTEC code is the European code, joint property of Institut de Radioprotection et de Sûreté Nucléaire (IRSN, France) and Gesellschaft für Anlagen und Reaktorsicherheit (GRS, Germany). Codes’ overall characteristics, physics models, and the analysis results are compared herein. Despite considerable differences between the codes’ modelling features, the general trends of the NPP behaviour are found to be similar, although discrepancies related to simulation of the processes in the containment cavity are also observed and discussed in the paper.

1. Introduction

The accident at the Fukushima Dai-ichi nuclear power plant (NPP) made utilities review the plants’ safety systems and procedures. Significant changes in nuclear safety approach, operational procedures, introduction of new systems, and modifications of current ones took place after the accident at the Three Mile Island NPP. Activity in the severe accident area was intensified again after the accident at the Chernobyl NPP due to public opinion, even though there is no direct connection between the LWR and RBMK reactor types. After the accident at the Fukushima NPP more countries paid additional attention to containment and containment safety systems’ performance.

Preserving the containment integrity limits the radioactive material release even in the case of a core meltdown and the reactor pressure vessel failure. Installation of a venting system coupled with appropriate filter devices may prevent damage of the containment wall due to overpressure and reduce release of fission products. Autocatalytic recombiners can be used to lower combustible gases, hydrogen, and carbon monoxide concentration by triggering chemical reaction with oxygen in the containment building.

Controlling containment conditions will be difficult if the plant relies solely on active systems while there is an interruption in electrical power supply. In that case, passive systems provide high level of protection during a long period of time. The passive containment filtered venting (PCFV) system and passive autocatalytic recombiners (PAR) are specially designed systems and equipment which operate in the most severe conditions without the need of an operator action. The reports and papers published in the open literature, some of them listed in [1–4], support the interest of nuclear facilities for such systems to restrict radioactive releases out of the NPPs.

Influence of safety systems on containment behaviour is examined by calculating a station blackout event without any recovery actions with three integral severe accident codes: ASTEC, MELCOR, and MAAP. ASTEC [5] and MELCOR [6] are predominantly mechanistic computer codes and the MAAP code [7], on the other hand, is a fast running parametric code. The code main features are given in a separate chapter.

The calculation is performed for a two-loop pressurized water reactor (PWR) nuclear power plant of a Westinghouse type. The plant model is based on the NPP Krško layout [8]. The input decks for the three codes are made as similar as possible and consistent initial and boundary conditions are placed on all of the analyses. The objectives of calculations are to analytically confirm the ability of passive containment safety systems to mitigate accident consequences and to justify the applicability of, currently, the most popular severe accident codes for such analyses. The analysis covers containment thermal hydraulics, releases of incondensable gases, basemat concrete melting, and all other phenomena relevant for the study of containment integrity. The fission product release and the transport of radio nuclides were calculated by all three codes but results are not presented. It is our opinion that thermal hydraulic differences should be addressed before paying attention to corresponding differences in prediction of radio nuclides behaviour. A comparison of codes’ calculation results is reported herein.

2. Description of Containment Passive Safety Systems

2.1. Passive Containment Filtered Venting System

The passive containment filtered venting system, a system specially designed for light water reactor NPPs to work without external power supply, is used to filter radioactive aerosols, gaseous iodine, and iodine organic compounds and to control the containment pressure.

The PCFV system consists of aerosol filters, the iodine filter, the rupture disc, valves, expansion orifices, instrumentation, and associated piping, Figure 1. The venting gas first passes aerosol filter modules, leaves the containment via piping through a containment penetration, passes the iodine filter, and discharges to the environment through the stack.

Figure 1: PCFV system layout, taken from [9].

Aerosol filters (five units on the right hand side in Figure 1) remove solid particles from the vented gases by mechanical filtering using a metal fiber filter. The iodine filter (green boxes in the middle of Figure 1) removes radioactive iodine and its organic compounds through chemical sorption with silver, active material inside the iodine filter.

System actuation occurs after the failure of the rupture disk at the burst pressure of 0.6 MPa. Vented gases flow then through the piping containing two isolation valves and the pressure relief valve which closes at the pressure 0.41 MPa. The line containing the rupture disc, isolation valves, and the relief valve lies between the outlet of aerosol filters and the inlet to the iodine filter.

2.2. Passive Autocatalytic Recombiners

The aim of passive autocatalytic recombiners is to prevent occurrence of flammable gas mixture in containment compartments by reducing concentration of combustible gases, hydrogen and carbon monoxide.

The recombination process is based on a self-actuated catalytic exothermic reaction between hydrogen/carbon monoxide and oxygen. Products of the reaction are steam and CO₂, respectively, accompanied with release of heat. The reaction is supported by natural circulation of gases. PAR consists of a stainless steel enclosure and plates with the catalyst material. The enclosure is opened at the bottom to allow the gases to enter the PAR unit and at the top for the gases to be discharged back in the containment. It extends above the catalyst elevation to provide a chimney to yield additional lift to enhance the effect of buoyancy. Catalyst material (palladium) is not consumed during the reaction.

3. Overview of the Codes and Their Containment Models

3.1. The ASTEC Code

The ASTEC is a modular computer code. It is being jointly developed since almost 20 years by IRSN, France, and GRS, Germany. The code version V2.0R3p3 was used in the calculation.

ASTEC consists of 13 coupled modules (CESAR, ICARE, CPA, MEDICIS, RUPUICUV, CORIUM, COVI, SYSINT, ELSA, SOPHAEROS, ISODOP, IODE, and DOSE) that model different phenomena, or different parts of a nuclear power plant. Some of the modules important for the presented calculation are briefly described.

The CESAR module [10] computes two-phase thermal hydraulics (TH) in the primary and secondary circuits. Modelling is based on a 1D, two-fluid, five-equation approach completed by a phase slip model. The finite volume method is used for the space discretization by means of the staggered mesh and the time discretization of the basic equations is done using the fully implicit first-order backward difference scheme ensuring the numerical stability. The ICARE module models in-vessel core degradation and vessel rupture [11]. The thermal hydraulics in the core is based on a 1D swollen water level approach completed by a 2D gas modelling. The corium behaviour in the lower plenum is based on a 0D modelling of corium layers (oxide, metallic, and debris layers) with a 2D meshing of the RPV lower head. The molten corium concrete interaction is simulated by the MEDICIS module [12].

The reactor containment thermal hydraulics, aerosol, and fission product behaviour is calculated by the CPA module [13]. The containment is divided into zones: thermal hydraulic volumes consisting of two parts, liquid and gaseous, also named the sump and the atmosphere. Gaseous zone part can contain liquid particles (fog), and liquid zone part can contain dissolved gaseous particles, including incondensable gases (nitrogen, oxygen, hydrogen, carbon monoxide, and carbon dioxide). The user can choose between equilibrium and nonequilibrium zone models. If the zone is in equilibrium thermodynamic (TD) state, temperatures of the liquid and gaseous parts are equal. If the nonequilibrium model is used, temperatures can differ. The latter model is usually used in order to correctly evaluate energy distribution between phases due to differences in heat capacities of water, steam, and incondensable gases.

Connections between the zones are realized by the two types of junctions, atmospheric and drainage junctions. Gaseous components which may carry liquid droplets are transported through atmospheric junctions. Water including dissolved gases is transported through drainage junctions. Water drainage may occur from one zone to another or across the wall in the case of steam condensation.

The thermodynamic behaviour of walls and internal structures in the containment is simulated by heat structures. One-dimensional heat conduction models are used to calculate temperature profiles and energy transfer. Heat transfer coefficients and condensation models are calculated using appropriate correlations built in the code.

3.2. The MELCOR Code

The MELCOR is a modular, engineering-level computer code whose primary purpose is to model the progression of severe accidents in light water reactor nuclear power plants. It is developed by the Sandia National Laboratories, USA, for the US Nuclear Regulatory Commission. The MELCOR version 1.8.6 was used in the calculation.

Initially, the MELCOR code was envisioned as being predominantly parametric but over the years, as phenomenological uncertainties have been reduced, the models implemented into MELCOR have become increasingly best estimate in nature. The use of models that are strictly parametric is limited, in general, to areas of high phenomenological uncertainties where there is no consensus concerning an acceptable mechanistic approach.

Unlike in ASTEC, thermal hydraulic calculation in MELCOR is performed by the same modules for all NPP systems. Control volume hydrodynamics (CVH), flow path (FL), and heat structure (HS) packages are used to model TH behaviour of the primary and secondary circuits, as well as of the containment. The CVH package is concerned with control volumes (CV) and their contents and the FL package with connections that allow transfer of these contents between control volumes. That approach is analogous to the ASTEC CPA module which uses concept of TH zone volumes connected by junctions. Furthermore, the contents of the MELCOR control volume may be divided between a pool containing water which may be subcooled (liquid) or saturated (two-phase), and an atmosphere containing water vapour, liquid water fog, and incondensable gases. Governing differential equations, converted to linearized implicit finite difference form, are solved for each phase separately. Also, like in the ASTEC code, two thermodynamic options are available: equilibrium and nonequilibrium. Equilibrium thermodynamics assumes that the pool and the atmosphere are in thermal and mechanical equilibrium, that is, that they have the same temperatures and pressures. Nonequilibrium thermodynamics, on the other hand, assumes mechanical equilibrium but not thermal equilibrium, so that pressures are equal but temperatures may be different and there may be a substantial driving force for condensation or evaporation.

The junction representation is different in the MELCOR code. In the ASTEC code, atmospheric and water junctions are separated. MELCOR junctions transport both phases at the same time. When preparing a NPP input database, it is easier for the user to use only one junction, instead of two, to connect two volumes that exchange water and gases. In that case, there is no need to define all possible drainage flow paths in complicated containment geometry.

The heat structure package calculates one-dimensional heat conduction within an intact, solid structure and energy transfer across its boundary surfaces into a TH control volume.

3.3. The MAAP Code

The modular accident analysis program (MAAP) is a computer code that can simulate response of light water reactor power plants, both current designs and advanced reactors, during severe accident sequences, including actions taken as part of the accident management. The code is developed for the Electric Power Research Institute (EPRI) by Fauske and Associates, LLC. The MAAP version 4.0.5 was used herein.

MAAP includes models for all of the important phenomena which might occur during accident transients involving degraded cores. It utilizes simplified and fast running models for thermal hydraulics description using a simple, fixed nodalization of the primary and secondary circuits in which the type and number of components and the geometry are predetermined. MAAP is a parametric code that includes combination of phenomenological and user defined parametric models necessary to describe the important trends in the behaviour of a NPP.

The code solves a set of lumped parameter, nonlinear, first-order, coupled ordinary differential equations in time. Most of differential equations express conservation of mass or energy. Momentum balances in MAAP are considered to be quasi-steady, which reduces them to algebraic expressions. Therefore, there are no differential equations in MAAP for the conservation of momentum. Combination of phenomenological, zero-dimensional models with nonexistence of momentum equations makes MAAP a fast running code, calculating 10–100 times faster than ASTEC and MELCOR.

Contrary to the primary and secondary circuit models which are prearranged in a way that the user cannot make its own nodalization but uses the default one, the containment model has a free format. The modelling follows the same pattern as in the ASTEC and MELCOR codes; the containment is divided into thermal hydraulic volumes connected by junctions with addition of heat structures which act as heat sinks.

The code distinguishes two types of heat sinks, distributed and lumped. Both can be used at the same time. Distributed heat sinks are one-dimensional structures with heat flow rate directed through the wall. They can be a wall or a floor between two compartments, or an internal wall within a compartment. Any significant masses of equipment, such as piping, piping supports, valves, pumps, and ironwork structures within the compartment (internal structures), are modelled as lumped heat sinks, represented only with their total masses and surface areas.

4. Computational Model of the Nuclear Power Plant

In order to conduct reliable comparison of the results, the NPP input models for ASTEC, MELCOR, and MAAP are prepared to be as similar as possible. Since processes in the containment depend on the mass and energy releases from the reactor coolant system (RCS), an integral analysis of the NPP behaviour is performed, although the emphasis is put on containment results. The primary, secondary, and containment systems are all modelled to support such analysis. Before the accident simulation, a steady state calculation was performed first to check model accuracy and to qualify it for transient simulations.

4.1. Models of the Primary and Secondary Systems

ASTEC and MELCOR codes allow the user to develop its own input deck almost without any constraints regarding the representation of the nuclear equipment, the reactor coolant system, steam generators, the containment, and all other systems important for the plant operation. Their nodalization scheme of the primary and secondary circuits is shown in Figure 2. The primary circuit is marked with the orange colour, the secondary circuit with the green colour, and the part of the reactor pressure vessel below the upper core plate with the blue colour. Each box in the scheme represents a single control volume. The reactor pressure vessel and steam generators are modelled with a fine mesh of control volumes in order to better predict heat transfer across the fuel rods in the reactor core and across the U-tubes in the steam generators. The core fuel elements are divided into five radial regions and twelve axial nodes. Pressurizer safety and relief valves, as well as proportional and backup heaters are also modelled because they are used for the primary pressure regulation. Many other systems, like pressurizer spray, RCS charging and letdown flows, SG valves, and RCP seal flows, are modelled as well but are not explicitly displayed in Figure 2 because it would make the scheme too complex. All those systems are important for the safe NPP operation and without them it would be difficult to obtain the steady state. However, they will not work in the case of a station blackout, so they are switched off at the beginning of the transient calculation.

Figure 2: ASTEC and MELCOR nodalization schemes of the primary and secondary systems.

The MAAP default nodalization is much simpler, Figure 3. The whole reactor coolant system is represented with 15 control volumes. The reactor pressure vessel is represented with four CVs, the SG U-tubes with two CVs and the hot leg, the cold leg, the intermediate leg, and the pressurizer, each with one CV, making it in total 15 control volumes, taking into account the fact that there are two loops in the RCS. The reactor core is modelled in the same way as in the ASTEC and MELCOR codes, with five radial rings and twelve axial nodes. A single control volume is used to calculate core thermal hydraulic behaviour. The 15-node representation is characteristic for the gas and fission product TH calculation. Water mass balance calculation is performed within the water pools. There are six water pools meaning that more gas volumes are lumped into one water volume. For example, hot legs are lumped together with the RPV upper plenum and the reactor core volumes and the cold legs with the reactor downcomer.

Figure 3: MAAP nodalization of the primary system, taken from [7].

4.2. The Containment Model

All the three codes use the same containment model. The containment nodalization scheme is shown in Figure 4.

Figure 4: Containment nodalization scheme (representation of control volumes and junction connections).

The containment building is represented with 10 control volumes:

(1)DOM (containment dome): cylindrical/spherical air space above the reactor pool, steam generators, and pressurizer compartments.(2)ANL (annulus): air space between the steel liner and the containment building.(3)SG1 (steam generator 1 compartment): air space in the SG1 compartment that contains components SG1 and RCP1.(4)SG2 (steam generator 2 compartment): air space in the SG2 compartment that contains components SG2 and RCP2.(5)PRZ (pressurizer compartment): air space in the compartment that contains pressurizer and primary system safety and relief valves.(6)BET (lower compartment): lower compartment below the containment dome placed between SG1, SG2, and PRZ compartments excluding the reactor pool and the reactor pressure vessel area.(7)RPO (reactor pool): air space above the reactor vessel filled with water during the shutdown, otherwise empty.(8)ARV (around reactor vessel): air space between the reactor vessel and the primary shield walls.(9)CAV (reactor cavity): air space below the reactor vessel including the instrumentation tunnel.(10)SMP (containment sump): the lowest control volume below the SG1 compartment and the lower compartment that contains the recirculation and drainage sumps.

An additional control volume with a large volume and fixed temperature (308 K) is used to represent the environment. This volume is necessary for the code to accurately calculate heat losses from the containment building. Heat transfer coefficient from the outside containment wall to the environment is calculated by the code.

The containment atmosphere is an air-vapour mixture initialized at the atmospheric pressure 101.3 kPa and the temperature 322 K, with 30% relative humidity.

Control volumes are connected with 24 junctions, Figure 4. More than one opening is used between the same volumes if they are located at different elevations to promote internal thermal mixing flow, what can be important for long term containment transients. For example, there are three connections between the lower compartment CV-BET and SG1 and SG2 compartments, respectively, at floor levels. There are also more connections between the containment dome CV and steam generator and pressurizer compartments. Pressurizer and steam generator compartments are open and junction areas between these compartments and the dome are large, between 6 m² and 35 m². Other connections, such as between ARV and SG1 and SG2 compartments which are through cold and hot leg openings in the primary shield walls, are smaller; their values are taken to be 1 m². Connections between the cavity and the ARV/BET volumes are established through small openings on the top of the cavity compartment. The connection between the sump and the cavity is based on cross section area of a 4-inch pipe. The largest connection area is between the reactor pool and the dome, 108.5 m². The reactor sump is just below the SG1 compartment with the connection area being 4.45 m².

The outer containment concrete wall, the steel liner, internal concrete walls and floors, the polar crane, fan coolers, platforms, the refuelling channel, embedment, and other miscellaneous stainless and carbon steel structures are modelled as heat structures. They act as heat sinks and exchange heat with water and gases inside and outside the containment.

Passive autocatalytic recombiners are modelled using correlations developed by German manufacturer NIS Ingenieurgesellschaft mbH (MAAP, MELCOR) and GRS (ASTEC), depending on the available correlations in the codes. Twenty-two PAR units are installed across the containment. There are 14 units in the containment dome, six units in the lower compartment, one unit in the SG1 compartment, and one unit in the SG2 compartment.

The PCFV system is modelled as a simple pipe between the containment upper dome and the environment. The junction connecting the dome and the pipe contains the rupture disc that breaks when the upstream pressure, on the containment side, reaches 0.6 MPa. The other junction, connecting the pipe with the environment, contains the pressure relief valve which is modelled taking into account its “hysteresis” characteristic; the size of the flow area alternates between being fully open and fully closed at the opening and closing set points: 0.49 MPa and 0.41 MPa, respectively. Thus, the PCFV system is not modelled explicitly and the operation of aerosol and iodine filters is not considered. Only the system’s function in controlling the containment pressure and temperature by releasing excessive containment inventory was taken into account.

4.2.1. The Cavity Model

In each code there is a package responsible for calculation of the molten corium concrete interaction (MCCI). The following concrete composition, including reinforcement, is used in the calculation: 35% CaO, 13% SiO₂, 4% H₂O, 21.5% CO₂, 2.5% Al₂O₃, 1% Na₂O, 0.5% MgO, 0.5% Fe₂O₃, and 22% Fe. Corium discharged from the reactor vessel will spread on the cavity floor. A corium spreading area is shown in Figure 5. Molten core materials discharged from the ruptured reactor vessel will react with concrete at the bottom of the cavity. The reaction results in melting of the cavity floor and is accompanied with releases of hydrogen, carbon monoxide, carbon dioxide, and steam. Accumulation of gases leads to containment pressurization and the release of fission products from the melt causes heating of the atmosphere.

Figure 5: Corium spreading area (a cross section and a floor plan).

Intensity of the containment pressurization and heating depends on the reactor cavity layout. In the performed analysis which is based on the NPP Krško containment design, the connection between the cavity and the containment dome is restricted to several small openings on the top of the cavity compartment. In such configuration it is hard to expect a significant dispersion of corium debris in the containment after the failure of the reactor pressure vessel lower head. Nevertheless, steam and gases will freely exit through these holes and cause containment pressure increase.

5. Analysis of the NPP Behaviour during the Accident

5.1. Accident Description

The analyzed transient is a station blackout (SBO) which includes the loss of both off-site and on-site AC power. The only systems available are passive safety systems: accumulators and the pressurizer and steam generator safety valves. Unavailability of electrical power means that reactor coolant pumps, main and auxiliary feedwater pumps, charging, high-pressure, and low pressure safety injection pumps are disabled. Containment safety systems, fan coolers and sprays, are also inoperable. Following the loss of power, RCP seals will overheat due to lack of cooling normally provided by charging pumps, a break will be formed, and coolant will be released from the reactor coolant system to the containment.

Reactor coolant pumps are equipped with staged shaft seals which are provided with cooling system designed to maintain seal integrity such that there is a low seal leakage rate at the nominal RCS pressure. For accident sequences in which there is no cooling of the RCP seals (e.g., SBO), the leakage rate through the seals will increase due to degradation of seal materials when exposed to the coolant at elevated RCS temperatures. The seal leakage rate is 0.03 m³/s, a value that corresponds to a scenario of a total seal rupture in pumps which use a high temperature o-ring RCP seal package [14], a typical arrangement in Westinghouse PWR plants. Leakage of the RCS fluid through the RCP seals, combined with unavailability of electrical power, is a small LOCA (loss of coolant accident) without makeup capability.

5.2. In-Vessel Severe Accident Progression

System thermal hydraulic behaviour and core damage progression are briefly described as the focus of the paper is on the containment analysis. ASTEC results of the calculation of the in-vessel phase of a station blackout accident are reported in [15].

Shutting off the reactor coolant pumps leads to decrease of the coolant mass flow rate. Shortly afterwards, the reactor and the turbine are tripped due to the low cold leg coolant flow. The turbine trip means the closure of the turbine stop valve and isolation of the steam line. Steam generator pressure rises instantly, as a consequence of the steam line isolation, forcing the opening of the SG safety valves and release of excess steam. Since there is no auxiliary feedwater supply, steam generators dry out after about one hour, deteriorating heat transfer from the primary to the secondary side across the SG U-tubes. The insufficient cooling of the RCS, in combination with generation of decay heat and the loss of coolant through the damaged RCP seals, leads to decrease of the core water level, production of steam, and increase of fuel elements’ temperatures. The core heat-up, additionally supported by oxidation of fuel rod cladding and other metallic materials, causes the core to melt. The melting process propagates to formation of an in-core molten pool and ends up with relocation of molten material to the lower head of the reactor pressure vessel. The RPV wall ultimately fails under thermal and mechanical stress and the corium is released in the containment cavity.

Time sequence of main events during the in-vessel phase is shown in Table 1. Calculated MELCOR events precede the other two by about 1000 s. The water mass flow rate out of the RCS through the break during initial 2500 s is 10–15% higher in MELCOR than in ASTEC or MAAP. The difference is not large but affects the ensuing accident progression. The core is, thus, uncovered earlier and the whole process of core degradation begins before that calculated by the other two codes. A larger release of liquid causes earlier transition to a two-phase flow. In the long term, the total coolant release in MELCOR is lower since the vapour flow rate is lower than the liquid flow rate. Table 2 summarizes mass and energy releases from the RCS into the containment as calculated by all three codes. Masses and temperatures of released coolant and molten material are rather well reproduced. Apart from the primary pressure, whose influence is described later, the biggest discrepancies between codes’ predictions are for the hydrogen mass generated in the reactor vessel and the time of the RPV failure. The hydrogen production depends on RCS thermal hydraulic conditions which, as noted before, differ between the code’s calculations. It should be emphasized here that such difference in hydrogen release is not a general trend; only in this specific scenario was a lower amount of hydrogen calculated by MELCOR. Regarding the total hydrogen releases in the containment, this behaviour only has a limited effect since the unoxidized corium in MELCOR will eventually oxidize in the cavity and the hydrogen production will continue after the start of the MCCI process.

Table 1: Time sequence of main events during the in-vessel phase.

Event
ASTEC
MELCOR
MAAP

Two-phase flow at the break
3800 s
2900 s
4000 s

Loss of the SG heat sink
4310 s
3540 s
4500 s

Core uncovery
5150 s
4050 s
4990 s

Onset of fuel rod cladding oxidation
5350 s
4350 s
5250 s

Start of the core melting
5960 s
4750 s
6330 s

Melt relocation to the lower head
6380 s
6100 s
9630 s

RPV failure
9440 s
9080 s
15020 s

Table 2: Main results of the in-vessel severe accident analysis important for the latter containment behaviour.

Parameter
ASTEC
MELCOR
MAAP

Mass of water released from the RCS before the RPV failure
128000 kg
105000 kg
130000 kg

Mean mass flow rate at the break
11 kg/s
9 kg/s
11 kg/s

Temperature of released water/steam
600–1100 K
600–1100 K
600–1000 K

Mass of hydrogen produced in the RPV
268 kg
211 kg
265 kg

RCS pressure at the time of the RPV failure
5.6 MPa
7.8 MPa
6.9 MPa

Mass of material released from the RPV
85700 kg
87500 kg
88000 kg

Temperature of released material
2400 K
2120 K
2330 K

Long term decay heat level in the material accumulated in the cavity
4–14 MW

The failure criteria employed in the MAAP code lead in medium and low pressure accident sequences (LOCAs) to later lower head failure times [16]. Containment conditions are considered at a larger time scale than RCS conditions during a severe accident. Accident progression in the RCS and the reactor core lasts for few hours, and, in the containment, the accident sequence lasts for days. Therefore, differences in the time of the reactor vessel failure are not very significant for the presented containment analysis.

5.3. Containment Behaviour and Comparison of Codes’ Results

5.3.1. Heat-Up and Pressurization

Discharge of reactor coolant in the containment is responsible for the initial containment pressure increase, Figure 6. (Results of ASTEC, MELCOR, and MAAP calculation are put together on the same graphs.) Mass and energy release from the RCS causes the containment to heat up, Figure 7. For the first 4000 s the released coolant is mainly water with a low void fraction of steam, but as the pressure continues to decrease, the steam fraction is increasing. The air heat-up also contributes to the pressure rise but not as effective as the release of steam at the RCP seal break, Figure 8.

Figure 6: Pressure in the containment upper dome.

Figure 7: Temperature in the containment upper dome.

Figure 8: Gas partial pressures in the containment (ASTEC calculation).

Discharge of hot molten corium (>2100 K) from the RPV to the containment cavity, followed by the blow-down of primary circuit gases, speeds up containment heating. Mass of corium released in the containment is about 90000 kg meaning that almost all fuel elements and a large portion of reactor internals have been melted and carried away out of the reactor vessel. Initial decay heat generation inside the melt is 14 MW and during the next seven days it gradually decreases to 4 MW. The reactor nominal power is 2000 MWt. The total coolant inventory in the primary system during normal plant operation is 133000 kg, out of which 105000–130000 kg is released in the containment before the vessel break. The small discrepancy between code simulation results, regarding the released coolant inventory, is mainly due to differences in predictions of thermal hydraulic conditions in the RCS and timing of the reactor pressure vessel failure. Coolant is released from the RCP breaks to steam generator compartments and from there it drains into the containment sump. Pipe connection between the sump and the cavity enables water to enter and to flood the cavity. Half of the released water accumulates in the cavity (~60000 kg) and the other half evaporates. Injection of corium leads to fast water evaporation and containment pressurization. Due to intensive evaporation, the reactor cavity dries out in less than one day, Figure 9. The effect of drying out is also visible on Figure 8 as a sharp drop in the steam partial pressure increase.

Figure 9: Mass of water in the reactor cavity.

Conditions in the RCS before the vessel rupture influence initial increase of pressure and temperature. The fastest early pressurization rate is calculated by MELCOR because the primary system pressure is the highest when the RPV failed. The RCS pressure is rapidly decreasing in the period between 9000 s and 10000 s due to uninterrupted loss of coolant through the break. The ASTEC calculates pressure that is 2 MPa lower than MELCOR’s result since ASTEC simulates the RPV failure 360 s later. The ASTEC’s containment pressure increase rate is thus smaller than in the MELCOR case. When the RCS pressure drops to 5 MPa, accumulators discharge water in the cold legs. The MAAP in-vessel sequence is long enough to account for accumulator actuation. Evaporation of the injected water increases the primary pressure and, after the vessel breach, causes containment pressurization rate to surpass the pressure calculated by MELCOR.

When the containment dome pressure reaches 0.6 MPa (the first pressure peak), the rupture disc in the PCFV line will break causing containment gases to be released in the environment. The pressure drops fast to 0.41 MPa prompting the relief valve in the PCFV line to close. Following the valve closure, the pressure rises once again. After reaching 0.49 MPa, the relief valve opens and again some containment inventory is released. Later, the pressure continues to cycle between 0.41 MPa and 0.49 MPa by the operation of the PCFV pressure relief valve, Figure 6. That kind of valve behaviour is important for preserving containment integrity and minimizing radioactive releases. Failure of the containment wall is assessed by using fragility curves which determine failure probabilities depending on the containment pressure. The containment fragility curve shows 5% failure probability at ~0.6 MPa and for pressures above 0.9 MPa, the probability for containment failure is about 90–95%. If there were no pressure relief systems inside the containment (e.g., PCFV), the pressure would reach critical value in less than a day (Figure 10).

Figure 10: Containment pressure behaviour and indication of failure criterion in the case without safety systems.

After each cycle of the relief valve operation, the new gas distribution is established. Concentrations of steam, nitrogen, and oxygen are being reduced, while those of hydrogen, CO, and CO₂, products of the MCCI, are going up (Figure 8). Apart from being released out of the containment, steam is also produced in the recombiners and by boiling of water bounded in the cavity concrete. Its concentration, therefore, tends to stabilize. Oxygen partial pressure drops to zero already during the first day because it reacts with hydrogen in PARs to produce steam. Nitrogen is neither produced nor consumed and its concentration decreases steadily.

5.3.2. Influence of the Molten Corium Concrete Interaction

Decay heat generated in corium dissolves concrete basemat at the bottom of the cavity.

Concrete is a mixture of calcium carbonate, water, and metal oxides, predominately silica. At temperatures 873–1173 K calcium carbonate is decomposed into calcium oxide and carbon dioxide [17]:

The reaction is endothermic; thus internal energy of the corium is used to dissolve CaCO₃. The released CO₂ and steam produced by evaporation of water from the concrete will react with free metals from the corium (Zr, Cr, and Fe) and iron from the concrete reinforcement (rebar).

Reactions between metals and steam are the following:

and reactions between metals and CO₂ are

Oxidation of zirconium and chromium in steam and CO₂ is an exothermic reaction, while iron oxidation is a slightly endothermic reaction. The amounts of Zr and Cr are limited as they are found only in the reactor core. Thus, the long term releases of H₂ and CO are due to oxidation of concrete rebar since there are no elementary metals in the concrete itself.

Intensity of incondensable gases production can be demonstrated by their partial pressures, shown in Figures 11–13. Differences are substantial, but a general trend can be identified. Considerable amounts of hydrogen and carbon monoxide are released during the first two days owing mostly to oxidation of metals inside the corium. At the same time hydrogen concentration decreases due to operation of recombiners and that is why the partial pressure of H₂ does not go up like that of CO. Carbon dioxide is largely consumed in that process, so its release becomes significant only after the reinforcement remains the only material that can oxidize. After approximately 4-5 days CO₂ is the most important incondensable gas, Figure 8, which causes containment pressurization.

Figure 11: Partial pressure of hydrogen in the containment.

Figure 12: Partial pressure of carbon monoxide in the containment.

Figure 13: Partial pressure of carbon dioxide in the containment.

Initial and final cavity temperature profiles as calculated by the ASTEC code, indicating concrete degradation during the MCCI, are shown in Figure 14. The initial corium thickness is about 10 cm because spreading area is relatively large (38.2 m²). Mass of eroded concrete is shown in Figure 15. As the core-concrete interaction progresses, concrete oxides are dissolved and the molten debris pool and the surface area grow in size. Hence, volumetric heat rate and the melt temperature decrease.

Figure 14: Initial and final cavity temperature profiles as calculated by the ASTEC code.

Figure 15: Mass of eroded concrete in the cavity during the process of the MCCI.

The surface of the concrete is ablated at a rate of 1-2 centimetres per hour. Gases released at the bottom of the pool are assumed to rise through it as bubbles. The rising bubbles also promote production of aerosols containing fission products stripped from the fuel debris. Removal of fission products leads to decrease of decay heat level in the pool. Heat losses from the surface are due to melt eruptions, radiation, and convection to containment gases or to an overlying water layer by means of water boiling. Melt eruptions and water evaporation are major mechanisms for corium cooling in the early phase of the accident. Later, as the corium surface stabilizes, convection from the melt to containment atmosphere gases prevails over the heat transfer caused by melt eruption.

Melt configuration is modelled to be homogenous; thus, there is no melt separation on oxide and metallic materials, although that is not completely fulfilled for the MELCOR calculation. MELCOR’s CORCON module, responsible for the cavity simulation, considers up to 15 possible debris configurations depending on the extent of oxides and metals entrainment into a molten corium mixture. ASTEC and MAAP codes also contain models for layer separation but not as detailed as the MELCOR models.

The cavity erosion progresses in axial and radial directions. The amount of liquefied concrete is calculated based on the data of the latent heat of fusion, liquids and solids temperatures for corium concrete mixtures, and the concrete composition.

The ablation rate of concrete is given by

where is the heat flux at the corium concrete interface, the density of concrete, and the latent heat for concrete melting.

Heat convection between the corium layer and concrete is enhanced by bubble formation at the corium concrete interface. Correlations [18–20] for the calculation of the heat transfer coefficient used by the ASTEC and MELCOR codes include superficial bubble transport velocities. For example, the Bali correlation that is used in the ASTEC calculation gives the following expression for the heat transfer coefficient:

where the Nusselt number is defined as

In the equations above is the heat transfer coefficient, , , are the thermal conductivity, density, and dynamic viscosity of the liquid debris, respectively, is the gas bubble radius, is the superficial gas rising velocity, is the gravity acceleration, and Pr is the Prandtl number.

The heat transfer coefficient in the MAAP code is not determined by experimental correlations, but it is directly entered by the user. It exponentially depends on the corium solid fraction where exponent is also a user defined value. Majority of MAAP models follow the similar approach; mechanistic models are replaced with simple algebraic equations whose parameters are selected by the user. Although the MAAP is relatively simple to use, broad knowledge about severe accident phenomena is necessary to correctly predict the NPP behaviour.

Mass of hydrogen removed by passive autocatalytic recombiners is shown in Figure 16. Hydrogen production during the oxidation in the core and the molten corium concrete interaction is shown in Figure 17. The PAR operation starts when hydrogen mole fraction reaches value of 0.02 and stops after oxygen mole fraction drops to 0.005. Despite being rather short, about 1.5 days, the process of recombination is very efficient, since 70–85% of hydrogen is removed. The time interval when the PARs are active coincides with the early phase of the MCCI process. This is very important for the severe accident management planning because, during that period, hydrogen production rate is the highest. Thereby, operation of passive safety systems provides crucial time for the members of a technical support centre and emergency response organizations in taking preventive and mitigating actions to restrict consequences of a severe accident.

Figure 16: Mass of hydrogen removed by PAR operation.

Figure 17: Total hydrogen production.

6. Discussion of Results

The most significant differences between the results important for the latter accident progression occur during the first two days. Figure 18 shows the temperature of the molten material. The initial cool-down of corium is followed by a temperature increase lasting from 3000 s (MELCOR) to 20000 s (ASTEC). Steam outflow to the neighbouring compartments is limited by the cavity design, and the temperature increases because of the reduced heat transfer rate and convection heat flux. The total temperature increase and duration of that time period depend on the water mass in the cavity. In the case with less water (MELCOR), the two-phase flow regime is established earlier and the higher vapour void fraction results in more efficient cavity ventilation. Unlike ASTEC and MELCOR predictions, there is a temperature rise in the MAAP simulation after water in the cavity dries out. The mass of molten material is low, Figure 15, and so is the heat capacity. Degradation of the heat transfer to the cavity atmosphere causes heat-up of the melt, and since the mass of the melt is low, there is a considerable temperature increase. The bulk of molten material, in the analyses with ASTEC and MELCOR, has a heat capacity large enough to prevent temperature increase following the change in heat transfer conditions on its upper surface. In general, MAAP calculates slower concrete erosion at the beginning of the MCCI when there is still water present in the cavity. The other two codes simulate almost unperturbed concrete degradation regardless of water layer floating on the debris. Despite the cooling of the upper debris surface, the downward and sideward heat flow remains strong enough to melt the concrete. In MAAP, the heat flow is concentrated to the upper surface, water is rapidly evaporating, Figure 9, and the cavity heats up slightly faster in that initial accident period, Figure 19.

Figure 18: Temperature of the molten material in the cavity.

Figure 19: Temperature in the cavity compartment.

Corium heats up water and gases in the cavity. Temperatures of steam and incondensable gases rise sharply, by 20–80 K, after the water evaporates, Figure 19. The lowest value is obtained by ASTEC which calculates that water remains in the cavity for the longest time. Both MELCOR and MAAP calculate higher heat-up. Temperature in the containment dome, Figure 7, follows the pattern of the cavity temperature profile. The ASTEC temperature trend is different from the other two. The reason is that the heat released from the corium surface to fluid is lower in the short term and larger in the long term (>2 days) than calculated by MELCOR and MAAP. Thus, the rate of temperature increase in ASTEC is steeper. This is evident from Figure 20 showing the heat losses to concrete. They begin to decrease after the first day, when, at the same time, cavity temperature starts to increase. (The MAAP curve is not shown because there is no such data in the MAAP plot files.) The upper debris surface is much more active in ASTEC than in the other codes. Melt eruptions cause fission products and melt particles to carry away fraction of the internal energy out from the debris. That energy output is heating up the cavity and the containment. Eruptions eventually subside after the first day, but less energy remains in the melt to be transferred to concrete. A larger heat transfer for the first two days corresponds to more intensive concrete melting (Figure 15). The melting rate later slows down, and masses calculated by MELCOR and MAAP exceed the mass predicted by the ASTEC after seven days of transient.

Figure 20: Heat flow from the corium to concrete.

The US codes (MELCOR and MAAP) show more similarities in results comparing to ASTEC. The difference arises from the fact that ASTEC is a new severe accident code; its physical models and correlations are mostly based on the latest modelling practices and experimental data, whereas MELCOR and MAAP initial MCCI models were prepared in the early 1990s and occasionally updated afterwards. The complexity of the MCCI process itself, the phenomenology and mathematical representation of concrete melting, corium growth, release of gases, and so forth, and differences in the methods of model implementation in the codes additionally contribute to modelling uncertainties and variations in results.

The largest differences in results are related to production of gases, hydrogen, carbon dioxide, and carbon monoxide (Figures 11–13). Even though the partial pressures of gases are not representative for the chemical reaction kinetics during the MCCI, the amount of generated H₂, CO, and CO₂ can be relatively correctly determined based on their values. The codes have extensive models of the MCCI when it comes to heat transfer modes, concrete erosion, molten front progression in radial and axial directions, layer formation, and stratification, but the chemistry data are rather scarce. It can be found that equilibrium chemistry is assumed based on minimization of the total Gibbs function for oxidic, metallic, and gaseous phases [6]. Corium and concrete are regarded as homogenous mixtures of metals and oxides which is not true for the cavity basemat reinforced concrete structure. The reinforcement arrangement affects the progress of concrete erosion and iron oxidation with H₂O and CO₂. The masses of unoxidized zirconium and other corium metals are limited and the only permanent source of combustible gases is the oxidation of the cavity rebar. The impact of codes’ differences in modelling chemical reactions between debris and concrete is particularly emphasized for the hydrogen release, Figure 17, where ASTEC calculates 80% more hydrogen production than MELCOR. Gaseous emission influences containment pressurization and potential formation of flammable mixture of combustible gases, air, and steam. Figure 21 shows the ternary diagram of air, combustible (H₂ and CO), and inert gases with respect to oxidation (steam and CO₂). Flammability limits for hydrogen [21] and CO [22] are also drawn. The oxygen is consumed by PARs early during the transient and, hence, the mixture of gases remains outside the limits of combustion. The accumulation of hydrogen and carbon monoxide in the latter phase of the accident does not pose any threat to containment safety because there is no oxygen to support the burning. The aforementioned differences in codes’ predictions of gases concentrations are noticeable on the bottom right end of the diagram. The closest point to the flammability limit curve is calculated by MAAP just before the first opening of the PCFV relief valve. Some PCFV systems are thus designed with additional line of manually operated valves that enable operator to control containment conditions and keep the atmosphere in an inert state by actively discharging containment inventory.

Figure 21: Flammability limits of the air – H₂ + CO – steam + CO₂ mixture.

The current code models of recombiners do not simulate recombination of carbon monoxide. The results are therefore conservative because the concentration of combustible gases in reality is lower than that represented on the ternary diagram. Nonetheless, the ternary diagram representation is very sensitive on the accident scenario and no conclusion could be drawn for a different course of the accident because small variations in steam and air concentrations can move the mixture properties closer to or further away from the flammability limits.

7. Conclusions

Nuclear power plant behaviour during a hypothetical severe accident event, with emphasis on containment conditions, was analyzed with the ASTEC, MELCOR, and MAAP computer codes. All three codes consist of models that enable integral simulation of NPP phenomena. The physics of an SA progression is interpreted through code routines with significant differences in the level of details of numerical modelling and the representation of the incorporated experimental findings. The MAAP code, for example, uses simpler phenomenological models to reduce the running times, while MELCOR and ASTEC rely on more complicated numerics at the expense of the duration of the analysis.

The analyzed accident showed that installation of the passive filtered venting system and autocatalytic recombiners preserves containment integrity and keeps its atmosphere in inert conditions. In the absence of the containment spray system and the active heat removal, the wall would be breached after one day following the loss of complete power supply if the relief valve of the PCFV system does not open. Recombination of hydrogen maintains the mixture of gases beyond flammability limits. The oxygen starvation in the early part of the transient compensates production of hydrogen and carbon monoxide during the MCCI. Results are rather conservative since there are no code models for CO recombination in PAR units. The releases of incondensable gases in the codes differ considerably which is especially pronounced for the release of hydrogen. Regarding the importance for the NPP safety, the chemistry data about oxidation processes taking place in the cavity with focus on the concrete reinforcement behaviour is an area that could be improved in the codes. The heat loss from the melt in the cavity atmosphere, intensity of concrete erosion, and melt eruptions also differ between the codes and they influence the release of gases, temperature, and pressure increase in the containment.

Despite considerable differences between the codes’ modelling features, calculation results showed that the thermal hydraulic phenomena are in good agreement. Overall trends were similar between all three used codes. Discrepancies between code predictions were due to a number of reasons: unequal core and containment degradation models, usage of parametric instead of mechanistic approach to simulate certain phenomena, differences in resolution of some thermal hydraulic issues, and so forth. The results were, nonetheless, comparable and the differences can be attributed to uncertainties associated with complex SA processes and the scaling of the experimental to the real plant data.

Comparison of calculations of different severe accident codes is the only possible practice for evaluating the codes’ accuracy in the situation when integral severe accident experimental results are missing. It also aids in improving existing input database and selecting the appropriate physical models and correlations for various severe accident phenomena. Continuous code application in analyses of a broad spectrum of accidents can reveal possible deficiencies in the representative models and, thus, help the code developers to improve code routines based on existing experimental findings.

Nomenclature

Latin Symbols

: | Gravity acceleration [m/s²]
: | Heat transfer coefficient [W/m²K]
: | Superficial gas rising velocity [m/s]
: | Latent heat for melting [J/kg]
Nu: | Nusselt number [—]
Pr: | Prandtl number [—]
: | Heat flux [W/m²]
: | Gas bubble radius [m]
: | Velocity [m/s].

Greek Symbols

: | Thermal conductivity [W/m·K]
: | Dynamic viscosity [kg/m·s]
: | Density [kg/m³].

Subscripts

: | Ablation
: | Concrete
: | Liquid debris.

Conflicts of Interest

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

Acknowledgments

The authors would like to express their gratitude to the IRSN for providing the ASTEC code and to the NPP Krško for providing the MAAP code and the plant data.

References

S.-W. Lee, T.-H. Hong, Y.-J. Choi, M.-R. Seo, and H.-T. Kim, “Containment depressurization capabilities of filtered venting system in 1000 MWe PWR with large dry containment,” Science and Technology of Nuclear Installations, vol. 2014, Article ID 841895, 10 pages, 2014.

Y. M. Song, H. S. Jeong, S. Y. Park, D. H. Kim, and J. H. Song, “Overview of containment filtered vent under severe accident conditions at Wolsong NPP unit 1,” Nuclear Engineering and Technology, vol. 45, no. 5, pp. 597–604, 2013.

Y. S. Na, K. S. Ha, R.-J. Park, J.-H. Park, and S.-W. Cho, “Thermal hydraulic issues of containment filtered venting system for a long operating time,” Nuclear Engineering and Technology, vol. 46, no. 6, pp. 797–802, 2014.

E.-A. Reinecke, I. M. Tragsdorf, and K. Gierling, “Studies on innovative hydrogen recombiners as safety devices in the containments of light water reactors,” Nuclear Engineering and Design, vol. 230, no. 1–3, pp. 49–59, 2004.

P. Chatelard, N. Reinke, S. Arndt et al., “ASTEC V2 severe accident integral code main features, current V2.0 modelling status, perspectives,” Nuclear Engineering and Design, vol. 272, pp. 119–135, 2014.

R. O. Gauntt, J. E. Cash, R. K. Cole et al., MELCOR Computer Code Manuals 1-2, Version 1.8.6, Sandia National Laboratories, U.S. Nuclear Regulatory Commission, 2005.

Fauske and Associates Inc, MAAP4—Modular Accident Analysis Program for LWR Power Plants, Vols. 1 to 4, Electric Power Research Institute, 1994.

D. Grgić, V. Benčik, and S. Šadek, “NEK RELAP5/MOD3.3 Post-RTDBE nodalization notebook,” Tech. Rep. NEK ESD-TR-02/13, FER-ZVNE/SA/DA-TR03/13-1, NPP Krško, University of Zagreb FER, 2013.

D. Grgić, V. Benčik, S. Šadek, and I. Bašić, “Independent review of NPP modifications and safety upgrades,” International Journal of Contemporary Energy, vol. 1, no. 1, pp. 41–51, 2015.

L. Piar, N. Trégourès, and A. Moal, “ASTEC V2 code: CESAR physical and numerical modelling,” Tech. Rep. ASTEC-V2/DOC/09-10, DPAM-SEMCA-2010-380, Institut de Radioprotection et de Sûreté Nucléaire, 2011.

P. Chatelard, N. Chikhi, L. Cloarec et al., “ASTEC V2 code: ICARE physical modelling,” Tech. Rep. ASTEC-V2/DOC/09-04, DPAM-SEMCA-2009-148, Revision 0, Institut de Radioprotection et de Sûreté Nucléaire, Gesellschaft für Anlagen- und Reaktorsicherheit, 2009.

F. Duval and M. Cranga, “ASTEC V2 MEDICIS MCCI module: theoretical manual,” Tech. Rep. DPAM/SEMIC 2008-102, Revision 1, Institut de Radioprotection et de Sûreté Nucléaire, 2008.

W. Klein-Heßling and B. Schwinges, “ASTEC V0—CPA module—program reference manual,” Tech. Rep. ASTEC-V0/DOC/01-34, Institut de Radioprotection et de Sûreté Nucléaire, Gesellschaft für Anlagen- und Reaktorsicherheit, 1998.

WOG, “WOG 2000 reactor coolant pump seal leakage model for Westinghouse PWRs,” Tech. Rep. WCAP-15603, Revision 1-A, Westinghouse Electric Company, 2003.

S. Šadek, M. Amižić, and D. Grgić, “Severe accident analysis in a two-loop PWR nuclear power plant with the ASTEC code,” Atw—International Journal for Nuclear Power, vol. 58, no. 12, pp. 694–700, 2013.

K.-I. Ahn and D.-H. Kim, “A state-of-the-art review of the reactor lower head models employed in three representative U.S. severe accident codes,” Progress in Nuclear Energy, vol. 42, no. 3, pp. 361–382, 2003.

Z. P. Bažant and M. F. Kaplan, Concrete at High Temperatures: Material Properties and Mathematical Models, Monograph and Reference Volume, Longman (Addison-Wesley), London, UK, 1996.

J. M. Bonnet, “Thermal hydraulic phenomena in corium pools for ex-vessel situations: the Bali experiment,” in Proceedings of the 8th International Conference on Nuclear Engineering, ICONE-8, ICONE-8177, Baltimore, Md, USA, 2000.

F. A. Kulacki and M. E. Nagle, “Natural convection in a horizontal fluid layer with volumetric energy sources,” Journal of Heat Transfer, vol. 97, no. 2, pp. 204–211, 1975.

M. T. Farmer, J. J. Sienicki, and B. W. Spencer, “CORQUENCH: a model for gas sparging-enhanced, melt-water, film-boiling heat transfer,” Transactions of the American Nuclear Society, vol. 62, pp. 646-647, 1990, Proceedings of the American Nuclear Society (ANS) Winter Meeting, Washington, DC, USA, 1990.

Z. M. Shapiro and T. R. Moffette, “Hydrogen flammability data and application to PWR loss-of-coolant accident,” Tech. Rep. WAPD-SC-545, U.S. Atomic Energy Commission, 1957.

N. Cohen, “Flammability and explosion limits of H₂ and H₂/CO: a literature review,” Aerospace Report TR-92(2534)-1, Space and Missile Systems Center, The Aerospace Corporation, El Segundo, Calif, USA, 1992.