Contents lists available at ScienceDirect

Progress in Nuclear Energy

journal homepage: www.elsevier.com/locate/pnucene

A new method to study the transient feasibility of IVR-ERVC strategy

Rui Guo, Wei Xu, Zhen Cao, Xiaojing Liu*, Xu Cheng

School of Nuclear Science and Engineering, Shanghai Jiao Tong University, Dong chuan Road 800, Shanghai, 200240, China

CrossMark

ARTICLE INFO

Article history: Received 16 April 2015 Received in revised form 13 November 2015 Accepted 14 November 2015 Available online 28 November 2015

Keywords: IVR ERVC CHF

Theoretical model Transient analysis

ABSTRACT

The traditional method to evaluate the feasibility of IVR-ERVC strategy is based on the steady state of the molten pool. But in the early stage, the transient behavior of the molten corium may impose a greater threat to the integrity of the reactor pressure vessel. A new method to study the transient feasibility is proposed in this paper. In order to calculate the critical heat flux in transient severe accident, a theoretical CHF model is developed suitable for the outer surface of the lower head. The effect of orientation on bubble movement is taken into consideration, and the method to deal with the non-uniform heat flux is also proposed. By comparing the prediction with the ULPU experimental data, the new model shows satisfying accuracy. Parametric analysis of the new model shows that an increased reactor pressure vessel diameter will lead to a decrease in critical heat flux at the lower head outer surface when the structure of the external flow channel keeps unchanged. A transient severe accident analysis of the large scale PWR shows that the transient behavior of the molten corium imposes a greater threat to the integrity of the reactor pressure vessel.

© 2015 The Authors. Published by Elsevier Ltd. This is an open access article under the CC BY-NC-ND

license (http://creativecommons.org/licenses/by-nc-nd/4.0/).

1. Introduction

In-vessel retention of molten core corium by external reactor vessel cooling (IVR-ERVC) is an important severe accident management strategy. The criterion of the IVR-ERVC is to assure the thermal load from the melt is lower than the coolability limit, everywhere on the lower head of the reactor pressure vessel. In that case, the decay heat will be removed by the reactor cavity flooding, to prevent the escape of radioactive material from the reactor vessel.

It's generally accepted that the IVR-ERVC strategy will succeed in AP600 and AP1000. But there is controversy that currently proposed strategy without additional measures could provide sufficient heat removal in higher power reactors. China is now developing higher power passive PWR with an operating power of 1400 MW and 1700 MW, in order to obtain the independent intellectual property rights which is very important for nuclear exports. In the design, the feasibility of applying the IVR-ERVC strategy to keep the integrity of the pressure vessel is one of the key problems. So it is necessary to investigate it carefully.

The traditional method to evaluate the feasibility of IVR-ERVC

* Corresponding author. E-mail address: xiaojingliu@sjtu.edu.cn (X. Liu).

strategy is based on the steady heat flux distribution in late stage of severe accident (Theofanous et al., 1997a). The researchers believe that the thermal load to the lower head is maximized when the debris pool has reached a steady thermal state. Heat transfer correlations available for steady molten pool are provided to calculate the thermal energy on the lower head. Critical heat flux as function of position on the lower head is also obtained based on the steady state, which will not change with inlet water temperature, mass velocity or heat flux distribution.

Although at steady state the total thermal load to lower head is maximized, it is not guaranteed everywhere because the heat flux is not distributed uniformly. So in the early stage, the transient behavior of the molten corium may impose a greater threat to the integrity of the reactor pressure vessel. Considering that the boundary condition at the outer surface of RPV lower head varies in the transient melting process, the existing critical heat flux expression is not enough. In order to evaluate the feasibility of IVR-ERVC strategy in a transient severe accident, we need to know the critical heat flux in different conditions. Experimental work is effective, but as we know, this kind of experimental work is very expensive and takes a long time. So it's useful to develop a theoretical model to predict CHF under this situation, in which the characteristics in ERVC condition such as inclined heating wall and non-uniform heat flux must be considered.

There are various mechanistic CHF models proposed so far.

http://dx.doi.org/10.1016/j.pnucene.2015.11.005

0149-1970/© 2015 The Authors. Published by Elsevier Ltd. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).

Among them, boundary layer separation model (Tong, 1975), bubble crowding model (Weisman and Pei, 1983), sublayer dryout model (Lee and Mudawwar, 1988) and interfacial lift-off model (Galloway and Mudawar, 1993) are receiving attention in the flow conditions. Guo et al. (2014) proposed a theoretical CHF model for subcooled flow boiling in a curved channel based on bubble crowding model. The model was verified in uniform heat flux condition, but it is unknown whether it is effective in non-uniform heat flux condition.

In this paper, a theoretical CHF model is developed suitable for the outer surface of the lower head, and the effect of various parameters on CHF are investigated. In order to evaluate the feasibility of IVR-ERVC strategy in a transient severe accident, the 1700 MW-class plant is simulated by the code MELCOR, to provide heat flux distribution on the lower head at different time, and the corresponding critical heat flux is calculated by the proposed CHF model.

2. The proposed CHF model

2.1. CHF mechanism

The CHF mechanism of bubble crowding model was proposed by Weisman and Pei (1983). They thought that under low quality condition the flow region could be divided into bubbly layer and bulk flow layer, and the limited turbulent interchange between them leads to the onset of CHF. The max value of void fraction in bubbly layer was postulated to be 0.82, which was determined by a balance between the outward flow of vapor and the inward flow of liquid at the bubbly layer and bulk flow interface. Guo et al. (2014) established an experiment apparatus to study the CHF phenomenon in the IVR-ERVC condition, with the width of 150 mm and depth of 156 mm. Fig. 1 presented the visual observation while the boiling crisis occurred. The red line showed the approximate location of heater surface. Bubble crowding and vapor blanketing appeared near the wall, and small bubbles were dispersed in the bulk flow region. It's rational to apply the bubble crowding model to the IVR-ERVC.

The equation of critical heat flux is expressed as:

Heater

surface ^ v

n, ^ hf - hid q = G(X2 - x1)hfgih—h;-d

where G is the mass flux due to turbulent interchange at the edge of the bubbly layer, x2 is the vapor quality of the bubbly layer, x1 is the vapor quality in the bulk flow, and hld is the liquid enthalpy at the point of bubble detachment, which is calculated from the Levy (1967) model. From the above equation, we know that the turbulent interchange at the interface and the vapor qualities of the bubbly layer and the bulk flow are precondition to calculate the critical heat flux.

It is assumed that only those velocity fluctuations which are larger than the vapor generation velocity could penetrate to the interface, the turbulent interchange was determined by:

G = Gji

where G is the mass flux, j is the velocity fluctuations that are effective in reaching the wall, and i is the turbulent intensity at the bubbly layer and bulk flow interface.

In the model of Weisman and Pei, vapor quality was calculated using homogeneous flow model. Guo et al. (2014) considered that the slip model was more suitable in IVR-ERVC condition, in which the slip ratio varied with the inclination of the heater surface. The vapor quality in bubbly layer is written as:

1+1=22 p 1

«2 Pg S

a2 is the void fraction in bubbly layer, which is assumed to be 0.82 at the CHF point. The slip ratio S is defined as the ratio ofvapor velocity to liquid velocity, and can be expressed as:

S = 1 +-

U2 - UtX2

u2 is the average velocity in the bubbly layer. It is assumed to be one half of the velocity at the interface, which can be calculated by the Karman velocity profile. ut is the bubble rise velocity, and will change with the inclination of the flow channel, given as:

ut = 1.41

sg sin q (pi - pg)

From Eq. (4) and Eq. (5), we can calculate the vapor quality and slip ratio of the bubbly layer.

The vapor quality in the bulk flow can be got by energy balance equation.

hgPgU2g(1 - h)«2 + hiPiU2i(1 - h)(1 - «2) + h1P1U1h = Gh (6)

where hl is the liquid enthalpy, h1 is the average enthalpy of bulk flow, h is the area fraction occupied by the bulk flow, u2i is the liquid velocity in bubbly layer and u2g is the vapor velocity in bubbly layer. From Eq. (6), the enthalpy of bulk flow is received. Then the vapor quality of the bulk flow is given as:

h1 - hi hg - hi

Fig. 1. Boiling photograph at CHF.

2.2. Non-uniform heat flux method

Currently, there are three approaches to account for the effect of

non-uniform power shape (Yang et al., 2006): overall power, local conditions, F-factor. The overall power approach assumes that the critical power for non-uniform heated tube is the same as that for uniform heated tube at the same cross-sectional geometry and heated length at given inlet conditions. The local conditions approach assumes that the CHF is independent of upstream history. The F-factor approach derived by Tong assumes that there is a superheated liquid layer between the bubbly layer and the heated wall. The enthalpy of the liquid layer at CHF is the same under uniform and non-uniform power shape. In the subcooled region or at low qualities, the upstream memory effect is small, and the local heat flux determines the CHF. At high qualities, however, the memory effect becomes strong, and the average heat flux determines the CHF.

In the IVR-ERVC conditions, given the water at the upper tank being saturated at atmospheric pressure, water subcooling at the heater inlet is about 10 K. Even at the CHF point, subcooled boiling is the dominant regime. So the local heat flux determines the CHF here.

The present method is to calculate CHF values at different angular degrees first. Then we get the ratio of CHF value to the local heat flux at different locations. If we gradually increase the local heat flux, boiling crisis will take place at the point with minimal ratio value as shown in Fig. 2.

2.3. Comparison of predictions with experimental data

ULPU experiments (Theofanous et al., 1997a,b; 2002a,b,c and Dinh et al., 2003) have been conducted to identify the coolability limit for AP600 and AP1000. The test facility was an effective full-scale simulation of the reactor axisymmetric geometry. Configuration I was focused on the bottom of the lower head in a saturated pool boiling condition. In configurations II and III, the CHF experiments of the overall inclination angle were conducted under natural convection conditions. The configurations IV and V studied the effect of the streamlined flow path. In the present study, ULPU IV with a streamlined flow path is selected to verify the proposed CHF model. The experimental facility was constructed as shown in Fig. 3. The height of the facility was about 6 m. The radius of the heater blocks was 1.76 m, as that of AP600 lower head. The heater blocks were made of 7.6 cm thick copper, with a width of 15 cm, and they were heated by imbedded cartridge heaters that were individually controlled to create any heat flux shape as will. Power shaping was used to simulate the axisymmetric geometry in the

Fig. 3. Schematic of ULPU IV (Theofanous et al., 2002c).

reactor.

The comparison of predictions with experimental data is as Fig. 4. It can be seen that with the increase of orientation, the critical heat flux on the heated wall first increases then decreases, which can be explained by our developed model. Bubble rise velocity increases with the orientation, so the slip ratio and steam quality in the bubbly layer tend to become bigger. Meanwhile, the thickness of the bubbly layer also increases. These factors result in increase of CHF. But at high orientation, the upstream heat length is longer, and the upstream overall power is bigger, which makes the vapor quality of the CHF point bigger and easier to reach boiling crisis. This results in decrease of CHF. At low orientation, the positive factor is dominant, but at high orientation, the negative factor is dominant.

Fig. 2. Schematic of calculation in non-uniform power shape.

Fig. 4. Comparison of prediction and experimental data in ULPU IV.

The critical heat fluxes predicted by the present model are compared with the experimental CHF data. The results are quantitatively evaluated by the quantity K, defined as:

CHFp CHFm

where subscripts p and m mean predicted and measured values respectively.

The max error is less than 25%, and the standard deviation of K is 9.3%. This shows a relatively good agreement of this developed model under IVR-ERVC condition.

2.4. Parametric effect

Fig. 5 shows the effect of inlet subcooling on CHF. The parameters in the calculation are the same as those in ULPU IV except the inlet subcooling. The CHF value increases with inlet subcooling, about 30% at 90° position of the lower head when the inlet temperature varies from 100 °C to 60 °C, which shows the potential of cooling ability at severe accident if enough cooling water is provided when accident happens.

Fig. 6 shows the effect of mass velocity on CHF. The parameters in the calculation are the same as those in ULPU IV except the mass velocity. The CHF value increases with mass velocity, about 41% at 90° position of the lower head when the mass velocity varies from 200 kg/m2s to 600 kg/m2s. If circulation mass flow rate is raised by reducing the resistance at the entrance and exit, the cooling limit of the IVR-ERVC strategy will be improved significantly.

Fig. 7 shows the effect of channel gap on CHF. The parameters in the calculation are the same as those in ULPU IV except the channel gap. The CHF value decreases with gap at low orientation, while increases with gap at high orientation, about 9% at 90° position of the lower head when the gap varies from 15 cm to 25 cm. The increased gap will bring decreased vapor quality and decreased turbulent interchange. The former will lead to increase of CHF and the latter just the opposite. At low orientation, the former factor is dominant, but at high orientation, the latter dominant.

In the ULPU IV experiment, the mass flow rate keeps nearly unchanged with the channel gap. So it is necessary to study the effect of channel gap on CHF with the same mass flow rate. Fig. 8 shows the results. The CHF value decreases with gap, about 10% at 90° position of the lower head when the gap varies from 15 cm to 25 cm. With the same mass flow rate, the mass velocity will

2400 -,

ä 1400-

CT 1200-

">—i—1—i—1—r

-200kg/m s 400kg/m2s 600kg/m2s

Angle(deg)

Fig. 6. Mass flux effect on CHF.

n—1—I—1—I—1—I—1—I—1—I—1—I—1—r

0 10 20 30 40 50 60 70 80 90 100

Angle(deg)

Fig. 7. Channel gap effect on CHF with the same mass flux.

Fig. 5. Inlet water temperature effect on CHF.

Fig. 8. Channel gap effect on CHF with the same mass flow.

decrease when the channel gap increases, which leads to the decrease of the CHF value.

Fig. 9 shows the effect of reactor vessel radius on CHF. The parameters in the calculation are the same as those in ULPU IV except the vessel radius. The CHF decreases with radius, about 8% at 90° position of the lower head when the radius varies from 1.76 m to 2.35 m. It means that the increased reactor pressure vessel volume caused by the increased power in the advanced plant will lead to decreased critical heat flux at the lower head outer surface when the structure of the external flow channel keeps unchanged. At the same time, heat load on the vessel wall increases with the power. So the thermal margin to keep the integrity of the lower head decreases, which will lead to failure of the lower head.

3. Transient feasibility of IVR-ERVC strategy

In this paper, the transient severe accident process is studied in a large scale passive PWR. Its nominal electric power is 1700 MW. The coolant system is composed by three loops, and each loop consists of one hot leg, two cold legs, one steam generator and two pumps. The passive safety systems are designed to avoid the loss of a heat sink and a core meltdown. When the core exit temperature is higher than 922.05 K, the reactor cavity will be flooded by water from in-containment refueling water storage tank (IRWST). The MELCOR nodalizationl of the 1700 MW passive PWR is shown in Fig. 10.

In the MELCOR simulation, the initial event is taken as a large break on cold leg together with station black-out (SBO) transients that lead to loss of coolant of the primary system, and both IRWST gravity injection and recirculation are assumed fail. At the onset of the accident, substantial amount of coolant is ejected to the containment, which will actuate the operation of the passive safety systems. After the accumulator (ACC) inventory is depleted, the liquid level in core keeps reducing, and the core begins to melt. The core materials fall into the lower plenum region and molten pools form. The configuration of MELCOR molten pool model is given in Fig. 11. MP2 represents the metallic molten pool, MP1 the oxide molten pool, and PD the solid particulate debris. Contiguous volumes containing molten pool components constitute coherent molten pools that are assumed to be uniformly mixed by convection so as to have uniform material composition and temperature.

Fig. 12 shows the oxide molten pool formation process in the lower plenum. Oxide molten pool appears at 8000 s, and grows to 12.6 m3 at 24,000 s. After that its volume keeps unchanged, which

0 10 20 30 40 50 60 70 80 90 100 Angle(deg)

Fig. 9. Lower head radius effect on CHF.

means a quasi-steady state achieves. Eventually, metallic molten pool of Fe—Zr is on the top, oxide molten pool of UO2—ZrO2 in the middle, and particulate debris at the bottom. It is worthwhile to note that not all the lower plenum volume is taken by the molten core material because the lower plenum volume is bigger than that of the previous designed PWR.

Fig. 13 shows the water temperature change at the inlet of the reactor cavity. At 2000 s, the water from the IRWST is 57 °C, which is subcooled. With increasingly absorbing decay heat, the water temperature reaches to 101 °C, which is slightly higher than the saturation temperature at atmospheric pressure.

Fig. 14 shows the mass flow rate change of the reactor cavity. Natural circulation is established in the cavity by the lower head surface heating. Initially, the fluid is single phase and the mass flow increases to 600 kg/s at 20,000 s. Later two phase flow is dominant in the cavity, and the mass flow rate increases greatly. At 24,000 s, the mass flow rate is 1217 kg/s, which is about twice of the single phase flow rate.

Fig. 15 shows the heat load on the lower head at different accident times. As indicated, at 12,000 s the peak heat flux locates at 37°, which is completely different from that of the steady molten pool. At the early stage of the molten pool formation, solid debris occupies most parts of the lower head, and some small oxide molten pools are scattered among them. Since there is an oxide molten pool close to the lower head wall at low angle, and the oxide molten pool imposes higher heat load than the solid debris, the heat flux there is high.

After knowing the detailed condition of the reactor cavity at different times, the corresponding critical heat flux can be calculated by the present critical heat flux model. Fig. 16 shows the critical heat flux distribution on the lower head at different times. Different from the fixed width in the ULPU experiments, the heating surface in the following calculation is hemispherical. So the mass flux changes along the flow direction, which leads to a different CHF distribution. At most of the angular positions, the critical heat flux at 24,000 s is greater than that at 20,000 s while their heat flux distribution is similar, which means that the thermal threat to the lower head at 20000 s is greater than that at 24,000 s. Comparing the parameters at the two moments, mass flow is the main difference. In the previous parametric effect study, we know that the critical heat flux increases while the mass flux increases. So the critical heat flux at 24,000 s is higher.

Fig. 17 shows the ratio of heat flux to CHF values at different times. At 24,000 s, when the steady molten pool is formed, the ratio everywhere on the lower head is less than one, which means the heat load is lower than the critical heat flux on the lower head, and the IVR-ERVC strategy is effective at this time. But at 20000 s, before the steady molten pool formed, situation is worse. At 68°and 72°position, the heat flux and CHF ratio is greater than one, which means the heat load is higher than the critical heat flux, and the IVR-ERVC strategy fails. The result indicates clearly that the traditional method to evaluate IVR feasibility based on the steady molten pool is not conservative always.

4. Conclusion

A new method to study the transient feasibility of IVR-ERVC strategy is proposed. Results are summarized as follows:

(1) A theoretical model based on bubble crowding has been developed to predict the CHF on the outer surface of the RPV lower head.

(2) The max error between the predicted and measured CHF in ULPU IV experiment is less than 25%, which shows the availability of the proposed model.

Fig. 10. MELCOR nodalizationl of the 1700 MW passive PWR.

Fig. 11. Molten pools in lower plenum (Gauntt et al., 2005).

50 i-1-1-1-1-'-1-'-1-'-1-'-r

0.0 0.4 0.8 1.2 1.6 2.0 2.4 2.i

Time(104s)

Fig. 13. Inlet water temperature of the reactor cavity.

14 1 1 1 1 1 1 1 i 1 i 1 i 1 i 1

12 - -

/ ■ .

10- / _

8- / -

O 4- /

> ■

0- ■ ■ ■ ■

-2 - ■ i ' i ' i 1 1 1 1 1 1 1 1

ra 800

0.0 0.4 0.8 1.2 1.6 2.0 2.4 2.8 3.2

~i-1-1-1-1-1-1-1-r

1.2 1.6 2.0 2.4 2.8 3.2

Time(104s)

Time(104s)

0.4 0.8

Fig. 12. Oxide molten pool volume.

Fig. 14. Mass flow rate of the reactor cavity.

Angle(deg)

Fig. 15. Heat flux distribution on the lower head.

1800 160014001200 1000 800600400 200 0

—1—

—1—

—1—

Angle(deg)

Fig. 16. Critical heat flux distribution on the lower head.

X 0.6-

i-1-1-1—i-1-1-1-1-1—i-1-1-1-1—i-1-1-r

0 10 20 30 40 50 60 70 80 90 100

Angle(deg)

Fig. 17. Heat flux and CHF ratio.

(3) CHF decreases with reactor vessel radius, which means IVR-ERVC method may lose effectiveness in high power reactor plant.

(4) The traditional method to evaluate IVR feasibility based on the steady molten pool is not conservative always.

References

Dinh, T.N., Tu, J.P., Salmassi, T., Theofanous, T.G., 2003. Limits of coolability in the AP1000-related ULPU-2400 configuration V facility. In: The 10th International Topical Meeting on Nuclear Reactor Thermal Hydraulics (NURETH-10). Korean Nuclear Society, Seoul, Korea. October 5—9 (paper G00407). Galloway, J.E., Mudawar, I., 1993. CHF mechanism in flow boiling from a short heated wall—II. Theoretical CHF model. Int. J. Heat Mass Transf. 36 (10), 2527—2540.

Gauntt, R.O., Cole, R.K., Erichson, C.M., et al., 2005. MELCOR Computer Code Manuals. In: Reference Manuals, vol. 2. Sandia National Laboratories, Albuquerque. NM 87185—073.

Guo, R., Kuang, B., Cheng, X., 2014. A theoretical CHF model for subcooled flow

boiling in curved a channel at low pressure. Ann. Nucl. Energy 69, 196—202. Lee, C.H., Mudawwar, I., 1988. A mechanistic critical heat flux model for subcooled flow boiling based on local bulk flow conditions. Int. J. Multiph. Flow 14 (6), 711—728.

Levy, S., 1967. Forced convection subcooled boiling—prediction of vapor volumetric

fraction. Int. J. Heat Mass Transf. 10 (7), 951—965. Theofanous, T.G., Liu, C., Additon, S., Angelini, S., Kymalainen, O., Salmassi, T., 1997a.

In-vessel coolability and retention of a core melt. Nucl. Eng. Des. 169,1—48. Theofanous, T.G., Maguire, M., Angelini, S., Salmassi, T., 1997b. The first results from

the ACOPO experiment. Nucl. Eng. Des. 169, 49—57. Theofanous, T.G., Tu, J.P., Dinh, A.T., Dinh, T.N., 2002a. The boiling crisis phenome-non:Part I. Nucleation and nucleate boiling heat transfer. Exp. Therm. Fluid Sci. 26, 775—792.

Theofanous, T.G., Tu, J.P., Dinh, A.T., Dinh, T.N., 2002b. The boiling crisis phenom-enon:Part II. Dryout dynamics and burnout. Exp. Therm. Fluid Sci. 26, 793—810. Theofanous, T.G., Tu, J.P., Salmassi, T., et al., 2002c. Quantification of Limits to

Coolability in ULPU-2000 Configuration IV. CRSS-02.05, 3. Tong, L.S., 1975. A Phenomenological Study of Critical Heat Flux. ASME Paper, 1975. Weisman, J., Pei, B.S., 1983. Prediction of critical heat flux in flow boiling at low

qualities. Int. J. Heat Mass Transf. 26 (10), 1463—1477. Yang, J., Groeneveld, D., Leung, L., 2006. An experimental and analytical study of the effect of axial power profile on CHF. Nucl. Eng. Des. 236 (13), 1384—1395.

Appendix nomenclature

General symbols

G: mass flux (kg/m2s)

G: lateral mass flux (kg/m2s)

h: enthalpy (J/kg)

hfg: latent heat (J/kg)

i: turbulent intensity

K: ratio of predicted to measured CHF

q: heat flux (W/m2)

S: slip ratio

u: velocity (m/s)

ut: bubble rise velocity (m/s)

x: steam quality

Greek symbols

a: void fraction h: portion of bulk flow region 6: angular position p: density (kg/m3)

J: effective portion of velocity fluctuation Subscripts

1: bulk flow 2: bubble layer d: bubble departure f: saturated liquid g: vapor l: liquid m: measured p: predicted