Science and Technology of Nuclear Installations
Volume 2020 (2020), Article ID 6798738, 12 pages
https://doi.org/10.1155/2020/6798738

Research Article

Numerical Simulation and Validation for Early Core Degradation Phase under Severe Accidents

Dekui Zhan,¹ Xinhai Zhao,¹ Shaoxiong Xia,¹ Peng Chen,¹ and Huandong Chen²

¹China Nuclear Power Technology Research Institute, Shenzhen 518000, China
²Sino-French Institute of Nuclear Engineering and Technology, Sun Yat-sen University, Zhuhai 519000, China

Correspondence should be addressed to Xinhai Zhao; zhaoxinhai2@cgnpc.com.cn

Received 20 December 2019; Revised 26 March 2020; Accepted 24 June 2020; Published 3 August 2020

Academic Editor: Kostadin Ivanov

Copyright © 2020 Dekui Zhan 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.

Early core degradation determines the amount of hydrogen generated by cladding oxidation as well as the temperature, the mass, and the composition of corium that further relocates into the lower head of reactor pressure vessel (RPV), which is essential for the effectiveness analysis of in-vessel retention (IVR) and hydrogen recombiners. In this paper, the mechanisms of controlling phenomena in the early phase of core degradation are analysed at first. Then, numerical models adopted to calculate (1) core heating up, (2) cladding oxidation, (3) dissolution between molten zirconium and fuel pellets, and (4) formation of a molten pool in the core active section are presented. Compared with integral codes for severe accident analysis (such as MAAP and MELCOR), the models in this paper are established at the fuel pin level and the calculation is performed in 3D, which can capture the detail local phenomena during the core degradation and eliminate the average effect due to equivalent rings used in integral codes. In addition, most of the control equations in this paper are calculated by implicit schemes, which can improve the accuracy and stability of the calculation. In the simulation, the calculation oxidation is calculated by using the oxygen diffusion model, while the dissolution is calculated with Kim, Hayward, Hofmann, and IBRAE models to perform uncertainty analysis. For the validation, the cladding oxidation model is verified by Olander theoretical cases in the conditions of both steam-rich and steam-starved. The dissolution models are validated by the RIAR experiment. The code is overall verified by Phebus FPT0 on the integral phase of core early degradation. According to the simulation results, it can be inferred that the dissolution reaction between the molten zirconium and fuel pellets is the main reason for the melting of UO₂ at low temperature. In the case of starved steam, part of the fuel pellets can melt down even at 2248 K and relocate to the bottom of the core, which is much lower than the melting point of UO₂ (3113 K).

1. Introduction

In the Chernobyl and Fukushima accidents, the melted reactor cores caused a large amount of radioactive material to be released into the environment. It is essential to study the severe accident phenomena for the design of dedicated mitigation measures to cease the accident process and reduce its consequences.

Accidents that can cause the core to melt down mainly include (1) the station blackout accident (SBO), (2) the large break loss of coolant accident (LBLOCA) superimposed with the failure of safety injection pumps, and (3) the anticipated transient without scram (ATWS) [1]. Although the initial events are different, the phenomena in the core early degradation are approximately identical. The controlling phenomena include core uncovery, core heating up, failure of control rods, dissolution of UO₂ and ZrO₂ by molten zirconium, failure of the cladding oxide layer, corium relocation and solidification, and formation of a molten pool in core active section.

For the cladding oxidation, under the condition that the solid zirconium alloy is thick enough and the steam is sufficient, the oxidation of zirconium alloy is dominated by the oxygen diffusion, in which the chemical reaction rate follows the parabolic law [2]. Correlations of Urbanic and Heidrick [3], Prater and Courtright [4], and Leistikow-Schanz [5] are widely used in the calculation of cladding oxidation.

De Luze et al. [6] summarized the failure temperatures of Ag-In-Cd control rods in different experiments, including CORA-9, QUENCH-13, and Phebus FPT1. The failure temperatures of control rods in most experiments are between 1400 and 1700 K. It can be concluded that the eutectic reaction between the solid stainless steel and the zirconium alloy is the reason for the partial melting of the cladding and guide tube at about 1400 K, while the dissolution of zirconium by the molten stainless steel and Ag-In-Cd absorber causes the guide tube to be broken in a wide range at around 1700 K.

Crucible experiments found that the dissolution process of UO₂ can be divided into two phases, including the saturation phase and the precipitation phase. Hofmann [7], Hayward and George [8], and Kim and Olander [9] established separately dissolution models based on their experiments, respectively. However, significantly different results were discovered during the precipitation period. The dissolution of UO₂ reaches saturation in Kim and Hayward experiments while Hofmann found the dissolved UO₂ mass continues to increase and follows the parabolic law. The ZrO₂ dissolution process is also separated into two stages, including the erosion stage and the corrosion stage [10], in which the correlation of Hofmann is widely used.

Based on the experimental data and mechanism research studies, many integral codes for severe accident analysis were developed, including MELCOR, SCDAP/RELAP5, and ASTEC. The temperature criteria and the oxide thickness criteria are widely used in most codes to predict the failure moment of the cladding oxide layer. For the relocation process of melt after the failure of the cladding oxide layer, candling model is established in MELCOR by analogy with the flow of molten droplets of the candle along the wall [11]. LIQSOL model is developed in SCDAP/RELAP5 based on the experimental results of CORA, PBF-SFD, and FLHT [12]. Decanting model and porous media model are adopted in ASTEC, which can simulate the radial and axial relocation of corium simultaneously [13].

Due to the consideration of computing efficiency, equivalent rings are widely used in integral codes. However, radial power factors of assemblies in the same ring are different, especially between assemblies with black rods and without black rods. When equivalent rings are adopted in the modelling, the power factors of high-power assemblies are actually averaged, causing the delay of the calculated degradation process. It means that integral codes may underestimate the time of some important events in the reactor, such as the start time of cladding oxidation, the failure of the zirconium cladding (namely, the fission product begins to be released in large quantities), and the downward relocation of the molten materials, which are important for accident analysis. In this paper, the models are established at the fuel pin level and the calculation is performed in 3D which can capture more detailed local phenomena during the core degradation phase and eliminate the average effect due to equivalent rings used in integral codes to achieve more accurate simulation for the early phase of core degradation.

2. Mechanism of Early Core Degradation

In this section, the mechanisms of controlling phenomena in the early core degradation phase are presented based on the literature review.

2.1. Oxidation of Zirconium Cladding

The reaction of zirconium cladding with steam is the most important phenomenon in the core degradation, which generates not only a large amount of heat but also a great quantity of hydrogen. With neglecting considering the hydrogen absorption in the cladding, the zirconium oxidation can be divided into the following stages [7]: (1) diffusion of steam to the outer surface of cladding; (2) chemical absorption of steam; (3) dissociation of steam molecules into oxygen and hydrogen; (4) diffusion of oxygen from the outer surface of cladding to the inside, forming ZrO₂, α-Zr(O), and β-Zr three layers in sequence; (5) recombination of hydrogen; and (6) diffusion of hydrogen molecules into the bulk of the flow channel, causing the hydrogen content in the flow channel to increase.

2.2. Dissolution of UO₂ and ZrO₂

Based on the analysis of the experimental results of Hayward, Kim, and Hofmann, Olander [14] proposed that during the saturation period, the UO₂ dissolution is dominated by natural convection mass transfer. Liquid U with a density of 1.7  10⁴ kg/m³ is produced by a chemical reaction at the interface between the molten zirconium alloy and the solid UO₂. However, the density of the liquid Zr is only 6  10³ kg/m³. The difference in density between liquids U and Zr causes significant convection in the melt. During this period, the U content in the melt increases continuously and gradually reaches the saturation. Dissolved U and O are quickly transferred to the central region of the melt and mixed well with the molten Zr. Therefore, it is found in the experiments that the ratios of U, Zr, and O are uniformly distributed in the molten zone. In the precipitation stage, Olander does not suppose that UO₂ will continue to be dissolved. At the end of the saturation period, the dissolution reaction of UO₂ is ceased, and the mass ratio of U and Zr in the melting zone is kept unchanged. During this period, O in the UO₂ crucible continues to diffuse into the molten zone, causing the previously saturated (U, Zr, O) melt to be supersaturated to precipitate the (U, Zr)O2-x ceramic phase.

The dissolution of ZrO₂ is similar to that of UO₂. Based on the ZrO₂ crucible experiments, Hayward and George [10] speculated that in the erosion period, the dissolution of ZrO₂ is caused by a significant convective effect. The newly dissolved ZrO₂ is quickly brought into the melt central zone. Therefore, no transition layer was created in the experiment in this period. Thereafter, the oxygen concentration gradient at the interface gradually decreases causing the convective effect to be weakened. Once the mass flux of oxygen entering the melt by diffusion is greater than that by dissolution, the erosion period is ceased. In the following period, oxygen in the ZrO₂ crucible continues to diffuse into the melt driven by the difference of oxygen potential. Near the interface region, ZrO₂ in the crucible is reduced to ZrO2-x and the metallic zirconium in the melt is oxidized to ZrO2−x. Therefore, it was shown in the experiment that the thickness of the ZrO2−x transition layer increased not only towards the crucible region but also towards the melt central region.

2.3. Cladding Oxide Layer Failure

The cladding of fuel rods is oxidized by steam at high temperature. Then, a ZrO₂ oxide layer is formed on the outside of the cladding. The oxide layer with a high melting point (2963 K) will hold up the molten (U, Zr, O) in the cladding, preventing the downward relocation of the melt. However, the ZrO₂ oxide layer is dissolved by the metallic Zr in the molten (U, Zr, O), causing the oxide layer to become thinner. In addition, the creep effect at high temperatures reduces the ability of the oxide layer to hold up the melt and its yield limit decreases with increasing temperatures. Therefore, when the thickness of the oxide layer is reduced to a certain threshold, the oxide layer can no longer hold up the (U, Zr, O) melt. Then, the (U, Zr, O) melt is released to relocate downwards.

2.4. Corium Relocation

According to CORA W1 and W2 experimental results, Veshchunov and Palagin [15] infer that in the first period, the molten cladding materials flow down rapidly in the form of droplets or rivulets at about 0.5 m/s. Owing to a large temperature gradient below the front, the molten materials refreeze in the interrod space forming local blockages (debris) and accumulating the melt in the form of molten pool. The melt mass in the molten pool gradually increases. Once the solidified crust at the bottom cannot support the weight of the melt, the molten pool will relocate downwards. During the relocation process, the molten zirconium alloy in the molten pool will dissolve adjacent fuel pellets and be oxidized by steam. The dissolution reaction and the oxidation reaction cause U in the melt to reach a supersaturated state to produce the (U, Zr)O2−x ceramic phase precipitate. The presence of the (U, Zr)O2−x ceramic phase precipitate increases the viscosity of the melt, thereby preventing the molten pool from rapidly relocating towards the bottom of the bundle. At the same time, the oxidation reaction generates a large amount of heat raising the temperature of the melt. Therefore, the solubility of U in the melt is increased, resulting in the reduction of the (U, Zr)O2−x ceramic phase and the viscosity of the melt. Consequently, the molten pool moves downwards with a slow characteristic velocity 1∼2  10⁻³ m/s. As the molten pool relocates downwards, fresh melt continues to accumulate in the molten pool, causing the molten pool to expand in the radial and axial directions to form a massive slug flow.

3. Numerical Simulation Method

In this section, the physical models adopted in the simulation are presented for the phenomena mentioned above.

3.1. Bundle Heat Transfer

The core heating up is accompanied by heat conduction, convection heat transfer, and radiation heat transfer. For the convective heat transfer, the Dittus–Boelter forced convection correlation is adopted when there are drive pumps and the coolant flow rate is large, while the natural convection correlations are used in the case that the flow of coolant is driven by the buoyancy. The models introduced in the following sections are used to calculate the transient heat conduction in the fuel rods and the radiation heat transfer between rods.

3.2. Transient Heat Conduction

The two-dimensional axisymmetric unsteady heat conduction equation is used to solve the transient temperature distribution in the fuel rods and control rods. The governing equation is given by

in which is the volumetric heat generation rate. For the absorbers in control rods, , while for the fuel pellets,

in which is the total power of core decay heat. is the axial power factor. is the radial power factor of fuel assemblies. is the number of fuel rods in an assembly. is the volume of a fuel rod in one axial node. The following equation is obtained by integrating equation (1):

The coefficients in equation (3) are shown in (4)∼(6), where the geometric parameters are defined in Figure 1.

Figure 1: Schematic of mesh generation for fuel pellets or absorbers.

3.3. Radiation Heat Transfer

The radiative heat transfer is calculated at the fuel pin level. Due to the shielding effect, central rod P only radiates heat directly to the surrounding 16 rods as shown in Figure 2. According to the geometric position, the surrounding rods can be divided into three types: 1, 2, and 3. The radiative view factors are , , and , respectively. Then, the effective radiation for central rod P is calculated as follows:

where , , and are the effective radiation at the outer surface of rods 1, 2, and 3, respectively. The radiative heat flux at the outer surface of rod P is given by

Figure 2: Schematic of the radiation model between rods.

3.4. Cladding Oxidation

The oxygen diffusion model proposed by Olander [16] is adopted to calculate the oxidation of cladding, in which the cladding is divided into the oxide layer (ZrO₂) and the metal layer (including a-Zr(O) and β-Zr). The oxidation rate depends on the amount of steam available, the thickness of the metal layer as well as the diffusion rate of oxygen in the cladding (namely temperature). According to the amount of available steam, the steam state in the flow channel can be divided into two cases: steam-rich and steam-starved. The amount of steam available is calculated by the molar steam balance:

in which is the molar density of gas. is the steam mole fraction in the bulk gas in node i. The total gas molar flow rate is denoted by and is the external surface area of the cladding. The steam absorption rate is given by

where is the steam molar fraction at the cladding outer surface. The oxygen diffusion rate in the cladding is given by

in which the molar density of Zr in ZrO₂ is denoted as and is the diffusion coefficient of oxygen in the oxide layer. , , and are the O/Zr ratio at the outer surface of the oxide layer, at the interface in the oxide layer side, and at the interface in the metal layer side, respectively, which are determined by Zr-O phase diagram. The oxidation rate is calculated as the minimum of and . In the case of oxidation when both oxide and metal phases are present, the thickness of the oxide layer is calculated as follows:

where is the diffusion coefficient of oxygen in α-Zr(O). The average oxygen content in the metal layer is determined by the following equation:

For the case of oxidation in all-metal cladding or in all-oxide cladding, other controlling equations are used to calculate the oxide thickness and the average oxygen content in the metal layer.

3.5. Dissolution of UO₂ and ZrO₂

Hofmann [7], Hayward and George [8], and Kim and Olander [9] established dissolution models based on their experiments, respectively. However, the results of each experiment in the precipitate period are significantly different. In the calculation, all the three models as well as the IBRAE model are used to simulate the dissolution of UO₂ and ZrO₂ and will be verified by the RIAR test in the next section.

3.5.1. Kim UO₂ Dissolution Model

When the Ra number in the (U, Zr, O) melt is greater than 1000, the natural convection model is adopted [9]:

Otherwise, the diffusion model is used [9]:

in which is the mass fraction of U in the bulk melt. , , and are the UO₂, (U, Zr, O) contact surface, the volume, and viscosity of (U, Zr, O), respectively. is the diameter corresponding to and is the (U, Zr, O) mixture thickness.

3.5.2. Hayward UO₂ Dissolution Model

Based on crucible experiments, Hayward found that the mass fraction of U in the melt can be calculated as follows [8, 17]:

in which is a constant related to the contact area and the initial oxygen content.

3.5.3. Hofmann UO₂ and ZrO₂ Dissolution Model

The crucible experiments of Hofmann show that UO₂ is rapidly dissolved in the initial period, and the mass fraction of UO₂ in the melt rapidly increases to 35.8%. Subsequently, UO₂ is continuously dissolved, and the UO₂ mass fraction in the melt follows the parabolic law [7]:

The ZrO₂ dissolution model is calculated in a similar form of equation (7):

3.5.4. Simultaneous Dissolution of UO₂ and ZrO₂

For the calculation of the simultaneous dissolution of UO₂ and ZrO₂, the first method is to assume it as two separate processes. Within one time step, the dissolution of UO₂ is calculated firstly according to the models mentioned above, and then the dissolution of ZrO₂ with Hofmann correlation is calculated until the mixture composition reaches the liquidus in the (U, Zr, O) ternary phase diagram. The second method is to adopt the IBRAE model proposed by Veshchunov et al. [18, 19], which can be used to calculate the simultaneous dissolution of UO₂ and ZrO₂ as well as the oxidation of (U, Zr, O) mixture. In the IBRAE model, both intact fuel rods and the molten pool can be considered as UO₂, (U, Zr, O) mixture, and ZrO₂ three-layer structure. The interface movement rate can be derived from the flux match. For the interface of UO₂ and (U, Zr, O),

For the interface of (U, Zr, O) and ZrO₂,

where , , and are positions of boundaries of the solid UO₂ and ZrO₂ layers. , , and are the molar densities of different components (, U, and Zr) in UO₂, (U, Zr, O) mixture, and ZrO₂, respectively. The molar density of the melt is denoted by . and are temperature-dependent parameters which determine the position of liquidus line in the ternary phase diagram.

The oxygen concentration in the melt in (21) and (22) can be calculated by the mass conservation equations, while the oxygen concentration gradient at interface is calculated by the diffusion equation of oxygen in solid UO₂ with moving boundary conditions. For the oxygen concentration gradient at the interface ,

3.6. Criteria of Cladding Oxidation Layer Failure

In the steam-rich case, the oxidation reaction in the cladding dominates. Oxidation of cladding continues until it is completely oxidized. In the calculation, it is assumed that the fuel pellets and cladding will melt together once the temperature reaches the eutectic point (2673 K).

In the steam-starved condition, the dissolution reaction in the cladding dominates. In the calculation, it is assumed that the cladding oxide layer will fail once the (U, Zr, O) mixture temperature exceeds the α-Zr(O) melting point (2248 K) and the oxide layer thickness is less than 300 μm.

3.7. Corium Relocation

Based on the explanation of Veshchunov and Palagin on CORA W1 and W2 experimental results [15], the corium relocation downward in the bundle is considered as the slug flow with fixed velocity 1  10−3 m/s. During the relocation, convective heat exchange with the steam, radiant heat exchange with the surrounding nonadjacent rods, and the heat generated by corium oxidation are considered. The energy equation during relocation is given by

4. Separate Effect Validation

In this section, the physical models for both the oxidation reaction and dissolution reaction which are the two of the most important phenomena in the early phase of core degradation will be verified.

4.1. Oxidation of Zirconium Cladding

The oxidation model is validated with theoretical cases of a single fuel rod oxidation proposed by Olander [16]. Both steam-rich and steam-starved conditions are simulated. In order to simplify the model and constrain other variables, the effect of heat generation due to oxidation is not considered. Instead, the fuel rod temperature of each node is determined by the following equations:

in which is the initial temperature. The temperature at which the microstructure of the oxide layer begins to change from tetragonal to cubic is denoted by with a value of 1773 K. and are the temperature ramp rates below and above the transition temperature.

In the steam-starved case, the inlet steam flow rate is 2  10⁻² m/s. The temperature ramp rates and are 1.0 K/s and 2.5 K/s, respectively. The steam mole fraction in the flow channel at different positions is illustrated in Figure 3. The cladding oxygen absorption rate increases as the temperature of the fuel rod rises. Therefore, the steam mole fraction in the flow channel gradually decreases. When the steam-starved condition occurs, the oxygen absorption from the flow channel is limited by the amount of steam available. Once the oxygen absorption rate of the cladding from the flow channel is less than the diffusion rate of oxygen from the oxide layer to the metal layer, the oxide layer will become thinner, as shown in Figure 4.

Figure 3: Comparison of steam mole fractions in the flow channel.

Figure 4: Comparison of fractional oxide layer thickness in cladding.

In the steam-rich case, the inlet steam flow rate is 1 m/s. The temperature ramp rates and are 0.7 K/s and 1.0 K/s, respectively. The steam mole fraction in the flow channel at different positions is presented in Figure 5. It can be seen that the steam mole fraction in the flow channel decreases steeply at around 600 s. It is explained by the fact that the temperature of the rod increases to around 1773 K at which the microstructure of the oxide layer changes from tetragonal to cubic, causing the diffusion of oxygen in the oxide layer increasing rapidly and large amount of steam consumed. After 750 s, the steam mole fraction in the 210 cm and 390 cm flow channels increase since the cladding below 210 cm has been sufficiently oxidized. Then, the consumption of steam is reduced so that the steam mole fraction arriving above 210 cm increases. After 880 s, the steam mole fraction in the flow channel drops again due to the increase in the temperature of the fuel rod causing the increase of the oxygen diffusion rate in the oxide layer and the metal layer. The variation of fractional oxide layer thickness in cladding is shown in Figure 6.

Figure 5: Comparison of steam mole fractions in the flow channel.

Figure 6: Comparison of fractional oxide layer thickness in cladding.

It can be seen from Figures 3–6 that the calculation results in this paper are consistent with the results of Olander, whether it is steam-rich or steam-starved. Then the rationality of the oxidation model in this paper can be verified.

4.2. Dissolution of UO₂ and ZrO₂

RIAR test [20] is used to verify physical models for the simultaneous dissolution of UO₂ and ZrO₂ and perform the uncertainty analysis.

Figures 7 and 8 show the U mass fraction in the (U, Zr, O) mixture calculated by different models compared with experimental values and Muller’s simulation results [20]. Since there is no limit on the solubility, the solid UO₂ is continuously dissolved by the molten mixture, which causes the U mass fraction in the melt to increase approximately linearly with time in the Hofmann model. But in Kim and Hayward models, the dissolution of UO₂ is limited by its solubility in the (U, Zr, O) mixture. It can be inferred from the ternary phase diagram that there is a competition between UO₂ dissolution and ZrO₂ dissolution. It means that the solubility of UO₂ is reduced due to the dissolution of ZrO₂. In the simulation with Kim and Hayward models at 2200°C, the dissolution rate of ZrO₂ calculated by Hofmann parabolic correlation is overpredicted. Excessive and rapid dissolution of ZrO₂ makes the U in the mixture quickly reach saturation, thereby reducing the volume of solid UO₂ dissolved. In the IBRAE model, the simultaneous dissolution of ZrO₂ and UO₂ is characterized by a complete mathematical model instead of the experimental correlations. The results calculated by the IBRAE model in this paper are in agreement with Muller’s results and experimental values.

Figure 7: U mass fraction in the melt changes with time at 2100°C.

Figure 8: U mass fraction in the melt changes with time at 2200°C.

Muller’s results are considered as a reference to perform uncertainty analysis. Dissolved volumes of UO₂ and ZrO₂ at 2200°C are presented in Figures 9 and 10, respectively. Within 100 s after the occurrence of UO₂ dissolution, dissolved volumes of UO₂ calculated by each model are approximately consistent, during which the maximum uncertainty range does not exceed 20%. Then, due to the premature saturation of U in the melt, the dissolution of UO₂ in Hayward and Kim models is ceased with 35% uncertainty at 200 s. However, because of no limit on solubility and high dissolution rate, the volume of UO₂ dissolved using the Hofmann model is maximum and with the uncertainty of 60% at 200 s. As shown in Figure 10, the uncertainty of ZrO₂ dissolved volume is close to 300%. The calculated volumes of ZrO₂ dissolved by Kim, Hayward, and Hofmann models are considerably different from Muller’s results, which can be explained by the overpredicted dissolution rate of ZrO₂ caused by the Hofmann correlation. Therefore, it can be inferred that since the dissolution rate of ZrO₂ is high but the dissolution rate of UO₂ is low for Hayward and Kim models, it may cause the failure time of the cladding oxide layer to be advanced while the mass of UO₂ dissolved at this stage is underestimated. On the contrary, it will overestimate the volume of UO₂ dissolved by using the Hofmann model, which is suitable for conservative calculations.

Figure 9: Comparison of dissolved UO₂ volume for different models at 2200°C.

Figure 10: Comparison of dissolved ZrO₂ volume for different models at 2200°C.

In severe accidents, the dissolution of UO₂ by molten zirconium is usually accompanied by cladding oxidation. In the steam-rich case, the zirconium cladding may be completely oxidized to ZrO₂. The increase of oxygen content in the (U, Zr, O) mixture will greatly reduce the solubility of U. Therefore, only when the temperatures of UO₂ and ZrO₂ exceed the melting point or reach the eutectic temperature, UO₂ and ZrO₂ can relocate to the bottom of the core as corium. But in the steam-starved condition, the zirconium water reaction occurs but is not intense, forming a thin ZrO₂ oxide layer on the outside of cladding. In this case, the molten zirconium dissolves UO₂ inward and dissolves ZrO₂ outward. As a result, some of the fuel pellets melt at the temperature far below its melting point. The dissolved UO₂ is relocated downward in the form of (U, Zr, O) mixture, while the undissolved UO₂ in the same pellet collapses in the form of debris. Therefore, the temperature corresponding to the meltdown of core is greatly reduced.

5. Integral Effect Validation

In-pile experiments Phebus FPT provide a reliable benchmark for numerical simulation of the core degradation process. In this paper, the integral effect of the core early degradation phase is verified by Phebus FPT0 experiment [21, 22]. The input bundle power and steam flow rate are shown in Figure 11. Considering the pressure tube is immersed in circulating cooling water with temperature 438 K, the temperature of the outer surface of the pressure tube is approximately supposed as 438 K.

Figure 11: Bundle power and steam flow rate in Phebus FPT0.

Figure 12 shows the comparison of temperature at 0.6 m in the calibration phase (0∼10000 s) between the simulation results of this paper, ATHLET CD, and ICARE2 as well as experimental values. In this period, due to the increase in power and the decrease in steam flow, the temperature variation of fuel rods can be divided into three rising stages and three platform stages. It can be seen from Figure 13 that the temperature of the fuel rod at 0.7 m increases sharply at around 12000 s. This period is so-called oxidation escalation phase (10000∼13000 s), during which the power of the bundle increases causing the temperature of the bundle rises again. After the cladding temperature exceeds 1500 K, cladding oxidation becomes significant, generating a large amount of heat and hydrogen. Figures 14 and 15 show the hydrogen generation rate and the cumulative hydrogen generation mass, respectively. The hydrogen generated mass calculated in this paper during this period is 90.83 g, which is close to the experimental value 85 g. Figure 16 presents the variation of temperature at 0.6 m in the corium relocation phase (13000 s∼18000 s). The distribution of temperature along the height throughout the simulation is shown in Figure 17.

Figure 12: Variation of outer ring fuel rod temperature at 0.6 m.

Figure 13: Variation of inner ring fuel rod temperature at 0.7 m.

Figure 14: Hydrogen generation rate.

Figure 15: Cumulative hydrogen generation mass.

Figure 16: Variation of outer ring fuel rod temperature at 0.6 m.

Figure 17: Temperature distribution along the height.

Figures 13 and 16 show that the maximum temperatures at 0.6 m and 0.7 m exceed 2500 K in the oxidation escalation phase, while the maximum temperature at 0.8 m is only 2250 K, as shown in Figure 18. It can be explained that in this phase, massive steam consumed in the middle and lower regions of the bundle causes a short period of steam starved in the upper region, as shown in Figure 19. In the steam-starved case, it is assumed in the calculation that once the temperature of fuel rods exceeds the α-Zr(O) melting point (2248 K) and oxide layer thickness is less than 300 μm, the cladding oxide layer will fail and the previous hold-up (U, Zr, O) melt will relocate downwards. As presented in Figure 17, some dissolved fuel pellets begin to relocate downward after 11970 s, which is quite in agreement with the experimental value 12050 s. Based on the simulation results, it can be inferred that the reason for rapid accident process of LBLOCA may be not resulted from the fuel rapidly rising to its melting point due to lack of cooling, but a large scale of cladding melts, forming a molten metal pool in the flow channel, and dissolving the adjacent fuel pellets, causing the fuel to melt down at much lower temperatures compared with its melting point 3113 K.

Figure 18: Variation of inner ring fuel rod temperature at 0.8 m.

Figure 19: Variation of steam mole fraction at different heights.

6. Conclusions

In this paper, based on the detailed analysis of the mechanism of controlling phenomena in the early phase of core degradation, fine models at the fuel pin level are established in 3D and programed by Fortran to calculate (1) core heating up, (2) zirconium cladding oxidation, (3) dissolution reaction between molten zirconium and fuel pellets, and (4) formation of molten pool in the core active section.

Analysis of dissolution models shows that there is a considerable uncertainty in the Hofmann ZrO₂ dissolution correlation, which may cause the failure time of the cladding oxide layer to be advanced while the mass of UO₂ dissolved is underestimated by using Kim or Hayward models. According to the simulation results of Phebus FPT0, it can be inferred that the dissolution reaction between the molten zirconium and fuel pellets is the main reason for the melting of UO₂ at low temperatures. In the steam-starved case, part of the fuel pellets can melt down even at 2248 K and relocate to the bottom of the core, which is much lower than the melting point of UO₂.

Nomenclature

: | Area
: | Bulk
: | Average oxygen content in the metal layer
: | Mass fraction
: | Specific heat
: | Diffusion coefficient
: | Parameter to solve oxygen content
: | Pilling–Bedworth ratio of ZrO₂
: | Convective heat transfer coefficient
: | Interface
: | Effective radiation
: | Mass transfer coefficient
: | Pressure
: | Heat flux density
: | Heat generation rate
: | Radius
: | Universal gas constant
: | Volumetric heat generation rate
: | Time
: | Temperature
: | Volume
: | Mole fraction of steam in mixed gas
: | Axial height
: | Thickness of oxide layer
: | Surface emissivity
: | Thermal conductivity
: | Thickness of remaining metal
: | Density
: | Stefan–Boltzmann constant.

Subscripts

: | Flow channel
: | Cladding
: | Diffusion
: | Gas
: | Node i
: | Molten pool
: | Oxygen
: | Oxidation
: | Radial direction
: | Radiation
: | Total
: | Axial direction.

Superscripts

: | Time step n
: | Saturation.

Data Availability

The data used to support the findings of this study are included within the article.

Conflicts of Interest

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

Acknowledgments

This work was supported by the project National Key R&D Program of China (Grant no. 2018YFB1900100).

References

L. U. O. Bang-Qi and J. Lin, “Severe accident mitigation measure and severe accident sequence of CPRl000 nuclear power plant,” Nuclear Power Engineering, vol. 31, no. s1, pp. 1–3, 2010.

G. Schanz, B. Adroguer, and A. Volchek, “Advanced treatment of zircaloy cladding high-temperature oxidation in severe accident code calculations,” Nuclear Engineering & Design, vol. 232, no. 1, pp. 85–96, 2004.

V. F. Urbanic and T. R. Heidrick, “High-temperature oxidation of zircaloy-2 and zircaloy-4 in steam,” Journal of Nuclear Materials, vol. 75, no. 2, pp. 251–261, 1978.

J. T. Prater and E. L. Courtright, Properties of Reactor Fuel Rod Materials at High Temperatures, NUREG/CR-4891, PNL-6164, Pacific Northwest Laboratory, Richland, WA, USA, 1987.

S. Leistkow and G. Schanz, “Oxidation kinetics and related phenomena of Zircaloy-4 fuel cladding exposed to high temperature steam and hydrogen-steam mixtures under PWR accident conditions,” Nuclear Engineering and Design, vol. 103, no. 1, pp. 65–84, 1987.

O. De Luze, T. Haste, M. Barrachin, and G. Repetto, “Early phase fuel degradation in Phébus FP: initiating phenomena of degradation in fuel bundle tests,” Annals of Nuclear Energy, vol. 61, pp. 23–35, 2013.

P. Hofmann, G. Schanz, and A. Skokan, Chemical Interactions of Reactor Core Materials up to Very High Temperatures (KfK-4485), Forschungszentrum Karlsruhe, Karlsruhe, Germany, 1989.

P. J. Hayward and I. M. George, “Dissolution of UO₂ in molten zircaloy-4 part 3: solubility from 2000 to 2500°C,” Journal of Nuclear Materials, vol. 232, no. 1, pp. 1–12, 1996.

K. T. Kim and D. R. Olander, “Dissolution of uranium dioxide by molten zircaloy,” Journal of Nuclear Materials, vol. 154, no. 1, pp. 102–115, 1988.

P. J. Hayward and I. M. George, “Dissolution of ZrO₂ in molten zircaloy-4,” Journal of Nuclear Materials, vol. 265, no. 1-2, pp. 69–77, 1999.

R. O. Gauntt, R. K. Cole, C. M. Erickson et al., MELCOR Code Manual-Version 1.8.6, USNRC NUREC/CR 6119 Rev. 3, SAND2005-5713, Sandia National Laboratoires, Albuquerque, CA, USA, 2001.

NRC, RELAP5/RELAP5/MOD3 Code Manual, Volume V: Assessment of Modeling of Reactor Core Behavior during Severe Accidents (NUREG/CR-6150), Idaho National Engineering and Environment Laboratory, Idaho Falls, ID, USA, 2001.

O. Coindreau, “ASTEC V2.1: physical modelling of the ICARE module,” Tech. Rep., IRSN, Fontenay-aux-Roses, France, 2017, PSN-RES/SAG/2016-00422.

D. R. Olander, “Interpretation of laboratory crucible experiments on UO₂ dissolution by liquid zirconium,” Journal of Nuclear Materials, vol. 224, no. 3, pp. 254–265, 1995.

M. S. Veshchunov and A. V. Palagin, “Modeling of chemical interactions of fuel rod materials at high temperatures II. Investigation of downward relocation of molten materials,” Journal of Nuclear Materials, vol. 252, no. 1-2, pp. 110–120, 1998.

D. R. Olander, “Materials chemistry and transport modeling for severe accident analyses in light-water reactors I: external cladding oxidation,” Nuclear Engineering and Design, vol. 148, no. 2-3, pp. 253–271, 1994.

P. J. Hayward and I. M. George, “Dissolution of UO₂ in molten zircaloy-4 part 1: solubility from 2000 to 2200°C,” Journal of Nuclear Materials, vol. 208, no. 1-2, pp. 35–42, 1994.

M. S. Veshchunov, A. V. Boldyrev, V. E. Shestak, and K. Mueller, “Analysis of molten pool physico-chemical interactions and interpretation of the Phebus FP tests observations,” Nuclear Engineering and Design, vol. 238, no. 7, pp. 1728–1742, 2008.

M. S. Veshchunov, K. Mueller, and A. V. Berdyshev, “Molten corium oxidation model,” Nuclear Engineering and Design, vol. 235, no. 22, pp. 2431–2450, 2005.

K. Müller, A. Goryachev, V. Smirnov et al., Simultaneous Dissolution of UO₂ and ZrO₂ by Molten Zircaloy. New Experiments and Modelling (FZKA-6947), Forschungszentrum Karlsruhe, Karlsruhe, Germany, 2004.

K. Müller, B. Toth, M. S. Veshunov et al., Final Interpretation Report of the PHEBUS Test FPT0: Bundle Aspects, European Commission 6th Euratom Framework Programme, JRC, Brussels, Belgium, 2007.

B. Clément, N. Hanniet-Girault, G. Repetto et al., “LWR severe accident simulation: synthesis of the results and interpretation of the first Phebus FP experiment FPT0,” Nuclear Engineering and Design, vol. 226, no. 1, pp. 5–82, 2003.